18 08 2010 chemoinformatics Python network igraph Tweet
MCSを求めるには比較したい二つの化合物のModular Productを求めて、それの最大クリーク探索をすればいい。
Modular Productとは
- 二つのグラフの属性の対応が同じノードが両方共隣接していてエッジの属性が同じ場合
- 二つのグラフの属性の対応が同じノードが両方共隣接していない場合
にエッジを張ったグラフ。
これを作ってあとはigraphの最大クリーク探索メソッドで
import openbabel as ob
from igraph import *
sm1 = 'OCC(N)=O'
sm2 = 'O=CC(N)=O'
obc = ob.OBConversion()
obc.SetInAndOutFormats('smi','smi')
source = ob.OBMol()
obc.ReadString(source,sm1)
source_atoms = [a for a in ob.OBMolAtomIter(source)]
target = ob.OBMol()
obc.ReadString(target,sm2)
target_atoms = [a for a in ob.OBMolAtomIter(target)]
pairs = []
for i in range(len(source_atoms)):
for j in range(len(target_atoms)):
if source_atoms[i].GetAtomicNum() == target_atoms[j].GetAtomicNum():
pairs.append((source_atoms[i].GetIdx(),target_atoms[j].GetIdx()))
p = []
for x in range(len(pairs)-1):
for y in range(x+1,len(pairs)):
u = pairs[x]
v = pairs[y]
if u[0] != v[0] and u[1] != v[1]:
sb = source.GetBond(u[0],v[0])
tb = target.GetBond(u[1],v[1])
if sb != None: sbo = sb.GetBondOrder()
if tb != None: tbo = tb.GetBondOrder()
if (sb == None and tb == None) or (sb != None and tb != None and sbo == tbo):
p.append((pairs.index(u),pairs.index(v)))
g = Graph(p)
mc = g.largest_cliques()
for c in mc:
print [pairs[i] for i in c]
実際にはline graph(エッジとノードを入れ替えたグラフ)のModular Productをもとめてクリーク探索したほうが計算量がかなり減るらしい。
が、今回はやっていない。
Handbook of Chemoinformatics Algorithms (Chapman & Hall/Crc Mathematical and Computational Biology)
Chapman & Hall / ¥ 8,895 ()
通常1~3週間以内に発送
Chapman & Hall / ¥ 8,895 ()
通常1~3週間以内に発送
Chemoinformatics and Computational Chemical Biology (Methods in Molecular Biology)
Humana Press / ¥ 14,386 ()
近日発売 予約可
Humana Press / ¥ 14,386 ()
近日発売 予約可