上一節(jié)我們對(duì)蛋白配體復(fù)合物體系進(jìn)行一個(gè)能量最小化與動(dòng)力學(xué)模擬的操作，這一節(jié)我們將對(duì)其結(jié)果進(jìn)行分析。

重新定位和重新包裝坐標(biāo)

與任何使用周期性邊界條件進(jìn)行的模擬一樣，分子可能會(huì)出現(xiàn) "斷裂 "或在單元格中來(lái)回 "跳躍"。要重新調(diào)整蛋白質(zhì)的位置并將分子重新包覆在單元格內(nèi)以恢復(fù)所需的菱形十二面體形狀，請(qǐng)調(diào)用 trjconv：

gmx trjconv -s md_0_10.tpr -f md_0_10.xtc -o md_0_10_center.xtc -center -pbc mol -ur compact

選擇 "Protein "為中心，"System "為輸出。需要注意的是，如果模擬時(shí)間較長(zhǎng)，涉及許多跨周期邊界的躍遷，則可能難以對(duì)復(fù)合物（protein-ligand, protein-protein）進(jìn)行居中。在這種情況下（尤其是在protein-protein復(fù)合物中），可能需要?jiǎng)?chuàng)建一個(gè)自定義索引組來(lái)進(jìn)行定 心，該索引組與復(fù)合物中一個(gè)蛋白質(zhì)的活性位點(diǎn)或一個(gè)單體的界面殘基相對(duì)應(yīng)。

要提取軌跡的第一幀（t = 0 毫微秒），請(qǐng)使用 trjconv -dump 和重新定心的軌跡：

gmx trjconv -s md_0_10.tpr -f md_0_10_center.xtc -o start.pdb -dump 0

為了獲得更流暢的視覺(jué)效果，可能需要進(jìn)行旋轉(zhuǎn)和平移擬合。執(zhí)行 trjconv 如下：

gmx trjconv -s md_0_10.tpr -f md_0_10_center.xtc -o md_0_10_fit.xtc -fit rot+trans

選擇 "Backbone（骨架）"對(duì)蛋白質(zhì)骨架進(jìn)行最小二乘擬合，選擇 "System（系統(tǒng)）"進(jìn)行輸出。請(qǐng)注意，同時(shí)進(jìn)行 PBC 重包和坐標(biāo)擬合在數(shù)學(xué)上是不相容的。如果您希望進(jìn)行擬合（這對(duì)可視化很有用，但對(duì)大多數(shù)分析程序來(lái)說(shuō)并非必要），請(qǐng)按照此處的說(shuō)明分別進(jìn)行這些坐標(biāo)操作。

分析蛋白質(zhì)配體相互作用和配體動(dòng)力學(xué)

本教程不可能涵蓋您希望執(zhí)行的所有分析方法。這里將對(duì)一些基本操作進(jìn)行說(shuō)明

2-propylphenol 配體可以與 Gln102 側(cè)鏈形成氫鍵。我們將計(jì)算 JZ4 的羥基與 Gln102 的羰基 O 原子間的距離。形成氫鍵的典型標(biāo)準(zhǔn)是供體原子和受體原子之間的距離≤ 3.5 ?（0.35 nm）。使用命令行選擇語(yǔ)法（示例和更多語(yǔ)法請(qǐng)參見(jiàn) gmx 幫助選擇），使用距離模塊計(jì)算軌跡上的距離。

gmx distance -s md_0_10.tpr -f md_0_10_center.xtc -select 'resname "JZ4" and name OAB plus resid 102 and name OE1' -oall

平均距離為 0.31 ± 0.05 nm，與氫鍵的形成一致。

確定是否存在氫鍵的第二個(gè)標(biāo)準(zhǔn)通常是供體原子、氫原子和受體原子之間形成的角度。計(jì)算角度有不同的慣例。在 GROMACS hbond 模塊中，角度定義為氫-供體-受體，該角度應(yīng)≤ 30°。要進(jìn)行此分析，首先要為供體原子（必須包括供體重 原子和鍵合氫原子）和受體原子創(chuàng)建索引組：

gmx make_ndx -f em.gro -o index.ndx

使用角度模塊分析這三個(gè)原子形成的角度：

gmx angle -f md_0_10_center.xtc -n index.ndx -ov angle.xvg

請(qǐng)注意，與其他大多數(shù) GROMACS 分析模塊不同，角度計(jì)算不使用 -s 參數(shù)。角度計(jì)算不需要質(zhì)量或周期信息、原子名稱等，因此在處理軌跡時(shí)不需要 .tpr 或坐標(biāo)文件。該命令返回的值大致為 23 ± 9°。這一結(jié)果可能有些出乎意料，因?yàn)槲覀儤?gòu)建的索引組是按照 OAB、H12、OE1 的順序排列的，這與供體-氫-受體距離相對(duì)應(yīng)，而我們預(yù)計(jì)供體-氫-受體距離接近線性（~180°）。檢查索引文件的內(nèi)容，您會(huì)發(fā)現(xiàn)

[ JZ4_&_OAB_H12_Protein_&_r_102_&_OE1 ]?

1616 2624 2636

make_ndx 自動(dòng)將原子編號(hào)從低到高排序，因此計(jì)算結(jié)果是受體-供體-氫角，與 hbond 模塊計(jì)算的角度相同。因此計(jì)算結(jié)果與氫鍵的形成是一致的，因?yàn)闅滏I的角度≤ 30°。為了得到所需的供體-氫-受體角度，我們必須手動(dòng)編輯文本文件中的索引組，重新排列原子編號(hào)（2624 2636 1616）。使用該索引組重新進(jìn)行角度計(jì)算得出的平均值為 147 ± 11°

最后，我們可能有興趣量化配體結(jié)合姿態(tài)在模擬過(guò)程中的變化程度。要計(jì)算 JZ4 的重原子 RMSD，請(qǐng)為其創(chuàng)建一個(gè)新的索引組：

gmx make_ndx -f em.gro -n index.ndx

執(zhí)行 rms 模塊，選擇 "Backbone "進(jìn)行最小二乘擬合，選擇 "JZ4_Heavy "進(jìn)行 RMSD 計(jì)算。這樣做可以通過(guò)擬合去除蛋白質(zhì)的整體旋轉(zhuǎn)和平移，所報(bào)告的 RMSD 就是 JZ4 位置相對(duì)于蛋白質(zhì)的變化程度，是模擬過(guò)程中結(jié)合姿態(tài)保持程度的良好指標(biāo)。

gmx rms -s em.tpr -f md_0_10_center.xtc -n index.ndx -tu ns -o rmsd_jz4.xvg

計(jì)算得出的 RMSD 應(yīng)該在 0.1 nm（1 ?）左右，這表明配體的位置變化非常小。

蛋白質(zhì)配體相互作用能量

為了量化 JZ4 和 T4 溶菌酶之間相互作用的強(qiáng)度，計(jì)算這兩種物質(zhì)之間的非鍵相互作用能可能是有用的。GROMACS 能夠分解任意數(shù)量的定義基團(tuán)之間的短程非鍵能。需要注意的是，這個(gè)量并不是自由能或結(jié)合能。事實(shí)上，大多數(shù)力場(chǎng)的參數(shù)化方式并不能使這個(gè)量具有實(shí)際的物理意義。CHARMM 的參數(shù)化專門(mén)針對(duì)與水的量子力學(xué)相互作用能，因此它在本質(zhì)上平衡了有意義的量，因此相互作用能是有用的。

相互作用能的計(jì)算是通過(guò) .mdp 文件中的 energygrps 關(guān)鍵字進(jìn)行的。盡管是 .mdp 關(guān)鍵字，但相互作用能計(jì)算不應(yīng)被視為正常模擬的一部分。短程能量的分解與在 GPU 上運(yùn)行不兼容，而且會(huì)不必要地減慢計(jì)算速度。mdrun 模塊不需要做這些額外的工作來(lái)執(zhí)行有效的模擬。因此，計(jì)算相互作用能只能作為分析的一部分，而不能作為動(dòng)力學(xué)的一部分。從定義了 energygrps = Protein JZ4 的 .mdp 文件創(chuàng)建一個(gè)新的 .tpr 文件，就像下面這樣：

gmx grompp -f ie.mdp -c npt.gro -t npt.cpt -p topol.top -n index.ndx -o ie.tpr

接下來(lái)，使用 -rerun 選項(xiàng)調(diào)用 mdrun，根據(jù)現(xiàn)有模擬軌跡重新計(jì)算能量：

gmx mdrun -deffnm ie -rerun md_0_10.xtc -nb cpu

注意使用 -deffnm 讀取 ie.tpr，并將所有輸出文件寫(xiě)入 ie.* 作為文件名。-rerun 選項(xiàng)使用了要重新計(jì)算能量的軌跡名稱，而 -nb cpu 則告訴 mdrun 只嘗試在 CPU 硬件上運(yùn)行，而忽略任何可能可用的 GPU。如上所述，這種計(jì)算無(wú)法在 GPU 上執(zhí)行。重新運(yùn)行的速度應(yīng)該很快，幾分鐘就能完成。

通過(guò)能量模塊提取感興趣的能量項(xiàng)。我們感興趣的項(xiàng)是 Coul-SR:Protein-JZ4 和 LJ-SR:Protein-JZ4

gmx energy -f ie.edr -o interaction_energy.xvg

平均短程庫(kù)侖相互作用能為 -20.5 ± 7.4 kJ mol-1，短程倫納德-瓊斯能為 -99.1 ± 7.2 kJ mol-1。從這些量的相對(duì)大小得出結(jié)論可能很有誘惑力，但盡管 CHARMM 是根據(jù)明確的水相互作用能進(jìn)行參數(shù)化的，將相互作用能進(jìn)一步分解為這些分量并不一定是真實(shí)的。沒(méi)有辦法通過(guò)實(shí)驗(yàn)驗(yàn)證這些量，因此無(wú)法知道它們是否有意義。不過(guò)，在這種情況下，總相互作用能是有用的。該值（根據(jù)兩個(gè)量相加的標(biāo)準(zhǔn)公式傳播誤差后）為 -119.6 ± 10.3 kJ mol-1。

至此，分析部分結(jié)束。