对分子模拟轨迹数据的分析绘图

article/2025/7/18 10:20:21

简介 

建立模型进行分子动力学模拟后,对得到的轨迹进行主成分分析(PCA),绘制相关性矩阵(Correlation Matrix)和自由能井图(Free Energy Landscape)和dssp图(Definition of Secondary Structure of Proteins)

绘制相关性矩阵图:

使用cpptraj计算轨迹的相关性矩阵

#进入amber的cpptraj程序:
cpptraj##第一步,计算轨迹的平均结构:
parm tyr.prmtop 
trajin tyr_500.crd#计算非H原子的rmsd,保存为数据集rmsd-dataset,选择第一帧为参考文件:
rms rmsd-dataset first !@H=#根据轨迹计算结构的平均构象,保存为数据集average-dataset:
average crdset average-dataset
run

计算Cα原子的相关性矩阵,输出为ca_matrix.gnu文件:

#利用1~30帧,计算所有CA的相关性矩阵,
#输出文件为ca_matrix.gnu,保存数据集为ca_matrix:
matrix covar out ca_matrix.gnu name ca_matrix start 1 stop 30 @CA

给数据前后补上脚本,具体例子如下:

#设置保存为pdf文件
set terminal pdfcairo enhanced color size 8,6 linewidth 3 font "Time New Roman,20"
set output "Correlation Matrix.pdf"set size square
set pm3d map
#设置色条颜色,rainbow色条:
set palette rgb 30,31,32
#差值法解决马赛克问题:0,0是插入的值,设置越大,马赛克消除效果越好。
#0,0表示让系统自动优化:
set pm3d interpolate 0,0
#设置标题和轴标签:
set title "Correlation matrix" font "Time New Roman,14"
set xlabel "atom CA" font "Time New Roman,10"
set ylabel "atom CA" font "Time New Roman,10"
#隐藏图例:
unset key
#设置x,y轴范围:
set yrange [1:265]
set xrange [1:265]
#设置轴的小刻度线:
set mxtics 3
set mytics 3
#刻度线向外:
set xtics out
set ytics out
#设置色条的刻度朝内:
set cbtics in
#隐藏对称的刻度轴:
set tics nomirror
#设置色条的刻度长度:
set cbtics scale 3
#绘图数据:
splot "-" with pm3d title "pca-out.gnu"1.000    1.000  1.0001.000    2.000  0.9041.000    3.000  0.7031.000    4.000  0.5231.000    5.000  0.7091.000    6.000  0.7531.000    7.000  0.6381.000    8.000  0.6491.000    9.000  0.6941.000   10.000  0.6651.000   11.000  0.6211.000   12.000  0.5981.000   13.000  0.5861.000   14.000  0.5491.000   15.000  0.5091.000   16.000  0.5051.000   17.000  0.5111.000   18.000  0.5491.000   19.000  0.6081.000   20.000  0.5721.000   21.000  0.582
[过多的数据省略...]set output
end
pause -1

在gnuplot里拖入该gun文件。

查看绘制的效果:

相关性矩阵的图像呈对角线对称,色条两端颜色代表的区域对应的残基有着越强的运动相关性,颜色越浅的区域对应的残基运动相关性越弱。

用cpptraj进行蛋白质二级结构分析(DSSP)

编写cpptraj脚本dssp.in

parm tyr.parm7 [top]   #导入拓扑文件
trajin tyr.nc parm [top]     #导入100frame的轨迹文件center origin :236&!@H= mass   #将残基236的质心作为盒子的中心对体系进行移动
image origin center familiar
reference tyr.pdb         #导入坐标参考文件
secstruct :15-68 out dssp15_68.gnu    #计算15-68残基的dssp
run
quit

注意:这里reference tyr.pdb很重要,不用pdb文件,它识别不了提取的残基id 

运行脚本:

cpptraj dssp.in

运行gnuplot进行绘制:

gnuplot dssp15_68.gnu

调整gnu脚本:

#设置保存为pdf文件
set terminal pdfcairo enhanced color size 8,6 linewidth 3 font "Time New Roman,20"
set output "DefinitionofSecondaryStructureofProteins.pdf"
set size square
set pm3d map corners2color c1
#隐藏图例:
unset key
#范围:
set yr [1.000:54.000]
set xrange [1.000:100.000]
#刻度值:
set ytics ("ILE:19"    5.000,"GLY:29"   15.000,"ARG:39"   25.000,"PHE:49"  \35.000,"GLN:59"   45.000,"PRO:68"   54.000) font "Times New Roman,10"
set xtics ("0" 0,"20" 10,"40" 20,"60" 30,"80" 40,"100" 50,"120" 60,"140" 70,"160"\80,"180" 90,"200"100) font "Times New Roman,10"#设置整个图表的空间(面积)大小:
##数值表示距离窗口左边界的距离:
set lmargin at screen 0.15
set rmargin at screen 0.85
##数值表示距离窗口下边界的距离:
set tmargin at screen 0.78
set bmargin at screen 0.20   
#标题和刻度标签:
set xlabel "Time (ps)" 
set ylabel "residues"
set title "Definition of Secondary Structure of Proteins 15-68"
#刻度线条外翻:
set tics out
#色条的范围是 -0.5~7.5
set cbrange [-0.0500:0.7500]
#调整色条刻度值:
set cbtics("None"    -0.0500,"Extended"    0.05000,"Bridge"    0.1500,"3-10"  \ 
0.2500,"Alpha"    0.3500,"Pi"    0.4500,"Turn"    0.5500,"Bend"    0.6500)
set palette maxcolors 8
set palette defined (0 "#000000",1 "#0000FF",4 "#00FF00",7 "#FF0000")
#调整色条刻度线:
set cbtics in
set cbtics scale 2.8
set tics nomirror
splot "-" with pm3d title "dssp15_68.gnu"1.000    1.000             0.01.000    2.000             0.01.000    3.000             0.71.000    4.000             0.01.000    5.000             0.01.000    6.000             0.01.000    7.000             0.11.000    8.000             0.11.000    9.000             0.11.000   10.000             0.1[过多的数据省略]set output
end
pause -1

利用cpptraj进行PCA分析绘制自由能井图

编写脚本pca.in:

parm tyr.prmtop [p]
trajin tyr.crd parm [p]
rms rms-data :1-264&!@H= first      #生成rmsd数据集
run
average crdset average-dataset :1-264&!@H=   #生成平均构象数据集
run
createcrd tyr-trajectories-dataset   #建立一个轨迹数据集
runcrdaction tyr-trajectories-dataset \  #计算这个轨迹数据集的rmsd,参考平均构象数据集
rms ref average-dataset :1-264&!@H=crdaction tyr-trajectories-dataset \     #计算协方差矩阵
matrix covar name tyr-covar-dataset :1-264&!@H=
runrunanalysis diagmatrix tyr-covar-dataset out tyr-evecs.dat \  #获取前3个特征向量
vecs 3 name my-evecs \
nmwiz nmwizvecs 3 nmwizfile tyr.nmd nmwizmask :1-264&!@H=crdaction tyr-trajectories-dataset projection TYR modes my-evecs \
beg 1 end 3 :1-264&!@H=hist TYR:1 bins 25 out tyr-hist.agr norm name TYR-1
hist TYR:2 bins 25 out tyr-hist.agr norm name TYR-2
hist TYR:3 bins 25 out tyr-hist.agr norm name TYR-3
hist TYR:1 TYR:2 bins 25 out pca-out.gnu  norm name all_1_2  #将PC1和PC2输出为一个gnu文件
run

运行脚本:

cpptraj pca.in

得到的图像有点马赛克感觉;而且色图外围还有一圈白色的留白;还有颜色也不是我们熟悉的红蓝色条的热图。

进行脚本的润色:

#设置保存为pdf文件:
set terminal pdfcairo enhanced color size 8,6 linewidth 3 font "Time New Roman,20"
set output "Free Energy Landscape.pdf"
#必须要在splot末尾加上"set output"语句,否则保存的文件会损坏!
set size square
set pm3d map
#设置色条颜色,rainbow色条:
set palette rgb -31,13,22
#差值法解决马赛克问题:0,0是插入的值,设置越大,马赛克消除效果越好。
#0,0表示让系统自动优化:
set pm3d interpolate 0,0
#设置标题和轴标签:
set title "Free Energy Landscape" font "Time New Roman,20"
set xlabel "PC1" font "Time New Roman,20"
set ylabel "PC2" font "Time New Roman,20"
#隐藏图例:
unset key
#这里必须根据实际数据进行修改:
set yrange [ -18.955:  22.277]
set xrange [ -25.709:  22.280]
#设置轴的小刻度线:
set mxtics 3
set mytics 3
#刻度线向外:
set tics out
#隐藏对称的刻度轴:
set tics nomirror........(表示数据)#输出:
set output
end
pause -1

查看保存的pdf图片:

 

编写脚本将概率值转换为相对自由能,再绘制FEL图:

使用origin软件绘制自由能井图:

将数据文件pca-out.gnu改为pca-out.dat,将里面和gnuplot有关的命令语法删掉,只保留数据值,导入origin里:

将第3列数据表示强度的数据列改为Z轴:

 对3列数据全选,在【工作表】-->【转化为矩阵】-->【xyz网格化】生成一个虚拟矩阵。

选择虚拟矩阵,选择【绘图】-->【等高线图】,调整调整

关于origin软件保存图片带有demo的水印的问题(提醒自己):

发现用修复工具之后,一段时间demo由会像蘑菇一样冒出来(F**k!),修复操作如下:

水印修复工具(就是那个Origin.exe文件),复制(不要移动)到安装目录里,替换里面原有的Origin.exe

打开安装目录下的Origin64.exe就能去掉水印啦。


http://chatgpt.dhexx.cn/article/Rr4XS1l3.shtml

相关文章

【学习笔记】ICLR2022-GNNRefine

一、代码Run 1. 运行 python /home/huangjiehui/Project/AIProtein/StructuralReconstruction/RefineGNN/ab_train.py --cdr_type 3 --train_path /home/huangjiehui/Project/AIProtein/StructuralReconstruction/RefineGNN/data/sabdab/hcdr3_cluster/train_data.jsonl --va…

利用pymol批量对PDB文件三维结构比对并输出RMSD值

经验上、文献上大量的蛋白质或者核酸比对工作是从一级序列开始的,这是基于一级序列决定二级结构,二级结构决定三级结构,而且一级序列有30%的相似性,那么两者的结构就具有较高的相似性这样的共识理论而来,这些理论也是基…

AMBER:使用Cpptraj计算RMSD 以及使用中遇到的问题

记录笨比生活又一天 输入文件rms.in设置运行cpptraj遇到的问题1.cpptraj不输出结果2.空格的问题 Tofirst:[空格]1-249&!H firstTofirst[空格]:1-249&!H first 输入文件rms.in设置 parm XXXXX.prmtop #载入拓扑文件 trajin XXXX_prod.nc #载入轨迹文件 rms ToFirst …

分子动力学模拟Amber/Gromacs结合自由能计算 药效团模型构建RMSD、RMSF

文章来源:公众号“科研讨论圈” 以下是使用AMBER、GROMAVCS的教程,希望对开始学习分子动力学的同学有帮助。 分子动力学入门理/论 分子力学简介 分子…

RMSD:通过旋转计算两个分子间的最小rmsd

使用旋转计算两个分子的均方根偏差(RMSD) 使用Kabsch算法(1976)或Quaternion算法(1991)进行旋转,在两个笛卡尔坐标之间.xyz或者.pdb格式中计算均方根偏差(RMSD)&#xf…

PSP - TMScore(US-align)、RMSD、Sequence 源码

欢迎关注我的CSDN:https://spike.blog.csdn.net/ 本文地址:https://blog.csdn.net/caroline_wendy/article/details/129125467 参考文档:Nature Methods | 蛋白、RNA、DNA及其复合物结构的比对算法US-align 官网地址:https://zha…

基于Gromacs的蛋白分子动力学模拟(RMSD、RMSF及蛋白的回旋半径)

一、实验要求 实验对象:目标体系为modeller或其他方法建模的结果中评价最好的模型。 软件: Gromacs-5.1.2 二、实验步骤 加立场 gmx pdb2gmx –h 打开帮助菜单。 选力场的时候选择 Amber99sb…,溶剂类型选Tip3p。 2、加模拟盒子,溶剂层…

RMSD与PMSF 解释与区别

我不需要力量源泉 对我来说 这只是需要完成的一件事,我应该一直全力以赴的做下去 放弃不是我的天性,我甚至不在乎乐观还是悲观。 二者都是对位移的平方和再求平方根,最后求得均值。区别在于:RMSD为在同一帧情况下,对不…

pymol pymol-align两分子或蛋白距离误差计算RMSD;spyrmsd库计算RMSD

参考:https://zhuanlan.zhihu.com/p/347743101 https://www.codenong.com/cs106148400/ RMSD 单位是埃 RMSD,root-mean-square deviation,也就是均方根偏差。 原子位置的均方根偏差是叠加蛋白质的原子(通常是骨架原子)之间的平均…

RMSD和RMSF

RMSD,Root Mean Square Deviation,均方根偏差;RMSF,Root Mean Square Fluctuation,均方根波动。 在轨迹分析中,最经常用,最简单,也最有用的就是这两巨头,二者都是对位移的…

统计学常用指标

目录 标准差(SD)均方根误差(RMSE/RMSD)均方误差(MSE/MSD)平均绝对误差(MAE)决定系数/拟合优度(R^2)平均偏差(Bias) 标准差(SD) 标准差(StandardDeviation),在概率统计中最常使用作为统计分布程度(statisticaldispersion&#xff…

Excel如何查找批注

Excel查找的选项,查找范围选择批注就可以了

Excel 打印显示批注(亲自实践)

有时候需要将Excel中的批注,随同正文一起打印出来 方法如下: 1.右键有批注的单元格,选择"显示/隐藏批注" 2.选择"分页预览"模式显示工作表 3.右键任意单元格,选择"页面设置" 4.在"页面设置"窗口中,标签"工作表"的&quo…

excel中深入理解批注

excel中深入理解批注 系统要求一、单个批注二、所有批注三、形状修改四、插入图 系统要求 装有office2010以上版本功能:批注 一、单个批注 右击显示或隐藏批注 显示表示批注不退出(鼠标移开时) 二、所有批注 审阅显示或隐藏批注 三、…

Excel中批量添加批注图片

excel中想实现这种悬停时显示图片的效果 1、将图片与单元格命名一致,并将图片与excel文件放置在同一目录下 2、选中需要设置的单元格,点击【开发工具】 -查看代码(如果没有开发工具往下看) 3、在编辑框中输入以下代码并运行 Sub …

Excel技能之批注超链接,你竟然真的不会用

资源整合是新时代职场人的顶配。资源整合,可以快速查找想要的资源。一个目录,如果包含了所有相关的内容,那么,价值是宝贵的。点击目录打开对应的内容,离不开超链接。超链接,改变了互联网,也改变…

Excel表格中重要的数据如何隐藏不显示

Excel表格中重要的数据如何隐藏不显示 目录 Excel表格中重要数据值如何隐藏不显示 1、选中需要隐藏的单元格数值 2、鼠标右键 点击“设置单元格格式” 3、点击“自定义”在“类型”一栏中输入三个“;”号即可(输入法切换在英文状态输入分号)。 4、想…

关于poi/Npoi创建批注后,EXCEL不能显示,wps能显示的问题(2020-08-25)

一般百度到这个: https://blog.csdn.net/zyr2206328732/article/details/48341191 实测不是作者描述的原因。 我的需求上:在列名(一个集合)增加注释。 代码如下:Row headRow sheet.createRow(0);XSSFDrawing drawing…

esayExcel自定义注解导出表头批注

注解 package com.baidu.activitidemo.annotation;import com.baidu.activitidemo.handler.ExcelRemarkHandler;import java.lang.annotation.*;/*** 设置表头的批注, 需要配合{link ExcelRemarkHandler}使用** author li* date 2022/09/24*/ Target(ElementType.FIELD) Reten…

把Excel批注的“红三角”放在单元格左上角_Excel的批注功能,全部知道的不足10%,你会用的仅仅是冰山一角...

Excel【审阅】功能区中,显示了5个最基本的功能,如下图红色矩形框所示,这也是我们最常用的基础。其实,Excel中批注的相关操作远不止这些,一起来看看。 插入删除批注 【插入批注】 插入批注的方法常见的有这3种: ❶ 【审阅】→【新建批注】。 ❷ 单元格点击右键,插入批注。…