やること

この投稿では、ロジスティック回帰分析を契約履歴に応用して、見込顧客の「契約確率」を予測するまでの一連のステップをまとめています。

ロジスティック回帰は対象データが特定のカテゴリーに分類される確率を予測する分析手法です。例えば「契約する」「契約しない」という２つのカテゴリーがあったとき、ある顧客が「契約する」（のカテゴリーに分類される）確率は何％である、と答ることが出来ます。
 
※”ロジスティック”は「記号論理の」という意味の形容詞で、物流最適化の文脈で使われる「ロジスティクス」とは別の単語ですので混同しないようにご注意を。
※この投稿は以下の記事について、翻訳、修正、追記をしたものです。
 

DataScience+

 

175 Shares

67 Pockets

Building A Logistic Regression in Python, Step by Step

https://datascienceplus.com/building-a-logistic-regression-in-python-step-by-step/

Logistic Regression is a Machine Learning classification algorithm that is used to predict the probability of a categorical dependent variable. In logistic regression, the dependent variable is a binary variable that contains data coded as 1 (yes, success, etc.) or 0 (no, failure, etc.). In other words, the logistic regression model predicts P(Y=1) as a

ロジスティック回帰の仮定

まず、ロジスティック回帰が想定する仮定は以下の通り：

1. 二値分類問題を解くために用いる。
   Yes/No、成功/失敗、契約成立/契約非成立 等どちらになるか判別するような問題は「二値分類」と言われます。気温や株価など、変動範囲が青天井の数値データを予測する場面には用いません。

 

2. 目的変数（Y）が1のときに求める事象が発生する。
   「契約成立」という事象が発生する確率に関心があるため、求める結果は契約あり=1、契約なし=0となるように設定します。

 

3. モデルの説明変数xは多重共線性（multicollinearity）を持たない。
   例えば、見込顧客に対して「送信したメールの数」と、「そのメールに対して返信をもらった数」という2つの変数は互いに連動してしまうので、2つともモデルに組み込むのは得策ではありません^[1]それぞれの変数について相関関係をチェックし、相関係数が 0.8 とか 0.9 を越えるような場合には，「多重共線性」を疑う必要があります。。目的変数に対するそれぞれの説明力を正確に測定出来なくなってしまうため、相関関係があるものを同時に投入してしまわぬよう注意する必要があります。

 

4. 説明変数は意味のある（因果関係を説明しうる）ものを用いる。
   当たり前ですが、例えば自動車のディーラーが新車を買ってくれそうかどうかを、そのお客さんの「身長」で判断したりしないですよね。全然関係ないですから。それよりも「家族構成」とか「年齢」とか「職種」とかの方が関係ありそうですし納得もできます。何か新しい因果関係の可能性を探りたい、という状況でなければ、意味のない変数を計算に入れるのはやめておきましょう。

 

5. 被説明変数（Y）は対数オッズ（ロジット）と線形の関係を持つ。
   「対数オッズ」の考え方については、実際の計算結果を見ながらの方が理解しやすいので、後ほど詳しく説明します。

 

6. パラメーター推定には大きめのサンプルデータが必要。
   ロジスティック回帰は表題のような問いに答えるのに非常に強力な手法ですが、そのモデルを作成するための「最尤推定法」という計算アルゴリズムは、正しかろう結果を出すために大きめのデータ量を必要とします。解きたい問題やデータセットの状況によりますが、多くの場合、正しかろう計算をするためには最低でも1,000程度のデータ量が必要となるかと思います。（← 主観）

実際にビジネス課題に応用する際は、解こうとしている問題が上記6つの仮定に沿っているか意識しながら進めてください。

利用データ

今回のデータは UCI Machine Learning Repository にある こちら^[2]S. Moro, P. Cortez and P. Rita. A Data-Driven Approach to Predict the Success of Bank Telemarketing. Decision Support Systems, Elsevier, 62:22-31, June 2014 を使います。ポルトガルの金融機関で電話によるダイレクトマーケティングキャンペーンを実施した際のデータ (Moro et al., 2014) です。顧客属性に応じて定期預金に申し込んでくれたかどうか、がテーマです。今回のロジスティック回帰の目的は「見込顧客が定期預金を契約する確率（P（Y=1））」です。

Pythonでデータに接続して、まずは簡単に中身を確認してみましょう。
※実行環境は Python 3.6 + Jupyter

#必要なライブラリのインポート
import pandas as pd
import numpy as np
import math
import matplotlib.pyplot as plt
import statsmodels.api as sm
from sklearn.feature_selection import RFE
from sklearn import model_selection
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import train_test_split
from sklearn.model_selection import cross_val_score
from sklearn import metrics
from sklearn.metrics import confusion_matrix
from sklearn.metrics import roc_auc_score
from sklearn.metrics import roc_curve
from sklearn.metrics import classification_report
from matplotlib.font_manager import FontProperties
from sklearn import datasets

#描画用プロパティ設定(Jupyter)
plt.style.use('ggplot')
%config InlineBackend.figure_format = 'retina'

#データの取得（インターネット経由）
data = pd.read_csv('https://raw.githubusercontent.com/madmashup/targeted-marketing-predictive-engine/master/banking.csv', header = 0)

#データの概要
print(data.shape)
print(list(data.columns))

[Out:]
(41188, 21)
['age', 'job', 'marital', 'education', 'default', 'housing', 'loan', 'contact', 'month', 'day_of_week', 'duration', 'campaign', 'pdays', 'previous', 'poutcome', 'emp_var_rate', 'cons_price_idx', 'cons_conf_idx', 'euribor3m', 'nr_employed', 'y']

全部で41,188行（41,188顧客分）、21列（21種類の情報）あります。ここでのデータセットが持つ各列はその顧客（行）に関する次のような情報です。

顧客情報（説明変数20列）：

1. age: 年齢（数値）
2. job: 職種（カテゴリー）
3. marital: 結婚歴 (カテゴリー）
4. education: 学歴（カテゴリー ）
5. default: 債務不履行あり (カテゴリー)
6. housing: 住宅ローンあり（カテゴリー）
7. loan: 個人ローンあり（カテゴリー）
8. contact: 連絡手段（カテゴリー）
9. month: 最終連絡月（カテゴリー）
10. day_of_week: 最終連絡曜日（カテゴリー）
11. duration: 最終連絡時会話秒（数値）
    ※”この値が0であれば当然「契約なし」となるため、ベンチマークとしてのみ用い、予測モデルに個の変数を用いないのが現実的です”（データ提供者より）
12. campaign: キャンペーン中の総連絡回数（数値）
13. pdays: 最後の連絡から経過した日数（数値）
    ※”999は、過去に連絡を取ったことがないことを意味します”（データ提供者より）
14. previous: このキャンペーンの前に連絡をとっていた回数（数値）
15. poutcome: 前回のキャンペーンの結果（カテゴリー）
16. emp.var.rate: 就労率（数値）
17. cons.price.idx: 消費者物価指数（数値）
18. cons.conf.idx: 消費者信頼指数（数値）
19. euribor3m: 3カ月物EURIBOR（数値）
20. nr.employed: 従業員数（数値）

キャンペーンの結果（教師データ、目的変数1列）:

y: その顧客は定期預金を契約してくれたか（1: はい、0: いいえ）

クレンジング

欠損値の確認
まず初めに、データ内に欠損している（Null）の値が無いか確認しましょう。

data.isnull().sum()

[Out:]

欠損している値はなさそうです。

次にデータの特徴を確認していきます。education（学歴）列にどんな種類のカテゴリーがあるのか見てみると、

data['education'].unique()

[Out:]
array(['basic.4y', 'unknown', 'university.degree', 'high.school', 'basic.9y', 'professional.course', 'basic.6y', 'illiterate'], dtype=object)

分類が少し細かすぎるようです。「基礎教育（basic）」として集約して列数を減らします。

data['education'] = np.where(data['education'] == 'basic.9y', 'Basic', data['education'])
data['education'] = np.where(data['education'] == 'basic.6y', 'Basic', data['education'])
data['education'] = np.where(data['education'] == 'basic.4y', 'Basic', data['education'])

実行後、学歴カテゴリーは以下のように集約されました。

data['education'].unique()

[Out:]
array(['Basic', 'unknown', 'university.degree', 'high.school', 'professional.course', 'illiterate'], dtype=object)

職業カテゴリーも整理しておきましょう。見たところ顧客データには「学生（student）」も含まれており、定期預金の販促キャンペーン対象として現実的ではなさそうです。非雇用、職業不明の人たちと一緒に、「働いていない人（notworking）」というカテゴリーに集約します。^[3] … Continue reading

data['job'] = np.where(data['job'] == 'student', 'notworking', data['job'])
data['job'] = np.where(data['job'] == 'unemployed', 'notworking', data['job'])
data['job'] = np.where(data['job'] == 'unknown', 'notworking', data['job'])

実行後、職業カテゴリーは以下のように集約されました。

data['job'].unique()

[Out:]
array(['blue-collar', 'technician', 'management', 'services', 'retired', 'admin.', 'housemaid', 'notworking', 'entrepreneur', 'self-employed'], dtype=object)

pdaysについて、”999は、過去に連絡を取ったことがないことを意味”しており、数値としてこのまま計算に加えてしまうと意味合いが異なる、間違った計算結果が出てくることは明らかです。

pd.DataFrame(data.groupby('pdays').size()).plot(kind = 'bar', figsize = (10, 4))
plt.show()

[Out:]

なので、「過去に連絡したことがある/ない」のカテゴリーデータに変換しましょう。

data['pdays'] = np.where(data['pdays'] == 999, 'nocontact', data['pdays']) #連絡したことがない
data['pdays'] = np.where(data['pdays'] != 'nocontact', 'contact', data['pdays']) #連絡したことがある

print(data['pdays'].unique())
pd.DataFrame(data.groupby('pdays').size()).plot(kind = 'bar', figsize = (10, 4))
plt.show()

[Out:]
['nocontact' 'contact']

カテゴリーデータに変換できました。
 

データ探索

データの状況をもう少し詳細に見ていきます。キャンペーンを実施した顧客の中で実際に契約してくれる割合はどのくらいでしょうか？

fig, ax = plt.subplots(figsize=(10, 4))

y0 = data[data['y'] == 0]['y']
y1 = data[data['y'] == 1]['y']

ax.hist([y0, y1], histtype = "barstacked", bins = np.arange(0,3,1), align = 'left', rwidth = 0.8)
ax.set_xticks(np.arange(0,3,1)[:-1])
plt.show()

[Out:]

0が契約なし、1が契約ありの顧客です。全体に対して1割強の契約率であることが分かります。

数値データの状況

「年齢」や「キャンペーン中の総連絡回数」等、数値データの状況について理解を深めるために、それぞれ「契約あり/なし」の分布を可視化します。赤色部分が「契約なし（y = 0）」、青色部分が「契約あり（y = 1）」を表しています。また、変数ごとに左側のグラフが各数値毎の「出現頻度」、右側のグラフが各数値毎の契約あり/なしの「割合」を示すように可視化します。

#数値データを指定
num_col = ['age', 'campaign', 'previous', 'emp_var_rate', 'cons_price_idx', 'cons_conf_idx', 'euribor3m', 'nr_employed']

fig, axes = plt.subplots(nrows = len(num_col), ncols = 2, figsize = (18, 20))

for x in num_col:
    
  freqs = pd.crosstab(data[x], data.y) #出現頻度データを作成
  freqs.plot(ax = axes[num_col.index(x), 0], kind = 'area') #左カラムにプロット
    
  props = freqs.div(freqs.sum(1).astype(float), axis = 0) #出現頻度データから割合（100%）データを作成
  props.plot(ax = axes[num_col.index(x), 1], kind = 'area') #右カラムにプロット

  fig.tight_layout()
    
plt.show()

[Out:]

まず年齢（age）についてみてみると、キャンペーンの対象となっている顧客は30-40代に集中していることが分かります（左グラフ）。ただし、契約してくれる割合でみると、60代以降が多いように見えます。（右グラフ）
キャンペーン中の総連絡回数（campaign）では、「1度だけ連絡をした」場合の頻度が最も高く、契約割合も徐々に減っていっています。何度も連絡を取ってみたところで顧客の気持ちはあまり変わらない、ということを示唆しています。
ただし、キャンペーン前に連絡をとった回数（previous）を見てみると、5、6回のコンタクトをしたことのある顧客の契約率が高いようです。日常的にコミュニケーションしている見込顧客は契約確度が高めである、と言えそうです。

カテゴリーデータの状況

数値データ同様、カテゴリーデータについても理解を深めるために、それぞれ「契約あり/なし」の状況を可視化してみましょう。こちらも左グラフが頻度、右グラフが契約割合を示せるように一覧化します。

#対象となるカテゴリーデータを指定
category_col = ['job', 'marital', 'education', 'default', 'housing', 'loan', 'contact', 'month', 'day_of_week', 'poutcome', 'pdays']

fig, axes = plt.subplots(nrows = len(category_col), ncols = 2, figsize = (18, 50))

for x in category_col:
  freqs = pd.crosstab(data[x],data.y) #出現頻度データを作成
  freqs.plot(ax = axes[category_col.index(x), 0], kind = 'bar', stacked = True) #左カラムにプロット
  axes[category_col.index(x)][0].set_xticklabels(freqs.index, rotation=45, size=12)
    
  props = freqs.div(freqs.sum(1).astype(float), axis = 0) #出現頻度データから割合（100%）データを作成
  props.plot(ax = axes[category_col.index(x), 1], kind = 'bar', stacked = True) #右カラムにプロット
  axes[category_col.index(x)][1].set_xticklabels(props.index, rotation = 45, size = 12)

  fig.tight_layout()
    
plt.show()

[Out:]

職種（job）について、「働いていない人（notworking）」と「退職者（retied）」の契約割合が大きくなっている通り、職種別に契約率が異なっていることが分かります。
その他にも、学歴（education）、債務不履行あり（default）、連絡手段（contact）、最終連絡月（month）、前回のキャンペーンの結果（poutcome）、過去に連絡をしたことがある/ない（pdays）、といったデータが予測の要素として使えそうだということが分かります^[4] … Continue reading。

ダミー変数の作成

カテゴリーデータをダミー変数に変換していきます。

#対象となるカテゴリーデータを指定
category_col = ['job', 'marital', 'education', 'default', 'housing', 'loan', 'contact', 'month', 'day_of_week', 'poutcome', 'pdays']

data_modified = data.copy()
for x in category_col:
  categ_list = 'var' + '_' + x
  categ_list = pd.get_dummies(data_modified[x], prefix = x, drop_first = True)
  data_dmy = data_modified.join(categ_list)
  data_modified = data_dmy
    
data_vars = data_modified.columns.values.tolist()
to_keep = [i for i in data_vars if (i not in category_col)] #カテゴリーデータ（ダミー元）を除く

ここで、7行目の引数 ”drop_first = True” ですが、ダミー変換した最初の1列を削除しています。
冒頭、ロジスティック回帰の仮定であったポイントですが、いずれか１つの項目が、他の項目の組み合わせで表現できる時も、「多重共線性」が出てしまいます。ダミー変換して作った変数をすべて使って解析をしないように注意しましょう。

最終的に利用する列は以下の通り:

data_final = data_dmy[to_keep]
data_final.columns

[Out:]
Index([
  'age', 'duration', 'campaign', 'previous', 'emp_var_rate',
  'cons_price_idx', 'cons_conf_idx', 'euribor3m', 'nr_employed', 'y',
  'job_blue-collar', 'job_entrepreneur', 'job_housemaid',
  'job_management', 'job_notworking', 'job_retired', 'job_self-employed',
  'job_services', 'job_technician', 'marital_married', 'marital_single',
  'marital_unknown', 'education_high.school', 'education_illiterate',
  'education_professional.course', 'education_university.degree',
  'education_unknown', 'default_unknown', 'default_yes',
  'housing_unknown', 'housing_yes', 'loan_unknown', 'loan_yes',
  'contact_telephone', 'month_aug', 'month_dec', 'month_jul', 'month_jun',
  'month_mar', 'month_may', 'month_nov', 'month_oct', 'month_sep',
  'day_of_week_mon', 'day_of_week_thu', 'day_of_week_tue',
  'day_of_week_wed', 'poutcome_nonexistent', 'poutcome_success',
  'pdays_nocontact'], dtype='object')

説明変数の選択(RFE)

さて、データの状況はおおむね分かり、回帰分析用にダミー変数をセットして最終的な分析用データセットが完成しました。ではどの項目を予測モデルの説明変数として採用すれば良いでしょうか？^[5] … Continue reading
ここでは、変数の選択アルゴリズムとして「再帰的特徴消去（Recursive Feature Elimination: RFE）」という手法を用いて変数を選択していきます。
RFEはすべての説明変数を加えてモデルを作成 → 最も重要度が低いものを削除 → 残ったものでモデルを作成 → 最も重要度が低いものを削除 → … と選択と削除を繰り返しながら重要度の高い説明変数を特定していくアルゴリズムです。

data_final_vars = data_final.columns.values.tolist()

Y = 'y'
all_X = [i for i in data_final_vars if i not in Y]

logreg = LogisticRegression()

rfe = RFE(logreg, 17) #ランキングトップ17の説明変数を選択する
rfe = rfe.fit(data_final[all_X], data_final[Y])

#RFEの結果から有効な説明変数を絞り込む
col_name = pd.DataFrame(data_final[all_X].columns).rename(columns = {0:'columns'})
target = pd.DataFrame(rfe.support_).rename(columns = {0:'target'})
all_features = pd.concat([col_name, target], axis = 1)
selected_features = all_features[all_features['target'] == True]

print(selected_features)

[Out:]

※スクリプト8行目で、17個の変数を指定していますが、以下のモデル推定結果（主にP値）と汎化性能の評価を見ながら探索的に判断したものです。説明変数は少ない方が良いですが、最適なモデルかどうかの判断は試行錯誤が必要になります。

モデルの推定

選択した説明変数を使ってロジスティック回帰モデルのパラメータを推定します。モデル推定のためのPythonライブラリは複数ありますが、この段階では推定モデルのサマリーと係数の評価を見たいので、統計モデルに特化したライブラリ statsmodels (sm) を使います。

#RFEから選択された変数だけのデータに整理する
selected_X = data_final[selected_features['columns']]

#モデルを推定する
logit_model = sm.Logit(data_final[Y], selected_X)
result = logit_model.fit()
print(result.summary())

[Out:]

いずれのp値も問題なさそうです。

ロジスティック回帰 係数の解釈
ロジスティック回帰で推定された係数は「対数オッズ（またはロジット）」と呼ばれる、確率の表現方法のひとつです。
確率の表現が0%から100%以外にもあるのは意外に思われるかもしれませんが、30%を0.3と表せたりするのと同じように、単純に見せ方が違うだけだと思えば簡単です。
推定した係数についている「対数オッズ」を0%から100%の日常的な確率表現に変換するには次の式に当てはめます。

\[ \displaystyle p = \frac{exp(対数オッズ)}{1 + exp(対数オッズ)} \]

\[ \displaystyle p = 事象が発生する確率 \]

例えば、job_blue-collar の係数は -0.2280 という対数オッズになっています。これを式に当てはめると、

#ロジスティック変換の逆変換（ロジット表現から確率表現へ）
logit = -0.2280
p = math.exp(logit) / (1 + math.exp(logit))
print(p) #事象が発生する確率

[Out:]
0.4432456471070593

約44.32%という数値に変換できました。仮に説明変数が job_blue-collar の１つだけ（単変量ロジスティック回帰）だった場合は、「変数が1単位増加したとき、事象が発生する確率が44.32%に増加する」というような解釈になりますが、今回のように、説明変数を複数組み合わせている多変量ロジスティック回帰では、これよりも重要な変換ポイントがあります。それは、上の式で示している

\[ \displaystyle exp(対数オッズ) \]

の部分です。指数関数expで対数を取り除いた「オッズ」と呼ばれるもうひとつの確率表現で、「事象が発生する確率は、事象が発生しない確率の何倍か」を示します。

\[ \displaystyle exp(対数オッズ) = \frac{p}{(1 – p)} = オッズ \]

\[ \displaystyle p = 事象が発生する確率 \]

\[ \displaystyle 1 – p = 事象が発生しない確率 \]

同じく -0.2280 の対数オッズから計算したオッズの値は、

logit = -0.2280
odds = math.exp(logit)
print(odds)

[Out:]
0.7961242598354538

約0.7961になりました。
ここでの解釈は「見込顧客の職種がブルーカラー（現場作業員）であるならば、キャンペーンで定期預金を契約する確率は約0.7961倍になる」となります。
このように「対数オッズ」（推定された係数）の符号がマイナスの場合、「オッズ」の値は1より小さくなります。1よりも小さい値を元の事象発生確率に掛算するので、目的変数Y（契約確率）は説明変数Xの変動と反対方向に振れます。 逆に、job_notworking や month_dec のように係数符号がプラスの場合は、説明変数の変動と同じ方向に振れる、ということになります。

尚、ここでは推定結果を見て相対的に「どの係数がどのくらい影響力をもつのか」を確認しますが、ロジスティック回帰の推定結果を統計レポートとしてまとめる際には、それぞれ説明変数の「オッズ」を比較しながら事象原因とその影響力を説明する必要があるかと思います。

モデル性能の検証

先ほどはモデル推定結果の検証を行うために、statsmodels (sm) で推定しましたが、ここからは機械学習モデルの検証メソッドが充実している、scikit-learn で汎化性能の検証を進めていきます。

#各変数について、学習用データとテスト用データに分割する
X_train, X_test, Y_train, Y_test = train_test_split(selected_X, data_final[Y], test_size = 0.3, random_state = 0)

#scikit-learnでロジスティック回帰モデルを推定する
logreg = LogisticRegression(fit_intercept = False)
logreg.fit(X_train, Y_train)

[Out:]
LogisticRegression(C=1.0, class_weight=None, dual=False, fit_intercept=False,
          intercept_scaling=1, max_iter=100, multi_class='ovr', n_jobs=1,
          penalty='l2', random_state=None, solver='liblinear', tol=0.0001,
          verbose=0, warm_start=False)

交差検証（Cross Validation）
データセットを分割し、その一部をまず解析して、残る部分でその解析のテストを行い、解析自身の妥当性の検証・確認に当てる手法です。モデルが過学習（Overfitting）していないか確認することが出来ます。以下では10分割（k = 10）にしたテストデータで「Accuracy（正解率）」の平均値を計算します。

kfold = model_selection.KFold(n_splits = 10, shuffle = True, random_state = 42) #"Answer to the Ultimate Question of Life, the Universe, and Everything"
model_to_CV = LogisticRegression()

results = model_selection.cross_val_score(model_to_CV, selected_X, data_final[Y], cv = kfold, scoring = 'accuracy')
print("10-fold Cross Validation 平均 Accuracy: %.2f" % (results.mean()))

[Out:]
10-fold Cross Validation 平均 Accuracy: 0.90

平均 Accuracy が90%と十分な値が確認できたので、それなりのモデルが出来ていると言えそうです。

混同行列（Confusion Matrix）

混同行列（Confusion Matrix）は、サンプルのうち、Y = 0 と Y = 1 それぞれに対して何個が正しく判定され、何個が誤って判定されたかをクロス表で確認するものです。

print(confusion_matrix(Y_test, Y_pred, labels=[1, 0]))

[Out:]
[[  266  1110]
 [  108 10873]]

この結果は、0と予測したデータの数、1と予測したデータの数と、それぞれの実際の数を当てはめたもので、

正しく分類したのが 11,139( = 266 + 10,873) 件、
誤って分類したのが 1,218( = 108 + 1,110) 件、
全体の「正解率（Accuracy）」は 0.9014 ( = 11,139 / (11,139 + 1,218))

であったことを示しています。整理すると次のような表になります。
 

ここで出てくる概念は、識別型機械学習モデルの汎化性能を評価する一般的な要素となります。

• True Positive(TP)： 実際に正（1）であるものを、正しく正（1）と予測できた数 → OK
• False Positive(FP)：実際は負（0）であるものを、間違って正（1）と予測してしまった数 → NG
• Flase Negative(FN)：実際は正（1）であるものを、間違って負（0）と予測してしまった数 → NG
• True Negative(TN)：実際に負（0）であるものを、正しく負（0）と予測できた数 → OK

ちなみに、交差検証（Cross Validation）でデータセットを分割して検証した平均 Accuracy も、全体の Accuracy も

Accuracy = (TP+TN) / (TP+TN+FP+FN)

と定義されます。
 

Recall、Precision、F-measure
Accuracy だけでは、Y = 1 と Y = 0 のいずれかデータ数が多く偏っている方の正解率を過大評価し、少ない方の不正解率を過小評価する可能性があります。
そのため、一般に機械学習モデルの性能検証基準は Accuracy だけでは十分ではありません。混同行列（Confusion Matrix）の結果を使い、さらに以下3つの視点からも汎化性能を確認してみます。

Precision（適合率） = TP / (TP+FP)
→ Positiveと予測された内の正解率

Recall（再現率） = TP / (TP+FN)
→ 実際にPositiveだった内、Positiveと予測された割合

F-measure（F値） = 2 × Precision × Recall / (Precision + Recall)
→ PrecisionとRecallの調和平均。
※統計学のF分布に従うF値とは別物。
 

print(classification_report(Y_test, Y_pred))

[Out:]

「support」はテストデータセットにある各データの出現回数を表します。
この結果に対する解釈ですが、「1(契約あり)と予測された顧客のうち、実際に71％が契約をしていた（1のPrecisionが0.71）」ので、ある程度信頼できる予測モデルであると言えます。
ただし、「実際に1(契約あり）であった顧客に対して予測を的中させた割合は19%しかない（1のrecallが0.19）」ため、もしこの予測モデルを使う担当者が「本当は契約してくれるはずの顧客を取りこぼしたくない」と思っている場合はあまり良いモデルとは言えません。
この課題を解こうとするシーンでは、「キャンペーンをフォローする営業リソースが有限であり、契約見込みを数値化することで、その配分を「最適化」したい」という状況が多いと思うので、これで良しとします。

ROC曲線
最後に「受信者操作特性（Receiver Operating Characteristic: ROC）曲線を示します。
ROC曲線は二値分類問題で一般的に使われるモデル性能の可視化手法で、ロジスティック回帰モデルが予測した「契約あり」に分類される可能性が高い順に左から並べていったとき、正しく予測された割合がどのように積みあがって推移するか（青線）を示すグラフです。

推定結果モデルである識別器が優秀であればあるほど、青の曲線は赤の点線（完全にあてずっぽうの場合）から離れて、左上にの隅に近づいていきます。
また、青線の下の面積のことを Area Under the Curve (AUC) と呼び、0 から 1 の比率で示されます。
識別器がパーフェクトに分類できていれば、ROC曲線は左上の隅を通過し、AUCは 1 となります。

logit_roc_auc = roc_auc_score(Y_test, logreg.predict(X_test))
fpr, tpr, thresholds = roc_curve(Y_test, logreg.predict_proba(X_test)[:,1])
plt.figure(figsize=(12, 6))
plt.plot(fpr, tpr, 'b', label='Logistic Regression (AUC = %0.2f)' % logit_roc_auc)
plt.plot([0, 1], [0, 1], 'r--')
plt.xlim([0.0, 1.0])
plt.ylim([0.0, 1.05])
plt.xlabel('False Positive Rate')
plt.ylabel('True Positive Rate')
plt.title('ROC Curve')
plt.legend(loc = 'lower right')
plt.show()

[Out:]

ここでは AUC = 0.59 という結果になりました。

予測精度を確認する方法をいくつか紹介しましたが、これらの数字が「良い数字」かどうかの判断は、モデルを使う担当者の主観とプロジェクトの状況により判断されます。あくまで指標であり、数学的に絶対的な”解”がもとめられるものではありません。

予測してみる

ここでは予測精度はOKと判断したこととし、実際に予測に使ってみます。例えば以下ような状況の見込顧客がいたとき、

costomer1 = pd.DataFrame({
  'euribor3m': [-0.323], #銀行間金利の状況
  'job_blue-collar': [0],
  'job_notworking': [0],
  'job_retired': [1], #退職された方
  'job_services': [0],
  'default_unknown': [0],
  'contact_telephone': [0], #携帯電話での連絡
  'month_dec': [0],
  'month_jul': [0],
  'month_jun': [0],
  'month_mar': [1], #3月に連絡
  'month_may': [0],
  'month_sep': [0],
  'day_of_week_mon': [0],
  'poutcome_nonexistent': [0],
  'poutcome_success': [1], #前回のキャンペーンは成功している
  'pdays_nocontact': [0] #今までに連絡したことがある
})[X_train.columns] #説明変数の並び順を調整

次のような結果を得ます。

spot_predict = logreg.predict_proba(costomer1)[0][1] #logregが分類される確率を予測
print('「契約する」に分類される可能性は {:.2f} % です。'.format(spot_predict * 100))

[Out:]
「契約する」に分類される可能性は 92.42 % です。

さらに、例えば以下のような状況の見込顧客がいたとき、

costomer2 = pd.DataFrame({
  'euribor3m': [1.343], #銀行間金利の状況
  'job_blue-collar': [0],
  'job_notworking': [0],
  'job_retired': [0],
  'job_services': [1], #サービス業の方
  'default_unknown': [1], #債務不履行があるかは分からない
  'contact_telephone': [1], #固定電話での連絡
  'month_dec': [0],
  'month_jul': [0],
  'month_jun': [0],
  'month_mar': [0],
  'month_may': [1], #5月に連絡
  'month_sep': [0],
  'day_of_week_mon': [1], #月曜日に連絡
  'poutcome_nonexistent': [1], #前回のキャンペーンは対象外
  'poutcome_success': [0],
  'pdays_nocontact': [1] #今までに連絡したことがない
})[X_train.columns] #説明変数の並び順を調整

次のような結果を得ます。

spot_predict = logreg.predict_proba(costomer1)[0][1] #logregが分類される確率を予測
print('「契約する」に分類される可能性は {:.2f} % です。'.format(spot_predict * 100))

[Out:]
「契約する」に分類される可能性は 6.39 % です。

いかがでしたでしょうか。ロジスティック回帰分析以外にも「受注確度」を計算できるモデルはいくつもあるので、いろいろ試してみてください。