I have a model of a protein-protein complex, where several residues in a loop in the interface are undefined. I want to use xplor to add the loop residues, and put them in a location that does not clash with the rest of the protein or with the binding partner.

I used addAtoms2.py to generate the loop residues without disturbing the rest of the structure. The loop that was generated clashed with the binding partner, so I wanted to run some dynamics to move the loop into a more suitable position, again without disturbing the rest of the model.

I found the dock.py script has some features to treat the proteins as rigid bodies, and to prevent sidechain movement, while allowing some movement in a defined part of the protein(s). So I'm trying to use the dock.py script, with lines like:
dyn_fix.group( 'resid 1022:1237')
dyn_fix.group( 'resid 1238:1329' )
dyn_fix.group( 'resid 1338:1428' )
dyn_fix.fix( 'resid 1022:1237')
dyn_fix.fix( 'resid 1238:1329')
dyn_fix.fix( 'resid 1338:1428' )

xplor.requireVersion("2.17")

numberOfStructures=1

ensembleSize=2           # number of ensemble members

outFilename = "SCRIPT_%d_STRUCTURE_MEMBER.sa" % ensembleSize

import protocol

protocol.initRandomSeed()

command = xplor.command

protocol.initParams("protein")

#load pdb and covalently minimize the protein

protocol.loadPDB("model.pdb")

xplor.simulation.deleteAtoms("not known")

protocol.fixupCovalentGeom(maxIters=1000,useVDW=1)

#initilize ensemble simulations

from ensembleSimulation import EnsembleSimulation

esim = EnsembleSimulation("ensemble",ensembleSize)

#

# a PotList contains a list of potential terms. This is used to specify which

# terms are active during refinement.

#

from potList import PotList

potList = PotList()

# parameters to ramp up during the simulated annealing protocol

#

from simulationTools import MultRamp, StaticRamp, InitialParams

rampedParams=[]

highTempParams=[]

# orientation Tensor - used with the dipolar coupling term

#  one for each medium

#   For each medium, specify a name, and initial values of Da, Rh.

#

from varTensorTools import create_VarTensor

media={}

expdatadir = 'xplor5050/'

import os,re

mediaList = open(expdatadir+'mediaList.txt').readlines()

#media[medium] = create_VarTensor(medium)

for mediaName in mediaList:

    mediaName = mediaName.strip()

    if not mediaName in media.keys():

        media[mediaName] = create_VarTensor(mediaName)

            

from rdcPotTools import create_RDCPot, scale_toNH,correctGyromagneticSigns

correctGyromagneticSigns()

bondTypes = ['NH','CAC','CAHA','CN','CHN']

#xplorBondTypes = {'NH':'NH','CN':'NCO','HNC':'CHN','CAHA':'CAHA', 'CAC':'CACO'}

#bondTypeScales = {'NH':15, 'CAHA':3, 'CN':100, 'CAC':100, 'CHN':3}

bondTypeScales = {'NH':1, 'CAHA':1, 'CN':1, 'CAC':1, 'CHN':1}

rdcs = PotList('rdc') 

for medium in media.keys():

    #rdc = create_RDCPot("%s_%s"%(medium,expt),file,media[medium])

    for btype in bondTypes:

        if os.path.isfile(expdatadir+medium+'_'+btype+'.tbl'):

            rdc = create_RDCPot("%s_%s"%(medium,btype),expdatadir+medium+'_'+btype+'.tbl',media[medium])    

            rdc.setScale(bondTypeScales[btype])

            rdc.setAveType("sum")    # its set to be average by default         

            rdc.setShowAllRestraints(1)            

            rdcs.append(rdc)

potList.append(rdcs)

rampedParams.append( MultRamp(0.05,5.0, "rdcs.setScale( VALUE )") )

# calc. initial tensor orientation

# and setup tensor calculation during simulated annealing

#

from varTensorTools import calcTensorOrientation, calcTensor

for medium in media.keys():

    calcTensor(media[medium])

    rampedParams.append( StaticRamp("calcTensor(media['%s'])" % medium) )

    pass

#add other things about NOE later

from avePot import AvePot

from xplorPot import XplorPot

potList.append( AvePot(XplorPot,'VDW') )

rampedParams.append( StaticRamp("protocol.initNBond()") )

rampedParams.append( MultRamp(0.9,0.8,

                              "command('param nbonds repel VALUE end end')") )

rampedParams.append( MultRamp(.004,4,

                              "command('param nbonds rcon VALUE end end')") )

                              

potList.append( AvePot(XplorPot,"BOND") )

potList.append( AvePot(XplorPot,"ANGL") )

potList['ANGL'].setThreshold( 5 )

rampedParams.append( MultRamp(0.4,1,"potList['ANGL'].setScale(VALUE)") )

potList.append( AvePot(XplorPot,"IMPR") )

potList['IMPR'].setThreshold( 5 )

rampedParams.append( MultRamp(0.1,1,"potList['IMPR'].setScale(VALUE)") )

# Give atoms uniform weights, except for the anisotropy axis

#

protocol.massSetup()

# IVM setup

#   the IVM is used for performing dynamics and minimization in torsion-angle

#   space, and in Cartesian space.

#

from ivm import IVM

dyn = IVM()

dyn.reset()

for m in media.values():

    m.setFreedom("varyDa, varyRh")      #vary tensor Rh, Da, vary orientation

protocol.torsionTopology(dyn)

#for ending with cartesian minimzation

minc = IVM()

protocol.initMinimize(minc)

for m in media.values():

    m.setFreedom("varyDa, varyRh")    #allow all tensor parameters float here

    pass

    

protocol.cartesianTopology(minc)

# object which performs simulated annealing

#

from simulationTools import AnnealIVM

init_t  = 3000.     # Need high temp and slow annealing to converge

cool = AnnealIVM(initTemp =init_t,

                 finalTemp=25,

                 tempStep =12.5,

                 ivm=dyn,

                 rampedParams = rampedParams)

#criteria to accept a structure

def accept(potList):

    """

    return True if current structure meets acceptance criteria

    """

    #if potList['noe'].violations()>0:

    #    return False

    if potList['rdc'].rms()>1.2: #this might be tightened some

        return False

    #if potList['CDIH'].violations()>0:

    #    return False

    if potList['BOND'].violations()>0:

        return False

    if potList['ANGL'].violations()>0:

        return False

    if potList['IMPR'].violations()>1:

        return False

    

    return True

def calcOneStructure(loopInfo):

    """ this function calculates a single structure, performs analysis on the

    structure, and then writes out a pdb file, with remarks.

    """

    # initialize parameters for high temp dynamics.

    InitialParams( rampedParams )

    # high-temp dynamics setup - only need to specify parameters which

    #   differfrom initial values in rampedParams

    InitialParams( highTempParams )

    # high temp dynamics

    #

    protocol.initDynamics(dyn,

                          potList=potList, # potential terms to use

                          bathTemp=init_t,

                          initVelocities=1,

                          finalTime=10,    # stops at 10ps or 5000 steps

                          numSteps=5000,   # whichever comes first

                          printInterval=100)

    dyn.setETolerance( init_t/100 )  #used to det. stepsize. default: t/1000 

    dyn.run()

    # initialize parameters for cooling loop

    InitialParams( rampedParams )

    # initialize integrator for simulated annealing

    #

    protocol.initDynamics(dyn,

                          potList=potList,

                          numSteps=100,       #at each temp: 100 steps or

                          finalTime=.2 ,       # .2ps, whichever is less

                          printInterval=100)

    # perform simulated annealing

    #

    cool.run()

              

              

    # final torsion angle minimization

    #

    protocol.initMinimize(dyn,

                          printInterval=50)

    dyn.run()

    # final all- atom minimization

    #

    protocol.initMinimize(minc,

                          potList=potList,

                          dEPred=10)

    minc.run()

    #do analysis and write structure

    loopInfo.writeStructure(potList)

    pass

from simulationTools import StructureLoop, FinalParams

StructureLoop(numStructures=numberOfStructures,

              pdbTemplate=outFilename,

              structLoopAction=calcOneStructure,

              genViolationStats=1,

              averagePotList=potList,

              averageTopFraction=0.5, #report only on best 50% of structs

              averageAccept=accept,   #only use structures which pass accept()

              averageContext=FinalParams(rampedParams),

              averageFilename="SCRIPT_ave.pdb",    #generate regularized ave structure

              averageFitSel="name CA",

              averageCompSel="not resname ANI and not name H*"     ).run()

Thanks for your support,

Santhosh

--
Marie-Laurence Tremblay, PhD candidate
Department of Biochemistry and Molecular Biology
Dalhousie University
5850 College St.
Halifax NS, Canada
B3H 1X5
mltremblay <at> dal.ca
Phone:902-494-4812

ATOM    856  X   ANI   500     322.383 -47.025 154.522  1.00  3.06          

ATOM    857  Y   ANI   500     322.080 -50.957 156.088  1.00  4.37          

ATOM    858  Z   ANI   500     323.128 -50.368 152.019  1.00  4.35          

ATOM    859  OO  ANI   500     324.205 -49.405 154.647  1.00  4.04          

Would this be fine as per VDW calculations or anything.  I am hoping that they could overlap as I have many RDC datasets and it would be tricky to create many such random points. 

Thanks a lot,

Santhosh