今天仍然是分子动力学相关内容,假如我MD模拟了一个酶(做了不同溶液条件下的模拟或互调控蛋白的结合模拟),如果我想知道其催化口袋受不同因素影响的开放程度应该怎样做?那么这就需要用到CAVER软件。
应用实例大家请看这篇文献《A contribution to lipid digestion of Odobenidae family: Computational analysis of gastric and pancreatic lipases from walrus》海象是保护动物,作者为研究海象脂肪代谢相关的酶利用已知物种的脂肪酶做了进化分析并进行了同源建模,再利用MD做了稳定性分析并使用CAVER探究了不同盐浓度下酶催化口袋的开放程度。
一、准备输入文件
在tunnel计算中我们希望结果能体现口袋的动态开闭过程,因此显然我们不能拿单个pdb进行计算,这里教给大家抽帧生成pdb的指令,500ps一帧,因此0-100ns共201个pdb模型,共同输出到一个caver_input.pdb中
# 在GROMACS工作环境下输入 gmx trjconv -s md.tpr -f md_final.xtc -o caver_input.pdb -dt 500 -pbc mol -center -n index.ndx加index是为了仅保留目标蛋白的主体部分,请按蛋白质序列实际情况自行编写index文件,详情可参考我的《在云服务器AutoDL实现分子动力学全流程》文章
二、下载CAVER及相关准备工作
1.下载CAVER
自行搜索CAVER Analyst进入官网 --> 上方点击download --> 弹出的界面下载第一个All platforms即可
全部解压到你的文件夹,最好不用C盘
2.准备java
CAVER是一个基于java的语言,需要java 1.8以上。
# adoptium是下载java依赖很好的途径 # https://adoptium.net/zh-CN/temurin/releases/ # 选择8 -LTS或者11 - LTS # 选择jdk # 点击下载获得的msi程序 # msi程序不像exe程序使用管理员权限安装会很复杂,可以直接下在C盘用户AppData文件夹,这是个隐藏文件夹,请记一下下载路径 # 打开windows搜索环境变量,新建path,粘贴下载java的路径,可以把这个java拖拽到第一条 # 做生信的人大概率是没有装java的,但是也请注意区分之前安装的版本和x86的java另外我们找到下载CAVER的文件夹,找到 etc --> caver_analyst.conf 找到#jdkhome="/path/to/jdk"这一行,jdkhome="C:\Program Files\Eclipse Adoptium\jdk8u412-b08"根据你的实际下载路径去掉前面的“#”然后改成实际路径。
3.打开CAVER的方式及memory设置
找到你的caver_analyst2 --> bin --> 点击 x64.exe 文件即可,也可以为其创建快捷方式添加到桌面;初始memory大概是1000MB,会弹出形如这样的弹窗
点yes,因为1000MB对100-200ns的MD轨迹太少了,如果您的电脑是16GB的可以设置4000-6000MB,32GB则可以设置8000MB。
三、文件导入及计算tunnel操作
先在下载caver_analyst2的dir里创建一个caver_work文件夹,把caver_input.pdb放进去(为了保存workspace的时候可以找到原始pdb)找到下载的 caver_analyst2 --> bin --> caver_analyst64.exe,直接右键运行或者创建快捷方式。在软件内点击 file --> Open Molecular Dynamics --> PDB -->add file选择自己的文件
在下方选项栏structure dynamics里面可以拖动滑条看各帧构象,如果不喜欢默认显示模式可以在上方工具栏按我的设置展示cartoon模式
在下方选项栏Sturctur Squence通过点击选择待计算的res,尽量根据先验选择4-5个res
之后点击上方工具栏Tunnel --> surrounding --> from selection --> 输入你的选择 --> output directory改成你的caver_work地址 --> Compute Tunnels
等待即可,当出现“是否计算surrounding”的提示时点击“yes”
出现结果后先 File --> Save Workspace 保存操作存档,勾选包含input pdb的选项,此后该次计算cws文件可直接通过CAVER软件打开。
下方选项栏Tunnel Statistic会有两个子选项SET #1和SET #1 bottlenecks:
SET #1:会给出summary信息,Max_BR、Avg_L、Throughput是主要评估标准;点击其中一个cluster会出现该聚类tunnel出现的pdb位置,一般我们会锁定BR最大的分析;再次点击某个tunnel则会出现各部位通道尺寸数据(由于工具采用微分思想用小球模拟tunnel形状,这里的radius与summary中的BR有微小差异)
SET #1 bottlenecks:会按1-201的顺序逐帧展示每个pdb中出现的tunnel,后面的信息则是起主要贡献的res
所有以上信息会被输出到预定文件夹里的各个.csv文件中,各位可自行查看,另外右侧的Structures Overview可点击条各cluster后面的条形图符号进入Tunnel Graph工具,该工具可展示tunnel各部位radius变化的趋势。
四、ChimeraX及pymol的可视化操作
1.存储通道的方式
上方工具栏应该是没有存储通道的方式的,我看Guide文件也没有找到,可以点击右边栏目你想下载的tunnel,右键会有下载选项,存成pdb。
输出的文件夹会存储各帧的obj文件,但chimerax无法识别这个格式,因为我使用的是可视化CAVER版本,没有直接构建与pymol的管道,如果您有什么更便捷的方法也可以分享给我。
2.可视化方法
ChimeraX打开你的caver_input.pdb,会有1.1-1.201等很多子模型,选择你要的那帧,命令选择或手动点击选择按钮,在下方命令栏输入delete ~sel,之后保存成新文件即可
另一个窗口打开tunnel的pdb文件,会发现该文件保存了该cluster所有的组,如果使用split命令可以拆分但这样很可能无法知道目的tunnel是哪一根,所以我的建议是用任意文本编辑器打开pdb文件,照着set #1里面details的各部位的radius找到目的tunnel,然后把其他的手动删了。给大家贴一段pdb文件内容,简单来说pdb文件其实只是用一种特定的格式记录了原子的类型、残基位数、三维坐标等
HEADER TUNNEL COMPND caver_input EREST VAL EREST SER EREST HID EREST ILE EREST HIE EREST LYS EREST GLN EREST PHE EREST PRO EREST TYR EREST HIP EREST GLU EREST HIS EREST TRP EREST GLY EREST 20_AA EREST ALA EREST ARG EREST CYS EREST ASN EREST LEU EREST MET EREST ASP EREST THR ATOM 1 H FIL T 496 99.152 95.849 30.985 1.15 ATOM 2 H FIL T 496 99.249 95.713 30.514 1.37 CONECT 1 2 ATOM 3 H FIL T 496 99.018 95.610 30.138 1.57 CONECT 2 3 ATOM 4 H FIL T 496 98.646 95.746 29.981 1.70 CONECT 3 4 ATOM 5 H FIL T 496 98.479 95.948 29.555 1.56 CONECT 4 5 ATOM 6 H FIL T 496 98.312 96.149 29.128 1.51 CONECT 5 6前面的EREST是指这些tunnel线不是真实的残基
如果我们想看到tunnel各部位radius的大小,那么首先选定tunnel,然后 点击Tools --> General --> Shell,在Shell里输入如下命令,注意不是命令行(这步的原理是bfactor一般被记成温度,是非必须的列,因此约定俗成的tunnel一般在这一列写入radius参数,将这个参数赋值给sphere球体大小的参数即可)
from chimerax.atomic import selected_atoms for a in selected_atoms(session): a.radius = a.bfactor再在style里改成sphere即可,然后可以把结构与tunnel一起存成新pdb
combine #1,#2 name #33.Align携带HETATM原子的pdb模型的方法
大家用过ChimeraX的一定知道,Matchmaker工具可以直接把模型拟合到一起,但如果其中带了非残基原子就无法被选中作为模板,就算实现了也会有个问题——那个tunnel不会跟着移动,chimerax有如下命令但我试了所有参数无法实现需求,大家可以试一下在Shell中以屏幕面为坐标系整体移动的方式
align #2 & @CA toAtoms #1 & @CA move nothing reporMatrix true那么可以用pymol实现旋转,再用ChimeraX渲染(这个软件最大的优势在于美观)
load tunnel_merged_1.pdb, 1 load tunnel_merged_2.pdb, 2 super 1 and name ca, 2 and name ca贴一张我做的蛋白的图,融合了两条相近的tunnel
注:如果出现了merge后的pdb在pymol中不显示但chimerax里显示的问题,可能是pdb编码格式被tunnel的写法扰乱了,重新保存新的tunnel再合并一下