abinitio软件 Rosetta从头预测蛋白结构的方法学研究

近年来，随着计算机技术的发展和人们对生物学的深入探究，预测蛋白质结构的方法也得到了快速的发展。其中基于abinitio软件的Rosetta从头预测蛋白结构的方法备受研究者的关注。该方法利用物理化学原理和大量的统计学习方法，通过蛋白质的氨基酸序列预测蛋白质的三维结构，实现了从蛋白质序列到蛋白质结构的全覆盖预测。本文将围绕该方法展开研究，探索其具体的算法原理和应用价值，并对其未来的发展进行展望。

Rosetta从头预测蛋白结构的方法学研究

总的来说，AbinitioRelax应用程序包括两个主要步骤：第一步，使用基于知识的“质心”评分函数（Abinitio）对构象空间进行粗粒度片段搜索；第二步（可选），使用Rosetta全原子力场（Relax）进行全原子细化。

1、输入文件（所有文件【input_files/】路径均有提供）

2、参数介绍（通过 -help 查看所有可接受的参数）

-in:file:native ./input_files/1elw.pdb              Native structure (optional)
(or -in:file:fasta ./input_files/1elwA.fasta)       Protein sequence in fasta format (required if native structure is not provided)
-in:file:frag3 ./input_files/aa1elwA03_05.200_v1_3  3-residue fragments (fragments file)
-in:file:frag9 ./input_files/aa1elwA09_05.200_v1_3  9-residue fragments (fragments file)
-database path/to/rosetta/main/database             Path to rosetta database
-abinitio:relax                                     Do a relax after abinitio ("abrelax" protocol), default=false.

-nstruct 1                                          Number of output structures
-out:file:silent 1elwA_silent.out                   Use silent file output, use filename after this flag, default=default.out
(or -out:pdb)                                       Use PDB file output, default=false
-out:path /my/path                                  Path where PDB output files will be written to, default '.'

-abinitio::rg_reweight 0.5 以下翻译可能有误

-abinitio::rsd_wt_helix 0.5 以下翻译可能有误

-abinitio::rsd_wt_loop 0.5 以下翻译可能有误

-use_filters true                               启用回转半径（RG）、接触顺序和折叠片过滤器。这个参数通过在筛选器失败时停止优化来节省计算。对于某些序列，
                                                可能会有很大一部分模型无法通过筛选。无法通过筛选的模型的名称以F_开头。
-psipred_ss2 ./input_files/1elwA.psipred_ss2    psipred_ss2 二级结构预测文件，当使用了 -use_filters 时为必填。
-abinitio::increase_cycles 10                   从头计算时每阶段所增加的循环次数。
-abinitio::rg_reweight 0.5                      使用此值作为 RG 的权重重新计算。
-abinitio::rsd_wt_helix 0.5                     对于“helix”中的残基，使用此值作为 env,pair,cb 的权重。
-abinitio::rsd_wt_loop 0.5                      对于“loop”中的残基，使用此值作为 env,pair,cb 的权重。
-relax::fast                                    使用性能无显著损失的快速relax算法。
-kill_hairpins ./input_files/1elwA.psipred_ss2  由于β或α-β蛋白会频繁采样到发夹结构，通过这个参数可以 kill 发夹。

-constant_seed 以下翻译可能有误

-constant_seed                                  使用一个常量“1111111”作为随机种子
-jran 1234567                                   指定一个随机种子，对于所有job应该具有唯一性。(要求指定 -constant_seed)
-seed_offset 10                                 这个值将被加到随机种子上。当使用时间作为种子的集群计算时会非常有用。如果
                                                使用Condor作为调度器，则参数值也可以填写为“$(Process)”

$ROSETTA/bin/AbinitioRelax.linuxgccrelease \
        -database /path/to/rosetta/main/database \
        -in:file:fasta ./input_files/1elwA.fasta \
        -in:file:native ./input_files/1elw.pdb \
        -in:file:frag3 ./input_files/aa1elwA03_05.200_v1_3 \
        -in:file:frag9 ./input_files/aa1elwA09_05.200_v1_3 \
        -abinitio:relax \
        -relax:fast \
        -abinitio::increase_cycles 10 \
        -abinitio::rg_reweight 0.5 \
        -abinitio::rsd_wt_helix 0.5 \
        -abinitio::rsd_wt_loop 0.5 \
        -use_filters true \
        -psipred_ss2 ./input_files/1elwA.psipred_ss2 \
        -kill_hairpins ./input_files/1elwA.psipred_ss2 \
    -out:file:silent 1elwA_silent.out \
        -nstruct 10

3、输出结果的处理

生成的主要文件有以下几种

output_files/S_00000001.pdb,
output_files/score.fsc, 
output_files/default.out(silent output file).

3.1、从silent中解出pdb（若输入参数未指明 -out:pdb，则最终输出的是一个 silent 文件）

生成的silent文件后缀为“.out”，里面包含了生成的pdb模型，以及打分数据。使用打分应用 score.linuxgccrelease 可以从 silent 解析出pdb文件。

3.2、使用聚类程序对结果进行聚类

cluster.linuxgccrelease 可对silent或pdb格式的模型进行聚类。当结构集过大时，建议使用优化过的聚类程序——Calibur（参考 http://www.biomedcentral.com/1471-2105/11/25和http://sourceforge.net/projects/calibur/）。

1、AbinitioRelax的第二步（Relax）花费的时间略长，对于demo（117个残基），普通计算机运行大约需要8分钟。

2、由于全原子能量函数对原子间相互作用力非常敏感，加之理论上给定序列的蛋白质构象空间巨大。因此应尽可能多的提高采样点才更有可能覆盖到天然态构象，受算力的约束，一般正式预测时，需要生成50000到100000个模型，才能达到较为理想的结果。

3、AbinitioRelax能够很好地预测长度小于100个残基的小单体蛋白，也有一定几率精准预测达到150个残基的蛋白质，但更大的蛋白质因构象空间骤增，需要更多的算力。以下参考资料提供了与采样问题相关的信息：

Bradley P, Misura KM, Baker D (2005). Toward high-resolution de novo structure prediction for small proteins. Science 309, 1868-71.
Kim DE, Blum B, Bradley P, Baker D (2009). Sampling bottlenecks in de novo protein structure prediction. J Mol Biol 393, 249-60.

4、Rosetta可以很好地利用你的实验数据进行约束性预测，为此你需要提供参数 -cst_file 和 -cst_weight 所需的文件，文件说明可参考 contraint files。

5、使用同源序列运行AbinitioRelax是很有用的，例如可通过PSI-BLAST搜索（NCBI nr数据库）或Pfam获取同源序列。使用类似Jaliew的序列比对查看器可以帮助选择要运行的最佳同源序列集，还可以帮助选择模型。

6、为了增加采样到正确拓扑的机会，也可以运行一组非同源序列，因为选一组同源序列一定程度上意味着选了一个较窄的采样域，因此降低了收敛到天然态的概率。通常我们寻找一组不同的同源序列（最多10个），它们在保守位置和缺失上存在差异，这些序列可能代表一个截断环或无序区域（truncated loop or disordered region）。

7、使用诸如Disopred或metaPrDOS之类的公共程序（publicly available programs）识别和修剪无序的端序列（disordered termini）也很重要。

8、信号序列也应该使用像SignalP这样的公共程序来识别和裁剪。该方案不适用于膜蛋白。如果使用TMHMM等程序预测跨膜螺旋，请参阅Membrane ab initio 程序。

官方建议生成20000到30000个目标序列模型和多至10个同系物模型，然后使用聚类应用程序或Calibur来识别最常采样的构象（most frequently sampled conformations）。在一般情况下，按大小排列的前5-10个集群中，至少有一个存在RMSD最低的模型（相对本征结构）。

在理想的情况下，你的序列应该是由一些同源序列识别搜索工具得到，如PSI-BLAST。序列比对在模型选择中非常有用。例如，保守的疏水位置很可能代表蛋白质的核心，因此在这些位置暴露有侧链的模型可能被丢弃。同样的逻辑也适用于最有可能出现在表面的保守的极性位点（conserved polar positions）。此外，保守的半胱氨酸对可能代表二硫化物。具有序列比对功能的Jalview和模型可视化的PyMOL等工具对模型选择非常有用。

分数-RMSD散点图（plots）可能有助于识别目标序列和同源序列向本征结构的收敛性。例如，当使用score.linuxgccrelease打分程序重新排序模型时，最低评分模型可作为 -in:file:native 的值。从得到的score文件画出的score-RMSD图可能会呈现出一个朝向最低得分模型的收敛（能量漏斗）。如果存在能量漏斗，得分最低的模型有更大的机会接近本征结构。有关创建分数-RMSD散点图（Score vs RMSD plots）的讨论，请参见https://www.rosettacommons.org/node/3813 和 https://www.rosettacommons.org/content/how-make-benchmark。PyRosetta工具包GUI也可以用于Score-RMSD输出，并在运行期间发现最低能量结构集。

集群中的得分最低的模型以及在PDB中表示拓扑的模型也更有可能是正确的。可以用结构-结构比对工具，如Dali或Mammoth，来搜索PDB数据库。

　　综上所述，基于Abinitio软件Rosetta从头预测蛋白结构的方法学是一项非常重要的研究工作，可以帮助我们更深入地了解蛋白质的结构及其功能。通过不断地优化算法和模型，这一方法将在未来的生物医学研究中发挥越来越重要的作用。