実践的（？）線形重回帰
Python を使って機械学習をしてみたい、という実験系化学者の方も多いと思います。
しかし、ググってみても機械学習の方法論とかは出てくるけど、痒いところに手が届かない、という事例は多いはずです（私もそうです）。
そこで、線形重回帰分析による LogS 値を予測するコードを例に、解析の方法を解説してみたいと思います。
間違っているぞ、もっとこうした方がいい、などのご意見がありましたら、コメント欄にご記入頂けますと幸いです。
下に、データの読み込みからモデルの評価まで行うコードを示します。
import numpy as np import pandas as pd from rdkit import Chem from rdkit.Chem.Draw import IPythonConsole from sklearn.linear_model import LinearRegression from sklearn.preprocessing import StandardScaler from sklearn.metrics import r2_score # データの読み込み X_train = pd.read_csv("../data/LogS_train_descriptors.csv", index_col = 0) X_test = pd.read_csv("../data/LogS_test_descriptors.csv", index_col = 0) y_train = pd.read_csv("../data/LogS_train_value.csv", index_col = 0) y_test = pd.read_csv("../data/LogS_test_value.csv", index_col = 0) # 文字列を返していない記述子を選択 X_train = X_train.select_dtypes(include = "number") #全て同じ値の記述子を削除 X_train = X_train[X_train.columns[X_train.std() != 0]] # 記述子名のラベルを取り出す descriptor_label = X_train.columns # 記述子間の相関係数を計算 corr_matrix = np.absolute(np.triu(np.corrcoef(np.array(X_train).T, ddof = 1), k = 1)) # 相関係数の絶対値が 0.6 以下の記述子のみを選択 descriptor_label = descriptor_label[np.max(corr_matrix, 0) <= 0.6] # ラベルを基にトレーニングデータの記述子を選択 X_train = X_train[descriptor_label] # テストデータはトレーニングデータと同じ記述子を選択 X_test = X_test[descriptor_label] # オートスケーリング scaler = StandardScaler() X_train_scaled = scaler.fit_transform(X_train) X_test_scaled = scaler.transform(X_test) # 線形重回帰モデル構築 clf = LinearRegression() clf.fit(X_train_scaled, y_train) y_pred = clf.predict(X_test_scaled) # R2 値を算出する r2_score(y_test, y_pred)
最初に、データの読み込みを行います。pd.read_csv で csv ファイルを読み込めますが、index_col = 0 を加えた方がいい場合が多いです。csv ファイルの 1 列目が行の番号だったりする場合、それ自体は記述子ではないため、それを行の名前として扱う操作です。そのため、1 列目から記述子が書かれている場合は、index_col を指定してはいけません。
読み込んだデータの中には、欠損値が nan や文字（string 型）で含まれることも多いです。欠損値の補完、削除方法は幾つかありますが、今回は単純に、全ての化合物に対して計算できている記述子のみを用います。select_dtypes(include = "number") で、pandas の DataFrame から、（int、float 問わず）数字のみで構成されている列のみを選択できます。
全て同じ値の記述子はデータを区別するのに役に立たず不要なので、削除します。
記述子間の相関係数が高い記述子は除きます。これは、数理統計学的には「共線性」と呼ばれる問題を回避するために行います。
数学的背景には立ち寄らず、イメージだけ述べると、x1 = 1, 2, 3, 4, 5 とx2 = 2, 4, 6, 8, 10 というデータがあり、y = w1*x1 + w2*x2 という式があったとします。y を作るのに、x2 がなくても x1 の係数 w1 の調整を行えば、y を算出することができます。x2 は必要なさそうですし、w1 と w2 を一義的に決めることができなくなってしまいます。そのため、相関係数の「絶対値」が大きい記述子を除きます（正の相関、負の相関があるため）。
具体的な操作としては、numpy 形式にして相関係数を算出し、上三角行列で相関係数の重複を除きます（k = 1 を指定することで対角成分を除去できる。対角成分は必ず 1）。そこから最大値を取ってきて、これが閾値以下の値を示す記述子のみを選びます。
テストデータは、トレーニングデータと同じ記述子を選べば良いので、計算を行う必要はなく、記述子ラベルを元に列を選択する操作だけを行います。
化合物間におけるトレーニングデータの記述子の値を平均 0、 標準偏差を 1 に変換します。記述子の絶対値の大きさがモデルに影響することを防ぐために行います。テストデータは、トレーニングデータの平均と標準偏差に合わせて処理を行うので、通常、X_train_scaled の平均と標準偏差はそれぞれ 0 と 1 にはなりません。
いよいよ線形重回帰モデルの構築を行うのですが、scikit learn は非常に便利で、3 行でモデルの構築と予測が行えてしまいます。最初に、モデルのインスタンスを定義してあげて、次にトレーニングデータを入力してモデルを構築します。最後に、テストデータを入力することで、予測が行えてしまいます。Scikit learn の機械学習モデルは、どの手法も基本的にこの形式に準拠しているので、慣れると楽です。
最後に、モデルの性能を R2 値で見てみます。式は他の書籍やサイトに譲りますが、値は 1 から負の無限大まで取り得ます。たまに R2 値が -2 になった、とかで慌てることがあるかもしれませんが、それはバグなどではなく、単にモデルの性能が恐ろしく低いだけなので安心しましょう（？）。ちなみに私は -70 とかを叩き出したことがあります。今回のモデルでは R2 値が 0.893 と、中々良好な結果が得られました。
基本的な流れを解説してみましたが、トレーニングデータとテストデータの分割など、細かい部分で未解説の部分もあるので、後日追記するかもしれません。
