fantom_zona’s diary

合成生物学(と最近では神経変性疾患)

deepchem動かしてみました①

deepchemというライブラリを勉強してみようと思います。チュートリアルからまずやっていこうと思います。日本語の記事が少ないので、自分用のメモとして記録に残そうと思います。deepchemの溶解度予測のチュートリアルを一個一個読んでいこうと思います。

データの読み込み

%load_ext autoreload
%autoreload 2
%pdb off

しょっぱなのautoreloadがまず不明でしたが、このページに解説がありました。チュートリアルの方では%load_ext autoreloadはありませんでしたが、これを加えるとipython notebookではエラーが出ないようです。
%pdb offの方は不明です。pdbはデバッガらしいですが......誰かわかったら教えてください。

from deepchem.utils.save import load_from_disk

dataset_file= "./deepchem/datasets/delaney-processed.csv"
dataset = load_from_disk(dataset_file)
#ソースコードを見るとload_from_diskはcsvじゃなくても読み込める関数
dataset.head()

f:id:fantom_zona:20180322070209p:plain
読み込んだcsvの形式はこういう形のようです。今回行いたいのは、smiles形式から溶解度の予測をすることなので、それ以外のものは全て使いません。

描画

import tempfile
from rdkit import Chem
from rdkit.Chem import Draw
from itertools import islice
from IPython.display import Image, HTML, display

def display_images(filenames):
    """Helper to pretty-print images."""
    imagesList=''.join(
        ["<img style='width: 140px; margin: 0px; float: left; border: 1px solid black;' src='%s' />"
         % str(s) for s in sorted(filenames)])
    display(HTML(imagesList))
#空白にpngをjoinで足していってHTML表示する
    
def mols_to_pngs(mols, basename="test"):
    """Helper to write RDKit mols to png files."""
    filenames = []
    for i, mol in enumerate(mols):
        filename = "%s%d.png" % (basename, i)
        Draw.MolToFile(mol, filename)
        filenames.append(filename)
    return filenames
#Draw.MolToFile: mol object->pngにする
#ファイル名のリストにする

num_to_display = 8
molecules = []
for _, data in islice(dataset.iterrows(), num_to_display):
    molecules.append(Chem.MolFromSmiles(data["smiles"]))
display_images(mols_to_pngs(molecules))
#MolFromSmiles: smiles->mol object
#display_imagesはファイル名をしているすると表示してくれる関数

RDKitではmol objectというのを生成して、そこから描画ができるようです。
上の流れでは、
Chem.MolFromSmiles: smiles→mol object
Draw.MolToFile : mol object→png
として画像を出力し、display_images(filename)で画像を表示しています。
f:id:fantom_zona:20180324080224p:plain

特徴量の取得

このチュートリアルでは、ECFP4という手法で化合物の特徴量を取得しています。
化合物の局所的な構造を元にした特徴量をfingerprintというらしく、ECFPは主だった特徴量として有名なようです。

import deepchem as dc

featurizer = dc.feat.CircularFingerprint(size=1024)

loader = dc.data.CSVLoader(
      tasks=["measured log solubility in mols per litre"], smiles_field="smiles",
      featurizer=featurizer)
dataset = loader.featurize(dataset_file)

dc.data.CSVLoader()は、tasksとsmikes_filedを指定するとcsv形式データを取り込むclassですが、methodとしてfeaturizeが用意されています。featurizeは、事前に定義したfeaturizerを利用してsmiles形式から勝手に特徴量を取得してくれるため便利です。
CSVloaderで出力されるオブジェクトは、deepchemではDatasetオブジェクトとして扱われています。Datasetオブジェクトは、(X, y, w, ids) という値を持っています。Xは特徴量、yは出力、wは重み、idsはID(今回の場合はsmiles形式を返す)です。deepchemのDatasetに関するチュートリアルもあるので今後それも記事にできたらと思います。

前処理

データをtrain/splitに分け、前処理は正規化を行います。

splitter = dc.splits.ScaffoldSplitter(dataset_file)
train_dataset, valid_dataset, test_dataset = splitter.train_valid_test_split(dataset)

transformers = [
    dc.trans.NormalizationTransformer(transform_y=True, dataset=train_dataset)]

dc.splits.ScaffoldSplitterでは、分子構造の偏り?が生まれないようにうまくsplitしてくれているらしいです。他にもrandom splitなども選択できるようです。

学習部分

deepchemはsklearnをラップしています。sklearnの異なるモデルに対して、SklearnModelに突っ込めば、fit(dataset object)で学習してくれます。学習したモデルはsave()で保存してくれます。save()では、sklearnのjoblib.dumpが背後で動いてくれます。

from sklearn.ensemble import RandomForestRegressor

sklearn_model = RandomForestRegressor(n_estimators=100)
model = dc.models.SklearnModel(sklearn_model)
model.fit(train_dataset)

評価部分

from deepchem.utils.evaluate import Evaluator

metric = dc.metrics.Metric(dc.metrics.r2_score)
#r2_scoreはdc.metrics内でimportしてあるので呼び出せる

evaluator = Evaluator(model, valid_dataset, transformers)
#class Evaluator(model, dataset, transformers, verbose=False)

r2score = evaluator.compute_model_performance([metric])
print(r2score)

Evaluator.compute_model_performance()で、metricの計算をしてくれます。
適当に決めた木の数(=100)では、決定係数は0.2くらいでした。

チューニング

木の数と特徴量の数を適当に振り、最適なものを選び出します。deepchem.hyper.grid_search moduleにあるHyperparamOptでグリッドサーチしてくれます。

def rf_model_builder(model_params, model_dir):
    sklearn_model = RandomForestRegressor(**model_params)
    return dc.models.SklearnModel(sklearn_model, model_dir)
params_dict = {
    "n_estimators": [10, 20,40, 80, 160, 320],
    "max_features": ["auto", "sqrt", "log2", None],
}

metric = dc.metrics.Metric(dc.metrics.r2_score)
optimizer = dc.hyper.HyperparamOpt(rf_model_builder)
#ただのclass
#methodでhyperparam_searchが使用できる
#dictionaryでパラメータを指定すれば検索してくれる
best_rf, best_rf_hyperparams, all_rf_results = optimizer.hyperparam_search(
    params_dict, train_dataset, valid_dataset, transformers,
    metric=metric)

木は10個、特徴量の数はルートが最適でした。R2の値は0.3くらいでした。
結果を出力してみるとそこそこ学習しているかなといったところでしょうか。ただし、大小関係のみが保たれている傾向があり、実際の値の予測はだめそうです。
deepchemのチュートリアルではランダムフォレストに加えて、ニューラルネットも使っていましたが、deepじゃなさそうなのにdeepって書いてあり、しかも使っている関数がdeepchem 2.0.0では削除されているようなので割愛しました。

task = "measured log solubility in mols per litre"
predicted_test = best_rf.predict(test_dataset)
true_test = test_dataset.y
plt.scatter(predicted_test, true_test)
plt.xlabel('Predicted log-solubility in mols/liter')
plt.ylabel('True log-solubility in mols/liter')
plt.title(r'RF- predicted vs. true log-solubilities')
plt.show()

f:id:fantom_zona:20180324080257p:plain

感想

以上、csvから予測までのチュートリアルを動かしてみたという話でした。
少しだけdeepchemが使えるようになった気がしてきました。今回はグラフ構造の学習をしなかったので、グラフ構造の学習とかもやっていきたいと思いました。