Gromacs——教程学习(2)

5. 分子动力学模拟结果分析与绘图自动化
   基于GROMACS分子动力学结果自动分析并绘图的工具–HeroMDAnalysis
   HeroMDAnalysis的工作流程
   

我们需要确保使用的GROMACS为2019.4或者更高的版本，HeroMDAnalysis使用Xmgrace进行绘图，因此需要系统已经安装了Xmgrace，如果尚未安装，可使用如下命令进行安装：

sudo apt-get install grace

HeroMDAnalysis工具的新版本要求安装SHC，可输入以下命令进行安装：

sudo add-apt-repository ppa:neurobin/ppa
sudo apt-get update
sudo apt-get install shc

下面我们需要下载HeroMDAnalysis工具包，并提取至工作目录中。下载地址如下：

ruhttps://heromdanalysis.wordpress.com/customized-version-of-the-program-for-a-specific-ligand-protein-complex-system/

本例中使用的HeroMDAnalysis版本为module1，程序的更新与定制版本可联系作者获取。
下面我们需要在如下网站注册一个账号，根据提示填写相应的信息即可。

https://heromdanalysis.wordpress.com/contact-us/

在工作目录中打开终端，注意工作目录中应包含GROMACS运行动力学之后的轨迹文件（.xtc, .edr和.tpr文件），输入以下命令：

chmod +x launch.sh
./launch.sh

弹出对话框提示在工作目录中识别到md_0_10.tpr、md_0_10.edr以及md_0_10.xtc 或者md_0_10_fit.xtc三个文件可用，HeroMDAnalysis输出的文件也同样保存在该目录下。
下面我们选择要分析的系统，可选择蛋白-配体复合物体系（Protein-ligand Complex）或脱辅基蛋白（纯蛋白，Apoprotein），根据自己要分析的体系选择即可，我们这里就选择复合物。
下面根据弹出的窗口提示选择自己要分析的项目，如RMSD分析可选择Calpha-PL、Backbone-PL、Sidechain-PL、Protein-PL、RMSD-L、Lig fit Prot以及Prot fit Lig，根据自己的需求进行选择即可。

同样的选择RMSF、SASA（Solvent accessible surface area，溶液可及表面积）、溶剂化自由能（Free energy of solvation）、氢键（H-bond contact）等的分析。
HeroMDAnalysis会根据我们勾选的要分析的项目自动调用GROMACS中对应的模块进行分析，如gmx rms、gmx rmsf、gmx gyrate、gmx sasa、gmx hbond、gmx energy等。
在弹出的“Would you like to plot data to images”窗口可以选择是否自动绘制图像，如果我们选择“Yes, let’s plot data to HD quality images”，那么HeroMDAnalysis会使用上步生成的动力学结果分析数据自动绘制出高质量的可供发表的图像；如果我们选择“No, I don’t want it”，则不会生成图像，但我们仍然可以以生成的分析数据使用各类绘图软件如Excel、Origin等自行绘图。
我们这里就选择让HeroMDAnalysis自动绘制图像，在弹出的窗口根据提示填写相应的信息，如蛋白的PDB ID、配体的名称、动力学运行的时间（如果运行了100 ns的分子动力学，则输入“100”）等。
此外，我们也需要根据窗口提示输入相应的绘图参数，如X轴和Y轴的最小最大刻度值等，同时会在相邻的位置打开Xmgrace窗口以供参考。
全部设置完成后，就等待程序运行完成吧。下面来看一下HeroMDAnalysis自动绘制的图像。

参考资料：
1.Rawat R, Kant K, Kumar A, et al. HeroMDAnalysis: an automagical tool for GROMACS-based molecular dynamics simulation analysis[J]. Future Medicinal Chemistry, 2021, 13(05): 447-456.
2.https://heromdanalysis.wordpress.com/user-manual-2/
（文末附视频教程）

6. 分子动力学模拟结果分析之MDAnalysis实例1—RMSD计算

计算分子动力学轨迹中蛋白质的均方根偏差（Root Mean Square Deviation, RMSD）是一种常用的分析方法，用于衡量蛋白质结构随时间的变化。
MDAnalysis 是一个用于分析分子动力学（MD）模拟的 Python 库。它是面向对象的工具，使用户能够处理来自多种流行格式的 MD 模拟轨迹。
MDAnalysis (https://www.mdanalysis.org/)
https://docs.mdanalysis.org/stable/index.html

（1）安装Anaconda
如果您还没有安装Anaconda，您可以从Anaconda官网下载并安装它。
（2）打开Anaconda Prompt
安装完成后，在开始菜单中找到并打开Anaconda Prompt。
（3）创建虚拟环境
输入命令 conda create -n my_mdanalysis_env python=3.8 来创建一个名为 my_mdanalysis_env 的新环境。这里使用的Python版本是3.8，但您可以选择适合MDAnalysis的其他版本。
按回车键，然后根据提示完成环境的创建。
（4）激活虚拟环境
一旦虚拟环境创建完成，您可以通过命令 conda activate my_mdanalysis_env 来激活这个环境。
（5）安装MDAnalysis
在激活的虚拟环境中，输入 pip install MDAnalysis 来安装MDAnalysis库。
（6）验证安装
安装完成后，您可以通过输入 python 进入Python的交互模式，然后尝试导入MDAnalysis来检查是否正确安装：import MDAnalysis as mda。
如果没有报错，表示MDAnalysis已成功安装在您的虚拟环境中。
（7）退出虚拟环境
当您完成工作后，可以通过命令conda deactivate来退出虚拟环境。

编写 Python 脚本计算 RMSD

import MDAnalysis as mda

from MDAnalysis.analysis import rms

# 加载拓扑结构和轨迹文件
u = mda.Universe('topology_file.pdb', 'trajectory_file.xtc')

# 选择蛋白质原子
protein = u.select_atoms('protein')

# 选择参考帧（通常是第一帧）
ref = mda.Universe('topology_file.pdb')
ref_protein = ref.select_atoms('protein')

# 初始化 RMSD 分析
R = rms.RMSD(protein, ref_protein, select='backbone')
R.run()

# 打印结果
import numpy as np
time = u.trajectory.times
rmsd = R.rmsd.T  # transpose makes it easier for plotting
print("时间（ps）: ", time)
print("RMSD（Å）: ", rmsd[2])

# 可以选择将结果绘制为图表
import matplotlib.pyplot as plt
plt.plot(time, rmsd[2])
plt.xlabel("时间 (ps)")
plt.ylabel("RMSD (Å)")
plt.show()

说明
替换‘topology_file.pdb 和 ‘trajectory_file.xtc’ 为你的实际文件名。
‘protein’ 和 ‘backbone’ 选择子用于选择蛋白质和其主链原子。这些选择子可能需要根据你的系统进行调整。
RMSD 是针对时间的每一帧计算的，相对于初始帧（或你选择的任何参考帧）。