Pythonで動いて多彩な機能を持っておりながらフリーである。
ということで自分はRDKitが大好きです。
次のバージョンのβ版がリリースされているのですがどうもMMPAアルゴリズムが実装されるようです。
たくさんの情報に埋もれた何かを見るにはいいかもしれないというのと、コードに興味があるのとで早速使ってみます。
SVMはここで公開されているので
待てない人は使ってくれよ！という話です。
ちなみにリファレンスはこちら。
アルゴリズム自体は前から報告されていたのですがソースコードは今回が初です。
まだソースをちゃんと読んでないですが使ってみます。
取りあえずChemblにいってABL/Humanで1504化合物のSDFをとってきました。
入力は
smiles\tidのテキストなので
そのようにかえます。
from rdkit import Chem
f = open("chembl.smi","w")
mols =[m for m in Chem.SDMolSupplier("cmpd_download_27530.sdf")]
for mol in mols:
sarregno = mol.GetProp("sarregno")
smi = Chem.MolToSmiles(mol)
f.write("%s\t%s"%(smi,sarregno))
f.write("\n")
...:
f.close()
でファイルができたんでフラグメント化します
lion:mmpa user$ wc chembl.smi 1504 3008 86733 chembl.smi lion:mmpa user$ time python rfrag.py <chembl.smi >chemblfrag.txt real 5m13.693s user 5m12.636s sys 0m0.489s lion:mmpa user$ wc chemblfrag.txt 37177 37177 9807140 chemblfrag.txt
frag.pyでsmilesの分子をフラグメントに切断します。
1504分子から37177のフラグメントができました。
元分子 ID フラグメントみたいなテキストがchemblfrag.txtにできます。
この辺をpythonで考えると辞書型で作る感じですよね。
で最後に
ペアを作ります。
lion:mmpa user$ time python indexing.py <chemblfrag.txt >chemblmmpa.txt real 0m39.278s user 0m39.138s sys 0m0.103s lion:mmpa user$ wc chemblmmpa.txt 886 886 158799 chemblmmpa.txt lion:mmpa user$ head chemblmmpa.txt Cc1ccc(C(=O)Nc2ccc(CN3CCN(C)CC3)c(C(F)(F)F)c2)cc1C#Cc1cnc[nH]1,Cc1ccc(C(=O)Nc2ccc(CN3CCN(C)CC3)c(C(F)(F)F)c2)cc1C#Cc1cncnc1,1784631,1784626,[*:1]C#Cc1cnc[nH]1>>[*:1]C#Cc1cncnc1,[*:1]c1cc(C(=O)Nc2ccc(CN3CCN(C)CC3)c(C(F)(F)F)c2)ccc1C Cc1ccc(C(=O)Nc2ccc(CN3CCN(C)CC3)c(C(F)(F)F)c2)cc1C#Cc1cnc[nH]1,Cc1ccc(C(=O)Nc2ccc(CN3CCN(C)CC3)c(C(F)(F)F)c2)cc1C#Cc1cnc(N)s1,1784631,1784640,[*:1]C#Cc1cnc[nH]1>>[*:1]C#Cc1cnc(N)s1,[*:1]c1cc(C(=O)Nc2ccc(CN3CCN(C)CC3)c(C(F)(F)F)c2)ccc1C
結局886個のペアができました。
構造1,構造2,ID1,ID2,変換のsmirkes
みたいになっています。
Spotfireに突っ込めば全部構造になるので良さそうですね。
切断部分で1504化合物で5分かかっているので大規模にやるのは時間がかかりそうだ。
合成した化合物を逐次フラグメント化してDBに入れておけばマッチングは早そうだからきびきび行くかな。
web.pyでがんばってwebインタフェースができればみんながアップしたら的なものできそうだけど。
POSTで受けたファイルのハンドリングがよくわからんからweb.pyも理解せねば
こういったこと、社内で言ってもきっと理解されないなぁ
何か考えるかな〜
Is life worth living? 