import sys
sys.path.append('/usr/lib/modeller9.12/modlib/')
sys.path.append('/usr/lib/modeller9.12/lib/x86_64-intel8/python2.5/')
import modeller
import _modeller
import modeller.automodel
#impo modeller.automodel import *
env=modeller.environ()
env.io.hetatm=True
MODELLER 9.12, 2013/06/12, r9480 PROTEIN STRUCTURE MODELLING BY SATISFACTION OF SPATIAL RESTRAINTS Copyright(c) 1989-2013 Andrej Sali All Rights Reserved Written by A. Sali with help from B. Webb, M.S. Madhusudhan, M-Y. Shen, M.A. Marti-Renom, N. Eswar, F. Alber, M. Topf, B. Oliva, A. Fiser, R. Sanchez, B. Yerkovich, A. Badretdinov, F. Melo, J.P. Overington, E. Feyfant University of California, San Francisco, USA Rockefeller University, New York, USA Harvard University, Cambridge, USA Imperial Cancer Research Fund, London, UK Birkbeck College, University of London, London, UK Kind, OS, HostName, Kernel, Processor: 4, Linux srv 3.13.0-88-generic x86_64 Date and time of compilation : 2013/06/11 21:33:06 MODELLER executable type : x86_64-intel8 Job starting time (YY/MM/DD HH:MM:SS): 2016/11/29 19:02:54
cd /tmp/mod-teach/
/tmp/mod-teach
%%bash
wget http://www.pdb.org/pdb/files/1lmp.pdb -O 1lmp.pdb
--2016-11-29 18:36:23-- http://www.pdb.org/pdb/files/1lmp.pdb Resolving www.pdb.org (www.pdb.org)... 132.249.231.77 Connecting to www.pdb.org (www.pdb.org)|132.249.231.77|:80... connected. HTTP request sent, awaiting response... 301 Moved Permanently Location: http://www.rcsb.org/pdb/files/1lmp.pdb [following] --2016-11-29 18:36:23-- http://www.rcsb.org/pdb/files/1lmp.pdb Resolving www.rcsb.org (www.rcsb.org)... 132.249.231.10 Connecting to www.rcsb.org (www.rcsb.org)|132.249.231.10|:80... connected. HTTP request sent, awaiting response... 301 Moved Permanently Location: http://files.rcsb.org/view/1lmp.pdb [following] --2016-11-29 18:36:23-- http://files.rcsb.org/view/1lmp.pdb Resolving files.rcsb.org (files.rcsb.org)... 132.249.213.140 Connecting to files.rcsb.org (files.rcsb.org)|132.249.213.140|:80... connected. HTTP request sent, awaiting response... 200 OK Length: unspecified [text/plain] Saving to: ‘1lmp.pdb’ 0K .......... .......... .......... .......... .......... 127K 50K .......... .......... .......... .......... .......... 254K 100K .......... .......... ........ 279M=0.6s 2016-11-29 18:36:25 (217 KB/s) - ‘1lmp.pdb’ saved [131301]
%%bash
wget http://www.uniprot.org/uniprot/P61626.fasta
--2016-11-29 18:40:00-- http://www.uniprot.org/uniprot/P61626.fasta Resolving www.uniprot.org (www.uniprot.org)... 128.175.240.211, 193.62.192.81, 193.62.193.81 Connecting to www.uniprot.org (www.uniprot.org)|128.175.240.211|:80... connected. HTTP request sent, awaiting response... 200 OK Length: 217 [text/plain] Saving to: ‘P61626.fasta’ 0K 100% 42.5M=0s 2016-11-29 18:40:00 (42.5 MB/s) - ‘P61626.fasta’ saved [217/217]
ls
1lmp.pdb P61626.fasta
alignm=modeller.alignment(env)
alignm.append(file='P61626.fasta', align_codes='all',alignment_format='FASTA')
mdl = modeller.model(env, file='1lmp.pdb', model_segment=('FIRST:'+'A', 'LAST:'+'A'))
alignm.append_model(mdl, atom_files='1lmp.pdb', align_codes='1lmp_A')
alignm[0].code='LYSC_HUMAN'
read_pd_459W> Residue type NAG not recognized. 'automodel' model building will treat this residue as a rigid body. To use real parameters, add the residue type to ${LIB}/restyp.lib, its topology to ${LIB}/top_*.lib, and suitable forcefield parameters to ${LIB}/par.lib. rdpdb___459W> Residue type NDG not recognized. 'automodel' model building will treat this residue as a rigid body. To use real parameters, add the residue type to ${LIB}/restyp.lib, its topology to ${LIB}/top_*.lib, and suitable forcefield parameters to ${LIB}/par.lib.
alignm.salign()
s= alignm[0]
for i in range(len(alignm)):
print "%s identical to wt in %4.1f perc" %(alignm[i].code,s.get_sequence_identity(alignm[i]))
alignm.write(file='all_in_one.ali', alignment_format='PIR')
SALIGN_____> adding the next group to the alignment; iteration 1 LYSC_HUMAN identical to wt in 100.0 perc 1lmp_A identical to wt in 68.2 perc
s=alignm[0]
print s.code
a = modeller.automodel.automodel(env, alnfile='all_in_one.ali', knowns= ('1lmp_A'), sequence = s.code )
a.name='mod'+s.code
a.starting_model = 1
a.ending_model = 1
a.make()
LYSC_HUMAN fndatmi_285W> Only 129 residues out of 132 contain atoms of type CA (This is usually caused by non-standard residues, such as ligands, or by PDB files with missing atoms.) fndatmi_285W> Only 129 residues out of 132 contain atoms of type CA (This is usually caused by non-standard residues, such as ligands, or by PDB files with missing atoms.) check_ali___> Checking the sequence-structure alignment. Implied target CA(i)-CA(i+1) distances longer than 8.0 angstroms: ALN_POS TMPL RID1 RID2 NAM1 NAM2 DIST ---------------------------------------------- END OF TABLE read_to_681_> topology.submodel read from topology file: 3 getf_______W> RTF restraint not found in the atoms list: residue type, indices: 18 148 atom names : C +N atom indices : 1155 0 getf_______W> RTF restraint not found in the atoms list: residue type, indices: 18 148 atom names : C CA +N O atom indices : 1155 1151 0 1156 fndatmi_285W> Only 129 residues out of 132 contain atoms of type CA (This is usually caused by non-standard residues, such as ligands, or by PDB files with missing atoms.) patch_s_522_> Number of disulfides patched in MODEL: 4 mdtrsr__446W> A potential that relies on one protein is used, yet you have at least one known structure available. MDT, not library, potential is used. iup2crm_280W> No topology library in memory or assigning a BLK residue. Default CHARMM atom type assigned: C1 --> CT2 This message is written only for the first such atom. 0 atoms in HETATM/BLK residues constrained to protein atoms within 2.30 angstroms and protein CA atoms within 10.00 angstroms 0 atoms in residues without defined topology constrained to be rigid bodies condens_443_> Restraints marked for deletion were removed. Total number of restraints before, now: 13194 12190 >> ENERGY; Differences between the model's features and restraints: Number of all residues in MODEL : 148 Number of all, selected real atoms : 1157 1157 Number of all, selected pseudo atoms : 0 0 Number of all static, selected restraints : 12190 12190 COVALENT_CYS : F NONBONDED_SEL_ATOMS : 1 Number of non-bonded pairs (excluding 1-2,1-3,1-4): 2351 Dynamic pairs routine : 2, NATM x NATM cell sorting Atomic shift for contacts update (UPDATE_DYNAMIC) : 0.390 LENNARD_JONES_SWITCH : 6.500 7.500 COULOMB_JONES_SWITCH : 6.500 7.500 RESIDUE_SPAN_RANGE : 0 99999 NLOGN_USE : 15 CONTACT_SHELL : 4.000 DYNAMIC_PAIRS,_SPHERE,_COULOMB,_LENNARD,_MODELLER : T T F F F SPHERE_STDV : 0.050 RADII_FACTOR : 0.820 Current energy : 804.2236 Summary of the restraint violations: NUM ... number of restraints. NUMVI ... number of restraints with RVIOL > VIOL_REPORT_CUT[i]. RVIOL ... relative difference from the best value. NUMVP ... number of restraints with -Ln(pdf) > VIOL_REPORT_CUT2[i]. RMS_1 ... RMS(feature, minimally_violated_basis_restraint, NUMB). RMS_2 ... RMS(feature, best_value, NUMB). MOL.PDF ... scaled contribution to -Ln(Molecular pdf). # RESTRAINT_GROUP NUM NUMVI NUMVP RMS_1 RMS_2 MOL.PDF S_i ------------------------------------------------------------------------------------------------------ 1 Bond length potential : 1177 0 0 0.006 0.006 11.387 1.000 2 Bond angle potential : 1592 0 7 2.150 2.150 147.53 1.000 3 Stereochemical cosine torsion poten: 738 0 24 47.751 47.751 266.00 1.000 4 Stereochemical improper torsion pot: 481 0 0 1.402 1.402 20.918 1.000 5 Soft-sphere overlap restraints : 2351 0 0 0.002 0.002 0.65140 1.000 6 Lennard-Jones 6-12 potential : 0 0 0 0.000 0.000 0.0000 1.000 7 Coulomb point-point electrostatic p: 0 0 0 0.000 0.000 0.0000 1.000 8 H-bonding potential : 0 0 0 0.000 0.000 0.0000 1.000 9 Distance restraints 1 (CA-CA) : 2405 0 0 0.131 0.131 40.957 1.000 10 Distance restraints 2 (N-O) : 2564 0 0 0.155 0.155 68.747 1.000 11 Mainchain Phi dihedral restraints : 0 0 0 0.000 0.000 0.0000 1.000 12 Mainchain Psi dihedral restraints : 0 0 0 0.000 0.000 0.0000 1.000 13 Mainchain Omega dihedral restraints: 147 0 3 4.636 4.636 37.258 1.000 14 Sidechain Chi_1 dihedral restraints: 120 0 0 68.491 68.491 25.477 1.000 15 Sidechain Chi_2 dihedral restraints: 86 0 0 63.198 63.198 33.297 1.000 16 Sidechain Chi_3 dihedral restraints: 35 0 0 89.980 89.980 22.231 1.000 17 Sidechain Chi_4 dihedral restraints: 20 0 0 102.640 102.640 12.968 1.000 18 Disulfide distance restraints : 4 0 0 0.006 0.006 0.27561E-01 1.000 19 Disulfide angle restraints : 8 0 0 2.789 2.789 1.3737 1.000 20 Disulfide dihedral angle restraints: 4 0 0 23.141 23.141 2.0451 1.000 21 Lower bound distance restraints : 0 0 0 0.000 0.000 0.0000 1.000 22 Upper bound distance restraints : 0 0 0 0.000 0.000 0.0000 1.000 23 Distance restraints 3 (SDCH-MNCH) : 1691 0 0 0.425 0.425 34.838 1.000 24 Sidechain Chi_5 dihedral restraints: 0 0 0 0.000 0.000 0.0000 1.000 25 Phi/Psi pair of dihedral restraints: 146 16 13 23.998 62.053 26.471 1.000 26 Distance restraints 4 (SDCH-SDCH) : 972 0 0 0.637 0.637 52.045 1.000 27 Distance restraints 5 (X-Y) : 0 0 0 0.000 0.000 0.0000 1.000 28 NMR distance restraints 6 (X-Y) : 0 0 0 0.000 0.000 0.0000 1.000 29 NMR distance restraints 7 (X-Y) : 0 0 0 0.000 0.000 0.0000 1.000 30 Minimal distance restraints : 0 0 0 0.000 0.000 0.0000 1.000 31 Non-bonded restraints : 0 0 0 0.000 0.000 0.0000 1.000 32 Atomic accessibility restraints : 0 0 0 0.000 0.000 0.0000 1.000 33 Atomic density restraints : 0 0 0 0.000 0.000 0.0000 1.000 34 Absolute position restraints : 0 0 0 0.000 0.000 0.0000 1.000 35 Dihedral angle difference restraint: 0 0 0 0.000 0.000 0.0000 1.000 36 GBSA implicit solvent potential : 0 0 0 0.000 0.000 0.0000 1.000 37 EM density fitting potential : 0 0 0 0.000 0.000 0.0000 1.000 38 SAXS restraints : 0 0 0 0.000 0.000 0.0000 1.000 39 Symmetry restraints : 0 0 0 0.000 0.000 0.0000 1.000 # Heavy relative violation of each residue is written to: LYSC_HUMAN.V99990001 # The profile is NOT normalized by the number of restraints. # The profiles are smoothed over a window of residues: 1 # The sum of all numbers in the file: 14675.1133 List of the violated restraints: A restraint is violated when the relative difference from the best value (RVIOL) is larger than CUTOFF. ICSR ... index of a restraint in the current set. RESNO ... residue numbers of the first two atoms. ATM ... IUPAC atom names of the first two atoms. FEAT ... the value of the feature in the model. restr ... the mean of the basis restraint with the smallest difference from the model (local minimum). viol ... difference from the local minimum. rviol ... relative difference from the local minimum. RESTR ... the best value (global minimum). VIOL ... difference from the best value. RVIOL ... relative difference from the best value. ------------------------------------------------------------------------------------------------- Feature 25 : Phi/Psi pair of dihedral restraints List of the RVIOL violations larger than : 6.5000 # ICSR RESNO1/2 ATM1/2 INDATM1/2 FEAT restr viol rviol RESTR VIOL RVIOL 1 4005 1M 2K C N 7 9 -64.57 -70.20 5.95 0.39 -62.90 179.26 23.19 1 2K 2K N CA 9 10 138.46 140.40 -40.80 2 4006 2K 3A C N 16 18 -74.64 -68.20 9.27 0.58 -62.50 167.58 28.07 2 3A 3A N CA 18 19 151.96 145.30 -40.90 3 4007 3A 4L C N 21 23 -98.45 -108.50 13.80 0.77 -63.50 -179.75 27.15 3 4L 4L N CA 23 24 141.97 132.50 -41.20 4 4010 6V 7L C N 44 46 -103.37 -108.50 6.06 0.27 -63.50 175.07 22.35 4 7L 7L N CA 46 47 129.27 132.50 -41.20 5 4011 7L 8G C N 52 54 -90.16 -80.20 68.15 3.19 82.20 -161.64 8.57 5 8G 8G N CA 54 55 106.68 174.10 8.50 6 4012 8G 9L C N 56 58 -125.41 -108.50 41.06 2.01 -63.50 161.23 26.22 6 9L 9L N CA 58 59 169.93 132.50 -41.20 7 4013 9L 10V C N 64 66 -58.86 -62.40 4.38 0.54 -125.40 -175.85 10.37 7 10V 10V N CA 66 67 -44.98 -42.40 143.30 8 4014 10V 11L C N 71 73 -111.10 -108.50 18.88 0.99 -63.50 174.23 27.12 8 11L 11L N CA 73 74 151.20 132.50 -41.20 9 4015 11L 12L C N 79 81 -72.91 -70.70 13.00 1.07 -63.50 170.25 23.12 9 12L 12L N CA 81 82 128.79 141.60 -41.20 10 4016 12L 13S C N 87 89 163.67 -136.60 83.64 2.69 -64.10 175.41 22.09 10 13S 13S N CA 89 90 -150.26 151.20 -35.00 11 4018 14V 15T C N 100 102 -104.81 -124.80 23.58 0.88 -63.20 178.03 20.35 11 15T 15T N CA 102 103 131.00 143.50 -42.10 12 4020 16V 17Q C N 114 116 88.80 -121.10 150.10 5.65 -63.80 -124.45 42.00 12 17Q 17Q N CA 116 117 139.14 139.70 -40.30 13 4021 17Q 18G C N 123 125 -101.84 -80.20 44.39 1.04 82.20 -125.09 7.89 13 18G 18G N CA 125 126 -147.14 174.10 8.50 14 4022 18G 19K C N 127 129 -77.16 -70.20 13.93 1.14 -62.90 169.73 21.30 14 19K 19K N CA 129 130 128.33 140.40 -40.80 15 4057 53E 54S C N 415 417 -136.30 -64.10 78.23 8.47 -64.10 78.23 8.47 15 54S 54S N CA 417 418 -4.87 -35.00 -35.00 16 4069 65A 66G C N 508 510 -66.26 -62.40 4.87 0.92 82.20 157.53 12.01 16 66G 66G N CA 510 511 -44.18 -41.20 8.50 report______> Distribution of short non-bonded contacts: DISTANCE1: 0.00 2.10 2.20 2.30 2.40 2.50 2.60 2.70 2.80 2.90 3.00 3.10 3.20 3.30 3.40 DISTANCE2: 2.10 2.20 2.30 2.40 2.50 2.60 2.70 2.80 2.90 3.00 3.10 3.20 3.30 3.40 3.50 FREQUENCY: 0 0 0 0 0 4 10 37 67 104 113 140 158 187 183 << end of ENERGY. >> Summary of successfully produced models: Filename molpdf ---------------------------------------- LYSC_HUMAN.B99990001.pdb 804.22363
import nglview
import ipywidgets
w1 = nglview.show_structure_file('LYSC_HUMAN.B99990001.pdb')
w2 = nglview.show_pdbid('1lmp')
class mymodel(modeller.automodel.automodel):
def special_restraints(self, aln):
rsr = self.restraints
for ids in (('OD1:98:A', 'O6A:131:A'),
('N:65:A', 'O7B:132:A'),
('OD2:73:A', 'O1C:133:A')):
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))