RDKit 电拓扑状态指数(E-State)模块详解:从 Kier–Hall 算法到描述符与指纹的完整实战指南
科学计算科研机器学习【免费下载链接】rdkitThe official sources for the RDKit library项目地址https://gitcode.com/gh_mirrors/rd/rdkit点击查看免费下载导读本文聚焦 RDKit 的rdkit.Chem.EState模块系统讲解基于 Kier 与 Hall 电拓扑状态指数Electrotopological State, E-State的一整套分子描述工具从 EState 指数的数学定义与 RDKit 源码实现、四个极值型派生描述符到 MOE 风格的 EState-VSA 组合描述符与 E-State 原子类型指纹。读完本文你将掌握该模块每个 API 的参数语义、缓存机制与算法细节并能直接把这些描述符接入 QSAR/QSPR 建模流程。该模块的官方文档入口是 Docs/Book/source/rdkit.Chem.EState.EState.rst它是一个 Sphinxautomodule文档桩真正的文档内容由 rdkit/Chem/EState/EState.py 等源码模块的 docstring 与代码自动生成因此本文以源码为准展开。一、模块总览rdkit.Chem.EState包含什么rdkit.Chem.EState是 RDKit 中专门实现电拓扑状态E-State理论的 Python 子包位于 rdkit/Chem/EState/包内按职责分为四个模块模块文件职责核心参考论文EState.pyEState 指数本体计算与极值描述符Hall, Mohney and Kier,JCICS31, 76–81 (1991)EState_VSA.py将 EState 指数与 Labute 近似表面积VSA组合成直方图描述符仿 MOE VSA 描述符设计AtomTypes.py79 种 E-State 原子类型的 SMARTS 定义与匹配Hall and Kier,JCICS35, 1039–1045 (1995) Table 1Fingerprinter.py基于原子类型的 E-State 指纹生成同上包入口init.py 通过from rdkit.Chem.EState.AtomTypes import *和from rdkit.Chem.EState.EState import *将原子类型工具与 EState 指数函数直接暴露在rdkit.Chem.EState命名空间下因此实际使用时既可以直接from rdkit.Chem.EState import EStateIndices也可以from rdkit.Chem import EState后通过EState.EStateIndices访问。二、核心算法EStateIndices()的完整实现拆解EStateIndices(mol, forceTrue)返回分子中每个重原子非氢原子的 EState 指数按原子索引顺序组织为 numpy 数组。它是整个模块的算法核心其余描述符与指纹均建立在其之上。其计算流程可分两步理解原子内禀状态intrinsic state计算与基于距离矩阵的分子环境修正。2.1 第一步内禀状态I的计算对每个原子i源码先计算内禀状态tbl Chem.GetPeriodicTable() for i in range(nAtoms): at mol.GetAtomWithIdx(i) d at.GetDegree() if d 0: atNum at.GetAtomicNum() dv tbl.GetNOuterElecs(atNum) - at.GetTotalNumHs() N GetPrincipleQuantumNumber(atNum) Is[i] (4. / (N * N) * dv 1) / d对应 EState.py其数学本质是 Kier–Hall 定义的内禀值I ( (4 / N²) · dv 1 ) / δ其中N为原子的主量子数principal quantum number由GetPrincipleQuantumNumber(atNum)依据原子序数区间查表得到1–2 → 13–10 → 211–18 → 319–36 → 437–54 → 555–86 → 687 以上 → 7见 EState.pydv为该原子的价电子数PeriodicTable.GetNOuterElecs减去其键连的氢原子总数GetTotalNumHs()即非氢价电子数d为原子的拓扑度degree即非氢键数。因此高周期大 N、高电负性差大 dv的原子内禀状态高而高度分支大 d的原子内禀状态被摊薄。从GetPrincipleQuantumNumber的实现可以推断N²出现在分母正是为了将不同周期的原子归一到同一能量尺度。2.2 第二步基于图距离矩阵的修正内禀状态只反映原子自身的电子结构E-State 的拓扑属性体现在第二步的环境修正中EState.pydists Chem.GetDistanceMatrix(mol, useBOFalse, useAtomWtsFalse) 1 for i in range(nAtoms): for j in range(i 1, nAtoms): p dists[i, j] if p 1e6: tmp (Is[i] - Is[j]) / (p * p) accum[i] tmp accum[j] - tmp res accum Is这里的dists是 RDKit 计算的全对最短路径距离矩阵不加 1 之前的图距离useBOFalse, useAtomWtsFalse表示不按键级加权重距离加 1 后将键连原子对的距离置为 2避免除零并保证同一键两端原子也有修正项。对每一对原子(i, j)差值(Is[i] - Is[j]) / p²被加入i的累积项、反向减入j的累积项——这正是 E-State 指数的核心思想原子的电子状态沿拓扑路径向周围原子传递距离越近影响越大与距离平方成反比。最终指数EState I Σ修正项。代码中的p 1e6是一个数值保护条件对不连通分子片段GetDistanceMatrix会产生极大的无穷远距离值该判断用于跳过此类无效对使算法天然兼容多片段分子如带反离子的盐见 2.4 示例。2.3 缓存机制force参数的语义EStateIndices(mol, forceTrue)的force参数控制结果缓存if not force and hasattr(mol, _eStateIndices): return mol._eStateIndices计算完成后结果被写入mol._eStateIndices属性。当forceFalse且分子上已有该属性时直接返回缓存避免重复计算forceTrue默认则无条件重算。该缓存行为在 UnitTestEState.py 的test_cacheEstate中有专门验证测试先断言新分子没有_eStateIndices属性计算后属性出现再把缓存属性改为字符串cached并确认forceFalse时原样返回该字符串而默认调用会重新计算覆盖缓存。这说明缓存与分子对象绑定同一分子对象的描述符计算是惰性共享的——EState_VSA 等依赖方传入forceFalse即可复用同一次计算。2.4 运行示例与数值直觉模块自带的_exampleCode()EState.py给出了可直接运行的演示smis [CCCC, CCCCC, CCCCCC, CC(N)C(O)O, CC(N)C(O)[O-].[Na]] for smi in smis: m Chem.MolFromSmiles(smi) print(smi) inds EStateIndices(m) print(\t, inds)配合测试数据可以建立数值直觉来自 UnitTestEState.py 的验证值即论文原始数据正构烷烃CCCC→[2.18, 1.32, 1.32, 2.18]端点碳高于中间碳异构烷烃CCC(C)(C)CC→ 季碳中心降至0.54分支效应显著引入杂原子后电负性原子数值骤升如CC(C)CO的氧端达8.14、CC(C)CF的氟端达11.11、氨基酸骨架CC(N)C(O)O的羰基碳9.57与两个羧基氧7.86/4.84芳香体系如Fc1ccc(C)cc1的氟原子达12.09而邻位碳出现负值-0.17说明 EState 可以刻画供电子/吸电子效应在环上的分布。测试代码中的_validate用tol1e-2的容差将 RDKit 计算值与论文值比对MaxEStateIndex/MinEStateIndex/MaxAbsEStateIndex/MinAbsEStateIndex四个极值也一并验证UnitTestEState.py。三、四个极值型派生描述符EState.py在EStateIndices之上提供了四个聚合描述符源码实现极为精简def MaxEStateIndex(mol, force1): return max(EStateIndices(mol, force)) def MinEStateIndex(mol, force1): return min(EStateIndices(mol, force)) def MaxAbsEStateIndex(mol, force1): return max(abs(x) for x in EStateIndices(mol, force)) def MinAbsEStateIndex(mol, force1): return min(abs(x) for x in EStateIndices(mol, force))四个函数均接受与EStateIndices相同的force参数并向下透传各自带有version 1.0.0版本标记。它们分别刻画分子中最富电子原子Max、最缺电子原子Min常为负值、偏离中性程度最大/最小的原子MaxAbs/MinAbs。在 RDKit 的顶层描述符体系里这四个函数被直接导入并注册为全局描述符Descriptors.py 中执行from rdkit.Chem.EState.EState import (MaxAbsEStateIndex, MaxEStateIndex, MinAbsEStateIndex, MinEStateIndex)随后经_setupDescriptors收集进descList。因此用户无需手动计算 EState 数组即可通过Descriptors模块的_descList/ 计算函数清单直接批量获得这四个分子级特征用于建模特征矩阵。四、MOE 风格混合描述符EState_VSA与VSA_EStateEState_VSA.py 实现了仿 MOE VSA 的混合描述符将 EState 指数与 MolSurf.py 的 Labute 原子表面积贡献_LabuteHelper as VSAContribs_组合成直方图。模块顶部注释明确说明分箱边界是用 PP3K 水溶性数据集挑选、按每个箱子原子数大致相等原则确定的EState_VSA.py。4.1 两组默认分箱vsaBins [4.78, 5.00, 5.410, 5.740, 6.00, 6.07, 6.45, 7.00, 11.0] # VSA_EState 用 estateBins [-0.390, 0.290, 0.717, 1.165, 1.540, 1.807, 2.05, 4.69, 9.17, 15.0] # EState_VSA 用4.2 两种组合方向VSA_EState_(mol, binsNone, force1)对每个原子按它的VSA 值落入的箱号把该原子的EState 指数累加进去bisect.bisect_right(bins, volContribs[i1])。即按表面积分箱、按 EState 加权EState_VSA_(mol, binsNone, force1)反过来按原子的EState 值分箱把它的VSA 贡献累加进去。即按 EState 分箱、按表面积加权。两者都支持传入自定义bins覆盖默认边界也都有与EStateIndices一致的缓存属性mol._vsaEState、mol._eStateVSA和force语义。最终数组长度恒为len(bins)1最后一个箱子收纳大于最大边界的值。4.3 描述符命名与注册_InstallDescriptors()在模块导入时动态生成全部描述符EState_VSA.pyVSA_EState1…VSA_EState10基于 9 个 vsaBins 边界 → 10 个箱子版本1.0.0EState_VSA1…EState_VSA11基于 10 个 estateBins 边界 → 11 个箱子版本1.0.1源码变更日志注明 1.0.1 为优化改动、数值不受影响。每个描述符的 docstring 由_descriptorDocstring动态生成明确写出区间语义例如第一个箱子的区间是-inf x 4.78最后一个箱子是11.00 x infEState_VSA.py。这些函数经 Descriptors.py 的from rdkit.Chem.EState import EState_VSA与_setupDescriptors注册进_descList因而EState_VSA11、VSA_EState10等名字会出现在rdkit.Chem.Descriptors的描述符清单中。这一组描述符的数值正确性由 UnitTestVSA.py 保障它读取参考数据 EState_VSA.csv文件头即为描述符名逐分子、逐描述符以delta1e-4的精度断言计算值与参考值一致。五、E-State 原子类型与指纹5.1AtomTypes.py79 种原子类型的 SMARTS 定义AtomTypes.py 依据 Hall 与 Kier 1995 年论文 Table 1定义了 79 种 E-State 原子类型_rawD列表覆盖 Li、Be、B、C、N、O、F、Si、P、S、Cl、Ge、As、Se、Br、Sn、I、Pb 等元素的各类化学环境。类型命名遵循键型前缀 元素 氢数后缀规则例如sCH3→[CD1H3]-*连单键、带 3 个 H 的饱和碳甲基dCH2→[CD1H2]*连双键的亚甲基sssCH→CD3H(-*)-*连三个单键的次甲基dssC→CD3H0(-*)-*一个双键加两个单键的季碳aaaC→C,c;D3H0(:*):*芳香碳sOH→[OD1H]-*、dO→[OD1H0]*、ssO→OD2H0-*氧的不同键合态。代码注释中的# mod标记如ddsN、aasN、ddssS表示这些 SMARTS 相对论文原始定义做过微调。所有模式统一用[元素][D/键]...[H]的 SMARTS 语法编码原子价态D1–D5、键型单键-、双键、芳香键:与氢数H0–H3由BuildPatts()惰性编译esPatterns全局变量首次使用时构建解析失败的条目会打印 WARNING 并跳过见 AtomTypes.py。5.2TypeAtoms()逐原子类型标注def TypeAtoms(mol): ... for name, patt in esPatterns: matches mol.GetSubstructMatches(patt, uniquifyFalse) for match in matches: idx match[0] ...由于类型模式只匹配单原子模式中心是第一个原子TypeAtoms对每个原子收集其命中的所有类型名返回一个长度等于原子数的元组列表一个原子可能匹配多个类型如杂环氮。uniquifyFalse保证同一原子的多次匹配都被记录。5.3Fingerprinter.FingerprintMol()E-State 指纹def FingerprintMol(mol): esIndices EStateIndices(mol) ... counts numpy.zeros(nPatts, dtypenumpy.int64) sums numpy.zeros(nPatts, dtypenumpy.float64) for i, (_, pattern) in enumerate(AtomTypes.esPatterns): matches mol.GetSubstructMatches(pattern, uniquifyTrue) counts[i] len(matches) for match in matches: sums[i] esIndices[match[0]] return counts, sums对应 Fingerprinter.py该函数返回两个等长长度 原子类型数 79的 numpy 数组countsint64每种原子类型在分子中出现的次数sumsfloat64命中该类型的原子的 EState 指数之和。这样既保留了类型分布counts又用 EState 指数对同型原子按电子环境加权sums指纹向量天然定长、适合直接做相似度计算或作为机器学习特征。注意uniquifyTrue的语义对单原子模式而言它保证每个原子在每个类型下只计一次命中。六、测试与数据验证数值可信度的来源整个 EState 子包配备了完整的单元测试是理解各 API 预期行为的最好参考UnitTestEState.py核心指数验证测试数据直接取自 1991 年原始论文覆盖直链/支链烷烃test_simpleMolecules、test_isomers、杂原子取代test_heteroatoms1/2含 N、O、F、Cl、Br、I、S 与羰基/硫羰基、芳香卤代体系test_aromatics同时验证GetPrincipleQuantumNumber的原子序数区间映射、force缓存语义与_exampleCode可运行性UnitTestVSA.py以 EState_VSA.csv 为基准逐描述符断言delta1e-4精度UnitTestTypes.py 与 UnitTestFingerprints.py分别覆盖原子类型标注与指纹生成。值得一提的是 UnitTestEState.py 中有一条注释# NOTE: this doesnt match the values in the paper硫醚CCSCC的测试值这是 RDKit 对自己与论文原始数值差异的显式记录提示用户若追求与 1991 年论文逐位一致个别含硫体系可能存在已知偏差而EStateIndices.version、各描述符的version属性则为追踪实现版本提供了机制。七、实战建议如何把这些描述符接入建模流程综合以上源码分析实际使用时的推荐路径快速体验直接调用rdkit.Chem.EState.EStateIndices(mol)查看逐原子指数使用Chem.MolFromSmiles构造分子后无需加氢预处理EState 计算基于非氢原子骨架。分子级特征通过from rdkit.Chem import Descriptors直接取Descriptors.MaxEStateIndex、MinEStateIndex、MaxAbsEStateIndex、MinAbsEStateIndex以及Descriptors.EState_VSA1…EState_VSA11、Descriptors.VSA_EState1…VSA_EState10它们均已注册进descList无需手动写循环。批量计算效率多次调用同一分子的相关描述符时利用各函数默认的缓存与forceFalse复用机制避免重复计算若分子结构被修改如加氢、成键务必使用默认forceTrue强制刷新。指纹用途需要定长向量做化合物相似性/聚类时使用Fingerprinter.FingerprintMol(mol)得到的(counts, sums)二元组需要逐原子类型信息时用AtomTypes.TypeAtoms(mol)。自定义分箱EState_VSA 与 VSA_EState 接受bins参数可按自己的数据集分布重新设定边界——源码注释与_descriptorDocstring都表明默认箱来自 PP3K 数据集自定义前可参考其近似等原子数分箱的选取思路。八、总结rdkit.Chem.EState是 RDKit 中对 Kier–Hall 电拓扑状态理论最完整的实现EStateIndices以内禀状态 距离平方反比修正两步算法输出逐原子指数四个极值函数与 21 个 EState-VSA 混合描述符通过 Descriptors.py 无缝接入 RDKit 顶层描述符体系79 种原子类型 SMARTS 与FingerprintMol则提供了基于 E-State 的定长分子指纹。算法数值由取自原始论文的单元测试与 CSV 参考数据双重校验缓存机制保证了批量计算性能。无论你是构建 QSAR 模型、计算分子相似度还是研究电子效应分布这个模块都能直接提供经论文验证、开箱即用的特征。赞分享科学计算科研机器学习【免费下载链接】rdkitThe official sources for the RDKit library项目地址https://gitcode.com/gh_mirrors/rd/rdkit点击查看免费下载相关推荐OpenArkWindows内核级Anti-Rootkit排查实操笔记OpenArkWindows内核级Anti Rootkit排查实操笔记 OpenArk 是一款 Windows 平台上的开源反 RootkitARK工具网络安全逆向工程桌面应用微信聊天记录怎么免费导出并永久保存WeChatMsg 新手指南微信聊天记录怎么免费导出并永久保存WeChatMsg 新手指南 你手机里那些几百上千天的微信聊天记录真的在云端吗其实没有。它们锁在 PC 版微信生成的科学计算科研机器学习Cline PR 审查工作流规则实战基于 gh CLI 的 Pull Request 评审全流程指南Cline PR 审查工作流规则实战基于 gh CLI 的 Pull Request 评审全流程指南 本篇指南围绕 Cline 仓库内 .clinerules科学计算科研机器学习上一篇如何快速开始使用libgdbm-rust10分钟安装与基础使用指南下一篇InternVL-C诞生全解InternViT-6B如何蜕变为CLIP式图像-文本检索模型创作声明:本文部分内容由AI辅助生成(AIGC),仅供参考