import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
このノートで、モデルの評価を中心に行うために、サンプル数の多い乳がんのデータを使用する。このデータは以下のように 569 サンプルからなる。目的変数は良性と悪性のいずれかとなっている。
import sklearn.datasets
cancer = sklearn.datasets.load_breast_cancer()
X = cancer.data
y = cancer.target
print(X.shape)
print(X[1:10, :])
print(y[1:10])
scikit-learn が提供している前処理用の関数は sklearn.preprocessing モジュールに実装されています。このノートではほとんど使わないが、今後、自分のデータを使ってモデルを作る場合に役だったりしますので、ここで使い方を軽くだけ紹介する。
カテゴリデータをゼロから始まる整数値に変更したい場合は OrdinalEncoder 変換器を使用する。例えば、1 列目が米の品種名、2 列目が県名となっているデータを OrdinalEncoder 変換器処理すると、両方とも整数化された結果が得られる。ただし、このまま整数化されたデータを特徴量として機械学習に代入すると、コシヒカリが 1.0、アキタコマチが 0.0 というなんらかの量と見なされて、学習されてしまう。そのため、OrdinalEncoder を使用して整数化するときは、注意すること。
import sklearn.preprocessing
enc = sklearn.preprocessing.OrdinalEncoder()
x = [['koshihikari',  'niigata'],
     ['akitakomachi', 'akita'],
     ['koshihikari',  'chiba'],
     ['koshihikari',  'toyama'],
     ['nanatsuboshi', 'hokkaido']]
print('Original data:')
print(x)
enc.fit(x)
x_enc = enc.transform(x)
print('Transformed data:')
print(x_enc)
次のようにして、自分で整数化する際の順番を設定することもできる。
x1 = ['nanatsuboshi', 'akitakomachi', 'koshihikari']
x2 = ['hokkaido', 'akita', 'niigata', 'chiba', 'toyama']
# numpy 配列にするとエラーがでるようなので 2 次元リストを利用する
enc = sklearn.preprocessing.OrdinalEncoder(categories=[x1, x2])
print('Original data:')
print(x)
enc.fit(x)
x_enc = enc.transform(x)
print('Transformed data:')
print(x_enc)
カテゴリデータの場合は一般的に one hot encodeing で数値化する。scikit-learn では OneHotEncoder を使用して変換を行う。この変換を行うと、出力結果は scikit-learn のスパース行列の形で保存されます。そのままで理解しにくいので、このスパース行列を普通の行列に直して表示させてみる。
import sklearn.preprocessing
enc = sklearn.preprocessing.OneHotEncoder()
x = [['koshihikari',  'niigata'],
     ['akitakomachi', 'akita'],
     ['koshihikari',  'chiba'],
     ['koshihikari',  'toyama'],
     ['nanatsuboshi', 'hokkaido']]
print('Original data:')
print(x)
enc.fit(x)
x_enc = enc.transform(x)
print('Transformed data (sparse matrix):')
print(x_enc)
x_enc = x_enc.toarray()
print('Transformed data:')
print(x_enc)
カテゴリ名と one hot 表現(列数)の対応は、categories_ からアクセスして確認できる。
print(enc.categories_)
import sklearn.preprocessing
x1 = ['nanatsuboshi', 'akitakomachi', 'koshihikari']
x2 = ['hokkaido', 'akita', 'niigata', 'chiba', 'toyama']
enc = sklearn.preprocessing.OneHotEncoder(categories=[x1, x2])
# numpy 配列にするとエラーがでるようなので 2 次元リストを利用する
x = [['koshihikari',  'niigata'],
     ['akitakomachi', 'akita'],
     ['koshihikari',  'chiba'],
     ['koshihikari',  'toyama'],
     ['nanatsuboshi', 'hokkaido']]
enc.fit(x)
x_enc = enc.transform(x)
x_enc = x_enc.toarray()
print('Transformed data:')
print(x_enc)
One hot 表現に変換するとき、最初の特徴量を取り除くという dummy encoding を行う場合は、次のようにする。
import sklearn.preprocessing
x1 = ['nanatsuboshi', 'akitakomachi', 'koshihikari']
x2 = ['hokkaido', 'akita', 'niigata', 'chiba', 'toyama']
enc = sklearn.preprocessing.OneHotEncoder(categories=[x1, x2], drop='first')
# numpy 配列にするとエラーがでるようなので 2 次元リストを利用する
x = [['koshihikari',  'niigata'],
     ['akitakomachi', 'akita'],
     ['koshihikari',  'chiba'],
     ['koshihikari',  'toyama'],
     ['nanatsuboshi', 'hokkaido']]
enc.fit(x)
x_enc = enc.transform(x)
x_enc = x_enc.toarray()
print('Transformed data:')
print(x_enc)
特徴量全体にいくつかの区画を設けて離散化する場合は KBinsDiscretizer を使用する。
import sklearn.preprocessing
x = np.array([[21.1, 16.3],
              [20.6, 12.1],
              [21.4, 13.2],
              [20.7,  9.9],
              [22.3, 11.1],
              [23.5, 10.8],
              [20.3, 10.2]])
est = sklearn.preprocessing.KBinsDiscretizer(n_bins=[3, 2], encode='ordinal').fit(x)
est.fit(x)
x_est = est.transform(x)
print(x_est)
特徴量の 2 値化は Binarizer を使用する。このとき、閾値を設けて、閾値を超えたら 1、そうでない場合は 0 にするような 2 値化が行われる。
import sklearn.preprocessing
x = np.array([[21.1, 16.3],
              [20.6, 12.1],
              [21.4, 13.2],
              [20.7,  9.9],
              [22.3, 11.1],
              [23.5, 10.8],
              [20.3, 10.2]])
est = sklearn.preprocessing.Binarizer(threshold=[21.0, 11.1]).fit(x)
est.fit(x)
x_est = est.transform(x)
print(x_est)
解析者がテストデータセットの予測結果と教師ラベルを比べて、TP、TN、FP、FN を計算し、定義式に基づいて precision や recall などを算出できるが、scikit-learn が提供する sklearn.metrics モジュールを利用する便利である。関数 1 つで、precision や recall などを計算できるようになり、間違いも回避できる。
以下は、ロジスティック回帰で分類モデルを構築し、様々な評価指標を計算する例を示している。
import sklearn.model_selection
import sklearn.linear_model
X_subset = X[:, 0:2]
X_train, X_test, y_train, y_test =\
  sklearn.model_selection.train_test_split(X_subset, y, test_size=0.4, random_state=2020)
clf = sklearn.linear_model.LogisticRegression(penalty='l2', C=0.01, multi_class='ovr', solver='lbfgs')
clf.fit(X_train, y_train)
y_pred = clf.predict(X_test)
評価指標のほとんどは sklearn.metrics モジュールに実装されている。このモジュール中の関数を使って、色々な指標を計算してみる。
import sklearn.metrics
confusion_matrix = sklearn.metrics.confusion_matrix(y_test, y_pred)
print('confusion matrix')
print(confusion_matrix)
accuracy = sklearn.metrics.accuracy_score(y_test, y_pred)
print('accuracy:', accuracy)
f1 = sklearn.metrics.f1_score(y_test, y_pred)
print('f1:', f1)
precision = sklearn.metrics.precision_score(y_test, y_pred)
print('precision:', precision)
recall = sklearn.metrics.recall_score(y_test, y_pred)
print('recall:', recall)
ROC 曲線や PR 曲線を描いたり、曲線下の面積(AUC)を計算したりする場合は、予測結果としてスコアを使用する必要がある。ここで、ロジスティック回帰を利用した予測で、clf.predict 関数で予測ラベルを出力するのではなく、clf.predict_proba 関数で予測確率を出力させて、その確率をスコアと皆して ROC 曲線の座標を計算したり、AUC を計算したりしてみる。
# 確率を予測する
y_prob_pred = clf.predict_proba(X_test)
print(y_prob_pred[1:10, :])
# class-0 の予測確率は 1 列目、class-1 の予測確率は 2 列目なので、y_prob_pred[:, 1] をスコアとして使用
# 1 を positive と皆したいので、pos_label=1 を指定
fpr, tpr, thresholds = sklearn.metrics.roc_curve(y_test, y_prob_pred[:, 1], pos_label=1)
print('fpr', fpr)
print('tpr', tpr)
print('thresholds', thresholds)
auc = sklearn.metrics.auc(fpr, tpr)
print('auc:', auc)
# ROC-AUC を計算するだけならば class-1 になる確率を `roc_auc_score` 関数に与えるだけで十分 
auc = sklearn.metrics.roc_auc_score(y_test, y_prob_pred[:, 1])
print('auc:', auc)
上で計算した ROC の座標(TPR および FPR)を利用して、ROC 曲線を描いてみよう。
import matplotlib.pyplot as plt
# 
# ROC 曲線
#
ROC 曲線と同様に、PR 曲線を描いてみよう。PR 曲線の座標を計算する関数の名前がわからない場合は、インターネットなどを活用しましょう。
# PR 曲線をプロットし、AUC を計算する
#
#
#
モデルのハイパーパラメーターを探索するときによく使われる検証方法である。訓練データセットをさらに k 分割して、k 回の学習と検証を繰り返して、その検証結果の平均を評価指標としているので、より汎化性能の高いモデルが得られると期待できる。
scikit-learn で k 交差検証を行う時、データセットを k 分割してそれを for 文で回して解析者自身が学習と検証を制御する方法と分割数を指定して scikit-learn に自動的に交差検証を行う方法の 2 通りがある。解析者自身で学習と検証を行う場合は次のようにする。なお、ここではパイプラインを作成し、そのパイプラインに対して k 交差検証でハイパーパラメーターの探索を行なってみる。
import sklearn.model_selection
import sklearn.metrics
import sklearn.linear_model
import sklearn.pipeline
X_train, X_test, y_train, y_test =\
  sklearn.model_selection.train_test_split(X, y, test_size=0.4, random_state=2020)
# 訓練データセットを k 分割する
skf = sklearn.model_selection.StratifiedKFold(n_splits=10, random_state=2020, shuffle=True)
#kf = sklearn.model_selection.KFold(n_splits=3, random_state=2020, shuffle=True)
# for 文で各ハイパーパラメーターで評価する
for c in [0.001, 0.01, 0.1, 1.0, 10, 100]:
    scores = []
    
    # for 文で学習と評価を繰り返す
    for train_index, test_index in skf.split(X_train, y_train):
        pipe = sklearn.pipeline.Pipeline([
            ('scaler', sklearn.preprocessing.StandardScaler()),
            ('clf', sklearn.svm.SVC(kernel='linear', C=c))
        ])
        pipe.fit(X_train, y_train)
        y_pred = pipe.predict(X_test)
        _score = sklearn.metrics.accuracy_score(y_test, y_pred)
        scores.append(_score)
    print(c, sum(scores)/len(scores))
分割数を指定して、scikit-learn で(for 文を使わずに)自動的に交差検証を繰り返す場合は、cross_val_scores または cross_validate 関数を使用する。後者の方が新しく追加された関数であり、複数の評価指標でモデル検証を行えるなどの機能がある。ここでは後者の cross_validate 関数の使い方を紹介する。
import sklearn.model_selection
import sklearn.metrics
import sklearn.linear_model
import sklearn.pipeline
X_train, X_test, y_train, y_test =\
  sklearn.model_selection.train_test_split(X, y, test_size=0.4, random_state=2020)
# 訓練データセットを k 分割する
skf = sklearn.model_selection.StratifiedKFold(n_splits=10, random_state=2020, shuffle=True)
#kf = sklearn.model_selection.KFold(n_splits=3, random_state=2020, shuffle=True)
# 評価指標
scoring = ['f1_macro', 'precision_macro', 'recall_macro']
# for 文で各ハイパーパラメーターで評価する
for c in [0.001, 0.01, 0.1, 1.0, 10, 100]:
    scores = {
        'test_f1_macro': [],
        'test_precision_macro': [],
        'test_recall_macro': []
    }
    
    pipe = sklearn.pipeline.Pipeline([
            ('scaler', sklearn.preprocessing.StandardScaler()),
            ('clf', sklearn.svm.SVC(kernel='linear', C=c))
        ])
    
    # cv=10 のように整数値を代入してもよい
    score = sklearn.model_selection.cross_validate(pipe, X_train, y_train, scoring=scoring, cv=skf)
    
    # 10-cv の描く指標の平均値を計算する
    for score_key in ['test_f1_macro', 'test_precision_macro', 'test_recall_macro']:
        scores[score_key] = np.mean(score[score_key])
    
    
    print(c, scores)
データセットを分割する関数として StratifiedKFold のほかに KFold や ShuffleSplit などのデータセット分割関数が用意されている。これらの関数の動作の違いについてでは、scikit-learn のウェブサイトで、画像付きで解説されている。詳細に知りたい方は下記のウェブページを参照してください。
https://scikit-learn.org/stable/modules/cross_validation.html 
 関連する日本語記事: https://qiita.com/Hatomugi/items/620c1bc757266b00e87f
ここで練習用に、マクロ F1 を評価指標として、上で定義した pipe パイプラインに対して、最適なハイパーパラメーター C を探してみよう。最適な C を探す時、最初は粗めに C の候補を決めて、交差検証を行う。
## 以下例、必ずしもこの通りにする必要はない
##
## for c in [0.001, 0.01, 0.1, 1, 10, 100, 1000]
##   交差検証で評価する
##
## 0.1 のときの性能が最大であれば、次の for 文を次のようにする
##
## for c in [0.01, 0.05, 0.1, 0.15, 0.2, ...., 0.95, 1.0]
##
##
## このように範囲を徐々に縮めていく
##
学習曲線および検証曲線を用いることで、モデルの学習不足や過学習などを評価することができる。ここで、学習曲線および検証曲線を描く例を示す。
学習曲線は横軸に訓練データのサンプル数、縦軸に評価指標をプロットした折れ線グラフである。for 文を使用して、サンプル数を調整しながら、学習と評価を繰り返してもよいが、ここでは scikit-learn が提供している learning_curve 関数を使用する。
import sklearn.model_selection
import sklearn.metrics
import sklearn.svm
import sklearn.pipeline
# 最初の 5 特徴量で学習曲線をプロットする
X_subset = X[:, 0:5]
X_train, X_test, y_train, y_test =\
  sklearn.model_selection.train_test_split(X_subset, y, test_size=0.4, random_state=2020)
# 標準化を行ってから分類を行うパイプラインを構築
pipe = sklearn.pipeline.Pipeline([
         ('scaler', sklearn.preprocessing.StandardScaler()),
         ('clf', sklearn.svm.SVC(kernel='linear', C=1))
       ])
# 学習と検証を繰り返して、学習曲線をプロットするための座標を計算
train_sizes, train_scores, test_scores =\
        sklearn.model_selection.learning_curve(estimator=pipe,
                                               X=X_train,
                                               y=y_train,
                                               train_sizes=np.linspace(0.1, 1.0, 21),
                                               cv=10,
                                               scoring='accuracy')
# 10-cv の平均値と標準偏差を計算
train_scores_mean = np.mean(train_scores, axis=1)
train_scores_std = np.std(train_scores, axis=1)
test_scores_mean = np.mean(test_scores, axis=1)
test_scores_std = np.std(test_scores, axis=1)
# プロット
plt.plot(train_sizes, train_scores_mean, color='orange', marker='o', markersize=5, label='training accuracy')
plt.fill_between(train_sizes, train_scores_mean + train_scores_std, train_scores_mean - train_scores_std,
                 alpha=0.2, color='orange')
plt.plot(train_sizes, test_scores_mean, color='darkblue', marker='o', markersize=5, label='validation accuracy')
plt.fill_between(train_sizes, test_scores_mean + test_scores_std, test_scores_mean - test_scores_std,
                 alpha=0.2, color='darkblue')
plt.xlabel('#training samples')
plt.ylabel('accuracy')
plt.legend(loc='lower right')
plt.show()
特徴量を少しずつ増やしたとき、学習曲線がどのように変化していくのかを調べてみよう。特徴量が少ないときは、training accuracy と validation accuracy のどちらも低い。特徴量が多い時は、training accuracy はほぼ 1.0 のままで、validation accuracy も 1.0 に近い値を取るようになる。
# use for-loop to check the changes of accuracy curve
## n_features = X.shape[1]
## for i in range(2, n_features):
##   plot learning-curve with i features
##
検証曲線もモデルの学習不足や過学習に使用できる。検証曲線の場合は validation_curve 関数を使用すると便利である。
import sklearn.model_selection
import sklearn.metrics
import sklearn.svm
import sklearn.pipeline
# 最初の 5 特徴量で学習曲線をプロットする
X_subset = X[:, 0:5]
X_train, X_test, y_train, y_test =\
  sklearn.model_selection.train_test_split(X_subset, y, test_size=0.4, random_state=2020)
pipe = sklearn.pipeline.Pipeline([
         ('scaler', sklearn.preprocessing.StandardScaler()),
         ('clf', sklearn.svm.SVC(kernel='linear'))
       ])
# 学習と検証を繰り返して、検証曲線をプロットするための座標を計算
svc_c_range = 10 ** np.linspace(-5, 2, 10)
print(svc_c_range)
train_scores, test_scores =\
        sklearn.model_selection.validation_curve(estimator=pipe,
                                                 X=X_train,
                                                 y=y_train,
                                                 param_name='clf__C',
                                                 param_range=svc_c_range,
                                                 cv=10,
                                                 scoring='accuracy')
# 10-cv の平均値と標準偏差を計算
train_scores_mean = np.mean(train_scores, axis=1)
train_scores_std  = np.std(train_scores, axis=1)
test_scores_mean = np.mean(test_scores, axis=1)
test_scores_std  = np.std(test_scores, axis=1)
# プロット
plt.plot(svc_c_range, train_scores_mean, color='orange', marker='o', markersize=5, label='training accuracy')
plt.fill_between(svc_c_range, train_scores_mean + train_scores_std, train_scores_mean - train_scores_std,
                 alpha=0.2, color='orange')
plt.plot(svc_c_range, test_scores_mean, color='darkblue', marker='o', markersize=5, label='validation accuracy')
plt.fill_between(svc_c_range, test_scores_mean + test_scores_std, test_scores_mean - test_scores_std,
                 alpha=0.2, color='darkblue')
plt.xlabel('C')
plt.ylabel('accuracy')
plt.legend(loc='lower right')
plt.xscale('log')
plt.show()
特徴量を少しずつ増やしたとき、検証曲線がどのように変化していくのかを調べてみよう。特徴量が少ないときは、training accuracy と validation accuracy のどちらも低い。特徴量が多い時は、training accuracy はほぼ 1.0 のままで、validation accuracy も 1.0 に近い値を取るようになる。
# use for-loop to check the changes of accuracy curve
## n_features = X.shape[1]
## for i in range(2, n_features):
##   plot validation-curve with i features
##
余裕があれば、サンプル数が多い時と少ない時、特徴量が多い時と少ない時の 4 通りの組み合わせで検証曲線をプロットして、サンプルの数・特徴量の数がモデルの性能に及ぼす影響を考察してみましょう。
# use for-loop to check the changes of accuracy curve
## n_features = [2, 30] # 自由に決めてよい
## n_samples = [50, 569] # 自由に決めてよい
## for each sample_set
##   for each feature_set
##     plot validation-curve with i samples x j features
##
最適なハイパーパラメーターを探し出すために、ハイパーパラメーターの候補となる値を複数用意し、すべての候補に対して交差検証を行えばよい。ハイパーパラメーターが 1 つだけあるとき for 文を利用したループでも十分にパラメーターの探索を行うことが可能。しかし、ハイパーパラメーターが複数になると、for 文を複数入れ子構造にする必要が生じてくる。このままでは不便であるから、scikit-learn で用意されている GridSearchCV や RandomizedSearchCV 関数を使用すると便利。
まずは理解しやすいグリッドサーチを行う GridSearchCV を行う例を示す。
import sklearn.model_selection
import sklearn.metrics
import sklearn.svm
import sklearn.pipeline
X_train, X_test, y_train, y_test =\
  sklearn.model_selection.train_test_split(X, y, test_size=0.4, random_state=2020)
pipe = sklearn.pipeline.Pipeline([
         ('scaler', sklearn.preprocessing.StandardScaler()),
         ('clf', sklearn.svm.SVC(kernel='linear'))
       ])
parameters = {
    'clf__kernel': ('linear', 'rbf'),
    'clf__C': 10**np.linspace(-5, 2, 10)
}
gs = sklearn.model_selection.GridSearchCV(pipe, parameters, scoring='f1', cv=10)
gs.fit(X_train, y_train)
グリッドサーチの結果が gs インスタンスに保存されている。性能の最も高いモデルは best_estimator_ に保存され、そのパラメーターは best_params_ に保存され、そのときのスコアは best_score_ に保存される。
best_params = gs.best_params_
best_score = gs.best_score_
print(best_params)
print(best_score)
clf = gs.best_estimator_
# テストデータセットを使って最終評価
y_pred = clf.predict(X_test)
score = sklearn.metrics.f1_score(y_test, y_pred)
print(score)
ハイパーパラメーターのすべての組み合わせを評価した結果は cv_results_ にディクショナリの形で保存されている。そのまま表示させるとみづらいので、これをデータフレームに変換して表示してみる。最初の 4 列は、各ハイパーパラメーターの組み合わせ時における 10-cv の学習時間の平均と標準偏差、テスト時間の平均と標準偏差である。
print(pd.DataFrame(gs.cv_results_))
SVM のようにカーネル関数が異なるとハイパーパラメーターをチューニングする必要がある。例えば、多項式カーネルのハイパーパラメーターは degree, gamma, coef0 があり、RBF カーネルのハイパーパラメーターは gamma、シグモイドカーネルのハイパーパラメーターは gamma, coef0 がある。そのため、次のようなグリッドサーチで、総当たりでのアプローチは効率が悪い。つまり、次のように書くとカーネル関数が RBF のときでも、degree を一つずつ動かして検証することになるので、その分だけ無駄になる。
pipe = sklearn.pipeline.Pipeline([
         ('scaler', sklearn.preprocessing.StandardScaler()),
         ('clf', sklearn.svm.SVC(kernel='linear'))
       ])
parameters = {
    'clf__C': 10**np.linspace(-5, 2, 10),
    'clf__kernel': ('linear', 'rbf'),
    'clf__degree': (2, 3, 4, 5, 6, 7),
    'clf__gamma': 10**np.linspace(-5, 2, 10),
    'clf__coef0': np.linspace(-10, 10, 10)
}
gs = sklearn.model_selection.GridSearchCV(pipe, parameters, scoring='f1', cv=10)
ここで、次のように、ハイパーパラメーターの組み合わせをいくつかのパターンにまとめて、与えると効率よくグリッドサーチができるようになる。
import sklearn.model_selection
import sklearn.metrics
import sklearn.svm
import sklearn.pipeline
X_train, X_test, y_train, y_test =\
  sklearn.model_selection.train_test_split(X, y, test_size=0.4, random_state=2020)
pipe = sklearn.pipeline.Pipeline([
         ('scaler', sklearn.preprocessing.StandardScaler()),
         ('clf', sklearn.svm.SVC(kernel='linear'))
       ])
parameters = [
        {
            'clf__kernel': ['linear'],
            'clf__C': 10**np.linspace(-5, 2, 3),
        },
        {
            'clf__kernel': ['rbf'],
            'clf__C': 10 ** np.linspace(-5, 5, 3),
            'clf__gamma': 10**np.linspace(-5, 2, 3),      
            'clf__degree': 10 ** np.linspace(-5, 5, 3),
            'clf__coef0': np.linspace(-10, 10, 3)
        },
        {
            'clf__kernel': ['rbf'],
            'clf__C': 10**np.linspace(-5, 2, 3),
            'clf__gamma': 10**np.linspace(-5, 2, 3)
        },
        {
            'clf__kernel': ['sigmoid'],
            'clf__C': 10**np.linspace(-5, 2, 3),
            'clf__gamma': 10**np.linspace(-5, 2, 3),
            'clf__coef': np.linspace(-10, 10, 3)
        }
    ]
gs = sklearn.model_selection.GridSearchCV(pipe, parameters, scoring='f1', cv=10)
gs.fit(X_train, y_train)
# 最適なモデルを取得
clf = gs.best_estimator_
best_params = gs.best_params_
best_score = gs.best_score_
print(best_params)
print(best_score)
# テストデータセットを使って最終評価
y_pred = clf.predict(X_test)
score = sklearn.metrics.f1_score(y_test, y_pred)
print(score)
ランダムサーチもグリッドサーチと同じ要領で実行できる。また、ランダムサーチの出力結果もグリッドサーチと同様に扱える。ただし、ハイパーパラメーターの候補は、確率分布の形で与える。
import sklearn.model_selection
import sklearn.metrics
import sklearn.svm
import sklearn.pipeline
import scipy
X_train, X_test, y_train, y_test =\
  sklearn.model_selection.train_test_split(X, y, test_size=0.4, random_state=2020)
pipe = sklearn.pipeline.Pipeline([
         ('scaler', sklearn.preprocessing.StandardScaler()),
         ('clf', sklearn.svm.SVC(kernel='linear'))
       ])
parameters = {
    'clf__kernel': ('linear', 'rbf'),
    'clf__C': scipy.stats.uniform(loc=0, scale=4)
}
gs = sklearn.model_selection.RandomizedSearchCV(pipe, parameters, scoring='f1', cv=10)
gs.fit(X_train, y_train)
best_params = gs.best_params_
best_score = gs.best_score_
print(best_params)
print(best_score)
clf = gs.best_estimator_
# テストデータセットを使って最終評価
y_pred = clf.predict(X_test)
score = sklearn.metrics.f1_score(y_test, y_pred)
print(score)
# 各ハイパーパラメーターの組み合わせおよびその時のスコア
print(pd.DataFrame(gs.cv_results_))
2 種類の分類アルゴリズムを比較するにはネストされた交差検証を行う必要がある。子の交差検証ループで各アルゴリズムの最適なハイパーパラメーター(SVM ならば C、決定木ならば max_depth など)を決定し、親の交差検証ループではアルゴリズム(SVM や決定木など)の性能評価を行う。
import sklearn.model_selection
import sklearn.metrics
import sklearn.svm
import sklearn.tree
import sklearn.pipeline
X_train, X_test, y_train, y_test =\
  sklearn.model_selection.train_test_split(X, y, test_size=0.4, random_state=2020)
pipe_1 = sklearn.pipeline.Pipeline([
         ('scaler', sklearn.preprocessing.StandardScaler()),
         ('clf', sklearn.svm.SVC(kernel='linear'))
       ])
parameters_1 = {
    'clf__kernel': ('linear', 'rbf'),
    'clf__C': 10**np.linspace(-5, 2, 10)
}
pipe_2 = sklearn.pipeline.Pipeline([
               ('clf', sklearn.tree.DecisionTreeClassifier())
           ])
parameters_2 = [
    {'clf__max_depth': [2, 3, 4, 5, 6, 7, 8]}
]
# 子の交差検証ループで最適なハイパーパラメーターを探索
gs_1 = sklearn.model_selection.GridSearchCV(pipe_1, parameters_1, scoring='f1', cv=2)
gs_2 = sklearn.model_selection.GridSearchCV(pipe_2, parameters_2, scoring='f1', cv=2)
# 親の交差検証ループで最適なアルゴリズムを評価
scores_1 = sklearn.model_selection.cross_validate(gs_1, X_train, y_train, scoring='f1', cv=10)
scores_2 = sklearn.model_selection.cross_validate(gs_2, X_train, y_train, scoring='f1', cv=10)
# 評価結果
print('pipeline 1', scores_1['test_score'].mean())
print('pipeline 2', scores_2['test_score'].mean())
最適なパイプラインは pipe_1 であることがわかったので、このパイプラインにすべての教師データを代入してモデルを構築し、最終評価を行う。
gs_1.fit(X_train, y_train)
clf = gs_1.best_estimator_
# テストデータセットを使って最終評価
y_pred = clf.predict(X_test)
score = sklearn.metrics.f1_score(y_test, y_pred)
print(score)
上の例では子の交差検証ループを GridSearchCV 関数にセットアップし、これらを親の交差検証ループを cross_validate 関数に代入して実行した。このような入れ子構造でネストされた交差検証を可能にし、評価結果も得られる。この結果からパイプライン 1 の方が優れていることはわかった。その後、すべての教師データをパイプライン 1 に代入し、モデルを構築した。
scikit-learn が提供している関数を使用すると、ネストされた交差検証のような複雑な機能を簡易に実現できる。しかし、複雑なモデルを構築する場合、調整が不可能になったりします。では、練習として、上で用いたネストされた交差検証とほぼ同じような機能を果たす交差検証を(cross_validate 関数を使わずに)for 文、KFold、GridSearchCV などの関数で書いてみよう。
##
## for i in 10-cv:
##    GridSearchCV with 2-cv for pipe_1
##    GridSearchCV with 2-cv for pipe_2
## 
## 10-cv 結果から pipe_1 と pipe_2 の評価結果を計算
##
#!/usr/bin/env python
"""
For windows users, install the graphviz with the installer listed below
    https://graphviz.gitlab.io/_pages/Download/Download_windows.html
If you got problems, you may copy the graphviz folder into the following location.
     AppData\Local\Continuum\anaconda3\Lib\site-packages\graphviz 
"""
import os
import subprocess
import graphviz.backend as be
import platform as p
import shutil
os_sys = ["Linux", "Darwin", "Windows"]
cmd = ""
dot_cmd = ['dot','-V']
plf = p.system()
rel = p.release()
print("OS : {}, Ver: {}".format(plf, rel))
print("System path:")
if plf == os_sys[2]:
    cmd = 'dot'
    print( os.getenv('Path') )
elif plf == os_sys[1] or plf == os_sys[0]:
    print( os.getenv('PATH') )
print()
print('Graphviz version: for dot command')
found = shutil.which('dot')
if not found:
    print("Do not found {}".format(cmd))
else:
    print("Run dot command with subprocess")
    proc = subprocess.Popen(dot_cmd, stdout=subprocess.PIPE, stdin=subprocess.PIPE, stderr=subprocess.PIPE)
    (output, err) = proc.communicate()
    proc.stderr.close()
    print(err.decode().strip())
    print()
    print("Run dot command with graphviz backend")
    stdout, stderr = be.run(dot_cmd, capture_output=True, check=True, quiet=True)
    print(stderr.decode().strip())