来一场蛋白和小分子的风花雪月

分子对接(Molecular Docking)理论

所谓分子对接就是两个或多个分子之间通过几何匹配和能量匹配相互识别找到最佳匹配模式的过程。分子对接对酶学研究和药物设计中有重要的应用意义。

分子对接计算是在受体活性位点区域通过空间结构互补和能量最小化原则来搜寻配体与受体是否能产生相互作用以及它们之间的最佳结合模式。分子对接的思想起源于Fisher E的”钥匙和锁模型”，主要强调的是空间形状的匹配。但配体和受体的识别要比这个模型更加复杂。首先，配体和受体在对接过程中会由于相互适应而产生构象的变化。其次，分子对接还要求能量匹配，对接过程中结合自由能的变化决定了两个分子是否能够结合以及结合的强度。

1958年D.E.Koshland提出分子识别过程中的诱导契合概念，受体分子活性中心的结构原本并非与底物完全吻合，但其是柔软和可塑的。当配体与受体相遇时，可诱导受体构象发生相应的变化，从而便于他们的结合进而引起相应的反应。

分子对接方法根据不同的简化程度分为三类：刚性对接、半柔性对接和柔性对接。刚性对接指在对接过程中，受体和配体的构象不发生变化，适合研究比较大的体系如蛋白-蛋白之间以及蛋白-核酸之间，计算简单，主要考虑对象之间的契合程度。半柔性对接常用于小分子和大分子的对接，在对接过程中，小分子的构象可以在一定范围内变化，但大分子是刚性的。这样既可以在一定程度上考察柔性的影响，又能保持较高的计算效率。在药物设计和虚拟筛选过程中一般采用半柔性的分子对接方法。柔性对接方法一般用于精确研究分子之间的识别情况，由于允许对接体系的构象变化，可以提高对接准确性但耗时较长。

分子对接的目的是找到底物分子和受体分子最佳结合位置及其结合强度，最终可以获得配体和受体的结合构象，但这样的构象可以有很多，一般认为自由能最小的构象存在的概率最高。搜寻最佳构象就要用到构象搜索方法，常用的有系统搜索法和非系统搜索法。系统搜索法通过改变每个扭转角评估所有可能的结合构象，进而选取能量最低的。这一方法计算量非常大。因此通常使用非系统搜索法来寻找能量较低构象，常用方法有：分子动力学方法、随机搜索、遗传算法、距离几何算法等。随机搜索又包括完全随机算法、蒙特卡罗法和模拟退火法等。

AutoDock Vina是The Scripp Research Institute的Olson科研小组开发的分子对接软件包，使用拉马克遗传算法提高效率。软件把遗传算法和局部搜索结合在一起，遗传算法用于全局搜索，而局部搜索用于能量优化。为了加快计算速度，AutoDock Vina采用格点(grid)计算。首先在受体活性氨基酸附近划定一个长方体区域作为搜索空间，扫描不同类型的原子计算格点能量，在搜索空间内，调整配体的构象、位置和方向，进而评分、排序获得能量最低的构象作为输出结果。

对范德华相互作用的计算：每个格点上保存的范德华能量的值的数目与要对接的配体上的原子类型数目相同。如果一个配体中含有C、H、O三种原子类型，那么每个葛店需要用单个探针原子与来计算其与受体之间的范德华相互作用值。当配体与受体进行分子对接时，配体中某个原子和受体之间的相互作用能通过周围8个格点上的这种原子类型为探针的格点值用内插法得到。

静电相互作用的计算采用静电势格点。当配体与受体对接时，某个原子和受体之间的静电相互作用能通过周围格点上静电势以及原子上的部分电荷计算得到。

蛋白和小分子可视化

例子文件是一个分辨率为2艾的X-射线衍射晶体结构(PDB ID: 1HSG)，其为HIV-1蛋白酶与药物茚地那韦(indinavir)结合在一起的构象。软件PyMOL用来观察HIV-蛋白酶、结合位点和药物分子的结构。

显示蛋白结构并做样式处理

1. 下载HIV-1蛋白酶的PDB结构(https://files.rcsb.org/download/1HSG.pdb)，存储到一个不含中文和空格的目录下。

2. 启动PyMOL，依次点选File-Open-1hsg.pdb导入PDB文件，会看到如下界面

3. 首先在右侧的对象控制面板，依次点选行all的H列-Hide: everything(如左图所示)，然后浩瀚无际的没有月亮的夜空出现在我们面前。

4. 在右侧的对象控制面板，依次点选行1hsg的S列-Show: cartoon，然后点选C列-By chain显色，这时可以看到如右图所示的同源二聚体。

显示与蛋白结合的小分子化合物和水分子

1. 从蛋白结构的PDB文件(PDB文件格式解析见后面)或PDB官网的信息(如下图所示)可以看到，1hsg结构中包含配体药物indinavir，其残基的名字为MK1。

2. 在PyMOL的命令行处输入PyMOL> select indinavir, resn MK1，回车，会看到 如下画面变化。

3. 在右侧的对象控制面板，依次点选行indinavir的S列-Show stick，再点选C列选择一种不同的颜色。在屏幕无图处点击鼠标，取消小分子药物的选择状态。这时可以清晰的看到小分子的结构和空间位置(如下左图)，随意拖动鼠标旋转或放大查看药物分子与蛋白的结合方式。PyMOL鼠标操作：按住左键移动旋转，按住右键移动放大，按住中键移动，观察结合位点所在的位置；滚动滚轮调节景深，化学结构会以 溶解形式出现。

4. 显示水分子。水分子的残基名字为HOH，运行命令PyMOL> select H2O, resn HOH调出水分子。然后点选S-Show spheres, C-red。再运行PyMOL> set sphere_scale, 0.2设置水球的大小。

5. 如果要存储结果，则在命令行输入png E:/docking/1shg.png保存当前结果。

准备docking需要的受体(蛋白)和配体(化合物)

Docking算法需要每个原子带有电荷并且需要标记原子的属性。这些信息通常未包含在PDB文件中。我们需要在对蛋白和小分子的PDB文件预处理，生成PDBQT文件同时包含以上信息和PDB文件中的原子坐标信息。进一步地对于“柔性配体docking”，我们还需要定义配体的柔性部分和刚性部分。所有这些都可以通过软件AutoDock Tools (adt)来完成。

Docking algorithms require each atom to have a charge and an atom type that describes its properties. However, the PDB structure lacks these. So, we have to prep the protein and ligand files to include these values along with the atomic coordinates. Furthermore, for flexible ligand docking, we should also define ligand bonds that are rotatable. All this will be done in a tool called AutoDock Tools (adt).

准备受体蛋白

1. PDB文件(1hsg.pdb)中包含了蛋白、配体和水分子；首先提取出蛋白的坐标，即以关键字ATOM和TER开头的行 (具体解释和例子见后面PDB格式解析)存储到文件1hsg_prot.pdb。

• 在windows下，我们可以手动选择，或者利用Excel的筛选功能。

• 在Linux下，使用命令egrep "^(ATOM|TER)" 1hsg.pdb >1hsg_prot.pdb。

• windows直接双击图标就可

• Linux可以使用命令adt &

• ADT中按住左键拖动旋转分子结构；点击中键滚动缩放；按住右键移动晶体位置。

• 1hsg_prot.pdbqt为只加了极性氢的结果

• 1hsg_prot_all_h.pdbqt为加了所有氢的结果

  这两个文件只在原子电量部分有所不同，经测试发现这两种处理对docking的结果没有影响，最后输出的日志文件和结果文件相同。

准备配体

1. 与蛋白结构类似，配体的结构也缺少氢原子，我们需要添加氢原子并且定义哪些键是可以旋转的以用于柔性docking。

2. 从PDB结构中提取配体的原子位置。indinavir的配体残基名字为MK1，以HETATM开头的行表示非核心多聚体的成分 (heteroatoms)(具体见PDB文件格式解释)。

• Linux系统下，运行grep "^HETATM.*MK1" 1hsg.pdb >indinavir.pdb

• Windows系统下，直接拷贝到文件indinavir.pdb

• indinavir.pdbqt为只加了极性氢的结果

• indinavir_all_h.pdbqt为加了所以氢的结果

准备docking配置文件

docking配置文件包含了输入的受体(蛋白)、配体(化合物)和搜索参数的信息，为一个文本文件，名字任意，可以为conf.txt，内容如下

receptor = 1hsg_prot.pdbqt
ligand   = indinavir.pdbqt
num_modes = 50
out = dockingResult.pdbqt
log = docking.log
center_x = 16
center_y = 25
center_z = 4
size_x = 30
size_y = 30
size_z = 30
seed = 2009

• receptor和ligand为输入文件的名字，与conf.txt在同一目录下; out为输出文件的名字；log为输出日志文件的名字。

• centerhe和size定义搜索空间的位置和大小。

• num_modes设置最多显示的结合模型，鉴于只输出符合能量值要求的结果，最后输出的结合模型数量可能少于这一数值。

• seed设置随机数生成的种子，可以为任意整数。如果想重现之前的分析结果就需要使用相同的seed。

Docking 小分子化合物indinavir到HIV-1蛋白酶

1. 使用AutoDock Vina执行docking预测。

• 在windows命令行提示符或linux终端下运行命令 vina --config conf.txt，大概需要几分钟时间。

Detected 4 CPUs
Reading input ... done.
Setting up the scoring function ... done.
Analyzing the binding site ... done.
Using random seed: 2009
Performing search ... done.
Refining results ... done.

mode |   affinity | dist from best mode
     | (kcal/mol) | rmsd l.b.| rmsd u.b.
-----+------------+----------+----------
   1        -11.5      0.000      0.000
   2        -10.6      1.425      4.304
   3        -10.4      2.042     10.990
   4        -10.3      2.034     10.326
   5        -10.2      2.517      4.774
   6        -10.1      1.933     10.911
   7         -9.9      2.176     10.884
   8         -9.8      1.794      3.600
   9         -9.6      1.981     10.865
  10         -9.5      2.431     10.943
  11         -9.3      2.417     10.370
  12         -8.9      2.404     10.285
  13         -8.8      4.058     10.904
  14         -8.7      5.574     11.291
  15         -8.7      4.441      8.312
  16         -8.6      5.659      8.929
  17         -8.6      4.404      8.275
  18         -8.5      5.630      8.900
Writing output ... done.

• dockingResult.pdbqt: 包含所有docking的模式，通常第一个为结合最好的构象，但如果前几个能量值相差不大时也有例外。

• docking.log: 日志文件，包含结合能量值(第一列，越低越稳定，默认由低到高排序，所以第一个为最好的构象)、每个构象与第一个构象的距离、每个构象与第一个构象的差别。

用PyMOL可视化docking结果。

打开PyMOL，依次点选File-Open文件类型选择All Files-选取结果dockingResult.pdbqt文件、原始蛋白和配体的pdb文件、原教程的pdbqt文件。

• dockingResult.pdbqt: 增加非极性氢的docking结果

• dockingResultAllH.pdbqt: 增加所有氢的docking结果

• original_tutorial_result.pdbqt：原教程中的docking结果