|
|
【答案】应助回帖
★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★
感谢参与，应助指数 +1
西门吹雪170: MolEPI+1 2017-07-28 07:44:12
西门吹雪170: 金币+6, 鼓励热心回帖交流 2017-07-28 07:44:45
cxl8400501: 金币+20 2017-07-29 16:47:50
正好熟悉这个，给你一个方法：利用MODLLER同源建模。
前期的准备工作（这里介绍并没有结构的蛋白序列的从头建立结构的过程）
首先：把你需要建立结构的蛋白序列写成modeller能够识别的文件格式
>P1;TvLDH
sequence:::::::::
MSEAAHVLITGAAGQIGYILSHWIASGELYGDRQVYLHLLDIPPAMNRLTALTMELEDCAFPHLAGFVATTDPKAAFKDIDCAFLVASMPLKPGQVRADLISSNSVIFKNTGEYLSKWAKPSVKVLVIGNPDNTNCEIAMLHAKNLKPENFSSLSMLDQNRAYYEVASKLGVDVKDVHDIIVWGNHGESMVADLTQATFTKEGKTQKVVDVLDHDYVFDTFFKKIGHRAWDILEHRGFTSAASPTKAAIQHMKAWLFGTAPGEVLSMGIPVPEGNPYGIKPGVVFSFPCNVDKEGKIHVVEGFKVNDWLREKLDFTEKDLFHEKEIALNHLAQGG*
第二行后面的字段基本没什么作用,只要把第一个字段sequence写上即可,注意，这一行一共10个字段用9个分号隔开。
第三行为将你要建立结构的序列的氨基酸写成单字母的形式，最后的 " * “ 符号表示序列的结束.
其次：
将PDB库中的序列准备好（modeller会经常更新PDB库中的序列,这些序列会经常更新,你可以从modeller的官网上下载）下面是将PDB库中的序列转化为二进制的文件格式,可以加快搜索具有高度相似的序列的速度. 转换的脚本如下所示：
from modeller import *
log.verbose()
env = environ()
sdb = sequence_db(env)
sdb.convert(seq_database_file='pdb_95.pir', seq_database_format='PIR', chains_list='ALL', minmax_db_seq_len=(30, 4000), clean_sequences=True, outfile='pdb_95.hdf5')
第四行 sequence_db()对象用于创建一个存放大量的蛋白序列的库对象。
第五行调用对象的convert()方法，以PIR格式的方式读入pdb_95.pir,去除残基序列的数目小于30并且大于4000的序列，用minmx_db_seq_len参数指出，去除不标准的残基用clean_sequences指出，最终文件以二进制的形式写到pdb_95.hdf5文件中。
第三：利用上一步中生成的二进制文件进行搜索，找到和未知序列相似的结构序列。
from modeller import *
log.verbose()
env = environ()
# Read in the sequence database in binary format
sdb = sequence_db(env, seq_database_file='pdb_95.hdf5', seq_database_format='BINARY', chains_list='ALL')
# Read in the target sequence in PIR alignment format
aln = alignment(env)
aln.append(file='TvLDH.ali', alignment_format='PIR', align_codes='ALL')
# Convert the input sequence "alignment" into profile format
prf = aln.to_profile()
# Scan sequence database to pick up homologous sequences
prf.build(sdb, matrix_offset=-450, rr_file='${LIB}/blosum62.sim.mat', gap_penalties_1d=(-500, -50), n_prof_iterations=1, check_profile=False, max_aln_evalue=0.01)
# Write out the profile in text format
prf.write(file='build_profile.prf', profile_format='TEXT')
# Convert the profile back to alignment format
aln = prf.to_alignment()
#- Write out a PIR alignment file
aln.write(file='build_profile.ali', alignment_format='PIR‘）
*****在build_profile.prf文件中，最重要的列为第2,10,11，12列，第二列指出了PDB库中的晶体结构的代码，
第11列指出了待建结构与晶体结构的序列相似度的百分数，第12列指明了叠加之后的e值。
第四：从build_profile.prf这个输出文件中找到和待建序列的相似度高的模版，评价模版的结构和序列的相似性。
from modeller import *
env = environ()
env.io.atom_files_directory = ['.', '../atom_files']
# Make a simple 1:1 alignment of 7 template structures
aln = alignment(env)for (pdb, chain) in (('1b8p', 'A'), ('1y7t', 'A'), ('1civ', 'A'), ('5mdh', 'A'), ('7mdh', 'A'), ('3d5t', 'A'), ('1smk', 'A')):
m = model(env, file=pdb, model_segment=('FIRST:'+chain, 'LAST:'+chain))
aln.append_model(m, atom_files=pdb, align_codes=pdb+chain)
# Sequence alignment
aln.malign()
# Structure alignment
aln.malign3d()
# Report details of the sequence/structure alignment
aln.compare_structures()
aln.id_table(matrix_file='family.mat')
env.dendrogram(matrix_file='family.mat', cluster_cut=-1.0)
对角线中的内容为每个结构的残基数目，上三角为相同的残基的数目，下三角为两个结构的序列相似度的百分数。
我们选择了1y7t作为模版结构，因为这个结构的结晶分辨率为1.6埃以及它与待建结构的序列的相似度为45%。
第五：将待建结构与模版结构进行align, MODELLER利用这些align的结果来提取出构建结构时必要的限制。
最好的方法是利用align2d()的方法，尽管align2d()方法是基于一种动态的算法，但是它不同于标准的序列-序列的align
方法，因为这种方法在进行align的时候会考虑到模版的结构信息。
from modeller import *
env = environ()
env.io.atom_files_directory = ['.', '../atom_files']
aln = alignment(env)
# Read in the 1y7t template structure and add to alignment
mdl = model(env, file='1y7t', model_segment=('FIRST:A','LAST:A'))
aln.append_model(mdl, align_codes='1y7tA', atom_files='1y7t')
# Add in the TvLDH sequence
aln.append(file='../step1_search/TvLDH.ali', align_codes='TvLDH')
# Sequence/structure alignment
aln.align2d(max_gap_length=40)
# Write out resulting alignment in both PIR and PAP formats
aln.write(file='TvLDH-1y7tA.ali', alignment_format='PIR')
aln.write(file='TvLDH-1y7tA.pap', alignment_format='PAP')
第六：创建比较模型
可以使用MODELLER中的automodel类，只需要告诉automodel alignment文件，使用哪个代码以及要建立多少个模型。
from modeller import *
from modeller.automodel import *
env = environ()
env.io.atom_files_directory = ['.', '../atom_files']
a = automodel(env,alnfile='../step3_align/TvLDH1y7tA.ali', knowns='1y7tA',sequence='TvLDH',assess_methods=assess.DOPE)
a.starting_model = 1
a.ending_model = 5
a.make()
第七：评价所建的模型
从上面的步骤中找到最好的模型，评价的方法有很多种，最简单的方法是molecule PDF (molpdf)方法，标准的MODELLER评分函数。这个评分结果会在每个PDB开头的REMARK中加入。其他的评价方法包括DOPE和GA341.
从上面产生的log文件中可以发现，第一个PDB文件模型是最好的，因为DOPE的值最小。
第八：Fit models into cryo-EM maps
通过将TvLDH的 cryo-EM数据导入来使产生的蛋白结构与这些数据相吻合，通过使用MODELLER中的Mod-EM方法。
from modeller import *
log.verbose()
env = environ()
struct='../step4_model/TvLDH.B99990001.pdb'
map='TvLDH.10A.mrc'
resolution=10.0
box_size=48
apix=1.88
x=-26.742; y=-9.5205; z=-10.375 #origin
steps=20
# Read in cryo-EM density map
den = density(env, file=map, em_density_format='MRC',
voxel_size=apix, resolution=resolution, em_map_size=box_size,
density_type='GAUSS', px=x,py=y,pz=z)
# Fit the PDB file into the map by MC simulated annealing
den.grid_search(em_density_format='MRC', num_structures=1, em_pdb_name=struct, chains_num=[1],
start_type='CENTER', number_of_steps=steps, angular_step_size=30., temperature=100.,
best_docked_models=1, translate_type='RANDOM', em_fit_output_file='modem.log') |
|
回复此楼
cxl8400501
304材料求调剂
+4