Result table
| image | mrblock_id | pdb_id | bmrb_id | cing | stage | position | program | type |
|
|
516915 | 2k3b RC | 16843 | cing | 1-original | 4 | XPLOR/CNS | unknown |
xplor.requireVersion("2.19")
# refine random loop structures
# pv 1-24-08
# Changed fixed atoms 2-14-07
xplor.parseArguments() # check for typos on the command-line
pvfile=open('refine_wlooprest_peg.current')
structnumst=pvfile.readline()
pvfile.close()
structnum=int(structnumst)
structname = "refinedloop_"+str(structnum)
outFilename ="refined_wlooprest_peg/" + structname + "_STRUCTURE.pdb"
outFilestart = "refined_wlooprest_peg/start/" + structname + ".pdb"
inFilename = "random_loop_pdbs/rand_loop_"+str(structnum)+".pdb"
numberOfStructures=96
# Stop is structnum is more than 10
if ( structnum > 10):
import sys
sys.exit(12)
# protocol module has many high-level helper functions.
import protocol
protocol.initRandomSeed(5324+(structnum*652)) #explicitly set random seed
# annealing settings
command = xplor.command
protocol.initParams("protein")
# Read PSF from file
protocol.initStruct("abpsh3.psf")
#read PDB
protocol.initCoords(inFilename)
#protocol.initCoords("abpsh3_xplor_psfgenerate.pdb")
# a PotList conatins 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, LinRamp
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={}
# medium Da Rh (are allowed to vary)
for (medium,Da,Rh) in [ ('pf1', -10.0, 0.35),
('peg', -7.37, 0.605)]:
oTensor = create_VarTensor(medium)
oTensor.setDa(Da)
oTensor.setRh(Rh)
oTensor.setDaMax(20) # Set DaMax
oTensor.setFreedom('varyDa,varyRh')
media[medium] = oTensor
pass
# dipolar coupling restraints for protein.
from rdcPotTools import create_RDCPot, scale_toNH
rdcs = PotList('rdc')
for (medium,expt,file, scale) in \
[('pf1','NH' , 'pf1_NH_scalephage_edit_scaleNH.tbl' ,1),
('pf1','TRP' , 'pf1_TRP_scalephage_scaleNH.tbl' ,1),
('pf1','CAHA','pf1_CAHA_scalephage_scaleNH.tbl' ,1),
('pf1','HNC', 'pf1_HNCO_scalephage_scaleNH.tbl' ,1),
('peg', 'NH', 'cpmgpeg.tbl', 1),
('peg', 'TRP','cpmgpeg_TRP.tbl', 1),
]:
rdc = create_RDCPot("%s_%s"%(medium,expt),file,media[medium])
rdc.setScale(scale)
rdc.setShowAllRestraints(1)
rdc.setThreshold(2.0) # in Hz
rdc.setPotType("square")
rdcs.append(rdc)
pass
potList.append(rdcs)
rampedParams.append( MultRamp(0.002,0.4, "rdcs.setScale( VALUE )") )
# Csa file
from csaPotTools import create_CSAPot
csas = PotList("csa")
for (medium,expt,file, scale) in \
[('pf1','CO', 'pf1_COcsa_scalephage_18off.tbl' ,1),
]:
csa = create_CSAPot("%s_%s"%(medium,expt),file,media[medium])
csa.setScale(scale)
csa.setShowAllRestraints(1)
csa.setThreshold(10.0) # in ppb
csa.setPotType("square")
csas.append(csa)
pass
potList.append(csas)
# rampedParams.append( MultRamp(0.00000333333333333333, 0.00066666666666666666, "csas.setScale( VALUE )") )
rampedParams.append( MultRamp(0.00005, 0.01, "csas.setScale( VALUE )") )
# calc initial Tensor orientation
from varTensorTools import calcTensorOrientation
for medium in media.values():
calcTensorOrientation(medium)
pass
# Set up dihedral angles
from xplorPot import XplorPot
# protocol.initDihedrals("talos_xplor_extra_wval_nottight.tab",useDefaults=0)
protocol.initDihedrals(["loop1/loop1_restraints.tbl","loop2_no31/loop2_no31_restraints.tbl","loop3/loop3_restraints.tbl","Da_pos_dihed.tbl"],useDefaults=0)
potList.append( XplorPot('CDIH') )
highTempParams.append( StaticRamp("potList['CDIH'].setScale(5)") )
#rampedParams.append( StaticRamp("potList['CDIH'].setScale(400)") )
#rampedParams.append( LinRamp(10,400, "potList['CDIH'].setScale( VALUE )") )
rampedParams.append( MultRamp(5,400, "potList['CDIH'].setScale( VALUE )") )
# set custom values of threshold values for violation calculation
potList['CDIH'].setThreshold( 10 ) #5 degrees is the default value, though
#Rama torsion angle database
protocol.initRamaDatabase()
potList.append( XplorPot('RAMA') )
rampedParams.append( MultRamp(.002,1,"potList['RAMA'].setScale(VALUE)") )
potList.append( XplorPot("BOND") )
potList.append( XplorPot("ANGL") )
potList['ANGL'].setThreshold( 5 )
rampedParams.append( MultRamp(0.4,1,"potList['ANGL'].setScale(VALUE)") )
potList.append( XplorPot("IMPR") )
potList['IMPR'].setThreshold( 5 )
rampedParams.append( MultRamp(0.1,1,"potList['IMPR'].setScale(VALUE)") )
# Set up LJ
from xplorSimulation import getXplorSimulation
def initLJ(wmin=0.01,
tolerance=0.5,
nbxmod=3,
cutnb=10,
ctonnb=8,
ctofnb=9,
selStr="known",
onlyCA=0,
simulation=0):
"""standard initialization of the non-bonded repel potential. The XPLOR
nonbonded potential term is described <l http://nmr.cit.nih.gov/xplor-nih/xplorMan/node117.html
here>
note that cutnb should be greater than rmax + 2*tolerance, where rmax
is the largest vdw radius.
"""
noSelStr="resname ANI" #exclude these nonbonded interactions
if onlyCA: selStr="name CA"
xSim = getXplorSimulation(simulation)
xSim.command("""
constraints
interaction (%s and (not (%s))) (%s and (not (%s)))
weights * 1 vdw 1 end
interaction (not (%s) and (not (%s))) (not (%s))
weights * 1 vdw 0 end
end""" % (selStr,noSelStr,selStr,noSelStr,selStr,noSelStr,noSelStr) )
xSim.command("""
parameters
nbonds
atom
wmin = %f ! warning off
tolerance = %f
nbxmod = %d
cutnb = %f ! nonbonded cutoff
ctonnb = %f ! switching on
ctofnb = %f ! switching off
rdie
vswitch
switch
end
end
""" % (wmin,tolerance,nbxmod,cutnb,ctonnb,ctofnb) )
potList.append( XplorPot('VDW') )
rampedParams.append( StaticRamp("initLJ()") )
# Give atoms uniform weights, except for the anisotropy axis
from atomAction import SetProperty
import varTensorTools
AtomSel("not resname ANI").apply( SetProperty("mass",100.) )
varTensorTools.massSetup(media.values(),300)
AtomSel("all ").apply( SetProperty("fric",10.) )
# IVM setup
from ivm import IVM
dyn = IVM()
dyn.reset()
dyn.fix(AtomSel("resid 2:4 or resid 20:28 or resid 41:47 or resid 57:59"))
for m in media.values():
m.setFreedom("varyDa, varyRh") #vary tensor Rh, Da, vary orientation
protocol.torsionTopology(dyn,oTensors=media.values())
# object which performs simulated annealing
from simulationTools import AnnealIVM
init_t = 5000. # Need high temp and slow annealing to converge
cool = AnnealIVM(initTemp =init_t,
finalTemp=25,
tempStep =12.5,
ivm=dyn,rampedParams = rampedParams)
# minc used for final cartesian minimization
# minc = IVM()
# minc.fix(AtomSel("resid 2:4 or resid 20:28 or resid 41:47 or resid 57:59"))
# protocol.initMinimize(minc)
# for m in media.values():
# m.setFreedom("varyDa, varyRh") #allow all tensor parameters float here
# pass
# protocol.cartesianTopology(minc,oTensors=media.values())
def accept(potList):
"""
return True if current structure meets acceptance criteria
"""
# if potList['CDIH'].violations()>0:
# return False
return True
from pdbTool import PDBTool
def calcOneStructure(loopInfo):
# initialize parameters for high temp dynamics.
InitialParams( rampedParams )
# command("restraints dihedral scale=25. end")
# high temp dynamics
protocol.initDynamics(dyn,
potList=potList, # potential terms to use
bathTemp=init_t,
initVelocities=1,
finalTime=100, # stops at 100ps
printInterval=1000)
dyn.setETolerance( init_t/100 ) #used to det. stepsize. default: t/1000
dyn.run()
# initialize integrator for simulated annealing
protocol.initDynamics(dyn,
potList=potList,
finalTime=2.0 , # 2ps
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
# Write starting structure
PDBTool(outFilestart).write()
from simulationTools import StructureLoop, FinalParams
StructureLoop(numStructures=numberOfStructures,
pdbTemplate=outFilename,
structLoopAction=calcOneStructure,
genViolationStats=1).run()