Modeller 学习记录（二）
前面讲完了安装及运行，那基本上就到了实际操作，建立自己要做的蛋白质模型了。开始建模型前我做的准备工作：1. 到salilab网站的tutorial去下载basic modeling example， advanced modeling example。因为我不会python，所以我只能依葫芦画瓢用这些例子里的脚本来做自己的工作。
2. 到uniprot上去找自己要做的蛋白质序列，把需要建模的那段序列的ali文件建好。怎么建？很简单，以葫芦画瓢，下面的例子是example里的ali
>P1;TvLDH
sequence:TvLDH:::::::0.00: 0.00
MSEAAHVLITGAAGQIGYILSHWIASGELYGDRQVYLHLLDIPPAMNRLTALTMELEDCAFPHLAGFVATTDPKA
AFKDIDCAFLVASMPLKPGQVRADLISSNSVIFKNTGEYLSKWAKPSVKVLVIGNPDNTNCEIAMLHAKNLKPEN
FSSLSMLDQNRAYYEVASKLGVDVKDVHDIIVWGNHGESMVADLTQATFTKEGKTQKVVDVLDHDYVFDTFFKKI
GHRAWDILEHRGFTSAASPTKAAIQHMKAWLFGTAPGEVLSMGIPVPEGNPYGIKPGVVFSFPCNVDKEGKIHVV
EGFKVNDWLREKLDFTEKDLFHEKEIALNHLAQGG*
用你的蛋白质序列替换TvLDH的序列，当然也需要用你的蛋白质名称替换TvLDH。然后存成ali后缀格式文件就好。
3. 准备好最新的pdb数据库文件，pir和bin。前面讲过了，这两个文件怎么生成。
4. 建一个以蛋白质文件命名的目录，将所有相关的文件都放这个目录，然后所有的工作都在这个目录里进行。（*.ali, *.pir, *.bin, build_profile.py,及从basic modeling example目录里拷贝过来的各种脚本）
建模开始：
1. build_profile.py的修改，下面贴出原始的build_profile.py文件，从basic example里来的。改动几个地方，就可以用这个来查找相似蛋白质信息。
from modeller import *
log.verbose()
env = environ()
#-- Prepare the input files
#-- Read in the sequence database
sdb = sequence_db(env)
sdb.read(seq_database_file='pdb_95.pir', seq_database_format='PIR',
chains_list='ALL', minmax_db_seq_len=(30, 4000), clean_sequences=True)
－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
我用了pdball.pir来替换这个pdb_95.pir
－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
#-- Write the sequence database in binary form
sdb.write(seq_database_file='pdb_95.bin', seq_database_format='BINARY',
chains_list='ALL')
#-- Now, read in the binary database
sdb.read(seq_database_file='pdb_95.bin', seq_database_format='BINARY',
chains_list='ALL')
－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
同上，我用了pdball.bin来替换这个pdb_95.bin
－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
#-- Read in the target sequence/alignment
aln = alignment(env)
aln.append(file='TvLDH.ali', alignment_format='PIR', align_codes='ALL')
－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
用自己的蛋白质ali文件替代TvLDH.ali
－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－
#-- 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 the alignment file
aln.write(file='build_profile.ali', alignment_format='PIR')
最后在命令行下运行mod9.14 build_profile.py
完成之后会在目录里生成三个文件：
build_profile.ali
build_profile.log
build_profile.prf
我只用log文件去抓蛋白质相似性信息。其它两个文件不怎么用。抓出下面这部分的信息就可以，主要的目的是为了挑选出相似性高的蛋白质来做进一步的分析。
Dynamically allocated memory at amaxprofile : 7496268 7320.574 7.149
> 1a5z 1 190 6900 312 335 28.14 0.85E-08 2 164 75 242 63 229
> 1b8pA 1 951 29400 327 335 42.50 0.0 3 316 7 331 6 325
> 1bdmA 1 1058 30700 318 335 44.52 0.0 4 309 1 325 1 310
> 1t2dA 1 1717 5200 315 335 24.70 0.67E-04 5 238 5 256 4 250
> 1civA 1 1784 23950 374 335 34.66 0.0 6 325 6 334 33 358
> 2cmd 1 1858 5900 312 335 27.24 0.17E-05 7 289 7 320 3 303
> 1o6zA 1 2060 5800 303 335 26.32 0.27E-05 8 278 7 320 3 287
> 1ur5A 1 3954 4500 299 335 30.67 0.25E-02 9 158 13 191 9 171
> 1guzA 1 3955 7100 305 335 25.27 0.29E-08 10 265 13 301 8 280
> 1gv0A 1 3962 5350 301 335 25.53 0.29E-04 11 274 13 323 8 289
> 1hyeA 1 4405 6800 307 335 29.28 0.14E-07 12 173 7 191 3 183
> 1i0zA 1 4433 5650 332 335 24.64 0.67E-05 13 207 85 300 94 304
> 1i10A 1 4435 5600 331 335 26.21 0.87E-05 14 196 85 295 93 298
> 1ldnA 1 5824 5000 316 335 25.76 0.19E-03 15 214 78 298 73 301
> 6ldh 1 5831 4600 329 335 23.08 0.17E-02 16 244 47 301 56 302
> 2ldx 1 5833 5400 331 335 25.83 0.25E-04 17 227 66 306 67 306
> 5ldh 1 5835 5800 333 335 25.59 0.30E-05 18 207 85 300 94 304
> 9ldtA 1 5836 6000 331 335 25.94 0.11E-05 19 207 85 301 93 304
> 1llc 1 5896 5000 321 335 25.82 0.20E-03 20 164 64 239 53 234
> 1lldA 1 5897 6650 313 335 30.67 0.32E-07 21 216 13 242 9 233
> 5mdhA 1 6084 32900 333 335 44.41 0.0 22 328 2 332 1 331
> 7mdhA 1 6085 22600 351 335 33.74 0.0 23 325 6 334 14 339
> 1mldA 1 6188 6100 313 335 25.93 0.58E-06 24 183 5 198 1 189
> 1oc4A 1 6736 5450 315 335 27.87 0.18E-04 25 174 5 191 4 186
> 1ojuA 1 6799 5700 294 335 27.98 0.44E-05 26 218 78 320 68 285
> 1pzgA 1 7253 6350 327 335 30.00 0.16E-06 27 114 74 191 71 190
> 1smkA 1 7850 8750 313 335 33.85 0.0 28 188 7 202 4 198
> 1sovA 1 7980 4700 316 335 26.59 0.94E-03 29 160 81 256 76 248
> 1y6jA 1 10912 5750 289 335 32.73 0.33E-05 30 109 77 191 58 167
把这部分信息整理，选高相似度及E值尽可能小的蛋白质来做下一步的compare分析。
页:
[1]