上一期刚刚讲了一种处理部分模板结构质量差的解决方法，不知道各位看官有没有看过瘾。今天我们要将两个相对比较简单的方法进行模型质量的优化提升。
loop环精炼
loop环优化依赖于评分功能，并且优化算法仅适用于loop环，可以发现上期处理完以后蛋白273-283附近结果相对不理想，检查是否能够再次进行优化。（chimera显示的二级结构有一些不准确，我们在这里不做探讨，实际上273-283是完整的loop环）
在这里我们使用loopmodel，而不是automodel,具体如下：
<font color="red" >相信这么多天的推文，大家对同源建模已经有了一定的基础，所以后面我们就仅进行例子查看，并讲解主要要注意的地方，不再写一个伪代码，一个示例代码</font>
# Loop refinement of an existing model
from modeller import *
from modeller.automodel import *
log.verbose()
env = environ()
# directories for input atom files
env.io.atom_files_directory = './:../atom_files'
# Create a new class based on 'loopmodel' so that we can redefine
# select_loop_atoms (necessary)
class MyLoop(loopmodel):
# This routine picks the residues to be refined by loop modeling
def select_loop_atoms(self):
# 10 residue insertion
return selection(self.residue_range('273', '283'))
m = MyLoop(env,
inimodel='TvLDH-mult.pdb', # initial model of the target
sequence='TvLDH') # code of the target
m.loop.starting_model= 1 # index of the first loop model
m.loop.ending_model = 10 # index of the last loop model
m.loop.md_level = refine.very_fast # loop refinement method; this yields
# models quickly but of low quality;
# use refine.slow for better models
m.make()
File: loop_modeling/loop_refine.py
这个例子我们创建了一个Myloop的子类，用来继承loopmodel类，并返回要精炼的loop环范围，上面的示例为273至283的区域. refine.very_fast使得建模更快，但是相对质量不佳，建议使用注释后面的refine.slow。最后再进行DOPE打分。
基于配体建模
个人理解基于配体建模最佳的方式是目标蛋白与模板蛋白有相似或者相同的小分子配体，这样建出来的模型效果相对较好。如果我的理解不对，欢迎指出问题：
首先我们按照同源建模入门到精通（二）的方法对loop环优化后的配体进行模板比对，将带有配体的1emd比对上去，或者手动比对上去也可以。结果如下：
_aln.pos 10 20 30 40 50 60
TvLDH MSEAAHVLITGAAGQIGYILSHWIASGELYGDRQVYLHLLDIPPAMNRLTALTMELEDCAFPHLA
TvLDH_model MSEAAHVLITGAAGQIGYILSHWIASGELYGDRQVYLHLLDIPPAMNRLTALTMELEDCAFPHLA
1emd -----------------------------------------------------------------
_consrvd
_aln.pos 70 80 90 100 110 120 130
TvLDH GFVATTDPKAAFKDIDCAFLVASMPLKPGQVRADLISSNSVIFKNTGEYLSKWAKPSVKVLVIGN
TvLDH_model GFVATTDPKAAFKDIDCAFLVASMPLKPGQVRADLISSNSVIFKNTGEYLSKWAKPSVKVLVIGN
1emd ----------------------GVRRKPGMDRSDLFNVN--------------------------
_consrvd *** * ** *
_aln.pos 140 150 160 170 180 190
TvLDH PDNTNCEIAMLHAKNLKPENFSSLSMLDQNRAYYEVASKLGVDVKDVHDIIVWGNHGESMVADLT
TvLDH_model PDNTNCEIAMLHAKNLKPENFSSLSMLDQNRAYYEVASKLGVDVKDVHDIIVWGNHGESMVADLT
1emd -----------------------------------------------------------------
_consrvd
_aln.pos 200 210 220 230 240 250 260
TvLDH QATFTKEGKTQKVVDVLDHDYVFDTFFKKIGHRAWDILEHRGFTSAASPTKAAIQHMKAWLFGTA
TvLDH_model QATFTKEGKTQKVVDVLDHDYVFDTFFKKIGHRAWDILEHRGFTSAASPTKAAIQHMKAWLFGTA
1emd -----------------------------------------------------------------
_consrvd
_aln.pos 270 280 290 300 310 320
TvLDH PGEVLSMGIPVPEGNPYGIKPGVVFSFPCNVDKEGKIHVVEGFKVNDWLREKLDFTEKDLFHEKE
TvLDH_model PGEVLSMGIPVPEGNPYGIKPGVVFSFPCNVDKEGKIHVVEGFKVNDWLREKLDFTEKDLFHEKE
1emd -----------------------------------------------------------------
_consrvd
_aln.pos 330
TvLDH IALNHLAQGG/..
TvLDH_model IALNHLAQGG/--
1emd ----------/..
_consrvd
比对好后我们在目标序列与带有配体的序列后加"/"，表示序列结束，然后加"."表示匹配配体，该示例后面加了两个"."表示匹配两个配体。
最后进行建模:
from modeller import *
from modeller.automodel import *
class MyModel(automodel):
def special_restraints(self, aln):
rsr = self.restraints
for ids in (('NH1:161:A', 'O1A:336:B'),
('NH2:161:A', 'O1B:336:B'),
('NE2:186:A', 'O2:336:B')):
atoms = [self.atoms[i] for i in ids]
rsr.add(forms.upper_bound(group=physical.upper_distance,
feature=features.distance(*atoms),
mean=3.5, stdev=0.1))
env = environ()
env.io.hetatm = True
a = MyModel(env, alnfile='TvLDH-1emd_bs.ali',
knowns=('TvLDH_model','1emd'), sequence='TvLDH')
a.starting_model = 1
a.ending_model = 5
a.make()
File: ligand/model-multiple-hetero.py
建模脚本需要注意的是``env.io.hetatm = True`表示非标准”残基"也可以被读取，这样才能读取到配体。
分享到这里其实同源建模用的最普遍的内容已经讲完，接下来准备和大家一起学习迭代建模和建模更深的一些常用函数。这一系列就将先告一段落，等到时候研究rosetta的时候或者看到关于建模的一些好文章再继续撰写该类型。
在这里我们使用loopmodel，而不是automodel,具体如下：
<font color="red" >相信这么多天的推文，大家对同源建模已经有了一定的基础，所以后面我们就仅进行例子查看，并讲解主要要注意的地方，不再写一个伪代码，一个示例代码</font>
# Loop refinement of an existing model
from modeller import *
from modeller.automodel import *
log.verbose()
env = environ()
# directories for input atom files
env.io.atom_files_directory = './:../atom_files'
# Create a new class based on 'loopmodel' so that we can redefine
# select_loop_atoms (necessary)
class MyLoop(loopmodel):
# This routine picks the residues to be refined by loop modeling
def select_loop_atoms(self):
# 10 residue insertion
return selection(self.residue_range('273', '283'))
m = MyLoop(env,
inimodel='TvLDH-mult.pdb', # initial model of the target
sequence='TvLDH') # code of the target
m.loop.starting_model= 1 # index of the first loop model
m.loop.ending_model = 10 # index of the last loop model
m.loop.md_level = refine.very_fast # loop refinement method; this yields
# models quickly but of low quality;
# use refine.slow for better models
m.make()
File: loop_modeling/loop_refine.py
这个例子我们创建了一个Myloop的子类，用来继承loopmodel类，并返回要精炼的loop环范围，上面的示例为273至283的区域. refine.very_fast使得建模更快，但是相对质量不佳，建议使用注释后面的refine.slow。最后再进行DOPE打分。
基于配体建模
个人理解基于配体建模最佳的方式是目标蛋白与模板蛋白有相似或者相同的小分子配体，这样建出来的模型效果相对较好。如果我的理解不对，欢迎指出问题：
首先我们按照同源建模入门到精通（二）的方法对loop环优化后的配体进行模板比对，将带有配体的1emd比对上去，或者手动比对上去也可以。结果如下：
_aln.pos 10 20 30 40 50 60
TvLDH MSEAAHVLITGAAGQIGYILSHWIASGELYGDRQVYLHLLDIPPAMNRLTALTMELEDCAFPHLA
TvLDH_model MSEAAHVLITGAAGQIGYILSHWIASGELYGDRQVYLHLLDIPPAMNRLTALTMELEDCAFPHLA
1emd -----------------------------------------------------------------
_consrvd
_aln.pos 70 80 90 100 110 120 130
TvLDH GFVATTDPKAAFKDIDCAFLVASMPLKPGQVRADLISSNSVIFKNTGEYLSKWAKPSVKVLVIGN
TvLDH_model GFVATTDPKAAFKDIDCAFLVASMPLKPGQVRADLISSNSVIFKNTGEYLSKWAKPSVKVLVIGN
1emd ----------------------GVRRKPGMDRSDLFNVN--------------------------
_consrvd *** * ** *
_aln.pos 140 150 160 170 180 190
TvLDH PDNTNCEIAMLHAKNLKPENFSSLSMLDQNRAYYEVASKLGVDVKDVHDIIVWGNHGESMVADLT
TvLDH_model PDNTNCEIAMLHAKNLKPENFSSLSMLDQNRAYYEVASKLGVDVKDVHDIIVWGNHGESMVADLT
1emd -----------------------------------------------------------------
_consrvd
_aln.pos 200 210 220 230 240 250 260
TvLDH QATFTKEGKTQKVVDVLDHDYVFDTFFKKIGHRAWDILEHRGFTSAASPTKAAIQHMKAWLFGTA
TvLDH_model QATFTKEGKTQKVVDVLDHDYVFDTFFKKIGHRAWDILEHRGFTSAASPTKAAIQHMKAWLFGTA
1emd -----------------------------------------------------------------
_consrvd
_aln.pos 270 280 290 300 310 320
TvLDH PGEVLSMGIPVPEGNPYGIKPGVVFSFPCNVDKEGKIHVVEGFKVNDWLREKLDFTEKDLFHEKE
TvLDH_model PGEVLSMGIPVPEGNPYGIKPGVVFSFPCNVDKEGKIHVVEGFKVNDWLREKLDFTEKDLFHEKE
1emd -----------------------------------------------------------------
_consrvd
_aln.pos 330
TvLDH IALNHLAQGG/..
TvLDH_model IALNHLAQGG/--
1emd ----------/..
_consrvd
比对好后我们在目标序列与带有配体的序列后加"/"，表示序列结束，然后加"."表示匹配配体，该示例后面加了两个"."表示匹配两个配体。
最后进行建模:
from modeller import *
from modeller.automodel import *
class MyModel(automodel):
def special_restraints(self, aln):
rsr = self.restraints
for ids in (('NH1:161:A', 'O1A:336:B'),
('NH2:161:A', 'O1B:336:B'),
('NE2:186:A', 'O2:336:B')):
atoms = [self.atoms[i] for i in ids]
rsr.add(forms.upper_bound(group=physical.upper_distance,
feature=features.distance(*atoms),
mean=3.5, stdev=0.1))
env = environ()
env.io.hetatm = True
a = MyModel(env, alnfile='TvLDH-1emd_bs.ali',
knowns=('TvLDH_model','1emd'), sequence='TvLDH')
a.starting_model = 1
a.ending_model = 5
a.make()
File: ligand/model-multiple-hetero.py
建模脚本需要注意的是``env.io.hetatm = True`表示非标准”残基"也可以被读取，这样才能读取到配体。
分享到这里其实同源建模用的最普遍的内容已经讲完，接下来准备和大家一起学习迭代建模和建模更深的一些常用函数。这一系列就将先告一段落，等到时候研究rosetta的时候或者看到关于建模的一些好文章再继续撰写该类型。
往期教程：
#分子模拟#同源建模从入门到精通（一）
#分子模拟#同源建模从入门到精通（二）
#分子模拟#同源建模从入门到精通（三）
