通过分子建模方法探索半胱天冬酶突变和翻译后修饰

本协议使用生物分子模拟包，并描述了用于模拟野生型半胱天冬酶及其突变形式的分子动力学（MD）方法。MD方法允许评估半胱天冬酶结构的动态演变以及突变或翻译后修饰的潜在影响。

细胞凋亡是一种程序性细胞死亡，可消除受损细胞并控制多细胞生物的发育和组织稳态。半胱天冬酶是半胱氨酸蛋白酶家族，在细胞凋亡的启动和执行中起关键作用。半胱天冬酶的成熟及其活性通过翻译后修饰以高度动态的方式进行微调。为了评估翻译后变化的影响，潜在的位点经常发生突变，残基对任何修饰都持续存在。例如，丝氨酸残基被丙氨酸或天冬氨酸取代。然而，这种取代可能会改变半胱天冬酶活性位点的构象，导致催化活性和细胞功能的紊乱。此外，位于关键位置的其他氨基酸残基的突变也可能破坏半胱天冬酶的结构和功能，导致细胞凋亡扰动。为了避免使用突变残基的困难，可以很容易地应用分子建模方法来估计氨基酸取代对半胱天冬酶结构的潜在影响。该协议允许使用生物分子模拟包（琥珀）和超级计算机设施对野生型半胱天冬酶及其突变形式进行建模，以测试突变对蛋白质结构和功能的影响。

细胞凋亡是调节多细胞生物形态发生和组织稳态的最广泛研究的细胞过程之一。细胞凋亡可由广泛的外部或内部刺激引发，例如死亡受体的激活、细胞周期信号紊乱、DNA 损伤、内质网 （ER） 应激以及各种细菌和病毒感染¹.半胱天冬酶 - 关键的凋亡参与者 - 通常分为两组:引发剂（半胱天冬酶-2，半胱天冬酶-8，半胱天冬酶-9和半胱天冬酶-10）和效应物（半胱天冬酶-3，半胱天冬酶-6和半胱天冬酶-7），这取决于它们的结构域结构和半胱天冬酶级联中的位置2^，3。在细胞死亡信号上，引发剂半胱天冬酶与接头分子相互作用，促进邻近诱导的二聚化和自动处理以形成活性酶。效应半胱天冬酶通过引发剂半胱天冬酶的切割被激活，并通过切割多个细胞底物4来执行下游执行^步骤。

引发剂和效应半胱天冬酶的成熟和功能受大量不同细胞内机制的调节，其中翻译后修饰在细胞死亡调节中起着不可或缺^{的作用5。}添加修饰基团（磷酸化、亚硝基化、甲基化或乙酰化）或蛋白质（泛素化或 SUMO化）会改变半胱天冬酶的酶活性或调节细胞凋亡的蛋白质构象和稳定性。定点诱变被广泛用于研究潜在的翻译后修饰位点并识别其作用。推定的修饰位点通常被另一种氨基酸取代，该氨基酸不能进一步修饰。因此，潜在的磷酸化丝氨酸和苏氨酸突变为丙氨酸，赖氨酸泛素化位点被精氨酸取代。另一种策略包括替换特别模仿翻译后修饰的氨基酸（例如，谷氨酸和天冬氨酸已被用于模拟磷酸化的丝氨酸或苏氨酸）⁶。然而，位于活性位点附近或关键位置的一些取代物可能会改变半胱天冬酶结构，干扰催化活性并抑制凋亡细胞死亡⁷。在半胱天冬酶基因中肿瘤相关的错义突变的情况下也可以观察到类似的效果。例如，半胱天冬酶-6 - R259H的肿瘤相关突变导致底物结合口袋中环的构象变化，降低了底物⁸的有效催化周转。头颈部鳞状细胞癌中半胱天冬酶-8中的G325A氨基酸取代可能阻碍半胱天冬酶-8活性，导致核因子-kB（NF-kB）信号传导的调节，促进肿瘤发生⁹。

为了评估氨基酸取代对半胱天冬酶结构和功能的潜在影响，可以应用分子建模。这项工作描述了分子动力学（MD）方法，用于使用生物分子模拟包（Amber）模拟野生型半胱天冬酶及其突变形式。MD方法给出了引入突变后蛋白质结构动态演变的观点。琥珀包最初由Peter Kollman的团队开发，成为生物分子模拟^10，11，^12，13中最受欢迎的软件工具之一。该软件分为两部分:（1）AmberTools，常规用于系统准备（原子类型分配，添加氢和显式水分子等）和轨迹分析的程序集合;（2）琥珀色，以PMEMD模拟程序为中心。AmberTools 是一个免费软件包（也是安装 Amber 本身的先决条件），而 Amber 则以单独的许可证和费用结构分发。在超级计算机上和/或使用图形处理单元（GPU）进行并行模拟可以大大提高蛋白质结构动力学科学研究的性能¹⁴。最新的可用软件版本是AmberTools21和Amber20，但所描述的协议也可以与以前的版本一起使用。

1. 系统准备

注意:天然和突变蛋白质形式的分子模型是基于从蛋白质数据库^15，16获得的适当晶体结构构建的。

1. 要检索选定的 PDB 结构，请使用" 下载文件 "下拉列表，然后单击 "PDB 格式"。删除备注和连接数据，并在 PDB 文件中的单独蛋白质链之间插入 TER 卡。或者，通过将残基名称 HIS 替换为 HIE、HID 或 HIP（补充文件 1）来设置组氨酸残基的质子化状态。
   注意:不去除晶体分辨的水分子（如果存在于 PDB 文件中）是合理的。
2. 要准备起始模型，请启动 tleap 程序（AmberTools 软件包），然后输入操作 tleap 对象的命令。
   1. 启动程序: 跳跃。加载ff14SB力场以用分子力学描述蛋白质: >来源leaprc.protein.ff14SB。
   2. 水分子和原子离子（Na⁺，Cl^-）的载荷力场:>来源leaprc.water.tip3p。
   3. 加载 PDB 文件 并构建氢的坐标，创建一个名为 mol 的对象: > mol = loadpdb 蛋白.pdb。
   4. 检查可能导致问题的内部不一致: >检查 mol。
      注意:二硫键（如果存在）必须手动指定。将共价结合的半胱氨酸CYS的名称替换为CYX，并在 tleap 中键入以下命令以在SG原子之间创建键: 键 mol.X.SG mol.Y.SG （X和Y是半胱氨酸残基数）。
   5. 在蛋白质周围创建一个溶剂盒（TIP3P水，12 Å距离）: >溶剂化盒分子TIP3PBOX 12.0。
   6. 检查总电荷: >电荷摩尔，并添加反离子（Na⁺ 或Cl^-）以中和系统:
      >添加剂摩尔Na⁺ 0。
   7. 创建拓扑 prmtop 文件和坐标 inpcrd 文件（用于运行模拟的输入）:> saveamberparm mol protein.prmtop protein.inpcrd。退出 tleap 计划:>退出。
      注意:琥珀色坐标文件可以使用 ambpdb 工具轻松转换为 PDB 格式（请参阅用户手册，请参阅 材料表）。

2. 能源最小化

注意:能量最小化对于消除启动系统中原子之间的任何不良接触和重叠是必要的，这些接触和重叠会导致运行MD时不稳定。

1. 进行能量最小化的第一阶段（最陡下降算法的2，500步+共轭梯度算法的2，500步），以优化添加的氢原子和水分子的位置，同时通过重原子的位置约束保持蛋白质坐标固定。
   1. 按如下方式运行 pmemd 程序:
      pmemd -O -i min1.in -p protein.prmtop -c protein.inpcrd -o protein_min1.out
      -r protein_min1.rst -ref protein.inpcrd.
   2. 遵循所需的参数:-i 文件控件数据;-p 文件分子拓扑、力场参数、原子名称;-c 文件初始坐标;-o 文件用户可读日志输出;-r 文件最终坐标;-ref 位置约束的文件参考坐标。
   3. 在输入文件中指定选项 min1.in:
      &cntrl
      IMin=1， maxcyc=5000， ncyc=2500，
      削减=10.0， NTB=1，
      NTC=1， NTF=1，
      NTPR=10，
      ntr=1，
      约束面具=':1-517 & ！@H='，
      restraint_wt=2.0
      /
      注意:输入选项:imin=1 执行最小化;maxcyc=5000 最大最小化循环次数;NCYC=2500最小化方法在NCYC循环后从最陡下降切换到共轭梯度;切割=10.0 非保税截止值 （Å）;NTB=1 强加周期性边界，恒定体积;NTC=1 键长不施加约束（以获得更好的能量收敛）;NTF=1 计算所有交互作用;NTPR=10 每个 NTPR 步骤将能源信息打印到输出文件中;ntr=1 标志，用于使用谐波势约束指定的原子;约束掩码=':1-517 & ！@H=' 指定约束原子;restraint_wt=2.0 重量 （千卡/摩尔∙Ų） 用于位置约束。
2. 执行能量最小化的第二阶段（5，000 个最陡峭下降步骤+ 5，000 个共轭梯度步骤），不受限制地优化整个系统。
   1. 使用 min2.in 和 protein_min1.rst 作为输入: pmemd -O -i min2.in -p protein.prmtop
      -c protein_min1.rst -o protein_min2.out -r protein_min2.rst.
   2. 在输入文件中指定选项 min2.in:
      &cntrl
      imin=1， maxcyc=10000， ncyc=5000，
      削减=10.0， NTB=1，
      NTC=1， NTF=1，
      NTPR=10
      /

3. 加热

注意:此阶段旨在将系统从0 K加热到300 K.初始速度分配给原子，因为基于PDB文件的起始模型不包含速度信息。

1. 在蛋白质原子（50 ps，恒定体积）上具有位置约束的情况下进行加热过程。
   1. 使用 heat.in 和 protein_min2.rst 作为输入: pmemd -O -i heat.in -p protein.prmtop
      -c protein_min2.rst -o protein_heat.out
      -r protein_heat.rst -x protein_heat.mdcrd -ref protein_min2.rst
   2. 遵循所需的参数:保存在 MD 轨迹上的 -x FILE 坐标集。
   3. 在输入文件中指定选项 heat.in:
      &cntrl
      imin=0， irest=0， ntx=1，
      nstlim=25000， dt=0.002，
      NTC=2， NTF=2，
      削减=10.0， NTB=1，
      NTPR=500， NTWX=500，
      ntt=3， gamma_ln=2.0，
      温度=0.0， 温度0=300.0，
      ntr=1， 约束面罩=':1-517'，
      restraint_wt=1.0，
      核磁共振=1
      /
      &wt 类型='TEMP0'，
      istep1=0， istep2=25000，
      值 1=0.1， 值 2=300.0 /
      &wt 类型='结束' /
      注意:输入选项:imin=0 运行 MD;irest=0 运行新的模拟;ntx=1 没有从 RST 文件中读取初始速度;nstlim=25000 MD步数;dt=0.002 时间步长 （ps）;NTC=2 涉及氢的键使用 SHAKE 算法进行约束;不计算约束键的NTF=2力;NTWX=500 个坐标每 NTWX 步骤写入 MDCRD 文件;NTT=3朗格文温度调节;朗格文动力学的碰撞频率gamma_ln=2.0（ps⁻¹）;温度=0.0初始温度;温度0=300.0参考温度;读取 &WT 名称列表中指定的 nmropt=1 参数;类型='TEMP0' 改变目标温度。

4. 平衡

注意:此阶段对于调整水的密度并获得蛋白质的平衡状态是必要的。

1. 在300 K下进行平衡，没有任何约束（500 ps，恒定压力）。
   1. 使用 equil.in 和 protein_heat.rst 作为输入: pmemd -O -i equil.in -p protein.prmtop
      -c protein_heat.rst -o protein_equil.out -r protein_equil.rst -x protein_equil.mdcrd.
   2. 在输入文件中指定选项 equil.in:
      &cntrl
      imin=0， irest=1， ntx=5，
      nstlim=250000，
      dt=0.002，
      NTC=2， NTF=2，
      削减=10.0， NTB=2， NTP=1， 灰褐色=2.0，
      NTPR=1000， NTWX=1000， NTWR=50000，
      ntt=3， gamma_ln=2.0，
      温度0=300.0
      /
      注意:输入选项: irest=1 重新开始模拟; ntx=5 坐标和速度是从先前生成的 RST 文件中读取的; NTB = 2 施加周期边界，恒定压力; NTP=1 各向同性压力标度; 灰褐色=2.0 压力松弛时间（PS）; NTWR=50000 每个 NTWR 步骤都会写入 RST 文件，确保从崩溃中恢复。

5. 生产动态

1. 成功达到平衡后，在恒定压力下进行生产MD模拟（10 ns或更长），并生成轨迹文件，用于蛋白质结构的后续分析。
   1. 使用 prod.in 和 protein_equil.rst 作为输入: pmemd -O -i prod.in -p protein.prmtop
      -c protein_equil.rst -o protein_prod.out -r protein_prod.rst -x protein_prod.mdcrd.
   2. 在输入文件中指定选项 prod.in:
      &cntrl
      imin=0， irest=1， ntx=5，
      nstlim=5000000，
      dt=0.002，
      NTC=2， NTF=2，
      削减=10.0， NTB=2， NTP=1， 灰褐色=2.0，
      NTPR=1000， NTWX=1000， NTWR=50000，
      ntt=3， gamma_ln=2.0，
      温度0=300.0， ig=-1
      /
      注意:pmemd 的并行版本 （pmemd.MPI）或GPU加速版本（pmemd.cuda）可用于计算机集群和超级计算机。长时间的MD模拟可以分成几个部分并按顺序执行。强烈建议每次模拟重启时设置 ig=-1（随机种子选项）。ptraj或cppptraj程序可用于坐标轨迹的分析和处理。有关详细信息，请参阅软件用户手册。

本方案可以很容易地应用于半胱天冬酶或致病突变的翻译后修饰的研究。在本节中，说明了MD建模工作流程（图1），该工作流程已成功用于半胱天冬酶-2⁷的研究。使用潜在磷酸化位点的 体外 定点诱变（Ser/Thr至Ala）和生化方法，证明Ser384Ala突变阻止了半胱天冬酶-2加工并阻断了酶活性和诱导凋亡细胞死亡（补充图1）。然而，质谱分析和Phos-Tag技术都没有证实Ser384经历磷酸化⁷。为了了解Ser384如何调节半胱天冬酶-2活性，对三维蛋白质结构进行了MD建模（图1）。

使用1pyo晶体结构构建了半胱天冬酶-2野生型和突变型的分子模型（可用的半胱天冬酶-2结构列于 补充表1中）。Ser384Ala突变体是通过去除Ser384残基中的O^γ 原子而创建的（在PDB中，Ser和Ala残基通过具有O^γ 原子来区分）。晶体结构中通常没有氢坐标。因此，将氢原子添加到蛋白质结构中，然后由12 Å厚的TIP3P水层溶解，以便进一步进行显式溶剂模拟。加入钠离子以中和系统。利用 min1.in、 min2.in、 heat.in、 equil.in 和 prod.in 控制数据文件，根据MD建模协议对半胱天冬酶-2的起始模型进行能量最小化、平衡和随后的10 ns MD模拟。

在半胱天冬酶-2的晶体结构中，Ser384与Arg219和Arg378一起参与在活性位点形成腔表面。Arg219和Arg378与基底的羧基形成氢键，而Ser384不干预与底物的直接相互作用。MD模拟证实，Ser384Ala取代不影响催化残基Cys320（亲核试剂）和His277（一般碱基），但诱导Arg378发生重大构象变化。它的胍基转向本体溶剂，因此N^ε 原子不能与底物的羧基形成必需的氢键（图1）。在天然酶中，Ser384（部分负电荷）的O^γ 原子和Arg378胍基团（正电荷）之间发生静电相互作用。然而，丝氨酸替代显然破坏了这种相互作用。因此，表明Ser384Ala取代影响了半胱天冬酶-2活性位点中精氨酸残基的底物识别，损害了酶活性和触发凋亡细胞死亡的能力。发现的机制似乎在进化上是保守的，并且对于其他半胱天冬酶家族成员来说很常见。因此，MD模拟的实施允许确认生化结果并获得对半胱天冬酶活性中心的分子结构的新见解。

图 1:用于研究半胱天冬酶结构的 MD 建模工作流程。 按照协议部分中的说明研究了野生型半胱天冬酶-2及其Ser384Ala突变体。结果表明，Ser384Ala取代诱导活性位点残基Arg378的重要构象变化。 请点击此处查看此图的大图。

补充图1:半胱天冬酶-2在细胞死亡中的作用。 为了响应外在和内在刺激，野生型半胱天冬酶-2可以通过自蛋白水解机制激活。活性半胱天冬酶-2切割Bid，进而促进线粒体外膜透化和凋亡细胞死亡。Ser384Ala突变可防止半胱天冬酶-2活化和诱导细胞死亡。 请点击此处下载此文件。

补充表1:蛋白质数据库中可用的半胱天冬酶-2结构。 结构按发布日期排序。 请按此下载此表格。

补充文件 1:编辑 PDB 文件的不同步骤。请点击此处下载此文件。

所描述的MD方法允许使用生物分子模拟包对半胱天冬酶的野生型和突变形式进行建模。这里讨论了该方法的几个重要问题。首先，需要从蛋白质数据库中选择半胱天冬酶的代表性晶体结构。重要的是，半胱天冬酶的单体和二聚体形式都是可以接受的。选择缺失残基数量最少的高分辨率结构是一个好主意。某些残基的质子化状态可以在输入PDB文件中手动设置。例如，在 图1中，氢原子连接到催化残基His277的N^ε2 原子（HIE形式）上。其次，必须在蛋白质周围创建一个溶剂盒，以便进一步进行显式溶剂模拟。需要能量最小化来优化添加的氢原子和水分子的坐标。第三，加热和平衡必须在生产动态模拟之前进行。在加热阶段，应确保蛋白质表面上的空腔和结合袋充满水分子（如有必要，可以增加MD步骤的数量）。在平衡阶段，蛋白质平衡构象的获得需要通过分析其主链原子与初始位置的均方根偏差来确认¹²。第四，需要在生产动力学模拟过程中监测蛋白质结构的构象变化。为此，应通过拟合主链原子将轨迹框架叠加到起始结构上，然后使用分子可视化工具（VMD，PyMOL等）分析构象变化。¹⁷. 如有必要，可以增加MD步骤的数量和/或将模拟分解为按顺序执行的几个部分。第五，建议使用GPU加速版本的MD仿真程序（pmemd.cuda），其性能明显超过传统CPU实现所能达到的性能^14，18。

所描述方法的一个可能的局限性是蛋白质数据库中某些半胱天冬酶结构的可用性。然而，在"半胱天冬酶"请求下发现了900多个PDB条目，表明已经详细研究了半胱天冬酶结构。当以低分辨率（即，在粗略的细节水平）获得晶体结构或缺少一些对半胱天冬酶功能很重要的片段时，可能会出现一些困难。

总之，MD建模方法已被证明可以有效地预测蛋白质基因突变，修饰或分子间相互作用后的结构变化。MD模拟在细胞死亡领域的应用为研究翻译后修饰半胱天冬酶调控的分子机制提供了巨大的机会。在这方面，长期MD模拟可以评估功能区域的动态行为和氨基酸取代的潜在影响，以阐明与突变或修饰的半胱天冬酶-2以及其他半胱天冬酶相关的分子原因。例如，在胃肠道癌症以及中枢和周围神经系统肿瘤中检测到引发性半胱天冬酶基因（半胱天冬酶-2/半胱天冬酶-8/半胱天冬酶-9/半胱天冬酶-10）的错义突变¹⁹。总之，半胱天冬酶的MD建模是一种强大的 计算机 方法，用于了解调节程序性细胞死亡和相关疾病的机制以及开发新的有效疗法。

作者没有利益冲突需要披露。

这项工作得到了俄罗斯科学基金会（17-75-20102，协议开发）的资助。代表性结果部分（磷酸化分析）中描述的实验得到了斯德哥尔摩（181301）和瑞典（190345）癌症协会的支持。