转载自:Computing Molecular Descriptors – Part 1 链接
一个由五部分组成的系列教程,内容包括实现python代码来计算在药物发现中高度使用的不同类型的二维分子描述符和指纹集。
什么是分子描述符或指纹
定义:“逻辑和数学过程的最终结果,该过程将分子符号表示中编码的化学信息转换为有用的数字或一些标准化实验的结果。”
它们可能呈现对数P、分子量、氢键供体、受体、可旋转键等分子的性质,这些性质可能与分子的实验证据相关联。或者它们可能是基于二维片段的指纹,由0s和1s的位数组显示,其中每个位位置指示是否存在结构片段。
对于每个级别的分子表示(见图2),如1D、2D、3D,甚至4D/MD(基于分子动力学时间序列的维度),我们可以计算数百或数千个结构特征。
二维(2D)描述符
二维表示解释了结构特征和原子的连接方式(例如邻接、连接)。它们通常是导出节点和边缘分子图表示的描述符,如大小、分支程度、原子在电子和立体效应方面的邻域、灵活性、整体形状等。
本项目python教程概述
我在二维描述符集上组合了四个python类,可以轻松导入和应用。这些python类如下,可以搜索以下:
MACCS描述符(第2部分)
ECFP6描述符(第3部分)
RDKit_2D(第4部分)
Mordred_MRC_描述符(第5部分)
从主python文件(名为main.py)调用这些python类。通过运行main.py脚本,您将计算上述所有描述符,并将其生成为单独的CSV文件。所以,最后,你会有如下图3所示的输出文件。
具体使用请看以下介绍:
MACCS Fingerprints in Python (第二部分)
MACCS介绍部分
MACCS键是关于化学结构的一组问题,是表征分子的二元(0或1)列表。以下是其中的一些key:
例:氧的种类是否少于3种?有S-S键吗?有四元环吗?是否至少存在一个F、Cl、Br或I?
最终得到true或false,也就是一行1/0的字节串。(共166个key的位置是不能颠倒的)
RDKit实现的166个MACCS公开的key(片段定义)可以在这里找到:点此链接
代码实现部分:
首先在环境中使用conda安装这两个必要的库。
conda install -c rdkit rdkitconda install pandas
MACCS的类 如下: 保存为MACCS.py
import pandas as pd
from rdkit import Chem
from rdkit.Chem import MACCSkeys
class MACCS:
def __init__(self, smiles):
self.smiles = smiles
self.mols = [Chem.MolFromSmiles(i) for i in smiles]
def compute_MACCS(self, name):
MACCS_list = []
header = ['bit' + str(i) for i in range(167)]
for i in range(len(self.mols)):
ds = list(MACCSkeys.GenMACCSKeys(self.mols[i]).ToBitString())
MACCS_list.append(ds)
df = pd.DataFrame(MACCS_list,columns=header)
df.to_csv(name[:-4]+'_MACCS.csv', index=False)
以上这个MACCS指纹的类可以被保存为单独的python文件:MACCS.py,以便我们轻松的加载它、import它。
例如,假设您要计算包含SMILES的CSV文件的MACC描述符,以下是您可以使用的代码:
(保证MACCS.py在相同工作路径下。)
import pandas as pd
from molvs import standardize_smiles
from MACCS import *
def main():
filename = 'data/macrolides_smiles.csv' # path to your csv file
df = pd.read_csv(filename) # read the csv file as pandas data frame
smiles = [standardize_smiles(i) for i in df['smiles'].values]
## Compute MACCS Fingerprints and export file.
maccs_descriptor = MACCS(smiles) # create your MACCS object and provide smiles
maccs_descriptor.compute_MACCS(filename) # 计算 MACCS 并且给文件命名。 命名可以直接使用初始名称,因为MACCS的类会在其后添加 "_MACCS.csv" 作为输出文件。
if __name__ == '__main__':
main()
MACCS代码逐行解释
class MACCS:
def __init__(self, smiles):
self.smiles = smiles
创建了MACCS这个类,在__init __()方法中启动类的属性。创建了一个所需的参数“smiles”,当我们调用这个类的时候,我们需要给出一串smiles。
将smiles变量保存为这个类的属性: self.smiles = smiles
使用‘Chem.MolFromSmiles()‘将smiles转换为mols文件。
self.mols = [Chem.MolFromSmiles(i) for i in smiles]
计算指纹并保存为CSV格式
def compute_MACCS(self, name):
MACCS_list = []
header = ['bit' + str(i) for i in range(167)]
for i in range(len(self.mols)):
ds = list(MACCSkeys.GenMACCSKeys(self.mols[i]).ToBitString())
MACCS_list.append(ds)
df = pd.DataFrame(MACCS_list,columns=header)
df.insert(loc=0, column='smiles', value=self.smiles)
df.to_csv(name[:-4]+'_MACCS.csv', index=False)
运行得到表格文件。
其中“header”如下:
header = ['bit0', 'bit1', 'bit2', 'bit3', 'bit4', 'bit5', ... , 'bit164', 'bit165', 'bit166']
每次ds都是一个列表,里面是166个“0、1”的指纹。
这样经过for in 循环,MACCS_list就有了n个列表。如下面是两个分子的MACCS指纹:
MACCS_list = [
['0', '0', '0', '0', '0', '0', '0', '0', '0', '0', '0', '0', '0',………… '1', '0'],
['0', '0', '0', '0', '0', '0', '0', '0', '0', '0', '0', '0', '0', ………… '0', '0']
]
代码使用(给出一个smiles列表,得出MACCS指纹)
from MACCS import *
smiles = ['Nc1nc(NC2CC2)c2ncn(C3C=CC(CO)C3)c2n1',
'CC(=O)NCCCS(=O)(=O)O',
'CCCC(=O)Nc1ccc(OCC(O)CNC(C)C)c(C(C)=O)c1',
'CC(=O)Nc1ccc(O)cc1',
'CC(=O)Nc1nnc(S(N)(=O)=O)s1',
'CC(=O)NO'
]
maccs_descriptor = MACCS(smiles)
ECFP6 Fingerprints in Python (第三部分)
Extended Connectivity FingerPrinting 对分子的每个原子进行迭代,并根据指定半径从该原子检索所有可能的分子路径。
1.Iteration 0
给每个原子添加一个标识符(相邻的非氢原子数、连接的键数、原子质量、相邻氢的数量、是否在环内,最终计算出一个整数哈希值)
2.Iteraltion 1、Iteration 2、Iteration……
Iteration1 根据每个原子最近相邻的原子更新标识符(为每个原子创建一个n行2列的数组,其中n为“相邻原子数+1”,第二列为各个原子的标识符。
例如上述的4号原子,它自身为(1,-1100000244),其他三个相邻原子的数列分别为(与4号相连键数,标识符):(1,1559650422),(1,1572579716),(2,-1074141656)),得到4号原子的最终数组:((1, -1100000244),(1,1559650422),(1,1572579716),(2,-1074141656)),再将其转换为哈希值:-1708545601。
其他5个原子重复上述操作,得到各自的“最近相邻原子”哈希值。
Iteration2 根根据两个原子半径内的原子 更新标识符
方法类似上述,最后得到哈希值列表。此时3号原子的Iteration2已经探索完整个分子,就不会再进行Iteration3了。
3.移除重复
在6个Iteration2中,由于1号、5号、6号位于边缘,探索的原子有限,他们完全可以被2号、3号、4号的Iteration2取代(2号、3号、4号的Iteration2 探索的原子包含了所有1、5、6能探索到的原子)所以只保留2、3、4号的Iteration2。
4.使用哈希算法将各个Iteration的标识符转为长度为2048位向量
实际应用中,根据实际需要,选择Iteration的级数得到对应的哈希值。如想探索毒性苯胺结构,就至少需要 Iteration 4。如果只是想预测疏水性,用Iteration1/2即可。
该指纹能提取分子各个子结构的信息,例如原子及其与相邻原子的连通性,然后是这些邻域的邻域。
代码使用 ECFP6 Fingerprints
配置环境
conda install -c rdkit rdkit
conda install pandas
conda install numpy
整个代码项目可以在这里找到:https://github.com/zinph/Cheminformatics/blob/master/compute_descriptors/ECFP6.py
以下是包含ECFP6的类的python文件:ECFP6.py
import numpy as np
import pandas as pd
from rdkit.Chem import AllChem
from rdkit import Chem, DataStructs
class ECFP6:
def __init__(self, smiles):
self.mols = [Chem.MolFromSmiles(i) for i in smiles]
self.smiles = smiles
def mol2fp(self, mol, radius = 3):
fp = AllChem.GetMorganFingerprintAsBitVect(mol, radius = radius)
array = np.zeros((1,))
DataStructs.ConvertToNumpyArray(fp, array)
return array
def compute_ECFP6(self, name):
bit_headers = ['bit' + str(i) for i in range(2048)]
arr = np.empty((0,2048), int).astype(int)
for i in self.mols:
fp = self.mol2fp(i)
arr = np.vstack((arr, fp))
df_ecfp6 = pd.DataFrame(np.asarray(arr).astype(int),columns=bit_headers)
df_ecfp6.insert(loc=0, column='smiles', value=self.smiles)
df_ecfp6.to_csv(name[:-4]+'_ECFP6.csv', index=False)
将上述文件存在工作目录文件夹中,再运行一段代码即可根据SMILES获得分子的指纹。以下是代码:
import pandas as pd
from molvs import standardize_smiles
from ECFP6 import *
def main():
filename = 'data/macrolides_smiles.csv' # path to your csv file
df = pd.read_csv(filename) # read the csv file as pandas data frame
smiles = [standardize_smiles(i) for i in df['smiles'].values]
## Compute ECFP6 Fingerprints and export a csv file.
ecfp6_descriptor = ECFP6(smiles) # create your ECFP6 object and provide smiles
ecfp6_descriptor.compute_ECFP6(filename) # compute ECFP6 and provide the name of your desired output file. you can use the same name as the input file because the ECFP6 class will ensure to add "_ECFP6.csv" as part of the output file.
if __name__ == '__main__':
main()