データ化学工学研究室(金子研)では、新しく配属になった学生にいろいろなトレーニングをしています。その１つがPython言語のトレーニングです。化学構造を扱ったりデータ解析・機械学習をしたりするときに、Python言語を使うわけです。
Python言語は無料で使えます。金子研では最初に Anaconda をインストールしています。Anaconda は Python 本体だけでなく、データ解析・機械学習をするときによく用いるライブラリ (ツール) がパッケージとして一緒に入っていて便利です。なお、いろいろな理由で Anaconda をインストールできない方は、以下で Anaconda を使わずにPython でデータ解析・機械学習する方法を紹介していますので、ご参照ください。
最初に Python やライブラリの説明を簡単にしたあと、学生たちは一人ずつ下に示す課題を順番に行います。課題をクリアしていくと、最終的に化学構造を扱ったりデータ解析をしたりするための Python によるプログラミング能力が身につく、というわけです。基本的に学生たちはwebや本で調べながら自分でクリアを目指します。ただし質問はいつでもOK。課題が解けたら、slackにあげてもらい、わたしがチェックします。問題や改善点があれば指摘し、クリアできたら次の課題へ、という流れです。
学生たちには、卒業までにデータサイエンティストとして中級レベルになりましょう、という話をしています。中級者は自分で考え学習しながら進められるレベルです。このプログラミングのトレーニングを経て、中級者に一気に近づくことでしょう。
応用化学科の学生ですので、化学とデータサイエンスの二刀流になるわけです。
ちなみに学生たちは、”Python”を聞いたことがありませんでした。プログラミングもほぼはじめての人たちです。そういったプログラミングの初学者が、今は課題に取り組みながらメキメキと力をつけているところです。
ちなみに、プログラミングをするときに最初に身につけるべき３つのテクニックはこちらです。
MATLAB を使ったことある人は、こちらにも目を通すと Python を始めやすいです。
データセット
課題を示す前に、課題に取り組むために必要なデータセットを載せます。プログラムと同じディレクトリ(フォルダ)に置いてください。
- あやめのデータ (speciesなし)
- あやめのデータ (speciesあり)
有名なFisher’s Iris Dataです。150個のあやめについて、がく片長 (Sepal Length)、がく片幅 (Sepal Width)、花びら長 (Petal Length)、花びら幅 (Petal Width) が計測されています。最初の方の課題ではこちらを使用します。最初は species なしのデータで教師なし学習、後に species ありのデータで教師あり学習のなかでもクラス分類のトレーニングをします。
E. Anderson, The species problem in Iris, Annals of the Missouri Botanical Garden, 23(3), 457–509, 1936. - 1290化合物の水溶解度データ(sdfファイル)
- 同じ化合物の化学構造に対するSMILES (Simplified Molecular Input Line Entry Specification syntax)
- 同じ化合物の化学構造に対してRDKitで196の構造記述子を計算したファイル
T. J. Hou らの論文にある 1290 の化合物に対する水溶解度の log 値のデータセットです。定量的構造活性相関 (Quantitative Structure-Property Relationship, QSPR) に関係する研究でよく用いられます。sdf ファイルと SMILES 形式のファイル、そして化学構造に対して RDKit で196の構造記述子を計算した csv ファイル (一番左に logS 値あり) です。
T. J. Hou, K. Xia, W. Zhang and X. J. Xu, ADME Evaluation in Drug Discovery. 4. Prediction of Aqueous Solubility Based on Atom Contribution Approach, J. Chem. Inf. Comput. Sci., 44(1), 266–275, 2004.
プログラミング課題と模範解答
以下に学生に取り組んでもらっている課題を載せます。課題は逐次更新します。模範解答(というには恐れ多い)プログラムも一緒に置いておきます。参考にしてください。できればコメントをいただけると嬉しいです。
課題1: iris.csv を読み込み各変数の最大値・最小値・平均値・中央値・分散・標準偏差を求め basic_statistics.csv に保存するプログラムを作成せよ。自分なりの方法で各統計量があっているか確認せよ。
課題2: iris.csv を読み込みすべての変数間の共分散と相関係数を求めそれぞれ covariance.csv, correlationcoefficient.csv に保存するプログラムを作成せよ。自分なりの方法で各統計量があっているか確認せよ。
課題3: まずは変数の標準化 (オートスケーリング) について調べよ。次に、iris.csv を読み込みオートスケーリングを行って autoscaled_iris.csv に保存するプログラムを作成せよ。自分なりの方法で結果があっているか確認せよ。
[補足] 模範解答では、rawdata は (pandasの) DataFrame 型なので、mean や std において axis=0 や ddof=1 はなくても問題ありませんが、仮に rawdata が (numpyの) array 型であったら、axis=0 や ddof=1 がないと標準化(オートスケーリング)されません。そのため、間違って array 型がきてもいいように、あえて axis=0 や ddof=1 を書いています。
課題4: scikit-learn を使って主成分分析 (Principal Component Analysis, PCA) を行うプログラムを作成せよ。寄与率と累積寄与率を図示できるようにしておくこと。
課題5: 乱数を使って変数間の相関係数の高い２変数のデータを適当に作成し、主成分分析を行い、ローディングの確認をせよ。データの前処理としてオートスケーリングを行うこと。主成分分析の前後でデータをプロットして考察せよ。
課題6: iris.csv を読み込みオートスケーリングをしてから主成分分析を行いデータの可視化をせよ。可視化の様子からいえることを考察せよ。
課題7: iris.csv を読み込みオートスケーリングをしてから、scipy を使ってウォード法による階層的クラスタリングを行え。クラスター数を変えて、結果を主成分分析 (PCA) により可視化をせよ。可視化の様子からいえることを考察せよ。iris_with_species.csv も参考にすること。
課題8: iris_with_species.csv を読み込みオートスケーリングをしてから、scikit-learn を使って線形判別分析 (Linear Discriminant Analysis, LDA) を行え。この際、setosa と versicolor とを合わせたものを１つのクラス、virginica をもう一つのクラスとして、２クラス分類をすること。さらに、実際のクラスと計算されたクラスとの間で混同行列を計算せよ。
課題8のプログラム (こちらでは、setosa, versicolor+virginica になっています)
課題9: logSdataset1290.csv を読み込み、オートスケーリングしてから logS と構造記述子との間で最小二乗法による重回帰分析を行え。計算できないもしくは計算結果がおかしい場合は、その原因について調べて対処法を実行せよ。r2 の計算、logS の実測値と計算値とのプロット、標準回帰係数の standard_regression_coefficients.csv への保存を実施せよ。
課題10: logSdataset1290.csv を読み込み、オートスケーリングしてから logS と構造記述子との間で PLS を行え。成分数は2とする。r2 の計算、logS の実測値と計算値とのプロット、標準回帰係数の standard_regression_coefficient.csv への保存を実施せよ。課題9の結果と比較して考察せよ
課題11: logSdataset1290.csv を読み込み、オートスケーリングしてから logS と構造記述子との間で PLS を行え。50 成分までの成分数ごとの r2・r2cv のプロット、最適成分数における logS の実測値と計算値とのプロットおよび logS の実測値とクロスバリデーション推定値とのプロットを作成せよ。クロスバリデーションは 2-fold クロスバリデーションとし、最適成分数は r2cv が最大のものとする。
課題12: 課題11で 50 成分までに r2cv の値が小さい(たとえば負の値になる)ことがあるときはその原因を調べて対応し、r2cv を向上させよ。ただし、化合物は削除しないこととする。
課題13: iris_with_species.csv を読み込み、ランダムに90サンプルを選択してトレーニングデータ (モデル構築用データ) とし、残りをテストデータ (モデル検証用データ) として、LDA モデルの作成およびモデルを用いた予測を行え。この際、setosa と versicolor とを合わせたものを１つのクラス、virginica をもう一つのクラスとして、２クラス分類をすること。またオートスケーリングを行うこと。さらに、トレーニングデータとテストデータそれぞれにおいて、実際のクラスと計算もしくは予測されたクラスとの間で混同行列を計算せよ。
課題13のプログラム (こちらでは、setosa, versicolor+virginica になっています)
課題14: logSdataset1290.csv を読み込み、1番目から878番目までのサンプルをトレーニングデータ (モデル構築用データ) とし、残りをテストデータ (モデル検証用データ) として、PLS モデルの作成およびモデルを用いた予測を行え。クロスバリデーションは 5-fold クロスバリデーションとし、最適成分数は50成分までで r2cv が最大のものとする。また説明変数も目的変数もオートスケーリングを行うこと。さらに、トレーニングデータとテストデータそれぞれにおいて、実測値と計算値もしくは予測値との間で r2, RMSE (Root-Mean-Square Error), MAE (Mean Absolute Error) を計算し、実測値と計算値もしくは予測値とのプロットを作成せよ。
課題15: 課題14と同様に PLS モデルを構築せよ。ただし、1番目から100番目までのサンプルをトレーニングデータ (モデル構築用データ) とし、残りをテストデータ (モデル検証用データ) とすること。k最近傍法 (k-nearest neighbor, kNN) によりモデルの適用範囲 (Applicability Domain, AD) を設定せよ。そして k の値を変えながら、AD の内側と外側とでテストデータの r2, RMSE, MAE を計算して比較し、考察せよ。また実測値と予測値とのプロットも確認すること。
課題16: logSdataset1290.csv を読み込み、1番目から878番目までのサンプルをトレーニングデータとし、残りをテストデータとして、PLS・リッジ回帰・LASSO・Elastic net・SVR (線形カーネル)・SVR (ガウシアンカーネル)・ランダムフォレストモデルをそれぞれ作成し、テストデータの予測を行え。クロスバリデーションは 5-fold クロスバリデーションとする。説明変数も目的変数もオートスケーリングを行うこと。すべてのモデルにおいて、トレーニングデータとテストデータそれぞれに対して、実測値と計算値もしくは予測値との間で r2, RMSE, MAEを計算し、実測値と計算値もしくは予測値とのプロットを作成せよ。モデル間で結果を比較し、考察せよ。
課題17: iris_with_species.csv を読み込み、ランダムに90サンプルを選択してトレーニングデータ (モデル構築用データ) とし、残りをテストデータ (モデル検証用データ) として、SVM の作成および SVM を用いた予測を行え。この際、setosa と versicolor とを合わせたものを１つのクラス、virginica をもう一つのクラスとして、２クラス分類をすること。またオートスケーリングを行うこと。さらに、トレーニングデータとテストデータそれぞれにおいて、実際のクラスと計算もしくは予測されたクラスとの間で混同行列を計算せよ。カーネル関数として線形カーネルとガウシアンカーネルを用いて結果を比較すること。ハイパーパラメータは 5-fold クロスバリデーション後の正解率が高くなるように設定すること。
課題17のプログラム (こちらでは、setosa, versicolor+virginica になっています)
課題18: RDKit を活用して、logSdataset1290_2d.sdf を読み込み、各化学構造において複数の構造記述子を計算して descriptors.csv に保存するプログラムを作成せよ。ただし各サンプル名を化学構造の SMILES 文字列にすること。
課題19: RDKit を活用して、logSdataset1290_2d.sdf を読み込み、各化学構造において fingerprint を計算して fingerprints.csv に保存するプログラムを作成せよ。ただし各サンプル名を化学構造の SMILES 文字列にすること。いろいろな種類の fingerprint を計算して確認すること。
課題20: RDKit を活用して、BRICSアルゴリズムにより仮想的な化学反応にもとづいて化学構造を生成して、生成した構造を generatedstructures.sdf に保存するプログラムを作成せよ。
課題21: logSdataset1290.csv を読み込み、1番目から878番目までをトレーニングデータとし、残りをテストデータとして、GPモデルを作成し、モデルを用いたテストデータの予測を行え。その際、化合物ごとに推定値だけでなく推定値の標準偏差も計算し、それに基づいて実際の logS の値が1.5から2の間に入る確率を計算すること。この確率が最大になる化合物と logS の推定値が最大になる化合物とを比較し、考察せよ。
発展課題1: scikit-learn を使わずに、固有値問題に基づく方法と NIPALS アルゴリズムによる方法の2通りで PCA を行うプログラムを作成せよ。scikit-learn の PCA の結果と等しくなることを確認せよ。
発展課題2: オートスケーリングをする場合としない場合とで主成分分析の結果を比較し、結果を考察せよ。
発展課題3: scikit-learn を使わずに、NIPALS アルゴリズムによる方法で PLS を行うプログラムを作成せよ。scikit-learn の PLS の結果と等しくなることを確認せよ。
発展課題4: scikit-learn を使わずに、SIMPLS アルゴリズムによる方法で PLS を行うプログラムを作成せよ。scikit-learn の PLS の結果と等しくなることを確認せよ。
以上です。
バグや不具合などありましたら、twitter・facebook・メールなどを通して教えていただけるとうれしいです。
