蛋白质聚集的分子动力学模拟（ Molecular dynamics simulations of protein aggregation）
蛋白质紊乱和聚集在许多神经退行性疾病（如阿尔茨海默病和帕金森病）的发病机制中起着重要作用。这些疾病中聚集过程的最终产物是高度结构化的淀粉样纤维。因此，合理设计靶向淀粉样蛋白寡聚体形成的药物，就需要全面了解驱动蛋白质聚集的物理化学力。原子级的分子动力学 (MD) 模拟提供了该过程最高的时间和空间分辨率，可以捕捉淀粉样蛋白寡聚体形成过程中的关键步骤。
1．MD 模拟和基本分析
以 Aβ16-22 为例，其两端都带有封端基团，序列为 ACE-KLVFFAE-NME。由于蛋白质数据库 (PDB) 不包含该肽的结构，因此可以从 Aβ42 的 PDB 结构的残基 16-21 的坐标中检索以下模拟的起始结构，例如 PDB ID 1Z0Q。在 PyMOL 的蛋白质模式下使用 Builder 工具，ACE 和 NME 封端组可以分别添加到 N 和 C 端。
1.1 六肽模拟体系建立
1. 第一步是为 Aβ16-22 单体产生一个松弛的构象。对于 Aβ16-22，建议模拟长度为 1 μs 或更长，使用构象聚类确定最稳定的单体结构（本文省略这步）。2. 使用 PACKMOL将六个单体随机放置到一个模拟盒子中。下面的示例脚本将六个 Aβ16-22 肽放置在一个大小为 ~10 nm × 10 nm × 10 nm 的模拟框中，它们之间的距离至少为 1.2 nm（或脚本中的 12 Å）。
packmol.inp#Six monomers of abeta16-22 peptide #minimum distance between two monomers tolerance 12.0 seed -1 #The file type of input and output files is PDB filetype pdb #The name of the output file output abeta16-22_hexamer.pdb #add TER to after each monomer chain add_amber_ter #distance from the edges of box add_box_sides 1.0 #path to input structure file #units for distance is measured in Angstrom #box size is 100 Å structure abeta16-22.pdb number 6 inside box 0. 0. 0. 100. 100. 100. end structure
运行packmol
packmol < packmol.inp
或者使用gmx insert-molecules命令
gmx insert-molecules -ci abeta16-22.pdb -nmol 6 -box 10 10 10 -o abeta16-22_hexamer.pdb
1.2 不同的模拟步骤创建目录
创建五个主要步骤的目录：拓扑文件构建、能量最小化、NVT 平衡、NPT 平衡和 MD 生产运行。每一步都需要一个 .mdp 文件。mdp 文件类型扩展名代表分子动力学参数，这些文件包含使用 GROMACS 设置 MD 模拟的所有关键参数。
mkdir 1-topol 2-em 3-nvt 4-npt 5-md mdp
1.3 拓扑构建
拓扑文件包含有关将被模拟的分子类型和分子数量的信息。上一步的 .pdb 文件作为输入，除了拓扑文件之外，还会生成一个 .gro 文件，该文件与 .pdb 文件一样，也包含模拟系统的坐标。它们之间的主要区别在于它们的格式。
(1)从http://macerell.umaryland.edu/download.php?filename=CHARMM_ff_params_files/charmm36-mar2019.ff.tgz下载Charmm36m力场，复制到1-topol目录下。切换到该目录：
cd 1-topol/
（2）运行GROMACS pdb2gmx 命令处理输入的结构文件并创建扩展名为.top 的拓扑文件、扩展名为.itp 的拓扑包含文件和扩展名为.itp 的位置约束文件。
gmx pdb2gmx -f ../abeta16-22_hexamer.pdb -o protein.gro -p topol.top -ignh -ter <<EOF
1
1
3
4
3
4
3
4
3
4
3
4
3
4
EOF
选项说明:
-f: 读取输入结构文件 abeta16-22_hexamer.pdb
-o and -p: 写入输出结构文件 protein.gro 和系统拓扑文件 topol.top
-ignh: 忽略输入文件中的氢原子，这是可取的，因为输入文件和力场中氢原子的命名约定不同。GROMACS 将使用所选力场的 H 原子名称添加新的氢原子。
-ter: N 端和 C 端的质子化状态可以-ter 标志交互选择
Option 1: choosing protein force field (charmm36-mar2019)
Option 1: choosing water force field (TIP3P)
Option 3: choosing N-terminus (None, as we use ACE capping)
Option 4: choosing C-terminus (None, as we use NME capping)
Options 3 and 4 are repeated for each peptide in the system, in this example six.
成功执行 GROMACS pdb2gmx 命令后，该目录还将包含六个拓扑包含文件和六个位置约束文件，每个肽一个：topol_Protein_chain_$chain.itp 和 posre_Protein_chain_$chain.itp，其中 $chain = A、B、C、D、E 或 F。后面的文件包含肽的所有非氢原子的位置限制条目，这是平衡 MD 步骤期间需要的。
（3） 接下来，创建一个模拟框。请注意，上面定义的框仅用于放置肽。和以前一样，选择了一个 10 × 10 × 10 nm3 的立方体盒子：
gmx editconf -f protein.gro -o box.gro -bt cubic -box 10 10 10
（4）在模拟盒里装入水分子
gmx solvate -cp box.gro -cs spc216.gro -o protein-solvated.gro -p topol.top溶剂分子的添加反映在 topol.top 文件中，除了 6 个肽段外，该文件现在还包括 32,323 个水分子。（5）为了使用周期性边界条件进行 MD 模拟，我们需要根据需要通过添加正 (Na+) 或负 (Cl-) 离子来中和系统的电荷。此外，我们可以为系统指定一个特定的离子浓度，通常约为 150 mM 以模拟生理条件：
gmx grompp -f ../mdp/ions.mdp -c protein-solvated.gro -p topol.top -o protein-ions.tpr
echo 13 | gmx genion -s protein-ions.tpr -o protein-ions.gro -p topol -neutral -conc 0.15
Explanation of flags and options: -neutral: 通过添加 Na+ 或 Cl- 离子将系统中和至零电荷 -conc: 将浓度更改为 150 mM NaClecho 13: 用离子替换溶剂分子，其中选项 13 代表 SOL（溶剂分子）。
在第一个命令中，grompp 模块，读取坐标文件、系统拓扑文件、ions.mdp 文件并将它们处理成 GROMACS 二进制格式生成 .tpr 文件。该文件包含模拟的起始结构、单个肽和水分子的拓扑信息以及所有模拟参数。第二个命令读取二进制输入文件 protein-ions.tpr 并中和净电荷，添加了 90 个 Na+ 和 90 个 Cl- 离子，相应地更新 topol.top 文件。
