Chem Sci｜分子形变能量计算：南方科技大学钟龙华团队推出D2AF软件

人们开发了许多理论和模型来了解化学过程。例如，过渡态理论 (TST)和马库斯理论，以及基于轨道的理论如前线分子轨道 (FMO)、伍德沃德-霍夫曼规则、价键 (VB)、自然键轨道 (NBO)从轨道角度提供了对反应更直观的见解。由Houk和Bickelhaupt等人开创的Distortion/Interaction-Activation Strain，简称为D/IAS模型能分析化学反应速率。该模型结合了分子几何形变（畸变）和反应物之间的相互作用，来解释和预测化学反应的活化能。

为了阐明不同阶段（来自内禀反应坐标IRC) 或分子动力学模拟的结构的形变分布，2024年12月23日，南方科技大学的钟龙华团队在Chemical Science上发表文章An efficient and flexible approach for local distortion: distortion distribution analysis enabled by fragmentation，提出了一种高效通用灵活的基于片段化的方法D2AF，来估计各种分子的局部形变能量。

概念介绍

首先介绍几个概念，畸变能（Distortion Energy）：指的是反应物从其平衡几何形态变为过渡态几何形态时所需要的能量。这部分能量通常来自反应物的结构调整和键的断裂等。相互作用能（Interaction Energy）：反应物在过渡态几何形态下的相互作用能。这与反应物的电子结构以及它们之间的相对定位密切相关，影响过渡态的稳定性。活化能（Activation Energy）：是化学反应过程中所需克服的能量屏障，通常包括畸变能和相互作用能的综合。

D/IAS模型的核心思想：强调反应物的畸变能和分子之间的相互作用如何共同决定反应的活化能。过渡态的稳定性由这两者的相互作用决定，而这种稳定性又直接影响反应速率。该模型将活化能（或沿反应坐标 (ζ) 的相对电子能；分别为ΔE或 ΔE (ζ)）分解为反应物片段的扭曲能（Edist）及其相互作用能（Eint）。

图1 D/IAS模型示意图

基于片段化的方法D2AF，来估计各种分子的局部形变能量，一来能可视化分子内的形变能量分布（即形变图）并识别关键的形变部分，还提供局部形变能量的值，这些值也可用作统计（例如多元线性回归，MLR）或机器学习（ML）分析/模型的新描述符。

该方法基于分子片段化策略，并结合ONIOM类型的边界和链接原子（Link Atom）方法，分析每个片段在目标状态和参考状态之间的形变能量差异。如图1所示，基于片段化的方法包括三个阶段：片段化、计算和可视化。片段化（Fragmentation）：首先选择目标分子（通常具有较大形变）和其参考形式，并将其分解成较小的片段。这些片段可以根据不同的目的和分子系统的复杂性进行调整。能量计算：对每个片段进行能量计算，通常使用量子力学（QM）或机器学习（ML）方法。形变能量计算与可视化：通过对各个片段的形变能量进行计算，生成形变图，从而直观地显示分子中不同部分的局部形变分布。

用户首先提供两个结构：目标分子（Tar；通常具有更多扭曲）及其参考（Ref）形式。分子内局部扭曲能量来自每个片段的参考形式和目标形式之间的能量差异。局部扭曲能量可以通过三种方法（M1-M3）中的一种来确定和分类，M1、M2、M3具有不同的片段化和坐标操作方案，具体取决于目的和系统的复杂性，用于不同的片段化策略。M1方法：适用于较简单的片段化，计算每个片段的形变能量。M2方法：对于涉及多个内部坐标（如键角扭转角等）的分子，进行更为细致的形变能量分解。M3方法：结合了M1和M2的优点，适用于复杂的共轭体系或金属配位体系。

图2 碎片化 (D2AF) 进行形变分布分析的工作流程。(a) 输入结构包含目标系统 ( Tar ) 及其参考 ( Ref ) 形式。(b) 示意图碎片化和坐标操作方法 ( M1-M3 )。(c) ONIOM 型连接原子 (LA) 处理。(d) 碎片的能量计算。(e) 畸变图可视化

连接原子处理

当涉及单键、双键或三键的边界断裂时，碎裂过程通常需要引入连接原子（LA）来覆盖所有生成的片段中的所有悬键。采用ONIOM型边界和LA方法。因此，连接原子的位置不是固定的，而是取决于两个边界原子的键距（图1c）。对不同连接原子的可靠性进行了系统的基准研究发现H-LA、C-LA和N-LA 通常分别适用于涉及单键、双键和三键的边界。

计算

使用用户指定的QM或机器学习势 (MLP) 方法和程序计算所有生成的片段的能量 (图1d)。当前的Python包与一些QM和MLP程序 (例如Gaussian、ORCA、xTB) 交互。用户也可以使用脚本为特定程序和方法准备输入文件，调用QM计算并提取每个片段的QM计算能量。

分析和可视化

计算完每个片段的能量后，可以根据M1-M3方案确定各个目标片段和参考片段的扭曲能量。此外，还生成了PyMOL脚本以可视化（相对）扭曲能量的分布（图1e）M1中的片段根据其相应的局部扭曲能量着色，而M2中的局部扭曲能量可以进一步分解为特定的键合项。

计算细节

所有这些计算都是使用开源Python包（D2AF，https://github.com/oscarchung-lab/D2AF）执行的，该包使用Open Babel包来处理片段的内部坐标。

接着将举几个例子说明。

逆Cope消除反应

图2a展示了一个分子内协同加成反应，其中质子和胺加到三键上，同时 O-H键断裂。Houk等地研究了该反应机理。DFT结果表明，主要的扭曲来自炔烃（∼17.6 kcal/mol）和羟胺（∼14.5 kcal/mol）部分，以及一些来自束缚部分（∼8.9 kcal/ mol）。当采用M1和M2方案对该反应的扭曲分布进行详细分析（图2b-g）。使用M1的扭曲分布支持显著的扭曲来自炔烃和羟胺部分（ΔE (I)：∼17.6 kcal/mol和 ΔE (II)：∼14.6 kcal/mol；图2c）。M2的扭曲分布通常显示出类似的分布（图 2e），即三键弯曲是主要扭曲贡献者（Δ（A1）：-41.9°，扭曲能量为∼9.6 kcal/mol；Δ（A2）：-19.6°，扭曲能量为∼3.6kcal/mol），羟胺部分的O-H键拉伸/O-N键压缩是显著贡献。

图3 逆 Cope消除反应的分析

Diels-Alder反应

第二个系统是环戊烯酮和环戊二烯之间的Diels-Alder反应（图3）。如图3c和e所示，主要的畸变主要由加成位点造成。此外，Houk等人开创性地发现在一系列 [4+2]加成反应中，分子畸变与活化能之间存在很强的线性相关性（R2=0.93）。进一步研究还观察到活化能和最大局部畸变分量（M1方案中的片段I；M2方案中的B1）之间的优异的线性相关性（图3h，R2∼0.96）。

正如 Bickelhaupt和Houk深刻指出的那样，两种反应物片段进入过渡态几何结构的畸变能与Marcus理论中的能量曲线相关（或类似于VB理论中的两条主态曲线）。因此，反应性被认为受热力学（影响过渡态的位置）和畸变能的影响。Houk团队提出，C-H键弯曲出反应烯烃C=C平面在新的C–C键形成过程中起着关键作用，这与最大的局部畸变分量（I）一致相关。这涉及关键键的延长，是关键坐标，从而控制着势垒。结果表明，主要的局部扭曲不仅有助于理解障碍的主要来源，而且还可以作为多元线性回归 (MLR) 和机器学习 (ML) 建模的新描述符/特征。

图4 Diels-Alder环加成反应概述

蛋白质-药物结合

以脾脏酪氨酸激酶蛋白结构中的细化伊马替尼药物结构为参考态，以其扭曲的X射线晶体结构为目标态（图4）。发现通过量化计算优化改进的药物结构表现出比其晶体结构（∼22.6 kcal/mol）低得多的扭曲能量。M1分析揭示出晶体结构的主要扭曲源自嘧啶（环1）和吡啶（环2）环（ΔE:分别为∼9.6和∼3.1kcal mol；图4c）。此外，M3结果进一步描绘出由吡啶和嘧啶环之间的 键距 (B1) 改变-0.07 Å而产生的畸变能∼2.7 kcal/mol。

图4 （a）伊马替尼 - 脾酪氨酸激酶结构概述。(b) 使用M1方案的畸变分布（kcal/mol）和（c）片段。(d)从量化计算优化结构到 X 射线衍射晶体结构的键长变化（Δ Bond in Å）。(e)使用M3方案的总畸变分布，(f) 键和 (g) 角度项的畸变贡献。

总结

本文例子很多 就不一 一介绍了，做下总结。研究团队提出并实施了一种基于片段化和ONIOM型边界策略的通用、高效且灵活的局部畸变分布分析方法。然后广泛比较了各种化学（有机、超分子、无机和金属配位系统）和生物（天然和人工酶以及含药物的蛋白质）系统的畸变分布。可以成功揭示目标结构与参考结构相比的定性（相对）畸变分布，其中可以采用和比较任何类型的结构（例如局部最小值、过渡态、MD或IRC衍生结构；甚至使用实验确定的晶体结构）。此外，可以对多个结构执行局部畸变分布分析，或使用不同的计算方法（如高级CCSD(T)、经济高效的xTB和MLP方法）来计算片段的能量，可以轻松应用和比较。可以采用快速xTB方法来获得初始畸变图并评估片段化设置。用户可以进一步采用高级计算方法来获得更高质量的形变图。

方法优势与挑战

优势：D2AF方法具有高效性和灵活性，能够对复杂分子系统进行形变能量分析，并通过不同的计算方法得到精确的局部形变分布。该方法的可扩展性使其能够应用于多种类型的分子，涵盖了有机、金属配位、超分子及生物分子等系统。

挑战：对于强去共轭体系（如某些金属配位系统或小环状分子），形变能量的分解可能会受到一定限制，尤其是在处理复杂的共轭系统时，可能需要结合不同的片段化方法以获得更好的结果。

课题组介绍

笔者第一次了解听说钟龙华（Lung Wa Chung）教授是源于2015年他在Chemical Reviews上发表的题为：The ONIOM Method and Its Applications的经典综述。

研究领域：

从事多尺度和机器学习模拟研究复杂均相催化化学和生物化学体系，包括生物与仿生催化反应机理与设计，可持续化学反应机理，光化学和量子隧道效应，和新模拟方法的发展及应用。

钟龙华博士(博导，独立课题组组长)，分别在00、03和06年获得香港科技大学一级荣誉学士、硕士和博士学位，师从吴云东院士，研究(贵金属和主族金属)金属催化有机反应机理。06-13年到京都大学福井謙一記念研究中心諸熊奎治 (Keiji Morokuma) 院士组，采用和发展多尺度模拟方法研究复杂生物体系。13年10月起加入南方科技大学化学系(预聘助理教授)，18年2月(破格)晋升预聘副教授，21年7月(破格)晋升长聘正教授。

从事多尺度和机器学习 (AI/ML) 模拟研究复杂均相催化化学和生物化学体系，包括生物与仿生催化反应机理与设计，可持续化学反应机理，光化学和量子隧道效应，和新模拟方法的发展及应用。部分提出新的反应机理被实验支持或证实，其中两个反应机理被Hartwig教授的《Organotransition Metal Chemistry》教科书和Roberts教授的《Encyclopedia of Biophysics》收录与讨论。团队24年获英国皇家化学学会的 Organic Chemistry Horizon Prize (Robert Robinson Prize in Synthetic Organic Chemistry)

学术成果

在国外期刊上发表98篇SCI论文(包括1篇Chem. Rev.，2篇Acc. Chem. Res. 21篇J. Am. Chem. Soc.，1篇Wiley跨学科综述，6篇Angew. Chem. Int. Ed.，5篇Nat. Commun.，7篇Chem. Sci.，5篇ACS Catal.，5篇J. Phys. Chem.，5篇J. Chem. Theory Comput.，5篇Org. Lett.，3篇Chem. Commun.，3篇Inorg. Chem.)和3篇Wiley书籍章节。Web of Science引用>5400次（他引>5100次），H-index: 41（23年6月）。入选了斯坦福大学发布的2021年，2022年和2023年的全球前2％顶尖科学家。

参考文献

1.L. W. Chung, W. M. C. Sameera, R. Ramozzi, A. J. Page, M. Hatanaka, G. P. Petrova, T. V. Harris, X. Li, Z. Ke, F. Liu, H.-B. Li, L. Ding and K. Morokuma, Chem. Rev., 2015, 115, 5678–5796.

2.F. M. Bickelhaupt and K. N. Houk, Angew. Chem., Int. Ed., 2017, 56, 10070–10086.

3.J.-X. Zhang, F. K. Sheong and Z. Lin, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2020, 10, e1469.