Hello,
I am having a great problem with modeller 9v1.
My aim is to model 2 chain protein from a 2 chain template with a
ligands. I want two ligands: Ca metal ion and SOX (Penicillin G
sulphoxide). My template pdb file briefly looks like this:
----------------------------------------------------------------------------------------------------
ATOM 1 N SER A 3 27.496 31.222 28.262 1.00
48.33 N
ATOM 2 CA SER A 3 27.785 31.042 26.803 1.00
43.97 C
... etc to the end of chain a
TER 1679 ALA A 209
ATOM 1684 N SER B 1 14.845 1.596 -1.736 1.00
6.93 N
ATOM 1685 CA SER B 1 15.548 2.917 -1.659 1.00
6.75 C
ATOM 1686 C SER B 1 16.489 3.059 -2.844 1.00
8.04 C
.....
ATOM 6113 NE ARG B 557 26.460 16.023 30.916 1.00
23.74 N
ATOM 6114 CZ ARG B 557 27.175 15.276 30.094 1.00
19.20 C
ATOM 6115 NH1 ARG B 557 28.016 15.921 29.310 1.00
21.83 N
ATOM 6116 NH2 ARG B 557 27.069 13.952 30.055 1.00
20.49 N
ATOM 6117 OXT ARG B 557 23.973 18.664 34.656 1.00
29.98 O
TER 6118 ARG B 557
.... Cutting useless HETATM
HETATM 6159 CA CA B1568 23.187 -16.230 6.242 1.00
10.86 CA << Ca ion
HETATM 6160 O8 SOX B1569 13.279 0.047 0.037 1.00
37.00 O << SOX ...
HETATM 6161 C7 SOX B1569 13.688 -1.062 0.060 1.00
38.52 C
..
HETATM 6166 O12 SOX B1569 13.790 -4.072 -2.049 1.00
54.27 O
...etc to the end of SOX (I guess resid for SOX = 569, chain B1, ?)
-----------------------------------------------------------------------------------------------------
I have made alignment in PIR using chain breaks "/" and block residues
"." for Ca and SOX.
My alignment looks like this:
-----------------------------------------------------------------------------------------------------
>P1;1gm9
structureX:1gm9.pdb:3:A:569:B1:::: #Starting from A resid =3(because
first two are missing) and finishing with SOX
--SSSEIKIVRDEYG-PHIYANDTWHLFYGYGYVVAQDRLFQMEMARRSTQGTVAEVLGK
DFVKFDKDIRRNYWPDAIRAQIAALSPEDMSILQGYADGMNAWIDKVNTNPETLLPKQFN
TFGFTPKRWEPFDVAMIFVGTMANRFSDSTSEIDNLALLTALKDKYGVSQGMAVFNQLKW
LVNPSAPTTIAVQESNYPLKFNQQNSQTA/ #Chain break = ending chain A
SNMWVIGKSKAQDAKAIMVNGPQFGWYAPAY
TYGIGLHGAGYDVTGNTPFAYPGLVFGHNGVISWGSTAGFGDDVDIFAERLSAEKPGYYL
HNGKWVKMLSREETITVKNGQAETFTVWRTVHGNILQTDQTTQTAYAKSRAWDGKEVASL
LAWTHQMKAKNWQEWTQQAAKQALTINWYYADVNGNIGYVHTGAYPDRQSGHDPRLPVPG
TGKWDWKGLLPFEMNPKVYNPQSGYIANWNNSPQKDYPASDLFAFLWGGADRVTEIDRLL
EQKPRLTADQAWDVIRQTSRQDLNLRLFLPTLQAATSGLTQSDPRRQLVETLTRWDGINL
LNDDGKTWQQPGSAILNVWLTSMLKRTVVAAVPMPFDKWYSASGYETTQDGPTGSLNISV
GAKILYEAVQGDKSPIPQAVDLFAGKPQQEVVLAALEDTWETLSKRYGNNVSNWKTPAMA
LTFRANNFFGVPQAAAEETRHQAEYQNRGTENDMIVFSPTTSDRPVLAWDVVAPGQSGFI
APDGTVDKHYEDQLKMYENFGRKSLWLTKQDVEAHKESQEVLHVQR/..*
#Chain break = ending chain B and ".." for two HETATM
>P1;1klc
sequence:1klc: : : : ::: 0.00: 0.00
ASPPTEVKIVRDEYGMPHIYADDTYRLFYGYGYVVAQDRLFQMEMARRSTQGTVSEVLGK
AFVSFDKDIRQNYWPDSIRAQIASLSAEDKSILQGYADGMNAWIDKVNASPDKLLPQQFS
TFGFKPKHWEPFDVAMIFVGTMANRFSDSTSEIDNLALLTAVKDKYGNDEGMAVFNQLKW
LVNPSAPTTIAARESSYPLKFDLQNTQTA/
SNMWVIGKNKAQDAKAIMVNGPQFGWYAPAY
TYGIGLHGAGYDVTGNTPFAYPGLVFGHNGTISWGSTAGFGDDVDIFAEKLSAEKPGYYQ
HNGEWVKMLSRKETIAVKDGQPETFTVWRTLDGNVIKTDTRTQTAYAKARAWAGKEVASL
LAWTHQMKAKNWPEWTQQAAKQALTINWYYADVNGNIGYVHTGAYPDRQPGHDPRLPVP-
DGKWDWKGLLSFDLNPKVYNPQSGYIANWNNSPQKDYPASDLFAFLWGGADRVTEIDTIL
DKQPRFTADQAWDVIRQTSLRDL-LRLFLPALKDATANLAENDPRRQLVDKLASWDGENL
VNDDGKTYQQPGSAILNAWLTSMLKRTVVAAVPAPFGKWYSASGYETTQDGPTGSLNISV
GAKILYEALQGDKSPIPQAVDLFGGKPEQEVILAALDDAWQTLSKRYGNDVTGWKTPAMA
LTFRANNFFGVPQAAAKEARHQAEYQNRGTENDMIVFSPTSGNRPVLAWDVVAPGQSGFI
APDGKADKHYDDQLKMYESFGRKSLWLTPQDVDEHKESQEVLQVQR/--*
---------------------------------------------------------------------------------------------------------
Then I used advanced example as a template to write my input script:
---------------------------------------------------------------------------------------------------------
from modeller import *
from modeller.automodel import *
class mymodel(automodel):
def special_restraints(self, aln):
rsr = self.restraints
for ids in (('NH2:471:B', 'O12:766:C'), # That is for SOX, no
restraints for Ca
('ND2:449:B', 'O16:766:C'),
('CE2:280:B', 'C5:766:C')):
atoms = [self.atoms[i] for i in ids]
rsr.add(forms.upper_bound(group=physical.upper_distance,
feature=features.distance(*atoms),
mean=3.04, stdev=0.1))
env = environ()
env.io.hetatm = True
a = mymodel(env, alnfile='1klc_1gm9.ali',
knowns='1gm9m, 1klc_model', sequence='1klc')
a.starting_model = 1
a.ending_model = 1
a.make()
------------------------------------------------------------------------------------------------------
Modeller starts with this error: "import site' failed; use -v for
traceback", (I guess it is not a big problem
because in other cases it does not cause any trouble) , then program
starts calculating, and then I always get :
-----------------------------------------------------------------------------------------------------
Traceback (most recent call last):
File "model-multiple-hetero.py", line 47, in ?
a.make()
File "C:\Program
Files\Modeller9v1\modlib\modeller\automodel\automodel.py", li
ne 108, in make
self.homcsr(exit_stage)
File "C:\Program
Files\Modeller9v1\modlib\modeller\automodel\automodel.py", li
ne 427, in homcsr
self.mkhomcsr(selection(self), aln)
File "C:\Program
Files\Modeller9v1\modlib\modeller\automodel\automodel.py", li
ne 537, in mkhomcsr
self.special_restraints(aln)
File "model-multiple-hetero.py", line 10, in special_restraints
atoms = [self.atoms[i] for i in ids]
File "C:\Program Files\Modeller9v1\modlib\modeller\coordinates.py",
line 127,
in __getitem__
(self.offset, self.length, self.suffix))
File "C:\Program Files\Modeller9v1\modlib\modeller\util\modutil.py",
line 76,
in handle_seq_indx
int_indx = lookup_func(*args)
File "C:\Program Files\Modeller9v1\modlib\modeller\coordinates.py",
line 42, i
n _indxatm
raise KeyError, ("No such atom: %s" % indx)
KeyError: 'No such atom: O12:766:C'
----------------------------------------------------------------------------------------------------------
That means, that it can not see atom O12:766:C.
I am sure, that my seq is 764 residues long (I mean residues in seq, not
position in alignment).
So Ca must be 765 and SOX - 766. I have tried different counts
(766,765,767...) and different chains
(B1, C) but it did not help. I have tried to model only chain B with
SOX, but it did not help.
If I remove user restraints and just use automodel, it output good
model, but without ligands at all!
And the latter is what I absolutely do no understand.
My log has this note about SOX:
-----------------------------------------------------------------------------------------------------------
read_pd_459W> Residue type SOX 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.
-----------------------------------------------------------------------------------------------------------
And nothing about CA (may be the program takes CA as Calpha?).
So my questions are:
1. What am I doing wrong?
2. Why modeller outputs model with no ligands if I use default automodel
with no restraints (env.io.hetatm =True)?.
What do you think? Thank you for your suggestions.
P.S. I am using windows compilation. May it be some platform specific
modeller bug?
