基于Martini 3力场的胶原纤维粗粒化模型开发与力学响应机制研究

本研究针对胶原蛋白力学特性研究中全原子模拟计算成本高、现有粗粒化模型分辨率不足的问题，开发了基于Martini 3力场的胶原纤维粗粒化模型。研究人员通过热力学计算和玻尔兹曼反演方法参数化二价HLKNL和三价PYD交联键，结合Gō模型维持三重螺旋结构，成功实现了从单分子到335 nm全长微纤维的多尺度模拟。该模型在平衡态和1000 pN拉力条件下与全原子模拟及实验数据高度吻合，揭示了交联构型对纤维力学性能的调控机制，为研究衰老相关组织力学变化提供了新工具。

　　

胶原蛋白作为动物界最普遍的蛋白质之一，在肌腱、骨骼等结缔组织中承担着机械力传递的关键角色，其独特的^{三重螺旋结构（triple helical structure）和共价交联键（crosslinks）网络共同赋予了组织卓越的力学性能。然而，理解交联键在胶原纤维稳定性中的作用一直面临计算方法的挑战：全原子分子动力学（all-atom MD）模拟虽能提供原子级精度，但受限于计算资源仅能模拟单D带（67 nm）片段；而现有粗粒化（coarse-grained, CG）模型又缺乏残基级分辨率和精确的三重螺旋结构表征。这种跨尺度建模的缺失，阻碍了人们对胶原纤维在生理载荷下力学响应机制的深入探索。}

为解决这一难题，来自德国马克斯·普朗克聚合物研究所和海德堡理论研究所的Matthias Brosz、Johanna Buck等研究人员在《Biophysical Journal》发表了突破性研究，成功开发出基于Martini 3力场的胶原纤维粗粒化模型。该模型首次实现了在保留残基级分辨率的同时，对全长335 nm胶原微纤维进行高效模拟，为揭示交联键构型与力学性能的关联提供了全新工具。

研究团队采用多尺度参数化策略，首先通过非平衡态热力学计算确定二价羟基赖氨酸酮正亮氨酸（HLKNL）和三价吡啶啉（PYD）交联键的珠类型，再通过玻尔兹曼反演拟合键长、键角和二面角参数。针对胶原特有的(Gly-X-Y)_n重复序列（X常为Pro，Y常为Hyp），创新性地采用双峰高斯分布描述主链键长的序列特异性。为维持螺旋结构，引入Gō模型施加虚拟位点间的Lennard-Jones势垒。模型验证涵盖三个层次：单三重螺旋、67 nm交联纤维和全长335 nm微纤维，在300-1000 pN拉力范围内与全原子模拟和实验数据对比。

关键技术方法包括：1) 基于CHARMM36m和AMBER99SB-ildn*力场的全原子模拟作为参数化基准；2) 非平衡态分子动力学计算交联键的油水分配系数（logP）；3) 通过ColBuilder工具自动生成不同尺度的胶原纤维拓扑结构；4) 在Fugaku超计算机上完成4,300万原子的全原子对照模拟；5) 采用修正的Crooks涨落定理计算自由能变化。

交联键映射与珠类型选择

研究人员通过对称性分析和极性论证，为HLKNL和PYD交联键设计了特异性的映射方案。HLKNL的中心仲胺基用SP1d珠表示，丙醇和丙酮基团分别用P1和N6a珠映射。PYD的芳香环采用四个微型珠（TQ2p、TP1q、TC6q和TC4）构成菱形结构维持环刚度，部分电荷分布遵循AMBER力场的离域特性。热力学计算显示CG模型与全原子模拟的logP值高度吻合（HLKNL：-19.06±0.07 kJ/mol；PYD：-60.63±0.24 kJ/mol），溶剂可及表面积（SASA）和Connolly表面分析进一步验证了分子尺寸的准确性。

三重螺旋主链键合参数化

通过分析全原子模拟的映射轨迹，发现主链珠距离呈双峰分布（μ_b,1=0.318 nm, μ_b,2=0.354 nm），首次揭示羟脯氨酸（HYP）和脯氨酸（PRO）贡献于较短键长峰值。据此开发了双电位方案：通用氨基酸对采用18,000 kJ/mol/nm2力常数，HYP/PRO-X组合采用34,000 kJ/mol/nm2力常数。键角和二面角参数通过单峰分布拟合获得。模型成功复现了胶原的关键结构参数：每残基上升高度（CG：0.304 nm vs 实验：0.290 nm）和每圈残基数（CG：3.15 vs 实验：3.28）。

力响应下的参数优化

在300-1000 pN拉力模拟中发现，标准Martini 3键合力常数（4,000 kJ/mol/nm2）严重低估胶原刚度，而优化后的参数与全原子模拟的力-延伸曲线高度一致。对于三价交联纤维，CG模型虽在键角分布上呈现更高刚性，但刻意保持此特性以维持螺旋结构稳定性。关键突破在于明确高拉力下（>1 nN）力学响应由键拉伸主导而非Gō模型势阱深度，确立了键合参数优化的优先性。

全长微纤维的力-延伸动力学

在335 nm微纤维的千皮牛拉力模拟中，CG模型与全原子模拟达到相同延伸量（图5A）。重叠区-间隙区应变比分析显示CG值（0.4）与全原子数据（0.5）接近，但略高估重叠区刚度，这源于交联区原子细节的简化。沿纤维轴的力分布分析揭示CG模型成功复现了重叠区低力、间隙区高力的交替模式，虽未完全捕获间隙-过渡区的力峰值，但整体力学响应高度一致。

该研究建立的Martini 3胶原模型突破了现有模拟方法的尺度限制，首次实现了在残基分辨率下对全长胶原微纤维的力学行为模拟。模型通过16倍计算加速（CG：5.0 ns/天/核心 vs AA：0.31 ns/天/核心），使研究交联密度与类型对组织力学的影响成为可能。值得注意的是，交联键的精确参数化是再现纤维力学响应的关键，特别是三价PYD交联的菱形映射方案有效维持了芳香环在高拉力下的稳定性。

研究启示我们，胶原纤维的力学性能由多层级结构共同决定：分子级别的键合刚度、介观级别的交联网络拓扑、以及宏观级别的纤维组装方式。该模型为探索衰老相关交联变化（如晚期糖基化终产物积累）对组织力学的影响提供了新范式，同时通过与ColBuilder工具的集成，为生物材料设计和组织工程应用提供了可预测的计算平台。未来工作可进一步拓展至其他交联类型，并结合反向映射技术实现化学反应性的原子级预测。