Cross validation + OptunaでMLPのモデル選択してみた話
やったこと
とあるテーブルコンペに出た時に, sklearnのMLPRegressorを使ってみました.
完全に素人がやってみた内容なので, 是非ともアドバイスください.
コード内容
optunaは最小化のためのパラメータ調整ライブラリらしいです.
sklearnのKFoldと合わせて使ってみました.
import pandas as pd import numpy as np from sklearn.model_selection import KFold from sklearn.neural_network import MLPRegressor from sklearn.metrics import r2_score import lightgbm as lgb import optuna
以下で, objectという関数を定義します. objectはパラメータを含む関数で, 関数内で調整するパラメータを明示的に記述する必要があります. optunaはobjectの出力を最小化します. 今回は決定係数でモデルの評価をしてみました. "X_"がinput, "y_"というのがoutputです.
kf = KFold(n_splits = 5, shuffle = True, random_state = 0) def objective(trial): ''' trial:set of hyperparameter ''' # hyper params max_layer_size = trial.suggest_int('hidden_layer_sizes', 50, 300) activation_func = trial.suggest_categorical('activation', ['relu', 'tanh']) alpha = trial.suggest_loguniform('alpha', 1e-5, 1e-3) # model model = MLPRegressor(hidden_layer_sizes=(max_layer_size, ), activation = activation_func, alpha=alpha) MLP_scores = [] for train_index, test_index in kf.split(X_): X_train, X_test = X_[train_index], X_[test_index] y_train, y_test = y_[train_index], y_[test_index] model = MLPRegressor(hidden_layer_sizes=(max_layer_size, ), activation = activation_func, alpha=alpha) model.fit(X_train, y_train) y_pred = model.predict(X_test) MLP_scores.append(r2_score(y_test, y_pred)) print("max_layer_size :", max_layer_size) print("activation_func :", activation_func) print("alpha :", alpha) print("result :", r2_score(y_test, y_pred)) score =-np.mean(MLP_scores) #optunaは最小化なので符号を反転させる. return score
関数を定義後, optunaを動かします.
study = optuna.create_study() study.optimize(func=objective, # 実行する関数 n_trials=20, timeout=None, #default: None n_jobs=-1 #multiprocess )
PythonでRNAlibを動かすときのメモ
RNAlibのdocumentなりhelpがあまりにも欠如しているのでちまちま調べていたことを軽くまとめる。
githubのissueを見ると割と勉強になることを学んだ。
温度の設定
import RNA print(RNA.__version__) #2.4.11 seq = RNA.random_string(20, "AUGC") md = RNA.md() print(md.temperature) #デフォルトで37.0 fc_37 = RNA.fold_compound(seq, md) md.temperature = 20 fc_20 = RNA.fold_compound(seq, md) md.temperature = 4 fc_4 = RNA.fold_compound(seq, md) print(fc_37.mfe()) print(fc_20.mfe()) print(fc_4.mfe())
制約条件下のfolding
制約条件の下でのfoldingの記載も見つけた。
constraint folding in python · Issue #31 · ViennaRNA/ViennaRNA · GitHub
constraint = ".......................xxxxxxx" #制約条件 sequence = RNA.random_string(len(constraint), "AUGC") #適当に生成 # create fold_compound object fc = RNA.fold_compound(sequence) # add hard constraints, specified by dot-bracket like string fc.hc_add_from_db(constraint) # compute partition function (prob, energy) = fc.pf() fc.hc_init() #初期化
duplexfoldの使い方について
fold_compoundを使うのが新しいRNAlibの使い方らしいが、duplexfoldみたいな古いもの未だに残っているらしい。
&を使ったstructureの出力が謎だったが、これを読んではっきり理解することが出来た。
Matching duplexfold structure to sequence · Issue #38 · ViennaRNA/ViennaRNA · GitHub
d = RNA.duplexfold('CGT', 'AAG') print(d.structure) d.energy print(d.i, d.j) #(.&.) #2 2
(.&.) の &の前後の番号を示しているらしい。なるほど。流石にdocumentに書け。
suboptの使い方
suboptの引数の数字deltaはx0.01 kcal/molの間にあるsuboptの構造を出力してくれる
a = RNA.fold_compound(sequence) solution = a.subopt(500) print "%d suboptimals" % len(solution) for s in solution: print "%s %6.2f" % (s.structure,s.energy)
deepchem.models.sequentialではグラフ構造を読み取れなかった話
注意
この記事を読んでも特に得るものはありません。ただdeepchem2.1.0ではsequentialにグラフ畳み込みが出来なかったというだけの話です。
概要
deepchem 2.1.0
python2.7
model=deepchem.models.Sequential() """ model.fit() → model._create_graph(self, feature_shape, label_shape) shapeは勝手に推論してくれる. → model._layer_listをfor文で回して各layer.in_layerに1個前のlayerとして定義していく. 最終的なlayerをoutputして保存. → TensorGraph.built()で各layerがbuild()される. →以下略 """
model._create_graphの中身
def _create_graph(self, feature_shape, label_shape): """This is called to create the full TensorGraph from the added layers.""" if self.built: return # The graph has already been created. # Add in features features = Feature(shape=feature_shape) # Add in labels labels = Label(shape=label_shape) # Add in all layers prev_layer = features if len(self._layer_list) == 0: raise ValueError("No layers have been added to model.") for ind, layer in enumerate(self._layer_list): if len(layer.in_layers) > 1: raise ValueError("Cannot specify more than one " "in_layer for Sequential.") layer.in_layers += [prev_layer] prev_layer = layer # The last layer is the output of the model self.outputs.append(prev_layer) if self._loss_function == "binary_crossentropy": smce = SoftMaxCrossEntropy(in_layers=[labels, prev_layer]) self.set_loss(ReduceMean(in_layers=[smce])) elif self._loss_function == "mse": mse = ReduceSquareDifference(in_layers=[prev_layer, labels]) self.set_loss(mse) else: # TODO(rbharath): Add in support for additional # losses. raise ValueError("Unsupported loss.") self.build()
debugの試み(まだ何も気づいていない時)
model=Sequential(loss="binary_crossentropy") model.add(GraphConv(out_channel=2)) model.add(GraphGather(32, activation_fn="sigmoid")) model.fit(train, nb_epoch=1)
以下がエラー
TypeError Traceback (most recent call last) <ipython-input-151-eccebc8e26f4> in <module>() 2 model.add(GraphConv(out_channel=2)) 3 model.add(GraphGather(32, activation_fn="sigmoid")) ----> 4 model.fit(train, nb_epoch=1) 5 ~/deepchem/models/tensorgraph/sequential.pyc in fit(self, dataset, **kwargs) 78 """ 79 X_shape, y_shape, _, _ = dataset.get_shape() ---> 80 self._create_graph((None,) + X_shape[1:], (None,) + y_shape[1:]) 81 super(Sequential, self).fit(dataset, **kwargs) 82 ~/deepchem/models/tensorgraph/sequential.pyc in _create_graph(self, feature_shape, label_shape) 114 raise ValueError("Unsupported loss.") 115 --> 116 self.build() 117 118 def make_estimator(self, ~/deepchem/models/tensorgraph/tensor_graph.pyc in build(self) 665 for layer in self.topsort(): 666 with tf.name_scope(layer.name): --> 667 layer.create_tensor(training=self._training_placeholder) 668 self.rnn_initial_states += layer.rnn_initial_states 669 self.rnn_final_states += layer.rnn_final_states ~/deepchem/models/tensorgraph/layers.pyc in create_tensor(self, in_layers, set_tensors, **kwargs) 2586 b_list = self.variables[self.num_deg:] 2587 else: -> 2588 W_list, b_list = self._create_variables(in_channels) 2589 2590 # Extract atom_features ~/deepchem/models/tensorgraph/layers.pyc in _create_variables(self, in_channels) 2564 initializations.glorot_uniform( 2565 [in_channels, self.out_channel], name='kernel') -> 2566 for k in range(self.num_deg) 2567 ] 2568 b_list = [ ~/deepchem/models/tensorgraph/initializations.pyc in glorot_uniform(shape, name) 59 def glorot_uniform(shape, name=None): 60 fan_in, fan_out = get_fans(shape) ---> 61 s = np.sqrt(6. / (fan_in + fan_out)) 62 return uniform(shape, s, name=name) 63 TypeError: unsupported operand type(s) for +: 'NoneType' and 'int'
どうやら_create_graphが出来てないみたいなので少し確認.
layerのリスト表示.
model._layer_list
[<deepchem.models.tensorgraph.layers.GraphConv at 0x1a395577d0>, <deepchem.models.tensorgraph.layers.GraphGather at 0x1a32d2f8d0>]
さらに、in_layer, outputsの表示. エラー的にはどこかのin_layerの定義がおかしそう?
print(model._layer_list[0].in_layers) print(model._layer_list[1].in_layers) print(model.outputs)
[<deepchem.models.tensorgraph.layers.Feature object at 0x1a3f89f450>] [<deepchem.models.tensorgraph.layers.GraphConv object at 0x1a3c72f3d0>] [<deepchem.models.tensorgraph.layers.GraphGather at 0x1a3c72f990>]
ソースコードを確認してみるとどうやらGraphConvのテンソル構築がうまくいっていないようです.
gc= model._layer_list[0] inputs = gc._get_input_tensors(gc.in_layers[0]()) #[<tf.Tensor 'Placeholder_20:0' shape=(?,) dtype=float32>] in_channels = inputs[0].get_shape() #TensorShape([Dimension(None)])
どうやらGraphConvに入れるサンプルの次元が無いよ、ということのようです。
そういえば、Sequentialのatom_featureの引数が無い......?
どうやらSequentialにはグラフ組めない......
ということで色々調べていたら、deepchem1.4まではSequentialGraphといってkearsっぽく書けたみたいなんですが、今となっては普通のnerural networkしか無理みたいですね. 以前触っていた時はdeepchem 1.4とかだったので気づいたら削除されていたようです. 悲しい. SequentialGraph絶対便利なのになんで消してしまったんだろうか。。。
deepchemのインストールについてのメモ
deepchemで遊んでいるうちに環境を壊してしまい、改めてインストールすることになりました。なのですが、これまでpython3系で動いていたはずなのにうまく行きません。deepchemのWebページでは3系でインストールできるみたいな雰囲気をかましているので、メモとして記録します。
2018年の4月くらいからこの問題は生じているみたいなのですが、ずっと解決していないようです。どうやっても3.5で動かすことができなかったので泣く泣く2.7で導入しました。
化合物-タンパク質相互作用の予測器を動かしてみた
colabでの環境構築
rdkitを入れなくてはいけないのでやむなくminicondaを入れました。
conda-forgeから入れるとこのバグは生じませんでした。
!wget -c https://repo.anaconda.com/miniconda/Miniconda3-latest-Linux-x86_64.sh !chmod +x Miniconda3-latest-Linux-x86_64.sh !bash ./Miniconda3-latest-Linux-x86_64.sh -b -f -p /usr/local conda install -c conda-forge rdkit pytorch scikit-learn
実行手順
あとはgithubに書いてある通りです。
wgetでフォルダを入れた後にCPI_prediction/codeまで降ります。ちなみに、colabはディレクトリの移動は%cdじゃなくてはいけないようですので注意してください。なぜかは知りません。
!mkdir cpi %cd cpi !wget -c https://github.com/masashitsubaki/CPI_prediction/archive/master.zip !unzip master.zip %cd CPI_prediction-master/code
!bash preprocess_data.sh !bash run_training.sh
結果
githubに書いてある通りです。
参考までに、epoch50でAUC=0.9596、Precision=0.9240、Recalll=0.8827です。
100 epoch回し終わるまでに1時間くらいだったかな?
びまん性正中グリオーマに対するCAR-T療法
この記事は、今年読んだ一番好きな論文2018の12月15日の記事です。
adventar.org
この一年は色々面白い話を見聞きしたのですが、いざ一番面白いのを選んで紹介するとなると、色々それぞれの良さがあるわけで選定するのって難しいですね*1。そういうわけで、色々悩んだ挙句、Potent antitumor efficacy of anti-GD2 CAR T cells in H3-K27M+ diffuse midline gliomasを紹介しようと思います。Bench-to-Bedside感あるやつで僕はCAR-T好きなので選定しました。
15秒で終わる紹介
・H3-K27M変異*2のびまん性正中グリオーマではGD2が高発現していることを発見した。
・GD2に対するCAR-T細胞を作成した。
・同所性に腫瘍細胞を移植されたマウスに抗GD2 CAR-T細胞を投与したところ、腫瘍細胞が除去された。
・一部マウスでは抗腫瘍効果による炎症反応で水頭症*3を発症し、致死的であった。
導入
まず、「キメラ抗原受容体T細胞(chimeric antigen receptor T cell, CAR-T細胞)」というワードと「びまん性正中グリオーマ(diffusive midline glioma, DMG)」というワードから説明したいと思います。
近年、B細胞性の急性リンパ球性白血病(B-ALL)に対する抗CD19 CAR-T細胞を始めとして、クローナル?*4な病変に対するChimeric antigen receptor T cell(CAR-T細胞)の有効性が注目されてきています。少し前にノバルティスからキムリアが出されて話題になっていたことを記憶している人もいるでしょう。CAR-T細胞のコンセプトは、細胞外領域に単鎖抗体、細胞内ドメインにCD3ζをくっ付けた受容体をT細胞に発現させることで、任意の細胞外抗原依存的にT細胞を活性化させるというものです。したがって、腫瘍細胞特異的な膜タンパクに対する人工受容体を作成・導入すれば、腫瘍細胞に対する炎症反応が惹起され、腫瘍細胞は除去されるだろうというのがCAR-T療法です*5。
(https://www.amed.go.jp/news/release_20171107.html より拝借)
分子標的薬など治療法の発展によって様々な腫瘍の生存期間が延長している一方、DMGという腫瘍ではそうではありません。こちらのページ*6を参照して、論文紹介に必要な項目だけ抜き出すと、DMGは、
・橋を中心とした脳幹部にびまん性に発生する小児腫瘍
・現時点の有効な治療法は、緩和的な放射線照射のみ
・予後は極めて不良
という腫瘍です。DMGは橋に発生することが多く、この場合、びまん性内在性橋グリオーマ(DIPG)という名前が付いています。橋には生命維持に重要な中枢が多く存在するために、切除することはできません。DMGは小児脳腫瘍の10%を占めているにも関わらず、DMGに対する治療法はほとんど進歩していないため、長期生存が見込める様な革新的な治療法が長らく求められています。
以上を踏まえて、DMGに対するCAR-T細胞を作成すれば新たな治療戦略になるのでは?という疑問に対する研究を今回紹介します。
結果
びまん性内在性橋グリオーマではGD2が高発現している
著者たちは、まず、DIPGの検体を用いて抗原を探しました。抗体のライブラリを用いて抗原になりそうな表面抗原を同定したところ、GD2が複数の検体に共通して最も高発現でした。GD2はdisialogangliosideというガングリオシドの仲間だそうです。また、免疫染色でもH3-K27M陽性の細胞にGD2が発現していること、H3-K27M陰性DIPG(VUMC-DIPG10)や高分化型グリオーマ(SU-pcGBM2)ではGD2の発現が弱いことが確認できました。
抗GD2 CAR-T細胞は、DIPGに対してGD2依存的に細胞障害性を示す
GD2が良い候補であったため、次にこれに対するCAR-T細胞がワークするのかを検討しました。作成したCAR-T細胞は、vitroでGD2陽性DIPGを除去することが分かりました。また、GD2依存的な活性化であることを、GD2陽性/陰性含む複数の検体やGD2合成酵素をノックアウトした細胞を作成して確認しています。人工受容体のネガティブコントロールとして抗CD19に対するCARを使用しています。GD2は当然正常組織にも発現しているため、正常組織に対しても細胞障害性を引き起こすように思われます。今回の人工受容体と同様のデザインによる抗GD2 CAR-T療法は他の疾患で先行研究があり、ヒトでの臨床研究が行われていますが、大きな神経毒性を持つことは報告されておらず、ヒトにおける安全性はある程度担保されていると考えられるようです。
CAR-T細胞は、DIPG同所性移植マウスで抗腫瘍効果を発揮する
次にin vivoでCAR-T細胞の効果を検証するため、ルシフェラーゼを発現させたDIPGを橋に移植したモデルマウスを作成しました(ちなみに、移植してもきちんとびまん性に浸潤する)。移植後7-8週間の後にCAR-T細胞を静注、静注後40日まで(40 DPT, 40 days post-treatment)追跡を行いました。すると、ルシフェラーゼの発光が時間とともに減少していき、実際にDPT50ではH3-K27M陽性細胞が組織はほとんど認められなくなりました。
実際の組織はこちら。(腫瘍はGFPを発現している)
CAR-T細胞は、DIPG同所性移植マウスの生存期間を延長する
最も重要な問いは、「実際にCAR-T細胞はマウスを救うのか?」というものですが、そのためにはモデルマウスが死んでくれないと困るわけです。ですが、先ほどのモデルマウスではDPT50を生きながらえてしまうようなので、著者たちは、腫瘍の増生が激しいタイプの患者由来細胞(SU-DIPG13P)を移植して、ひと月ほどで死ぬモデルマウスを作成し、これを用いて生存期間の解析を行いました。すると、CAR-T細胞によって腫瘍が除去され、生存率の改善を認めることができました。
しかし、CAR-T細胞によって治療したにも関わらず、死んでしまったマウスも存在します。このマウスの組織をみたところ、水頭症を示唆するような脳室の拡大が認められていました。移植部位以外での非特異的な脳実質の炎症などは認められませんでした。
CAR-T細胞は、他のタイプのDMG同所性移植マウスの生存期間も延長するか?
これまでは、DMGの中でも橋に発生するDIPGという腫瘍について検討してきましたが、同様の知見がDIPG以外でも成立するのかを最後に検討しました。
以下駆け足になりますが、視床と脊髄に発生したDMGでもGD2が発現していることを確認し、それぞれ同所性に移植しその後にCAR-Tを投与して反応をみました。脊髄の場合はDIPGのように上手く行ったのですが(略)、視床の場合、きちんと腫瘍が除去されているにも関わらず、マウスの生存率は下がってしまいました(図のk)。これは、炎症による腫瘍の腫脹によって*7、脳が押しやられ逸脱する現象(テント切痕ヘルニア)によって死亡するためのようです(図のj)。腫瘍の炎症による腫脹はDIPGの場合でも認められ、それは水頭症として現れましたが、視床は特にテント切痕ヘルニアが起こりやすい位置に存在するため、腫瘍よりも脳ヘルニアで死亡する割合が増えてしまうようです。モデルは免疫不全マウスに腫瘍を移植して作成したため、実際に移植した際の炎症反応よりも抑えられている可能性があると著者らは論文内で考察しています。
deepchemを動かしてみました?②
久しぶりの更新です。
deepchemのチュートリアルを動かそうとしたら、タンパク質の構造を表示してくれませんでした。ここで少しつまずいてしまったので、記事にまとめようと思います。deepchemと言うかdeepchem以前なのですが、、、。
はじめに
このチュートリアルをやろうと思いました。足りないパッケージは適宜追加しながらやれば問題ないだろうと思っていたのですが、途中、
ngltraj = visualize_complex(complex_mdtraj) ngltraj
と言う部分で、
Failed to display Jupyter Widget of type NGLWidget.
みたいな文字が表示され、複合体の構造が表示されませんでした。「あれ?NGLviewは途中で入れたはずなのにな?」と思い色々調べたことをまとめます。
結論
NGLViewを入れる際に、!pipで適当に付け足すだけではダメだったのです。
!pip install ipywidgets !jupyter nbextension enable --py widgetsnbextension --sys-prefix !pip install nglview !jupyter-nbextension enable nglview --py --sys-prefix
これをコピペしましょう。終わり。
試しに、
import nglview as nv view=nv.demo() view
としてタンパク質が出てくればオッケーです。クリクリして遊びましょう。