1 引言
在众多化学战剂和模拟剂中, 梭曼(soman, GD)、二异丙基氟磷酸酯(diisopropyl fluorophosphate, DIFP)、甲基膦酸二甲酯(dimethyl methyl phosphonate, DMMP)以及沙林(sarin, GB)等因其高毒性和广泛存在性, 成为化学战剂防护领域的重点关注对象[5]. 梭曼[6]作为一种毒性极强的神经性毒剂, 能迅速干扰人体神经系统正常功能, 极低剂量便可致人死亡. 在两伊战争期间, 就曾有使用梭曼等神经性毒剂造成大量人员伤亡的情况, 给当地人民带来了难以磨灭的伤痛. 沙林和DIFP[7]同样是毒性剧烈的神经性毒剂, 其对人体的危害极大, 曾在一些恐怖袭击事件中被使用, 造成了严重的人员伤亡和社会恐慌[8]. DMMP[9]作为常见的神经性毒剂模拟物, 虽毒性相对弱于典型化学战剂, 但在工业生产、军事科研场景广泛应用. 其若发生泄漏, 或遭不当处置, 仍会对周边环境产生污染, 接触人体时, 也会经呼吸道、皮肤等途径, 给健康带来不良影响. 因此, 开发安全高效的技术来捕获空气中微量化学战剂, 特别是在军事和消防领域, 已成为当前研究的重点[10].
金属-有机框架(Metal-organic Frameworks, MOFs)材料的诞生, 为攻克化学战剂捕获这一难题带来了新希望. MOFs是一类由金属离子或金属簇与有机配体通过自组装形成的具有周期性网络结构的多孔材料[14-18], 其独特的结构赋予了材料诸多令人瞩目的优异特性. 一方面, MOFs具有超高的比表面积[19], 可提供丰富的活性位点, 从而显著增强对各类分子的吸附能力. 另一方面, MOFs的孔道结构具有高度可调节性, 研究人员能够依据目标分子的尺寸、形状以及化学性质, 精心设计有机配体和金属节点, 精准调控孔道结构, 进而实现对特定化学战剂的特异性吸附[20]. 众多研究已充分证实, MOFs在膜分离[21]、药物传输[22]、催化[23]以及气体分离[24]等多个领域均展现出比传统材料更卓越的性能. 此外, MOFs在化学战剂吸附和降解领域同样展现出巨大潜力[25-26]. MOFs对神经性毒剂的降解主要是通过水解作用进行的, 通过水解作用使毒剂中的P—S、P—O或P—F键断裂, 将毒剂水解为无毒或低毒化合物[27]. Liu等[28]通过溶液吹纺聚苯乙烯/聚氨酯并结合水热原位生长法负载MOF-808, 制备出PS/PU/TTMA/MOF- 808纤维膜, 研究其对2-氯乙基乙基硫醚(2-chloroethyl ethyl sulfide, CEES)的吸附性能. 结果表明, 该纤维膜因具有蓬松结构、MOF-808的丰富微孔和大比表面积, 对油性CEES有良好吸附能力, 油滴接触表面后立即消失, 对柴油的最大吸附容量达41.7 mL•g-1, 且吸附为后续催化降解提供了充足反应物. Tao等[29]通过微波辅助溶剂热法在PP布料上垂直生长二维Zr-BTB金属有机框架纳米片, 制得的复合布料降解神经毒剂模拟物对硝基苯基磷酸二甲酯(Dimethyl p-nitrophenylphosphate, DMNP)的半衰期低至0.1~0.77 min, 负载率达50%, 单次过滤解毒率超88%, 凸显了其在化学防护和环境净化领域的应用潜力. Sun等[30]通过配体交换法, 在PCN-601框架中引入MnMo6-lg多金属氧酸盐簇, 构建了含原子级协同Ni-Mo双金属对位点的双簇MOF (MnMo6-lg-PCN- 601), 其对化学战剂模拟物CEES的降解半衰期小于1 min, 转化频率(TOF)达140 min-1, 且循环稳定性优异. 张等[31]采用TA-APTES涂层改性聚丙烯非织造布并通过原位生长法制备了高MOFs负载量的纤维基复合材料, 其对化学战剂模拟物DMNP的降解半衰期仅4.8 min, 30 min内即可实现完全转化. Wei等[32]研究了一种基于UiO-66-COOH/TiO2的一维光子晶体气体传感器, 通过氢键策略实现对神经毒剂模拟物DMMP的快速灵敏检测, DMMP在空气中的体积分数低至3.0×10⁻7时, 仍能被有效检测, 具有良好稳定性与选择性. Yost等[33]采用振动光谱和核磁共振等技术研究了不同粒径、缺陷类型和密度的Zr-UiO-66与Hf-UiO-66对梭曼模拟物亚磷酸二甲酯(DMP)的吸附及水解性能. 结果表明, DMP吸附受表面限制, 小粒径样品吸附量更高, 含大量缺失簇缺陷的Hf-UiO-66吸附与水解效果更优. Son等[34]研究10种Zr基 MOFs, 探讨其孔隙结构(比表面积和孔体积)、拓扑结构和二级结构单元(SBU)连接性对沙林(GB)和硫芥气模拟物(CEES)吸附及反应活性的影响, 发现高表面积、强氢键作用提升吸附量, 而且UiO-66等活性位点多的材料对GB反应性高, 开放金属位点在干燥环境有利, 潮湿时受水影响, 为MOFs靶向设计提供依据.
MOFs材料在化学战剂的防治领域展现出巨大潜力, 但MOFs材料家族庞大, 其结构与性能之间的关系错综复杂. 据统计, 目前已合成MOFs材料数量多达数万种, 若通过传统实验逐一筛选适用于捕获特定化学战剂的MOFs, 将耗费大量的时间、人力和物力资源, 效率极为低下[35-36]. 高通量计算技术的兴起, 为MOFs材料筛选工作注入了强大动力, 借助量子化学计算、分子动力学模拟等理论方法, 高通量计算能够在短时间内对海量MOFs结构进行系统计算与性能预测. 例如, Phadatare等[37]用巨正则蒙特卡罗(Grand Canonical Monte Carlo, GCMC)方法模拟了MOFs对DIFP和DMMP的吸附行为, 表明UiO-66-NH2对DMMP具有优异的吸附能力, 且DMMP和UiO-66-NH2的氨基之间的相互作用力为物理吸附类型. Yang等[38]利用GCMC模拟计算了606个疏水MOFs对芥子气及其模拟剂的吸附性能, 结果显示具有开放金属位点等特性的MOFs对化学战剂有较高的吸附性能. 不过, 单纯高通量计算在处理海量数据深度挖掘上有局限, 而MOFid作为一种系统化的MOF结构识别与索引工具, 通过将MOFs解构为金属节点、有机配体和拓扑结构等标准化单元, 为高通量计算中海量MOFs的结构解析、数据整合与关联分析提供了统一基础[39]. MOFid提供的结构化数据(如基于SMILES的构建块信息、拓扑代码等), 能为机器学习模型构建提供统一输入, 高效挖掘结构-性能规律. 例如, Cao等[40]开发的MOFormer模型以MOFid为输入, 实现了结构无关的MOFs性质预测, 无需依赖3D结构即可加速筛选, 其通过自监督学习在超40万MOF数据上预训练, 提升了预测精度, 且在训练数据有限时, 在量子化学性质预测上比基于结构的模型更具数据效率; Glasby等[41]则借助MOFid从1172个无序MOFs中识别出Cu paddlewheel MOFs, 为研究其活化过程中的结构坍塌现象提供了精准筛选基础, 体现了MOFid在特定MOF类型高通量筛选与性能关联研究中的价值. 通过这些方法, 能快速学习材料结构-性能规律, 精准预测众多材料特性、筛选潜力材料, 为高效吸附剂设计指明方向[42-43].
在本研究中, 采用GCMC模拟和机器学习方法探究了超过1.5万个计算就绪的实验金属有机框架(Computation-Ready Experimental MOFs, CoRE-MOFs)对4种毒剂或模拟剂(梭曼、沙林、DMMP和DIFP)的吸附性能. 首先, 分析了七个描述符即最大空腔直径(largest cavity diameter, LCD)、密度(density, ρ)、孔隙限制直径(pore limiting diameter, PLD)、孔隙率(φ)、体积比表面积(volumetric surface area, VSA)、亨利系数(Henry coefficient, K)和吸附热(heat of adsorption, Q₀st)对4种化学战剂/模拟剂的吸附量(N)、选择性(S)及两者权衡值(TSN)的构效关系; 然后, 利用机器学习训练和预测材料的吸附性能指标, 并结合SHAP可解释器(SHapley Additive exPlanation)来获取各个描述符对CoRE-MOFs吸附性能的贡献值; 最后, 基于权衡值分别筛选出了对每种毒剂具有优异吸附性能的MOFs, 用MOFid对其结构解析并进行共性分析, 提炼出关键的构筑要素, 包括理想的金属位点、拓扑构型及有机配体, 并基于上述要素利用Tobacco软件模拟组装新型MOFs, 从而实现了对新型吸附材料的理论设计与结构优化.
2 结果与讨论
2.1 单变量分析
为了深入理解MOFs吸附梭曼、沙林、DMMP和DIFP的结构-性能关系, 分析了MOFs的五种结构描述符(φ、VSA、LCD、PLD和ρ)和两种能量描述符(K和Q0st)与化学战剂吸附性能之间的关系. 本研究绘制了 7 种描述符与毒剂吸附量的相关性图谱, 下面将针对梭曼、沙林这两种典型毒剂展开深入讨论(其余两种毒剂表现出相似规律).
吸附量N是评估MOFs吸附性能的关键指标. 对于沙林和梭曼, 其吸附量与MOFs描述符(以φ、LCD、ρ、K为代表)的关联可从图中探究. 如图1a、1e显示, 当 φ<0.71时, Nsarin和Nsoman随φ增大而单调递增; 当φ>0.78时, 吸附量则出现衰减趋势, 最大值出现在φ在0.71~0.78范围之内. 当孔隙率逐渐增加时, 材料的比表面积相应增大, 为梭曼、沙林分子提供了更多的吸附活性位点, 气体分子能够更充分地与材料表面接触并发生吸附作用, 所以吸附量随之上升. 但当达到吸附饱和后, 无效孔空间的占比逐渐提升, 新吸附的CWAs分子与MOFs间相互作用力被削弱, 导致φ继续增大无法有效提升吸附量. 此趋势表明φ的最优区间较窄, 为筛选针对沙林、梭曼的最优MOFs提供了思路, 与Qiao等[44]的研究规律相似. 如图1b、1f显示, 沙林与梭曼吸附量随LCD变化呈现“先增后降”的特征. 当LCD<0.3 nm, 无论是沙林还是梭曼, 其吸附量都近乎为零, 这是因为过小的LCD限制了沙林和梭曼分子进入MOF孔道, 导致吸附无法有效发生[45]. 随着LCD的增大, 越来越多的气体分子可以进入到MOF内部, 吸附量随之增加. 然而, 当LCD继续增大, 吸附量呈衰减趋势, 这是由于过大的LCD使MOF孔道过度开放, N2、O2等空气组分更易进入孔道, 与沙林、梭曼分子形成竞争吸附, 分子-孔壁相互作用能降低, 导致吸附量下降. Nsarin和Nsoman与MOFs的ρ之间的结构性能关系如图1c、图1g所示. Nsarin和Nsoman随着ρ的增加先增加后减少, 从ρ=70 kg•m-3开始急剧上升, ρ=421~635 kg•m-3时Nsarin和Nsoman出现最大值, 达到ρ=635 kg•m-3时Nsarin和Nsoman逐渐降低. 从图中可以得到, 密度ρ在合适的区间内有着更合适的孔隙空间和VSA, 所以更有利于对沙林和梭曼的吸附. Nsarin和Nsoman与MOFs的K之间的结构性能关系如图1d、1h所示. 当K处于较低区间时, MOF框架与沙林、梭曼分子之间的相互作用较弱, 导致吸附量较低. 随着K逐渐增大, 吸附量开始迅速上升. 这表明随着K的增加, MOF框架与沙林、梭曼分子间的相互作用增强, 分子间的吸引力增大, 使得沙林和索曼分子更容易被吸附到MOF孔道中, 从而导致吸附量快速增加. 然而, 当K增大到一定程度后, 吸附量却开始逐渐下降. 这可能是由于过高的K改变了MOF与沙林、梭曼分子之间的相互作用方式, 或者导致其他竞争吸附过程加剧, 使得原本已吸附的分子更容易脱附, 进而导致吸附量降低.
图1 梭曼和沙林吸附行为的结构-性能关系Figure 1 Structure-performance relationships for the adsorption of sarin and soman (a) Nsarin-φ; (b) Nsarin-LCD; (c) Nsarin-ρ; (d) Nsarin-K; (e) Nsoman-φ; (f) Nsoman-LCD; (g) Nsoman-ρ; (h) Nsoman-K. The color represents the value of TSN. Each figure contains the data of MOFs |
图2 梭曼和沙林气体选择性的结构-性能关系Figure 2 Structure-performance relationships for the selectivity of sarin and soman (a) Ssarin-φ; (b) Ssarin-LCD; (c) Ssarin-ρ; (d) Ssarin-K; (e) Ssoman-φ; (f) Ssoman-LCD; (g) Ssoman-ρ; (h) Ssoman-K. The color represents the value of TSN. Each figure contains the data of MOFs |
而对于LCD(图2b、2f), Ssoman/N2+O2和Ssarin/N2+O2随着它们的增大呈现先增加后降低的趋势, 并在 0.317~0.761 nm范围内达到最大值, 此时吸附性能亦处于优良状态, 这是由于该尺寸范围恰好适配梭曼、沙林分子的动力学直径, 形成有效吸附空间, 这与之前的研究结论相似. 图2c、2g显示, Ssoman/N2+O2和Ssarin/N2+O2随着ρ的增加同样呈现先上升后降低的趋势, 与吸附量的关系变化图一致, 原理相同, 不再赘述. 值得注意的是, 图2d、2h揭示了K对选择性的调控效应: 当K过小时, 吸附驱动力不足导致目标分子吸附量极低, 选择性趋近于零; 而当K超过临界值, 过强的吸附能会加剧孔道内分子间的排斥作用, 导致吸附热力学平衡向脱附方向移动, 最终使选择性降低.
表1 MOFs吸附四种有毒物质时不同描述符的最佳范围Table 1 Optimal ranges of different descriptors for MOFs adsorbing four toxic agents |
| 气体 | 描述符 | 吸附量N/(mol•kg-1) | 选择性S |
|---|---|---|---|
| 梭曼、沙林 | φ | 0.71~0.78 | 0.22~0.90 |
| VSA/(m2•cm-3) | 2446~2647 | 328~1160 | |
| LCD/nm | 0.664~0.918 | 0.496~0.755 | |
| PLD/nm | 0.542~0.868 | 0.262~0.589 | |
| ρ/(kg•m-3) | 421~635 | 915~1953 | |
| K/(mol•kg-1•Pa-1) | 0.18×10-2~1.5×107 | 1.73×10-3~1.13×103 | |
| Q0st/(kJ•mol-1) | 19.37~86.02 | 24.52~68.22 | |
| DIFP | φ | 0.67~1.02 | 0.14~0.73 |
| VSA/(m2•cm-3) | 1065~3290 | 334.3~1099.7 | |
| LCD/nm | 0.665~1.336 | 0.504~1.164 | |
| PLD/nm | 0.493~1.362 | 0.219~1.164 | |
| ρ/(kg•m-3) | 392.7~690.9 | 1125.3~3260.5 | |
| K/(mol•kg-1•Pa-1) | 0.0001~5.33×106 | 0.12~5.8×1011 | |
| Q0st/(kJ•mol-1) | 24.19~85.41 | 45.20~105.73 | |
| DMMP | φ | 0.67~1.14 | 0.13~0.64 |
| VSA/(m2•cm-3) | 1065 ~3290 | 309.4~2383.7 | |
| LCD/nm | 0.664~0.930 | 0.433~0.787 | |
| PLD/nm | 0.493~1.092 | 0.242~0.659 | |
| ρ/(kg•m-3) | 392.72~808.97 | 728.76~3370.18 | |
| K/(mol•kg-1•Pa-1) | 0.62~9.41×106 | 0.003~2.28×109 | |
| Q0st/(kJ•mol-1) | 55.69~228.24 | 46.63~203.89 |
2.2 多变量分析
经过上述单变量分析, 已经初步得到各个描述符与吸附量及选择性之间的关系. 为进一步探究MOFs对四种毒剂(梭曼、沙林、DMMP和DIFP)共存时的吸附性能及挖掘高性能MOF的结构相似性, 采用机器学习算法进行了多变量分析. 另外, 考虑到毒剂分子自身的结构性能也会对MOFs的吸附性能起到一定的影响, 除了基础的MOFs结构及能量描述符, 同时引入了毒剂分子的物理化学参数(沸点Tb、临界温度Tc、临界压力Pc, 分子动力学直径Dia、偶极矩Dip、极化率Pol、熔点Tm、闪点Tf、分子密度ρgas、偏心因子ω)作为描述符, 采用权衡值TSN作为机器学习的目标预测值. 在此工作中, 采用了六种算法(决策树(DT)、随机森林(RF)、极端梯度提升(XGB)、反向传播神经网络(BPNN)、梯度提升回归(GBR)和轻量梯度提升机(LGBM))进行训练和预测. 在机器学习训练时, 将数据集分成了两部分, 70%的数据作为训练集和30%的数据作为测试集, 将每个机器学习模型训练五次, 以此确定模型的稳定性, 并使用得到的R2、MAE和RMSE作为评价标准. 相关的超参数设置详见支持信息表S9. 从图3中可以看到四种算法的预测R2值在0.846~0.923范围内, 六个机器学习算法按照R2值进行排序: XGB>LGBM>RF>GBR>BPNN>DT. 其中, XGB算法预测得到的R2值最高, 达0.923, 展现出最优预测能力. 结合图4数据可知, XGB算法训练集MAE仅为0.519、测试集为1.282, 训练集RMSE为1.104、测试集为2.219, 在6种算法中误差控制最优, 预测值与真实值偏差小且稳定性强. 相比之下, DT算法的MAE与RMSE均显著偏高, 训练集MAE达1.117、测试集1.604, RMSE训练集1.891、测试集1.750, 预测值与真实值偏差较大, 稳定性明显不足. 说明XGB算法在处理高维非线性数据时表现最优, 其预测精度显著优于DT模型. XGB算法表现优异的根源在于, 其依托梯度提升框架构建多棵决策树, 通过迭代优化残差精准捕捉复杂特征交互; 同时, 正则化机制有效抑制过拟合, 在保障高拟合能力的同时提升泛化性能. 此外, 并行化架构赋予其处理大规模数据集的高效性, 为实际工程应用奠定可行性基础. 反观DT模型, 因采用单变量分割策略, 难以捕捉特征间协同效应, 模型偏差与方差较高, 且树结构深度过深时易受训练数据噪声干扰, 在测试集呈现明显过拟合倾向; 加之对异常值敏感、适配高维数据能力不足, 进一步限制其预测可靠性, 与XGB算法形成鲜明性能差异.
图3 六种机器学习算法TSN预测性能对比Figure 3 Comparison of the simulated TSN and the predicted TSN obtained from six ML algorithms (a) DT; (b) RF; (c) XGB; (d) GBR; (e) BPNN; (f) LGBM. The dots represent the number of MOFs |
此外, 为了深入探究MOFs的各个描述符对MOFs吸附化学毒剂性能的影响机制, 采用了SHAP中的TreeExplainer模块[47]对MOF描述符与化学毒剂特性之间的关联进行解释和挖掘分析. 在SHAP模块中, 所有特征均被视为贡献因子, 通过计算SHAP值, 可量化特征对吸附性能的影响程度. 从图5左侧的蓝色柱状图看, 以“Percentage/%”量化各特征对模型输出的全局平均贡献, 计算全数据集样本SHAP绝对值的均值并转化为百分比, 直观呈现特征重要性排序. 其中特征K纵轴排序最靠前、柱状图占比最高, 是影响吸附性能的重要描述符, Tm、ρMOF、Tb等特征占比依次递减, 体现影响权重差异. 右侧蜂群图则展示了模型中每个特征在整个数据集上的贡献情况, 图中每一个点代表一个样本的特征 SHAP 值, 点的密集程度显示了特征值的集中趋势, 颜色表示特征值的高低, 红色表示高值, 蓝色表示低值; 横轴是SHAP值, 反映特征对预测结果的影响大小和方向(正值为促进作用, 负值为抑制作用), 纵轴为特征名称, 按重要性排序, 越靠上的特征对模型输出的影响越大. 以K为例, 其SHAP值分布跨度广、波动显著, 当K处于较低区间(显示为蓝色点)时, SHAP值多为负, 表明在较低的吸附驱动力下, 吸附过程需要依靠其他描述符的协同作用来维持. 而在K的有效高值区间(红色点向正SHAP值集中), K的提升增强了吸附驱动力, 促进了毒剂分子与MOFs吸附位点的结合, 从而导致吸附量上升. 孔隙率φ特征位于重要性排序第二, 仅次于K. 其蜂群图分布清晰呈现孔隙率对吸附性能的作用机制: 低孔隙率(蓝色点)的样本密集地分布在SHAP值–20至0的负值区域, 这意味着φ过低会导致吸附量减小——此时MOFs内部孔隙空间狭窄, 可供毒剂分子接触的有效吸附位点数量严重不足, 显著限制了吸附性能的提升. 随着孔隙率增大(特征值由蓝转红), 更多吸附位点暴露出来, 样本的SHAP值开始呈现出正向贡献, 吸附量也随之增加. 然而, 结合图1的结果, 吸附量在达到某个高水平后反而会下降. 图5右侧关于ρMOF的蜂群图很好地解释了这一现象: 低密度样本(蓝点)主要聚集在SHAP值0~40的正值区, 而高密度样本(红点)则集中在SHAP值–5至0的负值区(这一分布模式与φ形成了有趣的镜像). 这说明当孔隙率过大时, 必须同时调整ρMOF. 如果不加以控制, 过高的孔隙率可能会导致MOFs内部结构过于松散, 分子间作用力减弱, 最终引发已吸附的毒剂分子发生脱附[48]. 此外, 毒剂的物性参数对吸附性能也存在关键影响. 以沸点(Tb)为例, 从SHAP蜂群图中Tb的特征分布来看, 当毒剂具有较高沸点(红色点)时, 其SHAP值多为正值, 这表明高沸点毒剂分子与MOF吸附位点之间的范德华力等相互作用更强, 结合更牢固, 从而促进吸附过程[49]; 而低沸点(蓝色点)毒剂分子的热运动更剧烈, 易从吸附位点脱附, SHAP值多呈现负值, 抑制吸附效果. 再看偶极矩(Dip), 若毒剂偶极矩较高(红色点), 当MOFs自身具有适宜的极性(可通过Pol等描述符体现)时, 偶极-偶极相互作用会增强毒剂与MOFs的结合, 此时Dip的SHAP值多为正, 推动吸附量上升; 反之, 低偶极矩毒剂(蓝色点)因缺乏这类强相互作用, 对吸附的促进作用较弱, 甚至可能因与MOFs极性不匹配而产生抑制(SHAP值为负).
图5 全局特征重要性(左, 条形图)和局部解释摘要(右, 蜂群图), 其中点代表本研究中的单个数据点注意: x轴上的点显示每个数据点如何影响模型对每个特征的预测, 当多个点位于同一x位置时, 它们显示密度(颜色表示特征值, 其中红点表示高, 蓝点表示低; 正SHAP值表示对预测值的正特征影响, 负SHAP值表示对预测值的负特征影响) Figure 5 Global feature importance (left, bar charts) and local explanation summary (right, beeswarm plots), in which a dot represents an individual data point in this study Note: The dot on the x-axis shows how each data point impacts the model’s prediction for each feature, and when multiple dots are at the same x position, they show the density (Colors represent eigenvalues, where red dots represent high and blue dots represent low; positive SHAP values indicate a positive feature effect on the predicted value, and negative SHAP values indicate a negative feature effect on the predicted value) |
SHAP分析明确揭示MOFs吸附性能受多种描述符协同调控的机制, 结合SHAP蜂群图以及单变量分析图, 适当的K值, 较大的VSA和φ, 较小的ρ和PLD, 更加有利于MOFs对CWAs及模拟剂的吸附.
2.3 优异MOFs的共性
为精准挖掘MOFs对梭曼、沙林、DMMP及DIFP四种毒剂的吸附性能规律, 依据TSN递减排序确定了对四种毒剂具有最佳综合吸附性能的TOP1% MOFs. 为了进一步深入分析这些MOFs对吸附性能的影响机制, 借助Snurr团队[39]建立的MOFid软件对上述TOP1% MOFs进行结构解构, 通过其标准化的拆分算法精准提取出金属中心、有机配体及拓扑结构的关键信息.
首先, 研究了TOP1% MOFs所包含的金属中心类别及数量, 结果如图6所示. 图中, ratio代表TOP1% MOFs中含有该金属的MOFs数量与整个MOFs库中含有该金属的MOFs数量之间的比值, 颜色映射(即number)代表在TOP1% MOFs中包含该金属的MOFs数量. 蓝色划线将区域分成四个部分, 分别对应过渡金属、镧系金属、碱及碱土金属和贫金属. 从金属类别共性角度看, 梭曼、沙林、DMMP、DIFP 四个毒剂体系的分析结果高度一致: 过渡金属(如Cr、Nb)的“ratio” 与“number”均显著高于其他类别, 这个是因为过渡金属大多具有开放的金属位点, 经过前人的研究证明, 开放金属位点有利于MOFs的吸附性能[45,50]. 于是我们对TOP1%中三种具有开放金属位点的MOFs进行吸附密度分析(支持信息图S16、S17和S18). 这四种毒剂的最强吸附位点均紧密围绕在开放金属位点周围. 这一直观证据进一步证实了开放金属位点在吸附过程中的核心作用. 为量化这种相互作用的强度差异, 我们统计了四种毒剂在三种MOFs中的K(详见支持信息表S11). 如支持信息图S16分析表明: K值越高(如梭曼), 其吸附密度分布越紧密地包裹在金属位点周围; 而K值越低(如沙林), 其吸附位置的整体分布则相对偏离开放金属位点. 因此, K值的大小直接量化了毒剂分子与金属位点作用的强弱(支持信息图S17和S18也存在相似特征). 这一发现从微观机理层面印证了后文对宏观金属统计结果的解读, 不仅揭示了吸附强弱的本质, 也为过渡金属的主导地位提供了直接证据. 以梭曼体系(对应图6a)为例, Cr元素的“ratio”达4.2%, 对应“number”色条接近9, 表明含Cr-MOFs在高性能体系中既具较高富集度又有一定数量规模; Nb元素虽“number”绝对值低于5, 但其“ratio”高达10%, 结合MOFs库中含Nb材料基数小的背景, “低数量高比例”的特征凸显Nb-MOFs在高性能筛选中的极强富集性, 是提升梭曼吸附性能的关键金属之一. 沙林体系(对应图6b)中, Nb元素“ratio”为10%, Zn元素“number”对应深绿色区间(22~25), 进一步验证过渡金属在不同毒剂体系中对吸附性能的主导作用, 说明过渡金属的优势并非局限于单一毒剂, 而是具备跨体系的普适性. 综合来看, Cr、In、Al、Li、Nb等金属因所在MOFs多含开放金属位点, 在无水干扰时可借助位点与毒剂的强相互作用, 成为提升吸附性能的关键, 这也明确了含开放金属位点的金属中心在毒剂吸附中的核心地位, 为新吸附材料设计指明方向.
除金属中心外, 研究还对TOP1% MOFs开展标准结构特征搜索, 表2呈现了相关结果. 从有机配体(linker)看, 四个体系均有MOFkey[39]为QMKYBPDZA- NOJGF的配体; 梭曼、DIFP、DMMP体系共享MOFkey为GZJYHAYPVZNUJX的配体, 出现次数分别为16、8、15次, 说明特定配体在不同毒剂体系的高性能MOFs中均存在. 从配体结构特性看, 这些共性配体多含刚性芳香环与极性官能团(如羧基、羟基), 刚性芳香环可增强MOFs框架稳定性, 极性官能团则能通过氢键、偶极-偶极作用捕获毒剂分子, 其结构特性与吸附性能的关联值得深入探究. 以QMKYBPDZANOJGF配体为例, 其芳香环上的羧基可与毒剂分子(如沙林含磷氧键、梭曼含酯键)形成氢键, 强化分子-框架相互作用, 提升吸附选择性; 刚性骨架保障MOFs孔道结构稳定, 为毒剂分子提供持续可及的吸附空间, 这解释了其跨体系复用的内在原因.
表2 吸附四种毒剂对MOF吸附性能影响较大的结构Table 2 Critical structural properties governing MOF adsorption of four toxic agents |
| 毒剂 | Linker | 相应 次数 | Topology | 次数 | ||
|---|---|---|---|---|---|---|
| sarin |  QMKYBPDZANOJGF C9H6O6 |  OYFRNYNHAZOYNF C8H6O6 |  LSHYAWDVBASRBW PO4 | 8/5/4 |  srs | 13 |
| DIFP |  GZJYHAYPVZNUJX C8H6O4 |  SATWKVZGMWCXOJ C27H18O6 |  QMKYBPDZANOJGF C9H6O6 | 8/5/4 |  sql | 14 |
| soman |  GZJYHAYPVZNUJX C8H6O4 |  QMKYBPDZANOJGF C9H6O6 |  KKEYFWRCBNTPAC C8H6O4 | 16/5/3 |  kgm | 12 |
| DMMP |  GZJYHAYPVZNUJX C8H6O4 |  QMKYBPDZANOJGF C9H6O6 |  KKEYFWRCBNTPAC C8H6O4 | 15/8/3 |  sql | 15 |
拓扑结构(Topology)分析结果(表2)显示, 不同毒剂体系存在适配性优势拓扑: sql在DIFP、DMMP体系中最常见(分别有14种、15种MOFs具备该结构); 梭曼体系中cha、kgm、nbo拓扑出现频率高. 尤为关键的是, sql、kgm拓扑及QMKYBPDZANOJGF、KKEYFWRC- BNTPAC、GZJYHAYPVZNUJX配体, 在四个体系中均存在, 进一步印证这些结构要素与MOFs高效吸附性能的内在关联, 为新型MOFs设计明确了需重点关注的配体和拓扑结构类型.
本研究基于上述对金属中心、有机配体、拓扑结构的系统分析, 明确了提升MOFs毒剂吸附性能的核心要素: 选用含 Cr、In、Al、Li、Nb等具开放金属位点的过渡金属, 搭配QMKYBPDZANOJGF、KKEYFWRC- BNTPAC、GZJYHAYPVZNUJX等共性有机配体, 以及sql、kgm等优势拓扑结构. 以此为依据, 借助Tobacco软件组装新型MOF(见支持信息表S12), 通过模拟实践验证结构共性的应用价值, 为MOFs在毒剂吸附领域的研发, 从理论研究到实际材料设计搭建桥梁.
3 结论
在这项工作中, 通过系统筛选CoRE-MOFs数据库, 结合单变量分析与AI大数据算法, 揭示了MOFs吸附化学战剂和模拟剂之间的构效关系及关键影响因素. 单变量分析表明, φ和VSA的协同优化对提升吸附性能至关重要, 且当φ处于0.7~0.8和VSA为2000~2500 m2•cm-3时, MOFs对沙林、梭曼等毒剂的吸附量达到峰值. 最优XGB算法辅以可解释性工具SHAP分析, 进一步揭示了描述符对性能的复杂影响机制, 发现K是影响吸附性能的最关键描述符. 其正向贡献表明, 在有效的高K值区间内, 提升吸附驱动力能够显著促进毒剂分子的吸附. XGB算法与MOFid结合, 综合分析高吸附性能MOFs的结构相似性, 表现为含开放金属位点的过渡金属(如Cr、Nb、Al、In、Li)节点, 以及特定优势配体和拓扑结构(如sql、kgm). 研究发现, 开放金属位点与毒剂分子间的强配位作用是提升吸附性能的关键, 而某些过渡金属(如Cr、Nb)的电子特性可进一步增强对特定毒剂的选择性. 基于上述大数据挖掘与AI模型解析的发现, 我们设计了26种新型MOFs材料.
本研究为从现有MOFs库中高效筛选最优候选材料以及理性设计新型高性能MOF吸附剂提供了重要指导, 有助于推动高通量计算筛选与AI大数据技术在MOFs捕获化学战剂研究中的深度融合应用, 并为开发先进的化学战剂防护装备(如防护服)提供了更有效的新路径.
4 模拟与计算方法
4.1 模拟方法
采用GCMC方法, 对CoRE-MOFs捕获空气中微量有毒物质的性能展开模拟. 混合气体初始体积比设定为V(N2)∶V(O2)∶V(有毒气体)=0.789∶0.21∶0.001, 以模拟有毒气体体积分数为0.1%的大气环境. 模拟在温度298 K、压强101.325 kPa的条件下进行, 且假定MOFs结构为刚性. 模拟过程中, 维持化学势、体积和温度恒定. 通过Lorentz-Berthelot规则计算MOFs与吸附质分子间的相互作用, 利用Lennard-Jones (LJ)势描述吸附质与框架原子间的相互作用:
${{u}_{LJ+elec}}(r)\sum{4}{{\varepsilon }_{ij}}\left[ {{\left( \frac{{{\sigma }_{ij}}}{{{r}_{ij}}} \right)}^{12}}{{\left( \frac{{{\sigma }_{ij}}}{{{r}_{ij}}} \right)}^{6}} \right]\frac{{{q}_{i}}{{q}_{j}}}{4\pi {{\varepsilon }_{0}}{{r}_{ij}}}$
其中uLJ+elec(r)表示原子间势能, ɛij表示i原子和j原子间电势强度, σij表示原子间电势为零的距离, rij表示相互作用的原子之间的距离, q为原子电荷数, ɛ0=8.8542×10-12 C2•N-1, 表示真空介电常数. 长程校正的球形截止半径设为1.2 nm. 先通过前50000个周期平衡模拟系统, 随后50000个周期用于积累计算平均值. 每个GCMC模拟均独立开展, 整个模拟依托RASPA软件包[51]运行. 在每个循环中, 执行移位、旋转、交换、再生以及插入和删除等操作, 以促进系统平衡. 框架与气体分子之间以及气体分子间的静电相互作用, 运用Ewald求和方法计算. 通过MEPO-Qeq电荷算法[52-53]计算CoRE-MOFs的电荷. 所有CoRE-MOFs的LJ势参数源自通用力场(Universal force field, UFF)[54], 目标气体、N2和O2的力场参数则来自TraPPE力场[55-57]. 另外, 为了验证所选力场(UFF+TraPPE)及模拟参数对气体分子在MOF中吸附行为的预测精度, 我们进行了基准测试. 选取了四种具有不同孔结构和化学环境的典型MOF材料, 模拟了它们对O2的吸附等温线, 并将模拟结果与文献中的实验数据[58-60]进行了对比. 结果表明, 我们的模拟结果与实验数据吻合度良好(支持信息图S1), 验证了模型的可靠性. 为了评估真实大气环境(如湿度)对吸附性能的潜在影响, 我们考察了水蒸气(H2O)和二氧化碳(CO2)的竞争吸附效应. 在模拟体系中引入了H2O和CO2, 其含量比例为V(N2)∶V(O2)∶V(CO2)∶V(H2O)∶V(CWA)=0.7855∶0.210∶0.0004∶0.0031∶0.001. 模拟方法、力场参数和收敛标准与上述基础模拟保持一致. 本研究通过引入疏水性作为筛选指标(依据Snurr课题组[61]提出的疏水性评价方法), 从Top1%的材料中筛选出适用于高湿环境的候选MOFs, 并从中选取3种代表性MOFs进行混合气体竞争吸附模拟. 结果表明, 在模拟潮湿环境下, H2O和CO2的引入未对化学战剂吸附产生显著竞争抑制(详细对比数据参见支持信息S7). 基于此, 我们提出分场景应用策略: 高湿环境应优先采用Top1%中的疏水型MOF以抑制水分子竞争吸附(见支持信息表S8); 干燥环境则可选用对化学战剂吸附容量更高的材料, 以充分发挥其吸附优势. CoRE-MOFs的结构文件详见支持信息的晶体结构数据表S1.
在本研究中, 我们构建了15333个CoRE-MOFs[62]的晶体模型, 并计算了可代表其特性的7个不同描述符. 其中, LCD、ρ和PLD, 借助Zeo++[63]软件完成计算. φ和VSA, 分别利用RASPA[64]软件包, 以直径0.258 nm的He和直径为0.364 nm的N2作为探针进行计算. K和Q0st作为能量描述符, 通过RASPA软件的NVT-MC方法, 计算每个分子在无限稀释条件下于不同MOFs中的吸附能量得到. 上述七个描述符(LCD, ρ, PLD, φ, VSA, K, Q0st)的组合选择, 是基于其明确的物理化学意义: LCD和PLD表征孔道的空间大小, ρ和φ反映了框架的致密与开放程度, VSA提供了材料的比表面积, 而K和Q0st则量化了宿主-客体相互作用的强度. 该描述符组合已在众多MOF高通量计算筛选研究中被采用并验证[65-67], 能够全面捕捉影响MOF吸附性能的关键结构与能量特征. 此外, 运用MOFid[39]对MOFs的结构进行分割与分析.
4.2 机器学习
本研究依托AI大数据技术框架下的机器学习方法, 探索MOFs材料与对目标气体吸附性能间的潜在联系—通过从海量材料表征与吸附性能数据中挖掘规律, 并对材料性能(TSN值)[68]进行预测(见公式1). 变量TSN通过对选择性进行对数化处理, 平衡了选择性与吸附量的数值关系, 缩小了两者在数量级上的差距, 为后续大数据建模提供了更鲁棒的目标变量, 从而能够更好地综合考量这两种变量的共同作用. Ming等[69]正是利用TSN这一特征, 在MOFs材料大数据集中综合评估其对芥子气模拟物二乙基硫醚(DES)的吸附能力与选择性, 进而通过机器学习划分材料性能类别, 为高效MOFs设计提供了数据驱动的AI决策支持.
本研究针对四种目标气体的大规模综合数据集, 采用六种机器学习算法开展回归分析, 这些算法包括来自 Python3.8 中scikit-learn数据包[70]的DT、RF、XGB、BPNN、GBR和LGBM. DT[71]通过递归划分特征空间构建回归模型, 基于信息增益等准则选择最优分裂特征, 但其单棵树在高维大数据场景下易过拟合降低泛化能力; RF[72]基于Bagging框架集成多棵决策树, 通过自助抽样与随机特征子集选择生成基模型, 以投票或平均方式输出结果, 有效提升模型稳定性与抗噪性; GBR[73]属于Boosting算法, 逐阶段训练弱决策树拟合前序模型残差, 通过梯度下降迭代优化损失函数, 擅长捕捉数据非线性关系但需控制复杂度; XGBoost[74]在GBR基础上引入正则化项与二阶泰勒展开, 支持并行计算与缺失值处理, 显著提升大规模数据训练效率与泛化性能, 广泛应用于AI大数据体系中的大规模结构化数据分析; BPNN[75]通过多层感知机实现非线性映射, 基于误差梯度反向调整网络权重, 适合处理图像、语音等非结构化数据, 但其“黑箱”特性限制可解释性; LGBM[76]采用直方图算法离散化特征、梯度采样筛选样本、特征捆绑降维, 在保持精度的同时大幅提升训练速度, 尤其适用于亿级规模数据集与内存受限场景. 在回归分析中, 为在大数据环境下保障模型的鲁棒性, 采用决定系数(R2)、平均绝对误差(MAE)和均方根误差(RMSE)等AI模型评估体系中常用的量化标准评估模型性能(见公式2~4), 确保模型在海量的材料数据中仍能稳定输出可靠结果.
$TSN{{N}_{i}}\times \lg ({{S}_{i/{{N}_{2}}+{{O}_{2}}}})$
${{R}^{2}}1\frac{\sum\nolimits_{i=1}^{n}{{{({{y}_{i}}{{{\hat{y}}}_{i}})}^{2}}}}{\sum\nolimits_{i=1}^{n}{{{({{y}_{i}}{{{\bar{y}}}_{i}})}^{2}}}}$
$RMSE\sqrt{\frac{1}{n}\sum\limits_{i=1}^{n}{{{({{y}_{i}}{{{\hat{y}}}_{i}})}^{2}}}}$
$MAE\frac{1}{n}\sum\limits_{i=1}^{n}{\left| {{{\hat{y}}}_{i}}{{y}_{i}} \right|}$
其中, 公式(1)中的i分别代表四种目标气体: 梭曼、沙林、DMMP和DIFP. 公式(2)~(4)中的${{\hat{y}}_{i}}$和yi分别代表预测值和真实值, n代表数据量. 通常认为机器学习算法的评价标准中R2越高, 模型的拟合效果越好, 泛化能力越强, RMSE和MAE越小, 说明预测的误差越小, 这种机器学习算法更好. 有关机器学习算法的更多详细信息, 如支持信息图S10~S15所示, 有关超参数优化的细节见支持信息表S9.
(Cheng, B.)