需要对多个化合物批量进行特定片段或基团的拼接或替换时,可以参考以下代码。
一、片段的拼接
虽然rdkit中有一个Chem.CombineMols(mol_1,mol_2)的函数,仅仅是将两个分子放入同一个对象中,依然是未连接在一起的形式,不符合我们在指定位点拼接的需求。为此,我们需要自定义一个函数,并且通过marker原子指明拼接的位点来实现。
首先我们定义两个辅助函数,确认拼接位点后,连接两个片段
#函数一,获取marker邻居原子的index, 注意marker只能是一个单键连接核上的原子,否则邻居会多于一个
def get_neiid_bysymbol(mol,marker):
try:
for atom in mol.GetAtoms():
if atom.GetSymbol()==marker:
neighbors=atom.GetNeighbors()
if len(neighbors)>1:
print ('Cannot process more than one neighbor, will only return one of them')
atom_nb=neighbors[0]
return atom_nb.GetIdx()
except Exception as e:
print (e)
return None
#函数二,获取marker原子的index
def get_id_bysymbol(mol,marker):
for atom in mol.GetAtoms():
if atom.GetSymbol()==marker:
return atom.GetIdx()
然后是连接两个片段的主函数
def combine2frags(mol_a,mol_b,maker_b='Cs',maker_a='Fr'):
#将两个待连接分子置于同一个对象中
merged_mol = Chem.CombineMols(mol_a,mol_b)
bind_pos_a=get_neiid_bysymbol(merged_mol,maker_a)
bind_pos_b=get_neiid_bysymbol(merged_mol,maker_b)
#转换成可编辑分子,在两个待连接位点之间加入单键连接,特殊情形需要其他键类型的情况较少,需要时再修改
ed_merged_mol= Chem.EditableMol(merged_mol)
ed_merged_mol.AddBond(bind_pos_a,bind_pos_b,order=Chem.rdchem.BondType.SINGLE)
#将图中多余的marker原子逐个移除,先移除marker a
marker_a_idx=get_id_bysymbol(merged_mol,maker_a)
ed_merged_mol.RemoveAtom(marker_a_idx)
#marker a移除后原子序号变化了,所以又转换为普通分子后再次编辑,移除marker b
temp_mol = ed_merged_mol.GetMol()
marker_b_idx=get_id_bysymbol(temp_mol,maker_b)
ed_merged_mol=Chem.EditableMol(temp_mol)
ed_merged_mol.RemoveAtom(marker_b_idx)
final_mol = ed_merged_mol.GetMol()
return final_mol
示例,用Fr与Cs原子分别标注两个片段的连接位点后的连接过程
mol1=Chem.MolFromSmiles('c1ccccc1[Fr]')
mol2=Chem.MolFromSmiles('C1CC1C[Cs]')
combine2frags(mol1,mol2)
二、片段的替换
需要替换分子中特定基团为需求片段时候,可以简单的采用rdkit的ReplaceSubstructs函数, 示例如下,但该函数是封装得较好的高级函数,对于复杂的情况往往会出现难以debug的错误结果。
简单的片段替换可直接参考rdkit的示例,一个将乙酰胺的氨基替换为甲醚的案例。
from rdkit.Chem import AllChem
repl = Chem.MolFromSmiles('OC') # 期望替换为一个甲醚的片段
patt = Chem.MolFromSmarts('[$(NC(=O))]') #被替换的氨基部分,不熟悉SMARTS的朋友,也可以用Chem.MolFromSmiles('N')
m = Chem.MolFromSmiles('CC(=O)N') # 乙酰胺
rms = AllChem.ReplaceSubstructs(m,patt,repl)
rms
(<rdkit.Chem.rdchem.Mol object at 0x...>,)
Chem.MolToSmiles(rms[0])
'COC(C)=O'
但这样做有个很大的隐患,即我们新加入的甲醇片段,并没有指明应该从C原子还是O原子的部位连接到乙酰胺的核上。实际上ReplaceSubstructs这个函数有个可选argument为replacementConnectionPoint (default=0),用于指定连接的位点,当我们上面的代码略作修改,会得到另一个结果
rms = AllChem.ReplaceSubstructs(m,patt,repl,replacementConnectionPoint=1)
rms
Chem.MolToSmiles(rms[0])
因此在具体应用中,还是建议引入marker的方式,准确寻找连接位点的index,替换好了后再移除marker。
同样利用第一节中定义的辅助函数
def replace_sub(original_mol, sub_to_keep, new_sub, marker_a='Fr', marker_b='Cs'):
marked_sub=Chem.MolFromSmiles(Chem.MolToSmiles(tmp).replace('1*',marker_a))
new_mol=combine2frags(frag1,mol2,marker_a=marker_a, marker_b=marker_b)
return Chem.MolToSmiles(new_mol)
test_mol = Chem.MolFromSmiles('CNC(=O)c1cncnc1')
core = Chem.MolFromSmiles('c1cncnc1')
new_frag=Chem.MolFromSmiles('C1CC1C[Cs]')
replace_sub(test,core,new_frag)
'c1ncc(CC2CC2)cn1'
我们初始分子是一个带有甲基酰胺的嘧啶分子,希望保留嘧啶环,将侧链替换为一个甲基环丙烷,即可利用该函数实现目的。