[modeller_usage] Distance Restraints: are they working?

[João Rodrigues <>]

Hello all!

I'm new at Modeller and I'm trying my luck with Distance Restraints to model a small domain of a protein. I was running some tests on my modelling to see if I had those restraints working and, from my preliminary (and rookie) analysis, it seems it's not.

Attached are the files I'm using. The MyModel.py script has the MyModel class with my restraints defined. It is then called by the main modeller script (RNF11_1X4J_rest.pdb) (had to use this approach for parallel computation) which runs on modeller 9v5.

I've described each restraint like this:

rsr.make_distance(at['SG:1'],at['SG:4'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=5)

I just want both atoms to be, in a maximal distance of 5A. I made some further tests with that value of maximal_distance at 1.5 and 0.5, just to see if it ran ok, and I found out that the 3 best models (5A, 1.5A, 0.5A) are exactly the same :x (or nearly), with distances in the order of 2-3A, which would violate this restraint I've added in the 1.5A and 0.5A.

Is this behaviour normal, or is it something I'm doing that is not the most correct?

Thanks in advance for the answers and congrats for the software and for the whole team.

João Rodrigues
Bijvoet Center for Biomolecular Research
Utrecht, Netherlands
# Homology modeling by the automodel class
from modeller import *                       # Load standard Modeller classes
from modeller.automodel import *             # Load the automodel class

class MyModel(loopmodel):
    def special_restraints(self, aln):
        rsr = self.restraints
        at = self.atoms

        
        # Ring Finger: C3HC3H (Type 2)

        # First Zinc
        rsr.make_distance(at['SG:1'],at['SG:4'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['SG:1'],at['ND1:24'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['SG:1'],at['SG:27'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['SG:4'],at['ND1:24'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['SG:4'],at['SG:27'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['ND1:24'],at['SG:27'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['ZN:42'],at['SG:1'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['ZN:42'],at['SG:4'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['ZN:42'],at['ND1:24'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['ZN:42'],at['SG:27'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)

        # Second Zinc
        rsr.make_distance(at['SG:19'],at['ND1:21'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['SG:19'],at['SG:38'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['SG:19'],at['SG:41'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['ND1:21'],at['SG:38'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['ND1:21'],at['SG:41'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['SG:38'],at['SG:41'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)

        rsr.make_distance(at['ZN:43'],at['SG:19'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['ZN:43'],at['ND1:21'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['ZN:43'],at['SG:38'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)
        rsr.make_distance(at['ZN:43'],at['SG:41'],aln,spline_on_site=False,restraint_group=physical.xy_distance,maximal_distance=0.5)

# Homology modeling by the automodel class
from modeller import *                       # Load standard Modeller classes
from modeller.automodel import *             # Load the automodel class
from modeller.parallel import *              # Load the parallel class, to use multiple processors
from MyModel import *

# Use 8 CPUs in a parallel job on this machine

#j = job(modeller_path="/home/joao/Software/modeller9v5/bin/modslave.py") # Local
j = job(modeller_path="/home/software/bin/modeller9v5/bin/modslave.py") # Cluster

j.append(local_slave())
j.append(local_slave())
j.append(local_slave())
j.append(local_slave())
j.append(local_slave())
j.append(local_slave())
j.append(local_slave())
j.append(local_slave())

log.verbose()                                # request verbose output
env = environ()                              # create a new MODELLER environment to build this model in

# directories for input atom files
env.io.atom_files_directory = './'

# Read in HETATM records from template PDBs
env.io.hetatm = True

# Define my own automodel class to add optional restraints


# Loop Optimization
# Use dope_loopmodel or dopehr_loopmodel in place of loopmodel to obtain better quality loops (but slower)
a = MyModel(env,
            alnfile  = 'RNF11_1X4J.ali',      # alignment filename
            knowns   = ('1X4J'),           # codes of the templates
            sequence = 'RNF11',            # code of the target
            assess_methods = assess.DOPE)  # request assessment of the generated models by DOPE method
a.starting_model= 1                          # index of the first model
a.ending_model  = 100                        # index of the last model

# Optimization
  # CG
a.library_schedule = autosched.slow          # Very thorough VTFM optimization
a.max_var_iterations = 300                   # Select length of optimizations
a.max_molpdf = 1e6                           # do not stop unless obj.func. > 1E6

  #MD
a.md_level = refine.slow                     # model refinement level

  #All
#a.repeat_optimization = 2                    # Repeat the whole cycle 2 times

  #Loops
a.loop.starting_model = 1                    # First loop model
a.loop.ending_model   = 4                    # Last loop model
a.loop.md_level       = refine.slow          # Loop model refinement level

# Repeat the whole cycle 2 times and do not stop unless obj.func. > 1E6
#a.repeat_optimization = 2
a.max_molpdf = 1e6

a.use_parallel_job(j)                        # Use the job for model building
a.make()                                     # do the actual homology modeling


##################
# Postprocessing #
##################

# Cluster all built models
a.cluster(cluster_cut=1.5)

# Accessing output data after modeling is complete :
  # Get a list of all successfully built models from a.outputs
ok_models = filter(lambda x: x['failure'] is None, a.outputs)

  # Rank the models by DOPE score
key = 'DOPE score'
ok_models.sort(lambda a,b: cmp(a[key], b[key]))

  # Get top model
m = ok_models[0]
print "Top model: %s (DOPE score %.3f)" % (m['name'], m[key])