タグ: ruby

Call rdkit from Ruby

Recently I use web app made with ruby.
The app is very useful for me, but there is not chemoinformatics library in ruby.
So, I want to find way to use rdkit from ruby.
I found some web page that described how to run shell script from ruby. The simplest way is using back quotation.
It is means call “python xxx.py” from ruby.
I tested it.
At first, I wrote sample script named test.py. Just load SDF and print molecules as smiles.

 from rdkit import Chem
 sdf = Chem.SDMolSupplier( 'cdk2.sdf' )
  
 for mol in sdf:
     smi = Chem.MolToSmiles( mol )
     print( smi )
                               

Then, I wrote ruby script.
It just call test.py. That’s all. 😉

callpy.rb

         
res = `python test.py`
  print res

Run callpy.rb.

         
iwatobipen$ ruby callpy.rb 
...........
Cn1cnc2c(NCc3ccccc3)nc(NCCO)nc21
CCC(CO)Nc1nc(NCc2ccccc2)c2ncn(C(C)C)c2n1
COc1ccc(CNc2nc(N(CCO)CCO)nc3c2ncn3C(C)C)cc1
Nc1nc(N)c(N=O)c(OCC2CCCCC2)n1
COc1ccc2c(c1)C(=Cc1cnc[nH]1)C(=O)N2
COc1cc[nH]c1C=C1C(=O)Nc2ccc([N+](=O)[O-])cc21
CCc1cnc(CSc2cnc(NC(=O)C(C)C)s2)o1
CCCCOc1c(C(=O)c2nccs2)cnc2[nH]ncc12
COc1cc(-c2ccc[nH]2)c2c3c(ccc(F)c13)NC2=O
[NH3+]CCSc1cc(-c2ccc[nH]2)c2c3c(ccc(F)c13)NC2=O
NC(=O)Nc1cccc2c1C(=O)c1c-2n[nH]c1-c1cccs1
COc1ccc2c(c1)C(=Cc1[nH]cc3c1CCOC3=O)C(=O)N2
CNS(=O)(=O)c1ccc2c(c1)C(=Cc1[nH]cc3c1CCNC3=O)C(=O)N2
CC(=O)Nc1cccc2c1C(=O)c1c-2n[nH]c1-c1ccncc1
....

It worked. 😉

There are some ways to call shell script from ruby. I will try another method soon.

Lillyのフィルタ

J.Med.Chem. 2012, 55, 9763-9772
Lillyのこれまで積み重ねてきたフィルタに関しての報告です。
ソースコードがGitHubにあがっているので、ちょっと使ってみました。

pen$ git clone https://github.com/IanAWatson/Lilly-Medchem-Rules.git
pen$ cd Lilly-Medchem-Rules/
pen$ make

でインストールはよろしくやってくれます。簡単です。
続いてデータセットの準備をします。CHEMBLからDPP4阻害剤のデータを取ってきました。

pen$ head -2 ligand-13\ 12-44-32.txt
MOLREGNO	CHEMBLID	SYNONYMS	MOLWEIGHT	ALOGP	PSA	HBA	HBD	RO5	RTB	RO3	MED_CHEM_FRIENDLY	ACD_MOST_APKA	ACD_MOST_BPKA	ACD_LOGP	ACD_LOGD	AROMATIC_RINGS	HEAVY_ATOMS	NUM_ALERTS	QED_WEIGHTED	CANONICAL_SMILES
655929	CHEMBL1181979		339.43	.97	66.64	3	1	0	4	N	N		7.96	.04	-.29	1	25	0	.89	C[C@H]([C@H](N)C(=O)N1CCCC1)c2ccc(cc2)C3=CN(C)C(=O)C=C3
399135	CHEMBL235199		452.54	1.06	97.35	5	1	0	5	N	Y		6.8	.72	.67	1	32	0	.72	CN(C)C(=O)[C@H]([C@H](N)C(=O)N1CCC(F)(F)C1)[C@@H]2CC[C@H](CC2)C3CCc4ncnn4C3
162221	CHEMBL98869		309.83	1.58	58.36	3	2	0	6	N	Y		9.01	1.9	.04	1	21	1	.79	N[C@@H](CCNCc1ccc(Cl)cc1)C(=O)N2CCCCC2
521016	CHEMBL485123		392.42	2.61	72.11	4	1	0	4	N	Y		8.78	2.3	.91	2	28	1	.81	Cc1ncc(C)c(n1)C(=O)N2CCC(CC2)[C@H](N)Cc3cc(F)c(F)cc3F

入力に使うフォーマットはsmiles tab idにしないとだめなので
今度はPANDASで加工します。

import pandas as pd
table = pd.read_table( "ligand-13 12-44-32.txt",header=0 )
table_df = pd.DataFrame( table )
select_tab = pd.DataFrame( table_df, columns=["CANONICAL_SMILES","MOLREGNO"] )
select_tab.to_csv( "dpp4.smi", index=False, sep = "\t" )

で確認します。

pen$ head -4 dpp4.smi 
CANONICAL_SMILES	MOLREGNO
C[C@H]([C@H](N)C(=O)N1CCCC1)c2ccc(cc2)C3=CN(C)C(=O)C=C3	655929
CN(C)C(=O)[C@H]([C@H](N)C(=O)N1CCC(F)(F)C1)[C@@H]2CC[C@H](CC2)C3CCc4ncnn4C3	399135
N[C@@H](CCNCc1ccc(Cl)cc1)C(=O)N2CCCCC2	162221

ということでファイルができました。
実際に計算してみましょう

pen$ cd ~/Lilly-Medchem-Rules
pen$ ./Lilly_Medchem_Rules.rb dpp4.smi
~~~~~~~~~
CCN(CC)C1=NC2=C(C(=O)N(C)C(=O)N2C)C(=C1CN)C1=CC(F)=CC=C1Br 1425736 : D(60) too_many_atoms:bromine
N[C@H]1CN(CC2=CC=C(F)C(F)=C2)C[C@H]1C(=O)N1CCCC1 415800
CCN1C(=NC2=C1C=CC=C2)C1=CC=C(C=C1)C#CC1=CC=C2SCCC(C)(C)C2=C1 451778 : D(90) too_many_atoms:acetylene
NC(=O)[C@@H]1CCCN1C(=O)CCNC(=O)C1=CC=C(O)C=C1 1364243
CC(C)[C@H](N)C(=O)N1CCC[C@H]1C(=O)NC1=CC=C2N=CC=CC2=C1 577553
CN([C@@H]1CC[C@H](CC1)[C@H](N)CC1=CC(F)=CC=C1F)S(=O)(=O)C1=CC=CN=C1 1073174 : D(20) too_many_atoms
CCCCN1C(=NC2=C1C=CC=C2)C1=CC=C(C=C1)C#CC1=CC=CC=C1 451639 : D(73) too_many_atoms:C4:acetylene
FC(F)(F)C1=NSC(=N1)N1CCN(CC1)C(=O)[C@@H]1CN[C@@H](C1)C(=O)N1CCCC1 438446 : D(26) too_many_atoms
CCCC1=NSC(=N1)N1CCN(CC1)C(=O)[C@@H]1CN[C@@H](C1)C(=O)N1CCCC1 438495 : D(20) too_many_atoms
NCC(N1CCOCC1)C1=CC=CN=C1 341339
COC1=CC=C(C=C1)C1C=C(NC2=C1C(=O)CC(C)(C)C2)C1=CC=CC=C1 645647 : D(63) too_many_atoms:enamine_2
CN1C(=O)C2=C(N=C(N3CCC[C@@H](N)C3)N2CC=C(C)C)C2=CC(=CC=C12)C(=O)O 1449930 : D(33) too_many_atoms
N[C@@H]([C@@H]1CC[C@H](CC1)NS(=O)(=O)C1=CC=C(OC(F)(F)F)C=C1)C(=O)N1CCSC1 29718 : D(33) too_many_atoms
N[C@@H](CCNCC1=CC=C(OCC2=CC=CC=C2)C=C1)C(=O)N1CCCCC1 162236 : D(70) too_many_atoms:positive
CS(=O)(=O)NC1=CC2=C(C=C1)N1C[C@H](N)[C@H](CC1=N2)C1=CC(F)=C(F)C=C1F 671132 : D(20) too_many_atoms
N[C@@H]1CCCN(C1)C1=NC2=CC=C(Br)C=C2C(=O)N1CC1=CC=CC=C1C#N 377828 : D(54) too_many_atoms:bromine
C[C@H]([C@H](N)C(=O)N1CC[C@H](F)C1)C1=CC=C(C=C1)C1=CN(CC2CC2)C(=O)C=C1 658954 : D(26) too_many_atoms
COC1=CC2=C(C(=O)N(C(=N2)N2CCC[C@@H](N)C2)CC2=CC=CC=C2C#N)C=C1OC 377841 : D(40) too_many_atoms
O=C(C[C@@H]1C[C@H](NC1=O)C(=O)N1CCSC1)N1CC2=C(C=CC=C2)C1 702390
~~~~~~~~~~~~~

結果はsmiles / id / demerit / detailsのような出力です。
17クラスからなる275のルールから構成されるフィルターで
バイオレーションによってDemeritのスレッシュフォルドが決まります。
デフォルトは100でリラックス(ゆるめ)で160をカットオフにしています。
テスト用のデータが4万件くらいでしたが、手元のMBAで数十秒でした。
プロジェクトが進んでからのアプライは???ですが、
市販のセットやHTS用のセットに一度当てるのは面白そうですね。