多肽自组装计算 — 短肽β-sheet自组装的粗粒化MD模拟方法

为什么模拟短肽自组装比想象中难得多

短肽自组装在生物材料、药物递送和纳米技术中有着广泛的应用前景。一个经典的例子是二苯丙氨酸（FF）二肽，它能在水溶液中自组装成纳米管、纳米线等多种有序结构。从计算化学的角度看，模拟这类自组装过程在理论上似乎是直接的——把肽分子放在水里，设置好力场参数，跑足够长时间，”看”它们如何自发组装。

但现实残酷得多。我们首次尝试用全原子MD模拟FF二肽自组装时，在一个12×12×12 nm³的盒子里放了100个FF分子，跑了1微秒——结果只看到了一些不稳定的二聚体和三聚体，完全没有形成实验中观察到的高度有序纳米管。诊断结论是：全原子模拟能覆盖的时间尺度（微秒级）远小于自组装成核所需的实际时间（毫秒到秒级），而且全原子水溶剂极大地增加了计算负担。

这就是为什么多肽自组装计算必须依赖粗粒化模型。通过将多个重原子合并为一个”珠子”，粗粒化模型将模拟的时间尺度从微秒延伸到了毫秒甚至秒级——进入了自组装现象的真实时间窗口。本文详细复盘我们使用MARTINI粗粒化力场模拟短肽β-sheet自组装的全过程。

全原子到MARTINI：映射策略的正确打开方式

MARTINI力场采用4:1的映射策略——大约每4个重原子映射为1个粗粒化珠子。对于短肽，这是一个需要谨慎权衡的过程：映射太粗会丢失残基的化学特异性（比如无法区分天冬氨酸和谷氨酸），映射太细则计算收益有限。

对于FF二肽，每个苯丙氨酸残基映射为3个MARTINI珠子：骨架珠子（BB）、侧链珠子1（SC1，代表苯环的Cα-Cβ-Cγ部分）和侧链珠子2（SC2，代表苯环的其余部分）。二肽的总珠子数为6个——从原来的约70个原子减少到6个珠子，计算加速了约1000倍（因为时间步长从2飞秒增加到20飞秒，且需要计算的相互作用数量大幅减少）。

映射设计中最关键的决策是侧链珠子的划分方式。苯丙氨酸的苯环是一个整体电子结构，从化学直觉上似乎应该映射为1个珠子。但我们发现1珠映射无法再现苯环之间的π-π堆积几何——实际的苯环面对面堆积需要两个平面之间的特定间距和偏转角，而单个球形珠子无法编码这种方向性。因此我们选择2珠映射，并在珠子之间加入了角度约束来保留苯环的平面取向信息。这个”多余的”复杂度在后续的组装模拟中证明是完全必要的——1珠模型的组装产物是无定形聚集而非有序的β-sheet结构。

力场参数化：从实验和全原子数据构建粗粒化势

MARTINI力场的参数化哲学是”自下而上”与”自上而下”的结合。”自下而上”是指从高精度的全原子模拟中提取目标分布函数——键长分布、键角分布、二面角势——然后用这些分布来拟合粗粒化势。”自上而下”是指用实验热力学数据（如分配系数、水合自由能）校准非键相互作用强度。

在我们的FF项目中，非键参数的校准遇到了一个经典困境：MARTINI默认的非键相互作用强度对肽-肽堆积太弱，导致模拟中观察到的二聚化平衡常数比实验值低一个数量级。经过反复测试，我们将SC1-SC2珠子间的ε值从默认的3.5 kJ/mol提高到5.0 kJ/mol——这个调整被纯苯-水分配系数的实验数据所支持（苯在MARTINI中默认的分配系数偏亲水）。

参数化完成后必须进行验证：我们用100个FF分子的粗粒化模型跑了10微秒，观察到的组装行为与实验中的冷冻电镜图像定性一致——形成了宽度约8-12纳米、长度数十纳米的纳米带，β-sheet的层间距约0.47纳米与X射线衍射数据高度吻合。

组装动力学：从随机分布到有序结构的四个阶段

通过粗粒化轨迹的分析，我们识别出了FF自组装过程中四个清晰的动力学阶段。第一阶段（0-0.5 μs）：成核前的滞后期。分子在水中自由扩散，偶尔形成寿命极短（<10 ns）的二聚体。这个阶段的特征是高碰撞频率但极低的”成功”附着概率——99%以上的碰撞都以解离告终。

第二阶段（0.5-2 μs）：临界核形成。当一个三聚体或四聚体偶然形成了正确的β-sheet氢键排列——相邻分子之间的骨架珠子呈反平行排列——它就成为了一个稳定的”种子”。我们观察到，小于3个分子的聚集体几乎总是不稳定的，而一旦达到4-5个分子，解离的概率急剧下降。这个临界核大小的概念为我们理解自组装的触发条件提供了定量基础。

第三阶段（2-8 μs）：快速生长期。临界核形成后，单体以极高的效率附着到生长中的聚集体边缘。这个阶段的生长动力学近似线性——聚集体的分子数量随时间线性增长，说明生长是扩散控制的（单体的扩散是限速步骤），而非反应控制的。

第四阶段（8-10 μs）：成熟和结构调整期。组装体数量不再增加，但内部结构继续优化——不良堆积的分子从表面解离后被正确取向的分子取代。这是一个”奥斯特瓦尔德熟化”的分子版本。最终体系达到了一个稳态分布：约90%的分子纳入组装体，10%以单体形式存在于溶液中。

定量表征：如何评估你的组装模拟是否”正确”

除了视觉上”看起来像”实验结构外，多肽自组装计算需要严格的定量表征。我们使用了以下核心指标：

组装度——纳入最大聚集体的分子比例。在10微秒时这个值达到87%，与实验中动态光散射估算的聚集程度一致。β-sheet含量——通过相邻分子骨架珠子间距离和取向角来判别是否为β-sheet排列。在最大聚集体中，β-sheet氢键的占比为72%，其余为无序堆积或β-turn构象。

有序度参数——聚集体内部分子的取向有序性。我们定义了向列序参数S₂ = (3⟨cos²θ⟩ – 1)/2，其中θ是分子主轴方向与聚集体平均方向的夹角。S₂从0 ns时的0.05（几乎完全无序）升到10 μs时的0.82（高度取向有序），定量地展示了从无序到有序的转变。

溶剂可及表面积的变化——反映了疏水效应在驱动组装中的作用。在组装过程中，FF分子的总SASA从初始的约750 nm²减少到约250 nm²，减少的500 nm²对应于疏水苯环从水环境中被”埋藏”到组装体内部——这是组装的热力学驱动力在分子层面的直接体现。

——————————————————————————————

📐 科研计算代算服务

我们提供专业的第一性原理计算、分子动力学模拟、有限元仿真分析及材料表征测试服务。
团队由博士学历的计算化学、计算材料学专家组成，覆盖VASP、LAMMPS、GROMACS、ANSYS、COMSOL等主流工具链。

  🔗 立即咨询计算方案