面向NDDO后裔半经验方法的近黑箱参数优化器：MNDO与AM1的演示验证

《The Journal of Physical Chemistry A》：A Near Black Box Parameter Optimizer for NDDO-Descendant Semiempirical Methods: Demonstrations for MNDO and AM1

本文报道了一种增强版的几何校正参数优化算法，该算法能够直接整合关键几何变量（如键长、键角、二面角）的导数，克服了传统PARAM程序依赖衍生几何参考函数的局限性。通过对MNDO和AM1方法的重新参数化演示，研究揭示了衍生几何参考函数在参数优化中的不适用性，并可能无意中激励化学键力常数的低估。这项工作为提升半经验量子化学（SQC）模型性能提供了重要的算法基础和理论见解。

　　

理论背景与算法发展

半经验量子化学（SQC）方法，特别是NDDO（Neglect of Diatomic Differential Overlap）后裔模型，如PM6-D3H+和PM6-D3H4X，在蛋白质建模和其他生化系统研究中持续展现出与新型模型相媲美的性能。这些模型在过去半个世纪里精度的显著提升，主要源于对核心排斥函数的经验性修正、色散或氢键校正的加入以及参数化程序的改进。然而，参数优化过程中的一个关键挑战在于如何处理几何参考数据。直接包含键长（L_AB）、键角（θ_ABC）和二面角（?_ABCD）等几何变量十分困难，因为大分子的预测平衡几何结构对几何优化程序的精度极其敏感，且存在许多低能量构象体，使得直接包含二面角尤其困难。

在PM3开发过程中，有人提出使用在参考（实验）几何结构下评估的生成热相对于几何变量的导数，来代替几何变量本身。这些衍生的几何参考函数易于评估，且数值微分不会产生显著的浮点误差，它们的应用为PM6和PM7等通用半经验方法的发展铺平了道路。然而，从一开始就认识到，使用衍生几何参考函数也存在一些缺点：实际几何误差与评估导数之间的比例常数被确定为相关变形的力常数，其结果是，具有高力常数的化学键（即双键或三键）在参数优化过程中会被更优先考虑。

几何变量的导数计算

本研究在前人工作的基础上，详细描述了一种能够解释参考几何结构的增强版几何校正参数优化算法。该算法的核心在于计算预测分子性质相对于半经验参数的全导数。如公式（1a）和（1b）所示，该导数取决于冻结几何偏导数（在参数变化时保持几何结构不变）和一个几何校正项，该校正项描述了由于参数变化导致平衡几何结构改变如何影响分子性质。

对于键长、键角和二面角等几何变量，它们不显式地依赖于半经验参数，因此它们相对于半经验参数的全导数计算仅涉及几何校正项，具体公式见（3a）和（3b）。这些导数与Wilson的B矩阵中出现的项一致，该矩阵对于笛卡尔坐标和冗余内坐标之间的转换至关重要。

为了验证解析导数的准确性，研究以苯胺为例，展示了其C-N键长、H-N-C键角和H-N-C-C二面角对MNDO方法中氮原子2s轨道Slater指数（ζ_s^N）参数的依赖性。通过有限差分法（公式4a, 4b）进行比较，发现由于程序实现（如使用STO-6G展开而非STO基组进行分子积分计算）和几何优化 routines 的细微差异，导致基于MOPAC程序的验证不可行。因此，研究通过内部Java程序的有限差分比较来验证解析导数的准确性，并最终通过参数优化器成功识别参数曲面上的最小值证明了其可靠性。

CHNO分子的重新参数化

将参数优化算法扩展至包含几何变量后，研究尝试使用PM7训练集中几乎所有的CHNO分子，对MNDO和AM1进行完整的重新参数化。训练集最终包含了来自1,187个CHNO分子的1,616个参考数据。权重函数的选择力求与PARAM程序默认值匹配，其中键长权重（C_L）被视为可调参数，使用100 ?^-1的倍数进行了多组参数优化尝试。

参数优化的目标是最小化加权平方和误差（即误差函数S，公式5）。这通常通过非凸优化算法实现，其中第k个优化周期的搜索方向d_k是参数梯度g^{p_k}、参数Hessian矩阵H^{p_k}和信任半径r_k的函数（公式6a, 6b）。本研究采用基于修正Hessian矩阵（公式7a-c）的限制步长有理函数优化（RSRFO）程序（公式7d）来确定搜索方向，该算法在之前的几何校正优化器中表现出优异性能。

MNDO的重新参数化结果

从原始MNDO参数出发，采用键长权重C_L= 500 ?^-1进行RSRFO优化，经过193个优化周期后找到了一个局部最小值。随后尝试了键长权重为200至1000 ?^-1的优化。所有参数化运行均 uninterrupted，优化算法在优化MNDO参数时表现出合理的“黑箱”特性。

作为对比，使用相同的参考数据集在PARAM程序上进行了200个优化周期的参数化。由于训练集中各类实验数据量的不平衡，预测生成热和预测键长的误差贡献在优化过程中占据了总误差函数的绝大部分，这使得评估参数优化的有效性在很大程度上取决于这两类误差的改善情况。

对比分析表明，对于C_L= 400至800 ?^-1的五个参数集，其PARAM误差函数值均低于经过200个PARAM周期得到的参数集，这进一步证实了PARAM由于在导数评估程序中缺乏几何校正项而难以识别其自身误差函数上的最小值。总体而言，几何校正重新参数化带来了约5%的误差函数改进，这与之前排除几何参考函数时观察到的结果一致。

AM1的重新参数化结果

将已实现的MNDO参数优化算法扩展到AM1相对简单，因为AM1引入的额外参数仅影响核心排斥函数，无需对现有的矩阵和矩阵导数评估程序进行实质性更改。

由于AM1参数数量更多且核心排斥函数更复杂，收敛到梯度范数∥g^{p}∥ < 100且参数Hessian矩阵正定的点变得异常困难和耗时。从原始AM1参数开始，使用C_L= 500 ?^-1的优化程序经过975个优化周期后才识别出一个局部最小值。尝试使用键长权重为300、800和900 ?^-1的参数化未能成功识别参数曲面上的最小值。

与经过1000个PARAM优化周期得到的参数进行比较发现，尽管几何校正重新参数化在预测键长、键角和二面角方面取得了显著改进（平方和几何误差减少高达50%，无符号几何误差减少高达20%），同时也在预测生成热方面优于PARAM，但其参数集却导致了异常高的PARAM几何误差贡献。这一看似矛盾的现象促使研究进一步探究衍生几何参考函数的适用性。

衍生几何参考函数的不适用性分析

通过对H₂和CH₄等简单分子的模型计算，研究揭示了衍生几何参考函数存在的根本缺陷。以H₂为例，调整MNDO计算中的α^H参数会影响平衡键长。绘制PARAM几何误差贡献（E_geo^PARAM）随预测H₂键长变化的参数曲线发现，E_geo^{PARAM对较短的预测键长惩罚过重。}

这种对平衡键长高估和低估的不对称惩罚可以解释为：较小的α^H值导致较浅的势阱和较小的H-H键力常数，从而减少了相应衍生几何参考函数所施加的惩罚。因此，应用衍生几何参考函数不仅会在参数化过程中过度强调双键和三键，还会无意中激励更浅的势阱和更小的力常数。

更令人担忧的是，对于CH₄，改变ζ_s^C参数绘制E_geo^{PARAM随预测C-H键长变化的参数曲线表明，两个不同的ζ_sC值即使预测出完全相同的平衡键长，也可能由于预测力常数的差异而导致截然不同的E_geoPARAM值。这无可辩驳地证明了使用衍生几何参考函数存在重大缺陷。}

在AM1的PARAM参数化过程中，观察到许多分子（如季戊四醇四硝酸酯、甲烷、异丁烷）的E_geo^{PARAM值显著下降，但这并非源于预测几何精度的实际提高，而是由于相关力常数的减小所致。势能曲线分析证实，使用衍生几何变量进行参数化后，相关化学键的力常数均有所降低。}

非平衡几何结构上的评估

为了评估重新参数化对非平衡几何结构的描述能力，研究使用PM7-TS中部分CHNO反应能垒进行了测试。核密度估计（KDE）图显示，原始MNDO/AM1参数严重高估了大部分能垒，而经过PARAM或几何校正算法重新参数化后，这种趋势均得到显著改善。

对于MNDO，两种重新参数化方法预测的能垒误差分布几乎相同。然而，对于AM1，PARAM重新参数化的结果在预测能垒方面明显比几何校正重新参数化更准确。研究者认为，这可能是由于初始AM1参数系统性地高估了评估集中的能垒，而PARAM参数化过程中观察到的力常数强度降低可能无意中改善了预测能垒的准确性。尽管如此，几何校正AM1参数在预测能垒方面仍比原始AM1参数有约13%的改进，表明其能够描述参数化过程中未明确考虑的化学行为。

过拟合与参数优化

为了检验过拟合的可能性，研究将PM7的CHNO分子训练集随机划分为训练子集和验证子集，并在此基础上进行了额外的MNDO和AM1重新参数化。

对于MNDO，几何校正参数化在训练子集上经过177个优化周期后成功识别出局部最小值。评估结果显示，几何校正算法获得的参数在验证子集上的几何预测精度明显优于PARAM获得的参数，同时保持了相似的非几何误差贡献，表明其过拟合风险更低。

对于AM1，在训练子集上经过559个优化周期后获得了一个参数曲面上的低洼点（非严格最小值）。结果同样表明，几何校正算法获得的参数在验证子集上的几何误差贡献显著优于PARAM参数，且非几何误差贡献也更低，进一步证实了几何校正AM1重新参数化的结果不太可能是过拟合造成的。

结论与展望

本研究报道了一种能够包含几何参考数据的几何校正算法的扩展，表明直接几何变量相对于半经验参数的导数可以利用现有计算化学程序中已实现的例程轻松计算。将完整的几何校正参数优化程序应用于MNDO和AM1的结果表明，使用直接几何变量可能是实现对现有优化程序进一步改进的关键，并证实了PARAM甚至无法在其自身误差函数上识别最小值的断言。

几何校正重新参数化观察到的显著精度提升表明，仅通过更好的参数优化即可显著提高NDDO后裔半经验方法的精度。一个保留PMx家族公认的SCF收敛