镧系元素键参数拓扑指数的算法改进

杨寒彪，杨琳

（1.五邑大学 数学与计算科学学院，广东 江门 529020；2.江门职业技术学院 化学材料系，广东 江门 529000）

定量结构活性关系（quantitative structure activity relationship，QSAR）研究旨在通过合理的数理统计方法建立起一系列化合物的生理活性或某种性质（如药物的毒性、药效学性质、药物代谢动力学参数与生物利用度等）与其理化性质参数或者结构参数（包括二维分子结构参数、三维分子结构参数等）之间的定量关系. 其中拓扑指数能用于描述化学结构、预测化合物的物理化学性质，是现代计算化学、结构化学与量子化学相互结合交叉的产物[1]. 稀土元素及其化合物在现代科技各领域都有重要应用，有关镧系元素键参数拓扑指数与其物化性质的多元线性回归，已经取得了良好的线性相关性[2-6]. 但是，相关研究以镧系元素的离子半径作为拓扑指数参数，实际上，镧系元素的离子半径收缩比较快速，且是非等差收缩的；而镧系元素的原子半径是缓慢的等差收缩，相邻原子半径之差约为1 pm. 因此，本文考虑将元素键参数拓扑中镧系元素的离子半径改为了镧系元素的金属半径，以期得到相关性更好、更稳定的线性回归方程.

1 镧系元素的键参数拓扑指数

键参数拓扑指数的一般形式为

其中，pi,qi为分子结构图中第i条边两端顶点原子的度数；V i,Yi为对应于两顶点原子的键参数；i是分子结构图中所有边求和；t代表顶点原子键参数的类型.

把式（1）中的exp(V-Y) 项泰勒展开，取前4项得：

进一步简化式（2），对于某元素的原子，相当于分子图中i条边退化为一个点，因此可以用原子参数R替代(Vi-Yi)，原子参数S替代pi,qi，得到元素键参数拓扑指数HE：

在式（3）中，R为镧系元素的金属半径，S为镧系元素原子的价电子数，随着镧系元素原子序数从小至大，R严格变小，S严格变大：因此，式（3）并不是连续的满足拓扑指数的拓扑不变量. 但是镧系元素的数据是离散的，在(1 +R) 中，随着元素序数从小到大，(1 +R) 下降跨度最大为97.24%，上升跨度最小为103%. 即使在最极端的情况下，式（3）也整体上升了100.16%. 也就是说，在镧系元素的离散数据下，随着镧系元素原子序数从小到大，式（3）拓扑指数是严格的减函数，是拓扑不变量[7]. 该结论剔除了Eu 和Yb 两个离群元素的干扰.

另外，式（3）中仍包含有点价平方根的倒数，与应用广泛的拓扑指数分子连接性指数的主要含义一致.

2 本文的算法改进

2.1 本文算法

为了研究元素的物理化学性质，我们必须具体定义式（3）中R,S的参数内容. 文献[4]把R取为镧系元素 Ln3+离子半径，S取镧系元素原子的价电子数N. 例如，对于镧系元素， Ln3+离子半径为1.06 A，价电子数为3，代入式（3）计算得到其拓扑指数=0.28 .

但是本文针对元素键参数拓扑指数进行了算法的改进. 将R取为镧系元素金属半径. 例如，对于镧系 Ln3+，其离子半径为187.7 pm，价电子数为3，同样代入式（3）计算得到新的拓扑指数=798 .

特别需要指出的是，将R从镧系元素的离子半径改为金属半径，需要做一个剔除离群点的操作.镧系元素的离子半径绝大部分均呈现随着原子序数上升而均匀下降的趋势，但是Eu、Yb 两个元素大幅度背离了这一规律. 如果不剔除这两个元素的干扰，将会导致计算中相关系数会低于0.9. 因此本文的计算过程中首先剔除了Eu 和Yb 两个元素的金属半径数据，然后通过MATLAB中的fillmissing函数将Eu和Yb两个元素的金属半径数据进行“补全”.

水化性质是镧系元素极为重要的物理化学性质，它主要由镧系元素金属离子的水化能(H) 、水解常数(PK1) 、氢氧化物溶解度(PKm) 、与 EDTA 的络合稳定常数(lgK) 等参数衡量. 得到镧系元素的新的元素键参数拓扑指数H1后，我们可以得到其水化性质的线性回归方程：

其中，P为镧系元素的某种水化性质；a、b为线性回归方程的参数，镧系元素各水化性质线性回归方程的系数见表1.

表1 镧系元素各水化性质线性回归方程的系数

元素键参数拓扑指数算法改进后，我们得到了如表2所示的各项水化参数. 利用这些数据，我们得到了比文献[6-7]线性回归方程更高的相关系数. 也就是说以镧系金属的原子参数为标准，可以更准确地利用回归方程来预测镧系元素的物理化学性质.

表2 算法改进后镧系元素的各项水化参数与回归相关系数

2.2 算法改进的数学意义

镧系元素有镧系收缩的物理化学性质，其中镧系元素原子半径收缩的较为缓慢，相邻原子半径之差仅为1 pm 左右并基本是等差收缩的，但是镧系元素的离子半径收缩却比较快速，而且是非等差收缩的. 因此从线性相关性的角度，用镧系元素的金属半径作为拓扑指数的参数进行计算可以得到更高的相关系数. 如表 3 所示，除Eu 和Yb 两个离群元素外，镧系元素的金属半径和离子半径均随着原子序数的增大而单调下降. 但金属半径基本是稳定的等差下降，而离子半径是逐步放缓下降的.同时我们发现，lgK的数值虽然随着镧系元素原子序数增大而上升，但其也是等差增长的. 也就是说，lgK与镧系元素的金属半径更相关.

表3 镧系元素离子半径、金属半径之差和lg K 之差

这里谈一下为何前人更多使用镧系元素的离子半径作为拓扑指数参数计算的原因. 这可能是因为镧系元素的离子半径并不像金属半径那样有Eu 和Yb 两个离群元素，其几乎完全呈现了单调下降的特点；并且水化性质更多地跟离子性质有关而非元素金属半径有关. 我们通过 matlab 中的fillmissing 函数将Eu 和Yb 两个元素的金属半径数据进行“补全”，完全避免了这两个离群元素的干扰，得到了线性系数更好的回归方程.

3 新旧算法下，线性回归方程其他统计学指标的比较

本文算法下，计算值和文献值的平均偏差为2.85，而文献[6-7]旧算法的计算值和文献值的平均偏差为4.7. 故，本文算法的偏差明显优于文献[6-7]的算法.

但在回归方程各系数的显著性检验 t 检验中，本文算法与文献[6-7]算法就各有千秋. 从表 4 可以看到：1）对于水化能H和 EDTA 络稳定常数lgK而言，本文算法系数的 t 检验更强，均方差更小. 说明回归方程可以相对更加显著地反映其化学性质，而均方差更小说明本算法的稳定性更高. 2）对于水解常数PK1和氢氧化物溶解度PKm而言，本文算法系数的t 检验更弱，均方差更大. 这说明新算法的回归方程不能更显著地反映其化学性质，而且稳定性变更差了. 这也解释了镧系元素的离子半径对PK1、PKm有更显著和更稳定的影响.

表4 新旧算法中各水化性质的t检验和均方差

4 结论

本文得到了键参数拓扑指数算法改进下的镧系元素水化性质的回归方程. 新算法下回归方程的相关系数仅略有提高，但是线性回归方程平均偏差比之前降低，说明方程稳定性更强了，新算法可以为水化物性质的理论研究提供一定的铺垫，并可进一步推广到诸如电离性等其他化学物理性质的研究之中.