The following script for CNS simulated annealing was submitted by
trent. He recommends to use this script as a part of
RECOORD structure determination and water refinement protocol to prevent wrong geometry of aromatic rings (as reported
here and
here). You may find this script useful if you have a ligand in your model that makes CNS go crazy and produce protein structures with collapsed aromatic rings.
Thanks
trent!
Mark
************************************************** *****************
remarks file sa_cns.inp
remarks Simulated annealing protocol for NMR structure determination with cns.
remarks The starting structure for this protocol can be any structure with
remarks a reasonable geometry, such as randomly assigned torsion angles or
remarks extended strands.
remarks based on sa_new.inp from xplor
remarks high_steps and cool_steps increased 4 fold to improve convergence
remarks uses parallhdg5.3.pro and topallhdg5.3.pro from recoord without error
remarks works with additional ligand files
remarks Author: Trent Bjorndahl
{====>}
evaluate ($init_t = 1000 )
{====>}
evaluate ($high_steps= 24000 )
{====>}
evaluate ($cool_steps = 12000 )
parameter
{====>}
@parallhdg5.3.pro
!@yourparameterfilehere.par
end
{====>}
structure @yourstructurefile.mtf end {*Read the structure file.*}
{====>}
coordinates @yourcoordinatefile.pdb {*Read the coordinates.*}
noe
{====>}
nres=10000 {*Estimate greater than the actual number of NOEs.*}
class all
{====>}
@noes.tbl {*Read NOE distance ranges.*}
end
{====>}
restraints dihedral reset
nass = 1000
@dihedrals.tbl
end
flags exclude * include bonds angle impr vdw noe cdih end
{*Friction coefficient for MD heatbath, in 1/ps. *}
do (fbeta=10) (all)
{*Uniform heavy masses to speed molecular dynamics.*}
do (mass=100) (all)
noe {*Parameters for NOE effective energy term.*}
ceiling=1000
averaging * cent
potential * soft
scale * 50.
sqoffset * 0.0
sqconstant * 1.0
sqexponent * 2
soexponent * 1
asymptote * 0.1 {*Initial value--modified later.*}
rswitch * 0.5
end
parameter {*Parameters for the repulsive energy term.*}
nbonds
repel=1. {*Initial value for repel--modified later.*}
rexp=2 irexp=2 rcon=1.
nbxmod=3
wmin=0.01
cutnb=4.5 ctonnb=2.99 ctofnb=3.
tolerance=0.5
end
end
restraints dihedral
scale=5.
end
{====>}
evaluate ($end_count=10) {*Loop through a family of 10 structures.*}
coor copy end
evaluate ($count = 0)
while ($count < $end_count ) loop main
evaluate ($count=$count+1)
coor swap end
coor copy end
{* Initial minimization.*}
restraints dihedral scale=5. end
noe asymptote * 0.1 end
parameter nbonds repel=1. end end
igroup interaction
(all) (all) weights * 1 vdw 0.002 end end
minimize powell nstep=50 drop=10. nprint=25 end
{* High-temperature dynamics.*}
igroup interaction (all) (all)
weights * 1 angl 0.4 impr 0.1 vdw 0.002 end end
evaluate ($nstep1=int($high_steps * 2. / 3. ) )
evaluate ($nstep2=int($high_steps * 1. / 3. ) )
dynamics cartesian
nstep=$nstep1
timestep=0.003
temperature=$init_t
tcoup=true
nprint=50
ntrfrq=9999
cmremove=true
cmperiodic=0
end
{* Tilt the asymptote and increase weights on geometry.*}
noe asymptote * 1.0 end
igroup interaction
(all) (all) weights * 1 vdw 0.002 end end
{* Bring scaling factor for S-S bonds back *}
parameter
bonds ( name SG ) ( name SG ) 1000. TOKEN
angle ( name CB ) ( name SG ) ( name SG ) 500. TOKEN
end
dynamics cartesian
nstep=$nstep2
timestep=0.003
tcoup=true
temperature=$init_t
nprint=50
ntrfrq=9999
cmremove=true
cmperiodic=0
end
{* Cool the system.*}
restraints dihedral scale=200. end
evaluate ($final_t = 100) { K }
evaluate ($tempstep = 50) { K }
evaluate ($ncycle = ($init_t-$final_t)/$tempstep)
evaluate ($nstep = int($cool_steps/$ncycle))
evaluate ($ini_rad = 0.9) evaluate ($fin_rad = 0.75)
evaluate ($ini_con= 0.003) evaluate ($fin_con= 4.0)
evaluate ($bath = $init_t)
evaluate ($k_vdw = $ini_con)
evaluate ($k_vdwfact = ($fin_con/$ini_con)^(1/$ncycle))
evaluate ($radius= $ini_rad)
evaluate ($radfact = ($fin_rad/$ini_rad)^(1/$ncycle))
evaluate ($i_cool = 0)
while ($i_cool < $ncycle) loop cool
evaluate ($i_cool=$i_cool+1)
evaluate ($bath = $bath - $tempstep)
evaluate ($k_vdw=min($fin_con,$k_vdw*$k_vdwfact))
evaluate ($radius=max($fin_rad,$radius*$radfact))
parameter nbonds repel=$radius end end
igroup interaction (all) (all)
weights * 1. vdw $k_vdw end end
dynamics cartesian
nstep=$nstep
time=0.005
tcoup=true
temperature=$bath
nprint=$nstep
ntrfrq=9999
cmremove=true
cmperiodic=0
end
{====>} {*Abort condition.*}
evaluate ($critical=$temp/$bath)
if ($critical > 10. ) then
display ****&&&& rerun job with smaller timestep (i.e., 0.003)
stop
end if
end loop cool
{* ================================================= Final minimization.*}
igroup interaction (all) (all) weights * 1. vdw 1. end end
parameter
nbonds
repel=0.80
rexp=2 irexp=2 rcon=1.
nbxmod=3
wmin=0.01
cutnb=6.0 ctonnb=2.99 ctofnb=3.
tolerance=1.5
end
end
minimize powell nstep=1000 drop=10.0 nprint=25 end
{* =================================== Write out the final structure(s).*}
print threshold=0.5 noe
evaluate ($rms_noe=$result)
evaluate ($violations_noe=$violations)
print threshold=5. cdih
evaluate ($rms_cdih=$result)
evaluate ($violations_cdih=$violations)
print thres=0.05 bonds
evaluate ($rms_bonds=$result)
print thres=5. angles
evaluate ($rms_angles=$result)
print thres=5. impropers
evaluate ($rms_impropers=$result)
remarks ================================================== =============
remarks overall,bonds,angles,improper,vdw,noe,cdih
remarks energies: $ener, $bond, $angl, $impr, $vdw, $noe, $cdih
remarks ================================================== =============
remarks bonds,angles,impropers,noe,cdih
remarks rms-d: $rms_bonds,$rms_angles,$rms_impropers,$rms_noe,$rm s_cdih
remarks ================================================== =============
remarks noe, cdih
remarks violations.: $violations_noe, $violations_cdih
remarks ================================================== =============
{====>} {*Name(s) of the family of final structures.*}
evaluate ($filename="sa_cns_"+encode($count)+".pdb")
write coordinates output =$filename end
end loop main
stop