【番外編】AI（機械学習）やってみた。【第6回】乳がん診断AI　その1

ソムリAI（ワインソムリエAI）で使った手法で乳がん診断AIを作ってみます。まずは決定木とニューラルネットワークの2つのモデルを試してみます。かなり高性能なモデルができました。名医誕生かもしれません　笑

ソムリAIで使った手法で乳がん診断AIを作る

ソムリAI（ワインの評価を当てるソムリエAI）では、いろいろな手法を試しましたが66%程度の精度しか出ませんでした。イマイチの精度しか出ず悔しかったので、ソムリAIで使った決定木やニューラルネットワーク・交差検証・多数決などの手法を使って、別のデータでAIを作ってみようと思います。

UCI（University of California, Irvine・カリフォルニア大学アーバイン校）が様々なビッグデータを提供しています。個人的には医療×AIがアツいと考えており、今回はウィスコンシンで収集された乳がんの診断データを使ってみようと思います。

乳がんデータのダウンロード

乳がんの診断データは以下のリンクからダウンロードできます。
https://archive.ics.uci.edu/ml/machine-learning-databases/breast-cancer-wisconsin/wdbc.data

また、データの内容（項目の説明など）は以下のページに書かれています。
https://archive.ics.uci.edu/ml/machine-learning-databases/breast-cancer-wisconsin/wdbc.names

早速Python（Jupyter Notebook）で読み込んでみます。

#入力データセットを読み込み
import pandas as pd #データハンドリング用ライブラリ呼び出し
data = pd.read_csv('BreastCancer.csv', encoding='SHIFT-JIS')
data.head() #先頭5つを表示して読み込まれていることを確認

各列の項目はUCIの解説を基に日本語にしています。患者のID・診断結果（Mはmalignantで悪性、Bはbenignで良性)に始まり、数値化した腫瘍データが並んでいます。数値化されたデータには腫瘍の半径・きめ・凹みの個数など10個のデータがあり、それぞれに対して平均・標準誤差・最大値の3つが用意されています。つまり、1患者当たり10×3で30個のデータがあります。腫瘍の画像データをそのまま使うのではなく、腫瘍の画像を数値化したデータを使うイメージなのです。

機械学習用のデータを準備

ダウンロードしたデータを機械学習用に加工します。腫瘍の悪性・良性を目的変数に、腫瘍の半径などのデータを説明変数にします。説明変数に重みづけなどを行って目的変数（悪性・良性）を当てようということです。

また、目的変数である診断結果がMとBというアルファベットであり、コンピュータでは扱いにくいため、悪性（M）を1・良性（B）を0とするように加工しています。コンピュータは数値データの方が扱いやすいのです。

from sklearn.preprocessing import LabelEncoder

#機械学習で求める解である「評価」以外の項目をdata_Xに、「評価」をdata_Yに格納
data_X = data.copy()
del data_X['ID']
del data_X['診断結果']
data_Y = data['診断結果']
 
#目的変数data_YはM（malignant：悪性）・B（benign：良性)で表現されており、
#コンピュータで扱いやすいように悪性＝1・良性＝0のように数値変換
labelencode = LabelEncoder()
data_Y = labelencode.fit_transform(data_Y)
labelencode.transform(['M', 'B'])
#返還結果を表示
print(data_Y)

0・1変換された目的変数（診断結果）が表示されました。

続いて、機械学習用のデータを機械学習のモデルを作る訓練データと、モデルの精度（汎化性能）を測る検証データに分割します。ソムリAI同様、全量の70%を訓練用・30%を検証用に分割します。

#入力データセットの70%を訓練データ・30%を検証データに分割
from sklearn.model_selection import train_test_split
X_train, X_test, Y_train, Y_test = train_test_split(data_X, data_Y,random_state=0, test_size=0.3)
#それぞれの件数を表示
print(len(X_train))
print(len(X_test))

ここまでで事前準備は完了です。ソムリAIで行ったデータの正規化（平均0・標準偏差が1になるように各データのモノサシを統一すること）は後程行います。

機械学習①：決定木

ベストな階層数の探索

まずは決定木を使った機械学習から。決定木は木の階層数が変数です。とりあえず3階層から10階層まで試してみます。

#決定木の機械学習
from sklearn.tree import DecisionTreeClassifier
from sklearn import metrics
 
#決定木の深さを3階層から10階層の範囲で探索
for i in range(3,11):
    #モデル定義
    treeModel = DecisionTreeClassifier(max_depth=i, random_state=0)
    #機械学習実行
    treeModel.fit(X_train, Y_train)
    #精度表示
    Y_train_score = metrics.accuracy_score(Y_train, treeModel.predict(X_train))
    Y_test_score = metrics.accuracy_score(Y_test, treeModel.predict(X_test))
    print("%i階層　訓練精度: %0.4f 検証精度： %0.4f" % (i, Y_train_score, Y_test_score))

すごい！軒並み90%を超える精度です。訓練精度と検証精度のバランスを見ると、訓練精度99.75%・検証精度94.74%の6階層の決定木が一番バランスが良さそうです。

重要な変数の特定

次は決定木において、どのように条件分岐しているか（どの変数を重要とみなしているか）を見ていきます。まず6階層のモデルを作り直してから決定木における重要度を出力させます（3階層から10階層までループを回して探索しており、探索後最後のモデルは10階層で最適化されているモデルになっているので、6階層で再学習させる必要があります）。

#6階層で学習しなおし
treeModel = DecisionTreeClassifier(max_depth=6, random_state=0)
treeModel.fit(X_train, Y_train)
 
Y_train_score = metrics.accuracy_score(Y_train, treeModel.predict(X_train))
Y_test_score = metrics.accuracy_score(Y_test, treeModel.predict(X_test))
print("6階層　訓練精度: %0.4f 検証精度： %0.4f" % (Y_train_score, Y_test_score))

#モデルの決定係数を表示（入力データの各項目に対する重みを表示）
Importance=pd.DataFrame({'変数名':data_X.columns,'重要度':treeModel.feature_importances_})
#重要度を降順で表示
Importance.sort_values('重要度', ascending=False)

これを見ると「凹み個数」という変数の重要度が突出しています。続いて「面積」・「凹み度合い」が続きます。どうやら腫瘍の凹みに関連する情報が診断ポイントのようです。

決定木の可視化

最後に決定木を可視化してみます。日本語の項目名が文字化けしてしまいうまく表示出来なかったので、元の英語の項目名で可視化しました。

#決定木の可視化
from sklearn import tree
from graphviz import Source
from sklearn.externals.six import StringIO
import pydotplus as pdp

#決定木モデルの書き出し
dot_data = StringIO()
tree.export_graphviz(treeModel, out_file="tree.dot", feature_names=data_X.columns, filled=True, rounded=True)

6階層すべて表示するとかなり大きな図になるので部分的に切り抜いています。先ほど見た通り、決定木のモデルについても「凹み個数（concave points)」・「面積（area）」の順で条件分岐しているのが分かります。

機械学習②：ニューラルネットワーク

続いてニューラルネットワークで機械学習してみます。ニューラルネットワークは活性化関数・ソルバー・パーセプトロンの数と階層数がモデルの精度を決める主要なハイパーパラメータ（人間が手で決める機械学習の主要パラメータ）です。とりあえず適当に設定して精度を見てみます。活性化関数はReLU・ソルバーはLbfgs・パーセプトロンは500個ずつ3階層を指定してみます。

#ニューラルネットワークによる機械学習
from sklearn.neural_network import MLPClassifier
 
#モデル定義と機械学習実行
mlpModel = MLPClassifier(activation='relu',solver='lbfgs',
                         hidden_layer_sizes=(500,500,500), random_state=0)
mlpModel.fit(X_train, Y_train)
 
#精度を表示
Y_train_score = metrics.accuracy_score(Y_train, mlpModel.predict(X_train))
Y_test_score = metrics.accuracy_score(Y_test, mlpModel.predict(X_test))
print("訓練精度: %0.4f 検証精度： %0.4f" % (Y_train_score, Y_test_score))

訓練精度93.47%・検証精度95.32%です。いい感じなので、このパラメータをベースにベストなハイパーパラメータを探索します。

最適な活性化関数・ソルバーの探索

まずは、活性化関数とソルバーから。パーセプトロンの数と階層数を先ほどの500個ずつ3階層に固定して、活性化関数とソルバーを総当たりで試してみます。デフォルトで用意されている活性化関数4種類とソルバー3種類をグリッドサーチ関数（指定したハイパーパラメータを総当たりで試してベストなものを教えてくれる便利な関数）に放り込みます。結構時間かかるので所要時間も表示するようにしています。

#グリッドサーチ1回目（活性化関数・ソルバーの探索）
from sklearn.model_selection import GridSearchCV
import datetime
 
#時間表示（開始）
import datetime
start = datetime.datetime.now()
print("開始：",start)
 
#シミュレーションするパラメータを指定（活性化関数・ソルバーを指定）
parameters = {'activation': ['logistic', 'identity', 'relu', 'tanh'],
              'solver': ['lbfgs', 'sgd', 'adam']}
 
#グリッドサーチの実行
mlpModel2 = GridSearchCV(MLPClassifier(hidden_layer_sizes=(500,500,500),
                                       random_state=0, max_iter=500,
                                       early_stopping=True), parameters)
 
mlpModel2.fit(X_train, Y_train)
 
#ベストなパラメーターの表示
print(mlpModel2.best_params_)
 
#時間表示（終了）
end = datetime.datetime.now()
print("終了：",end)

CPUをマックスに使って15分くらいかかりました。最適な活性化関数はIdentity・ソルバーはLbfgsとのことです。このハイパーパラメータの精度を見てみます。

#1回目のグリッドサーチで得られたパラメータで機械学習
mlpModel = MLPClassifier(activation='identity',solver='lbfgs',
                         hidden_layer_sizes=(500,500,500), random_state=0,
                         max_iter=500, early_stopping=True)
mlpModel.fit(X_train, Y_train)
 
#精度を表示
Y_train_score = metrics.accuracy_score(Y_train, mlpModel.predict(X_train))
Y_test_score = metrics.accuracy_score(Y_test, mlpModel.predict(X_test))
print("訓練精度: %0.4f 検証精度： %0.4f" % (Y_train_score, Y_test_score))

訓練精度96.23%・検証精度95.32%と先ほどから改善された結果になりました。

最適な階層数の探索

続いて階層数の探索です。活性化関数とソルバー・パーセプトロン数を固定して階層数を変動させます。具体的には3階層・4階層・5階層で試してみました。

#グリッドサーチ2回目（階層数の探索）
 
#時間表示（開始）
start = datetime.datetime.now()
print("開始：",start)
 
#シミュレーションするパラメータを指定（階層数を変える）
parameters = {'hidden_layer_sizes': [(500,500,500),(500,500,500,500),(500,500,500,500,500)]}
 
#グリッドサーチの実行
mlpModel2 = GridSearchCV(MLPClassifier(activation='identity', solver='lbfgs',
                                       random_state=0, max_iter=500, early_stopping=True), parameters)
 
mlpModel2.fit(X_train, Y_train)
 
#ベストなパラメーターの表示
print(mlpModel2.best_params_)
 
#時間表示（終了）
end = datetime.datetime.now()
print("終了：",end)

階層を増やした方が精度が上がるかと思いきや、3階層がベストとのこと。

最適なパーセプトロン数の探索

最後に最適なパーセプトロン数の探索をしてみます。500個と1000個でグリッドサーチしてみます。

#グリッドサーチ3回目（パーセプトロン数の探索）
 
#時間表示（開始）
start = datetime.datetime.now()
print("開始：",start)
 
#シミュレーションするパラメータを指定（パーセプトロン数を変える）
parameters = {'hidden_layer_sizes': [(500,500,500),(1000,1000,1000)]}
 
#グリッドサーチの実行
mlpModel2 = GridSearchCV(MLPClassifier(activation='identity', solver='lbfgs',
                                       random_state=0, max_iter=500, early_stopping=True), parameters)
 
 
mlpModel2.fit(X_train, Y_train)
 
#ベストなパラメーターの表示
print(mlpModel2.best_params_)
 
#時間表示（終了）
end = datetime.datetime.now()
print("終了：",end)

1000個にした方が良いと思いきや、500個の方が良いとの結果に。人力ハイパーパラメータ探索の結果としては、活性化関数はIdentity・ソルバーはLbfgs・パーセプトロン数と階層は500個を3階層ということになりました。ニューラルネットワークの精度は訓練精度96.23%・検証精度95.32%。

とりあえず今回はここまでで、次回は交差検証と多数決モデルを使って更なる精度向上にトライしようと思います。

今回のまとめ

データに含まれる情報量が多いほど機械学習の精度が高くなるのかもしれない（ソムリAIで使ったワインデータは1データあたり11個の情報だったが、乳がん診断情報は1データあたり30個の情報）。
決定木はニューラルネットワークなど他のモデルと比較すると精度が低いことが多いが、重要度の高い変数や条件分岐など機械学習の内部プロセスが目に見える点が優れている（医療など意思決定のプロセスが透明化されていなければならない領域にAIを適用する場合に有効だと思われる）。
ニューラルネットワークのハイパーパラメータ探索は探索を分割して絞り込んでいくのが良い（多くのパラメータを網羅的にグリッドサーチに掛けることも可能だが、計算量が雪だるま式に増えるためかなり時間がかかる。メモリがオーバーフローする可能性もある）。
乳がん診断AIの性能はかなりいい感じ（決定木の検証精度は94.74%・ニューラルネットワークの検証精度は95.32%）。