RDKit|分子指纹提取、相似性比较及应用

  • 一、分子指纹提取
    • 1.Topological Fingerprints
    • 2.MACCS
    • 3.Atom Pairs and Topological Torsions
    • 4.Morgan Fingerprints (Circular Fingerprints)
  • 一、相似性
    • 1.相似性比较
    • 2.相似性地图
    • 3.相似性应用——构建多样化分子库

一、分子指纹提取

1.Topological Fingerprints

也叫RDKit topological fingerprint。该算法受Daylight fingerprint启发而产生,它会计算所有介于minPath和maxPath之间的分子路径(子图,subgraphs),对每个子图进行哈希运算,产生一个原始的bit ID。对bit ID取模(除数是fpSize),并在相应的位上进行设置。
对子图做哈希运算,其实就是对原子和键做哈希,它考虑了原子类型、芳香性、键的类型三种特征

对每个子图生成指纹的大致流程:

  • 1.对子图计算一个哈希值
  • 2.使用哈希值生成一个raw ID
  • 3.对raw ID取模
  • 4.如果nBitsPerHash大于1,根据raw ID再随机生成一个数,直到达到指定值
  • 5.组成分子指纹相应的位。

使用rdkit来生成拓扑分子指纹

  • 计算拓扑分子指纹:RDKFingerprint(mol, minPath,maxPath, fpSize, nBitsPerHash, useHs, tgtDensity, bitInfo, ...)
    mol:mol对象
    minPath:默认1。生成子图所需的最少键的数量
    maxPath:默认7。生成子图时最多键的数量
    fpSize:生成指纹的长度,默认2048
    nBitsPerHash:默认2。每个子图对应的位数。如果大于1,就将raw ID作为种子随机生成若干个数。
    useHs:默认True,将显式氢考虑进指纹中。
    tgtDensity:默认0。如果大于0,指纹长度会不断减半,直到大于等于设置的密度,或到达minSize的设置值。
    bitInfo:用于接收指纹结果的字典。键为位数,值为相同指纹的列表(列表的每个元素为键的序号)。
>>> ms = [Chem.MolFromSmiles('CCOC'), Chem.MolFromSmiles('CCO'), Chem.MolFromSmiles('COC')]
>>> fps = [Chem.RDKFingerprint(x) for x in ms]
>>> fps
[<rdkit.DataStructs.cDataStructs.ExplicitBitVect at 0x7538300>,
 <rdkit.DataStructs.cDataStructs.ExplicitBitVect at 0x75383f0>,
 <rdkit.DataStructs.cDataStructs.ExplicitBitVect at 0x7538670>]

先简单地进行一下相似性比较

  • 相似度比较:DataStructs.FingerprintSimilarity(fp1, fp2)
>>> DataStructs.FingerprintSimilarity(fps[0], fps[1])
0.6
  • 通过传入空字典rdkinfo获取非空分子指纹信息
    黄色表示带有芳香性质的原子
>>> mol = Chem.MolFromSmiles('c1ccccc1CC1CC1')
>>> rdkinfo = {}
>>> rdkfp = Chem.RDKFingerprint(mol, bitInfo=rdkinfo)
>>> Draw.DrawRDKitBit(mol, list(rdkinfo.keys())[0], rdkinfo)
1

2.MACCS

是一种基于SMARTS的,长度为167的分子指纹,每一位所代表的含义可以参考openbabel的介绍

>>> from rdkit.Chem import MACCSkeys
>>> fps = [MACCSkeys.GenMACCSKeys(x) for x in ms]
>>> DataStructs.FingerprintSimilarity(fps[0],fps[1])
0.5
>>> fps[0].GetNumBits()
167

3.Atom Pairs and Topological Torsions

这两种指纹算法的概念在30多年前提出,比较相似,指纹会包括原子序号、π电子数和相邻原子数三个维度的信息。此外还可以通过参数,加入手性的信息。两种指纹都可以用sparse form和explicit form形式来表示,并且可以查看指纹所代表的化学信息。

  • 计算atom pair指纹:GetAtomPairFingerprint(mol, minLength, maxLength, ...)
  • 同GetAtomPairFingerprintAsIntVect()
    这个方法生成的指纹非常庞大稀疏。
>>> from rdkit.Chem.AtomPairs import Pairs
>>> fps = [Pairs.GetAtomPairFingerprint(x) for x in ms]
  • 可以通过该方法获取atom pairs算法产生的非空元素:GetNonzeroElements()
    返回值是一个字典,键为指纹所在的位,值为出现的频数。
>>> fps[2].GetNonzeroElements()
{541730: 1, 1606689: 2}
  • 也可以查看位所代表的指纹:Pairs.ExplainPairScore()
>>> Pairs.ExplainPairScore(1606689)
(('C', 1, 0), 1, ('O', 2, 0))

该结果表示这样一种子结构:有1个相邻原子和0个π电子的碳原子,通过一个键,与有2个相邻原子和0个π电子的氧原子相连。

  • 也可以用一个定长向量来表示atom pair:GetAtomPairFingerprintAsBitVect()
    注意:该方法不会记录子结构的频数,只用0和1表示是否存在。
>>> pairFps = [Pairs.GetAtomPairFingerprintAsBitVect(x) for x in ms]
>>> pairFps[0].GetNumBits()
8388608
>>> DataStructs.DiceSimilarity(pairFps[0],pairFps[1])
0.2222222222222222
  • 生成Topological Torsion指纹:GetTopologicalTorsionFingerprint(mol, targetSize, ...)
  • 同GetTopologicalTorsionFingerprintAsIntVect()
    mol:mol对象
    targetSize:子图的长度
>>> from rdkit.Chem.AtomPairs import Torsions
>>> ttfp = Torsions.GetTopologicalTorsionFingerprintAsIntVect(ms[0])
  • 获取非空元素:GetNonzeroElements()
  • 查看位所代表的指纹:ExplainPairScore()
>>> ttfp.GetNonzeroElements()
{4320149536: 1}
>>> Torsions.ExplainPathScore(4320149536)
(('C', 1, 0), ('C', 2, 0), ('O', 2, 0), ('C', 1, 0))

4.Morgan Fingerprints (Circular Fingerprints)

摩根分子指纹,也成为圆形指纹,是采用摩根算法而产生。使用时,需要提供原子半径。这里只展示最基本的使用方法,更多关于指纹生成、提取与展示的操作可以参考这篇文章

  • 以SparseBitVects方式生成摩根指纹:GetMorganFingerprint(mol, radius)
    radius:考虑半径
>>> mfp = [AllChem.GetMorganFingerprint(x, 2) for x in ms]
>>> mfp[0].GetLength()
4294967295
  • 查看非空元素:GetNonzeroElements()
    这里记录的指纹是包括频数的
>>> mfp[2].GetNonzeroElements()
{864674487: 1, 2154640335: 1, 2246728737: 2, 3975275337: 2}
  • 以ExplicitBitVects方式生成摩根指纹:GetMorganFingerprintAsBitVect(mol, radius, nBits)
    radius:考虑半径
    nBits:指纹长度
    这种方法将不再记录频数
>>> mfp = [AllChem.GetMorganFingerprintAsBitVect(x, 2, nBits=10) for x in ms]
>>> print(mfp[0].GetNumBits())
10
>>> print(mfp[0].ToBitString())
0010100100
>>> print(DataStructs.DiceSimilarity(mfp[0], mfp[1]))
0.75

二、相似性

1.相似性比较

分子指纹的应用之一就是比较分子间相似性。

  • 相似度比较:DataStructs.FingerprintSimilarity(fp1, fp2, metric)
    fp1, fp2:待比较的分子
    metric:距离度量方法,默认为DataStructs.TanimotoSimilarity
>>> ms = [Chem.MolFromSmiles('CCOC'), Chem.MolFromSmiles('CCO')]
>>> fps = [Chem.RDKFingerprint(x) for x in ms]
>>> DataStructs.FingerprintSimilarity(mfp[0], mfp[1])
0.6
  • 度量方法非常多,随便举几例:Tanimoto, Dice, Cosine, Sokal, Russel, Kulczynski, McConnaughey等,按需选择。
>>> metic_list = ['DataStructs.TanimotoSimilarity',
>>>               'DataStructs.DiceSimilarity',
>>>               'DataStructs.CosineSimilarity',
>>>               'DataStructs.SokalSimilarity',
>>>               'DataStructs.RusselSimilarity',
>>>               'DataStructs.KulczynskiSimilarity',
>>>               'DataStructs.McConnaugheySimilarity']
>>> for i in metic_list:
>>>     print(DataStructs.FingerprintSimilarity(fps[0], fps[1], metric=eval(i)))
0.6
0.75
0.7745966692414834
0.42857142857142855
0.0029296875
0.8
0.6
  • 以上方法也可以直接调用,以Tanimoto为例:DataStructs.TanimotoSimilarity(bv1, bv2)
    bv1、bv2:两个分子指纹的bit vetor
>>> DataStructs.TanimotoSimilarity(fps[0], fps[1])
0.6

2.相似性地图

相似性地图是一种用来可视化分子中原子相似性程度的方法,具体可以参考rdkit中的rdkit.Chem.Draw.SimilarityMaps模块

  • 先来初始化两个分
>>> mol = Chem.MolFromSmiles('COc1cccc2cc(C(=O)NCCCCN3CCN(c4cccc5nccnc54)CC3)oc21')
>>> refmol = Chem.MolFromSmiles('CCCN(CCCCN1CCN(c2ccccc2OC)CC1)Cc1ccc2ccccc2c1')

相似性地图模块支持三种分子指纹的比较和可视化:atom pairs、topological torsions和Morgan fingerprints。

  • 使用atom pairs和topological torsions时,fpType可选值有:normal(默认值)、hashed、bit vector(bv)。
  • 使用morgan指纹时,fpType可选值为:bv(默认值)和count(count vector)
  • 该模块内的分子指纹方法也可以通过AllChem直接调用
>>> from rdkit.Chem.Draw import SimilarityMaps
>>> fp = SimilarityMaps.GetAPFingerprint(mol, fpType='normal')
>>> fp1 = AllChem.GetAtomPairFingerprint(mol)
>>> print(fp == fp1)
True
>>> fp = SimilarityMaps.GetTTFingerprint(mol, fpType='normal')
>>> fp1 = AllChem.GetTopologicalTorsionFingerprint(mol)
>>> print(fp == fp1)
True
>>> fp = SimilarityMaps.GetMorganFingerprint(mol, fpType='bv')
>>> fp1 = AllChem.GetMorganFingerprintAsBitVect(mol, 2)
>>> print(fp == fp1)
True
  • 生成相似性地图:GetSimilarityMapForFingerprint(refMol, probeMol, fpFunction, metric, ...)
    refMol:参照分子
    probeMol:探针分子
    fpFunction:指定所用的分子指纹方法
    metric:相似性度量方法(默认为DataStructs.DiceSimilarity)
>>> fig, maxweight = SimilarityMaps.GetSimilarityMapForFingerprint(refmol, mol, SimilarityMaps.GetMorganFingerprint)
2

3.相似性应用——构建多样化分子库

有时会遇到这种任务,需要从整个分子库中,挑出一个子集,使子集的多样性最高,或者能最大程度地代表原始分子库的化学空间。Rdkit在rdkit.SimDivFilters模块中提供了一系列的方法用来处理这种需求,最有效的是MaxMin算法,有文章报道该方法优于k-means。使用MaxMin算法从N个分子的库中挑选k个分子作为子集的大致流程如下:

  • 1.选择一个种子化合物用于初始化子集,并设置x=1。
  • 2.遍历分子库中剩余的N-x个分子,计算它与子集x中每个分子的相异性(可以认为是“1-相似度”),并取x个相异性值中最小的一个作为最终相异值(该值反应了库中分子与子集中最相似分子的相异性)。
  • 3.选择分子库中最终相异值最大的化合物(MaxMin名称的由来),将其加入到子集中,并更新x += 1,重复步骤2,直到x==k。
>>> from rdkit.SimDivFilters.rdSimDivPickers import >>> >>> MaxMinPicker
>>> ms = [Chem.MolFromSmiles('CCOC'),
>>>       Chem.MolFromSmiles('CCO'),
>>>       Chem.MolFromSmiles('COC'),
>>>       Chem.MolFromSmiles('CCC=C'),
>>>       Chem.MolFromSmiles('CCCO'),
>>>       Chem.MolFromSmiles('CCC=O'),
>>>       Chem.MolFromSmiles('CC(=O)C'),
>>>       Chem.MolFromSmiles('CCNC'),
>>>       Chem.MolFromSmiles('CCCCN'),
>>>       Chem.MolFromSmiles('CCC(N)C'),
>>>       Chem.MolFromSmiles('c1ccccc1'),
>>>       Chem.MolFromSmiles('c1cocc1'),]
>>> fps = [AllChem.GetMorganFingerprint(x, 2) for x in ms]
>>> Draw.MolsToGridImage(ms, molsPerRow=5)
3
  • 该算法需要先定义一种相异性度量方法,用于比较分子间的距离,以DiceSimilarity为例
>>> def distij(i, j, fps=fps):
>>>     return 1-DataStructs.TanimotoSimilarity(fps[i],fps[j])
  • 再创建一个picker对象:MaxMinPicker()
  • 挑选分子:picker.LazyPick(distFunc, poolSize, pickSize, firstPicks, seed, ...)
    distFunc:相异性度量函数,该函数可以接收两个索引,并返回相应的距离。
    poolSize:总分子数量
    pickSize:要挑出的分子数量
    firstPicks:第一个初始化的分子
    seed:随机种子
>>> picker = MaxMinPicker()
>>> pickIndices = picker.LazyPick(distij, len(ms), 4, seed=1)
>>> picks = [ms[x] for x in pickIndices]
>>> Draw.MolsToGridImage(picks, molsPerRow=4)
4
  • 顺便介绍下k-means基于聚类构建多样性分子库的方法,大致流程如下:
    • 1.初始化k个类别的代表(k为用户提供),并计算类的质心centroids
    • 2.遍历分子库中的每个分子,计算该分子与各个类的代表间的距离,并划分到距离最近的类中
    • 3.更新各个类别的代表,如果类别没有收敛,则重复步骤2
    • 4.从分好的各个类别中,选取若干分子组成最终的子集。

本文参考自rdkit官方文档
代码及源文件在这里

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 214,588评论 6 496
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 91,456评论 3 389
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 160,146评论 0 350
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 57,387评论 1 288
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 66,481评论 6 386
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 50,510评论 1 293
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 39,522评论 3 414
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 38,296评论 0 270
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 44,745评论 1 307
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 37,039评论 2 330
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 39,202评论 1 343
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 34,901评论 5 338
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 40,538评论 3 322
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 31,165评论 0 21
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,415评论 1 268
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 47,081评论 2 365
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 44,085评论 2 352