論文読んでたら、コンベンショナルなMMPの求め方が載ってたので、似たようなやり方を考えてみた。
要するに網羅的に結合を切断していって
{Scaffold; [fragment, fragment, fragment]}
っていうハッシュを作っていく。fragmentは共通のScaffoldを持つので配列中の任意のフラグメントのペアがMMPということになる。
僕の切断ルール
僕の切断ルール
- シングルボンド
- 環の一部ではない
- ボンドのどちらかの端がCであること
CHEMBL327743の結果
('CC[C@@H]([C@H](NC(=O)[C@@H]1CCCCN1CC(=O)*)/C=C/c1ccccc1)C', 'COc1cccc(c1*)OC', 'CHEMBL327743')
('CC[C@@H]([C@H](NC(=O)[C@@H]1CCCCN1C*)/C=C/c1ccccc1)C', 'COc1cccc(c1C(=O)*)OC', 'CHEMBL327743')
('CC[C@@H]([C@@H](/C=C/c1ccccc1)N*)C', 'COc1cccc(c1C(=O)CN1CCCC[C@H]1C(=O)*)OC', 'CHEMBL327743')
('CC[C@@H]([C@H](NC(=O)*)/C=C/c1ccccc1)C', 'COc1cccc(c1C(=O)CN1CCCCC1*)OC', 'CHEMBL327743')
('CC[C@@H](C(NC(=O)[C@@H]1CCCCN1CC(=O)c1c(OC)cccc1OC)*)C', '*/C=C/c1ccccc1', 'CHEMBL327743')
('COc1cccc(c1C(=O)CN1CCCC[C@H]1C(=O)N[C@H](/C=C/c1ccccc1)*)OC', 'CC[C@H](C)*', 'CHEMBL327743')
('*c1ccccc1', '*/C=C\\[C@H]([C@H](CC)C)NC(=O)[C@@H]1CCCCN1CC(=O)c1c(OC)cccc1OC', 'CHEMBL327743')
('*OC', 'CC[C@@H]([C@H](NC(=O)[C@@H]1CCCCN1CC(=O)c1c(*)cccc1OC)/C=C/c1ccccc1)C', 'CHEMBL327743')
('*OC', 'CC[C@@H]([C@H](NC(=O)[C@@H]1CCCCN1CC(=O)c1c(*)cccc1OC)/C=C/c1ccccc1)C', 'CHEMBL327743')
('*CC', 'COc1cccc(c1C(=O)CN1CCCC[C@H]1C(=O)N[C@@H](C(C)*)/C=C/c1ccccc1)OC', 'CHEMBL327743')
('CC[C@@H]([C@H](NC(=O)[C@@H]1CCCCN1CC(=O)c1c(OC)cccc1OC)/C=C/c1ccccc1)*', 'C*', 'CHEMBL327743')
('*C[C@@H]([C@H](NC(=O)[C@@H]1CCCCN1CC(=O)c1c(OC)cccc1OC)/C=C/c1ccccc1)C', 'C*', 'CHEMBL327743')
コード
どこで切断されたかどうかを明確にするために切断されたボンドにダミーアトム(atomicnumが0のアトム)をおいた。
canonicalSMILESで出力して.の文字でsplitした文字列をそのままキーに使っているが、openbabelのcanonicalizationアルゴリズムはMorgan Algorithmを改変したものを使っているのでおそらく大丈夫。
下のコードでは標準出力に書きだしたが、実際はデータベース化することになる。多分Mongo使う。
smi_string,id = clone.write(format="can")[:-1].split('\t')
frag1, frag2 = smi_string.split('.')
dbh.mmp.update({"core": frag1}, {"$push": {"frags": (frag2, id)}}, True)
dbh.mmp.update({"core": frag2}, {"$push": {"frags": (frag1, id)}}, True)
こんな感じだろうか。
以下MMPを求めるのに使ったコード
import pybel
import openbabel
import sys
input = sys.argv[1]
def mkclone(mol):
obc = openbabel.OBConversion()
obc.SetInAndOutFormats("mol", "mol")
molstring = obc.WriteString(mol)
new_mol = openbabel.OBMol()
obc.ReadString(new_mol, molstring)
return new_mol
def ring_member(rings, atom):
for i, ring in enumerate(rings):
if ring.IsMember(atom):
return i
return -1
def is_same_ring_member(rings, a1, a2):
a1_idx = ring_member(rings, a1)
a2_idx = ring_member(rings, a2)
if a1_idx == -1 or a2_idx == -1:
return False
return a1_idx == a2_idx
def detect_del_bonds(mol):
rings = mol.OBMol.GetSSSR()
del_bonds = []
for a in mol:
if a.atomicnum == 6:
for na in mol:
n_a = na.OBAtom
b = a.OBAtom.GetBond(n_a)
if b != None and b.GetBondOrder() == 1 and not is_same_ring_member(rings, a.OBAtom, n_a):
if a.OBAtom.GetIndex() < n_a.GetIndex():
del_bonds.append((a.OBAtom.GetIndex(), n_a.GetIndex()))
return del_bonds
for mol in pybel.readfile("sdf", input):
chembl_id = mol.data["CHEMBL ID"]
print chembl_id
del_bonds = detect_del_bonds(mol)
for a1, a2 in del_bonds:
clone = pybel.Molecule(mkclone(mol.OBMol))
# !!! index starts 0 but GetAtom starts 1 !!!
ob_a1 = clone.OBMol.GetAtom(a1 + 1)
ob_a2 = clone.OBMol.GetAtom(a2 + 1)
bond = ob_a1.GetBond(ob_a2)
clone.OBMol.DeleteBond(bond)
new_atom1 = clone.OBMol.NewAtom()
new_atom1.SetAtomicNum(0)
new_bond1 = clone.OBMol.NewBond()
new_bond1.SetBegin(ob_a1)
new_bond1.SetEnd(new_atom1)
new_bond1.SetBondOrder(1)
new_atom2 = clone.OBMol.NewAtom()
new_atom2.SetAtomicNum(0)
new_bond2 = clone.OBMol.NewBond()
new_bond2.SetBegin(ob_a2)
new_bond2.SetEnd(new_atom2)
new_bond2.SetBondOrder(1)
smi_string, id = clone.write(format="can")[:-1].split('\t')
frag1, frag2 = smi_string.split('.')
#print (frag1, chembl_id)
#print (frag2, chembl_id)
print (frag1, frag2, chembl_id)
Pythonプロフェッショナルプログラミング
Gitポケットリファレンス
エキスパートPythonプログラミング
Essential Sqlalchemy
Python Testing Cookbook
Python Testing: Beginner's Guide
みんなのPython 第3版
NHK 趣味の園芸 やさいの時間 2012年 09月号 [雑誌]