Dear all,
Many thanks for your previous help with my copper binding proteins and I will apologies in advance for the length of this query! I require some guidance with a modelling problem I am having. I am attempting to model the structure of a 258 residue segment of a decaheme cytochrome, I have template models covering the positions of the 10 hemes and the first 223 residues. I have also implemented several distance and bond angle restraints on the heme binding residues in an attempt to correctly orientate them around the hemes, however when I run modeller some of the restraints dont seem to take affect. A specific example is residue 258 that (if I have implemented the script correctly) should be within 2.0A of the Fe atom from heme 268, however, in the final model it is 50A away from the heme and has been pulled 50A apart.
The following is the complete script I have used for my modelling attempts:
--------------------------------------------------------------------------- # Addition of restraints to the default ones from modeller import * from modeller.automodel import * # Load the automodel class log.verbose() class MyModel(automodel):
def select_atoms(self): return selection(self.residue_range('1:', '258:')) def special_restraints(self, aln): rsr = self.restraints at = self.atoms
# Restrict secondary structure, Residues that should be an alpha helix: rsr.add(secondary_structure.alpha(self.residue_range('16:', '19:'))) rsr.add(secondary_structure.alpha(self.residue_range('42:', '43:'))) rsr.add(secondary_structure.alpha(self.residue_range('48:', '49:'))) rsr.add(secondary_structure.alpha(self.residue_range('70:', '73:'))) rsr.add(secondary_structure.alpha(self.residue_range('78:', '80:'))) rsr.add(secondary_structure.alpha(self.residue_range('84:', '87:'))) rsr.add(secondary_structure.alpha(self.residue_range('132:', '134:'))) rsr.add(secondary_structure.alpha(self.residue_range('150:', '151:'))) rsr.add(secondary_structure.alpha(self.residue_range('162:', '164:'))) rsr.add(secondary_structure.alpha(self.residue_range('170:', '178:'))) rsr.add(secondary_structure.alpha(self.residue_range('216:', '222:'))) rsr.add(secondary_structure.alpha(self.residue_range('254:', '256:')))
# Add distance restraints to heme binding CXXCH residues: rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:41'], at['CAB:260']), mean=1.825, stdev=0.01)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:44'], at['CAC:260']), mean=1.825, stdev=0.01)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:45'], at['FE:260']), mean=2.0, stdev=0.1)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:77'], at['CAB:261']), mean=1.825, stdev=0.01)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:80'], at['CAC:261']), mean=1.825, stdev=0.01)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:81'], at['FE:261']), mean=2.0, stdev=0.1)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:101'], at['CAB:262']), mean=1.825, stdev=0.01)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:104'], at['CAC:262']), mean=1.825, stdev=0.01)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:105'], at['FE:262']), mean=2.0, stdev=0.1)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:124'], at['CAB:263']), mean=1.825, stdev=0.01)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:127'], at['CAC:263']), mean=1.825, stdev=0.01)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:128'], at['FE:263']), mean=2.0, stdev=0.1)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:150'], at['CAB:264']), mean=1.825, stdev=0.01)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:153'], at['CAC:264']), mean=1.825, stdev=0.01)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:154'], at['FE:264']), mean=2.0, stdev=0.1)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:174'], at['CAB:265']), mean=1.825, stdev=0.01)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:177'], at['CAC:265']), mean=1.825, stdev=0.01)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:178'], at['FE:265']), mean=2.0, stdev=0.1)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:196'], at['CAB:266']), mean=1.825, stdev=0.01)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:199'], at['CAC:266']), mean=1.825, stdev=0.01)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:200'], at['FE:266']), mean=2.0, stdev=0.1)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:219'], at['CAB:267']), mean=1.825, stdev=0.01)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:222'], at['CAC:267']), mean=1.825, stdev=0.01)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:223'], at['FE:267']), mean=2.0, stdev=0.1)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:8'], at['CAB:259']), mean=1.825, stdev=0.01)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['SG:11'], at['CAC:259']), mean=1.825, stdev=0.01)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:12'], at['FE:259']), mean=2.0, stdev=0.1)) rsr.add(forms.upper_bound(group=physical.upper_distance, feature=features.distance(at['SG:254'], at['CAB:268']), mean=1.825, stdev=0.1)) rsr.add(forms.upper_bound(group=physical.upper_distance, feature=features.distance(at['SG:257'], at['CAC:268']), mean=1.825, stdev=0.1)) rsr.add(forms.upper_bound(group=physical.upper_distance, feature=features.distance(at['NE2:258'], at['FE:268']), mean=2.0, stdev=0.1))
# Add distance restraints to the other heme binding His residues: rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:51'], at['FE:259']), mean=2.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:108'], at['FE:260']), mean=2.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:26'], at['FE:261']), mean=2.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:94'], at['FE:262']), mean=2.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:93'], at['FE:263']), mean=2.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:203'], at['FE:264']), mean=2.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:141'], at['FE:265']), mean=2.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:189'], at['FE:267']), mean=2.0, stdev=0.1)) rsr.add(forms.upper_bound(group=physical.xy_distance, feature=features.distance(at['NE2:228'], at['FE:268']), mean=2.0, stdev=0.1))
# Set distance restraints for Heme binding His-His residues at 4.0A rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:12'], at['NE2:51']), mean=4.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:45'], at['NE2:108']), mean=4.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:81'], at['NE2:26']), mean=4.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:105'], at['NE2:94']), mean=4.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:128'], at['NE2:93']), mean=4.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:154'], at['NE2:203']), mean=4.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:178'], at['NE2:141']), mean=4.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:223'], at['NE2:189']), mean=4.0, stdev=0.1)) rsr.add(forms.gaussian(group=physical.xy_distance, feature=features.distance(at['NE2:258'], at['NE2:228']), mean=4.0, stdev=0.1))
# Set bond angles for iron coordinating HIS residues of CXXCH motifs rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:12'], at['FE:259'], at['CHA:259']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:45'], at['FE:260'], at['CHA:260']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:81'], at['FE:261'], at['CHA:261']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:105'], at['FE:262'], at['CHA:262']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:128'], at['FE:263'], at['CHA:263']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:154'], at['FE:264'], at['CHA:264']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:178'], at['FE:265'], at['CHA:265']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:200'], at['FE:266'], at['CHA:266']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:223'], at['FE:267'], at['CHA:267']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:258'], at['FE:268'], at['CHA:268']), mean=1.63, stdev=0.1))
# Set bond angles for other iron coordinating HIS residues rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:51'], at['FE:259'], at['CHA:259']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:108'], at['FE:260'], at['CHA:260']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:26'], at['FE:261'], at['CHA:261']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:94'], at['FE:262'], at['CHA:262']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:93'], at['FE:263'], at['CHA:263']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:203'], at['FE:264'], at['CHA:264']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:141'], at['FE:265'], at['CHA:265']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:189'], at['FE:267'], at['CHA:267']), mean=1.63, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['CG:228'], at['FE:268'], at['CHA:268']), mean=1.63, stdev=0.1))
# Restrain heme binding His-Fe-His bond angles rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['NE2:12'], at['FE:259'], at['NE2:51']), mean=3.07, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['NE2:45'], at['FE:260'], at['NE2:108']), mean=3.07, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['NE2:81'], at['FE:261'], at['NE2:26']), mean=3.07, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['NE2:105'], at['FE:262'], at['NE2:94']), mean=3.07, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['NE2:128'], at['FE:263'], at['NE2:93']), mean=3.07, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['NE2:154'], at['FE:264'], at['NE2:203']), mean=3.07, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['NE2:178'], at['FE:265'], at['NE2:141']), mean=3.07, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['NE2:223'], at['FE:267'], at['NE2:189']), mean=3.07, stdev=0.1)) rsr.add(forms.gaussian(group=physical.angle, feature=features.angle(at['NE2:258'], at['FE:268'], at['NE2:228']), mean=3.07, stdev=0.1))
env = environ() env.io.hetatm = True # directories for input atom files env.io.atom_files_directory = ['.', '../atom_files']
a = MyModel(env, alnfile = 'alignment.ali', # alignment filename knowns = ('2OZY','MtrAst'), # codes of the templates sequence = 'MtrA', # code of the target assess_methods=(assess.DOPE, assess.GA341))
# To keep the hemes from being moved def select_atoms(self): self.pick_atoms(selection_segment=('1:', '258:'), selection_search='segment', pick_atoms_set=1, res_types='all', atom_types='all', selection_from='all', selection_status='initialize') a.starting_model= 1 # index of the first model a.ending_model = 1 # index of the last model # (determines how many models to calculate) a.make() # do homology modeling
---------------------------------------------------------------------------
Many thanks for any help you can provide, and once again apologies for the length of this post.
Stephen Hearnshaw