remarks Solvent refinement protocol from ARIA1.2 (Nilges and Linge), modified for XPLOR-NIH remarks Site the following references: set message on echo on end {*==========================================================================*} {*====================== SET FILENAMES AND VARIABLES =======================*} {*==========================================================================*} {* type of non-bonded parameters *} {* choice: "PROLSQ" "PARMALLH6" "PARALLHDG" "OPLSX" *} {* The water refinement uses the OPLSX parameters *} evaluate ( $par_nonbonded = "OPLSX" ) {* parameter file(s) *} evaluate ( $par.1 = "TOPPAR:parallhdg5.3.pro" ) evaluate ( $par.2 = "TOPPAR:parallhdg5.3.sol" ) evaluate ( $par.3 = "" ) evaluate ( $par.4 = "" ) evaluate ( $par.5 = "" ) {* solvent topology file *} evaluate ( $solvent_topology = "TOPPAR:topallhdg5.3.sol" ) {* structure file(s) *} evaluate ( $struct.1 = "STRUCTURES:1xxx.psf" ) evaluate ( $struct.2 = "" ) evaluate ( $struct.3 = "" ) evaluate ( $struct.4 = "" ) evaluate ( $struct.5 = "" ) {* input coordinate file(s) *} evaluate ( $pdb.in.file.1 = "INPUTCOORDINATES:analyzed_3.pdb" ) {* output coordinate file(s) *} evaluate ( $pdb.out.file.1 = "OUTPUTCOORDINATES:refined_3.pdb" ) {* input distance restraint file(s) *} evaluate ( $noe.file.1 = "TABLES:1xxx_noe.tbl" ) {* Averaging for NOE restraints *} evaluate ( $noe.ave = sum ) {* input dihedral restraint file(s) *} evaluate ( $cdih.file.1 = "TABLES:1xxx_dihe.tbl" ) {* seed for velocity generation *} evaluate ( $seed = 12396 ) {* set number of md steps for the heating stage *} evaluate ( $mdsteps.heat = 200 ) {* set number of md steps for the hot md stage *} evaluate ( $mdsteps.hot = 2000 ) {* set number of md steps for the cooling stage *} evaluate ( $mdsteps.cool = 200 ) {* seed for velocity generation *} evaluate ( $seed = 12396 ) {* set acceptance criteria *} evaluate ( $accept.noe = 0.50 ) evaluate ( $accept.cdih = 5.00 ) evaluate ( $accept.coup = 1.00 ) evaluate ( $accept.sani = 0.00 ) evaluate ( $accept.vean = 5.00 ) {*==========================================================================*} {*=========== READ THE STRUCTURE, TOPOLOGY AND PARAMETER FILES =============*} {*==========================================================================*} structure @$struct.1 end topology @$solvent_topology end parameter @$par.1 @$par.2 end {*==========================================================================*} {*================== SET VALUES FOR NONBONDED PARAMETERS ===================*} {*==========================================================================*} parameter nbonds nbxmod=5 atom cdiel shift cutnb=9.5 ctofnb=8.5 ctonnb=6.5 eps=1.0 e14fac=0.4 inhibit 0.25 wmin=0.5 tolerance 0.5 end end {*==========================================================================*} {*============== READ THE COORDINATES AND COPY TO REFERENCE SET ============*} {*==========================================================================*} coor @@$pdb.in.file.1 vector do (refx = x) (all) vector do (refy = y) (all) vector do (refz = z) (all) {*==========================================================================*} {*========================= GENERATE THE WATER LAYER =======================*} {*==========================================================================*} vector do (segid = "PROT") (segid " ") @SOLVENT:generate_water.inp vector do (segid = " ") (segid "PROT") {*==========================================================================*} {*========================= READ THE EXPERIMENTAL DATA =====================*} {*==========================================================================*} noe reset nrestraints = 10000 ceiling = 100 @$noe.file.1 averaging * $noe.ave potential * soft scale * 50 sqconstant * 1.0 sqexponent * 2 soexponent * 1 rswitch * 1.0 sqoffset * 0.0 asymptote * 2.0 end restraints dihedral reset @$cdih.file.1 scale = 200 end #{ncs constraints for symmetric dimer} #evaluate ($kncs=0.1) #ncs restraints #initialize #group #equi (resid 6:15 or resid 24:59 or resid 46:61 or resid 66:72 or resid 78:88) #equi (resid 306:315 or resid 324:359 or resid 346:361 or resid 366:372 or resid 378:388) #weight = $kncs #end #? #end {*==========================================================================*} {*============================ SET FLAGS ===================================*} {*==========================================================================*} flags exclude * include bond angle dihe impr vdw elec noe cdih coup oneb carb ncs dani vean sani prot harm end {*==========================================================================*} {*========================= START THE REFINEMENT ===========================*} {*==========================================================================*} set seed $seed end ! We loop untill we have an accepted structure, maximum trials=3 evaluate ($end_count = 3) evaluate ($count = 0) while ($count < $end_count ) loop main evaluate ($accept = 0) ! since we do not use SHAKe, increase the water bond angle energy constant parameter angle (resn tip3) (resn tip3) (resn tip3) 500 TOKEN end ! reduce improper and angle force constant for some atoms evaluate ($kbonds = 1000) evaluate ($kangle = 50) evaluate ($kimpro = 5) evaluate ($kchira = 5) evaluate ($komega = 5) parameter angle (not resn tip3)(not resn tip3)(not resn tip3) $kangle TOKEN improper (all)(all)(all)(all) $kimpro TOKEN TOKEN end ! fix the protein for initial minimization constraints fix (not resn tip3) end minimize powell nstep=40 drop=100 end ! release protein and restrain harmonically constraints fix (not all) end vector do (refx=x) (all) vector do (refy=y) (all) vector do (refz=z) (all) restraints harmonic exponent = 2 end vector do (harmonic = 0) (all) vector do (harmonic = 10) (not name h*) vector do (harmonic = 20.0)(resname ANI and name OO) vector do (harmonic = 0.0) (resname ANI and name Z ) vector do (harmonic = 0.0) (resname ANI and name X ) vector do (harmonic = 0.0) (resname ANI and name Y ) constraints interaction (not resname ANI) (not resname ANI) interaction ( resname ANI) ( resname ANI) end minimize powell nstep=40 drop=10 end vector do (refx=x) (not resname ANI) vector do (refy=y) (not resname ANI) vector do (refz=z) (not resname ANI) minimize powell nstep=40 drop=10 end vector do (refx=x) (not resname ANI) vector do (refy=y) (not resname ANI) vector do (refz=z) (not resname ANI) vector do (mass =50) (all) vector do (mass=1000) (resname ani) vector do (fbeta = 0) (all) vector do (fbeta = 20. {1/ps} ) (not resn ani) evaluate ($kharm = 50) ! heat to 500 K for $bath in (100 200 300 400 500 ) loop heat vector do (harm = $kharm) (not name h* and not resname ANI) vector do (vx=maxwell($bath)) (all) vector do (vy=maxwell($bath)) (all) vector do (vz=maxwell($bath)) (all) dynamics verlet {time step was 3 fs} nstep=$mdsteps.heat timest=0.003{ps} tbath=$bath tcoupling = true iasvel=current nprint=50 end evaluate ($kharm = max(0, $kharm - 4)) vector do (refx=x) (not resname ANI) vector do (refy=y) (not resname ANI) vector do (refz=z) (not resname ANI) end loop heat ! refinement at high T: constraints interaction (not resname ANI) (not resname ANI) weights * 1 dihed 3 end interaction ( resname ANI) ( resname ANI) weights * 1 end end vector do (harm = 0) (not resname ANI) vector do (vx=maxwell($bath)) (all) vector do (vy=maxwell($bath)) (all) vector do (vz=maxwell($bath)) (all) dynamics verlet nstep=$mdsteps.hot timest=0.004 {ps} tbath=$bath tcoupling = true iasvel=current nprint=50 !trajectory=1xxx_hot.dat nsavc=5 end constraints interaction (not resname ANI) (not resname ANI) weights * 1 dihed 5 end interaction ( resname ANI) ( resname ANI) weights * 1 end end ! cool evaluate ($bath = 500) while ($bath >= 25) loop cool evaluate ($kbonds = max(250,$kbonds / 1.1)) evaluate ($kangle = min(250,$kangle * 1.1)) ! equi-ed after 17 steps evaluate ($kimpro = min(300,$kimpro * 1.4)) ! 12 evaluate ($kchira = min(1000,$kchira * 1.4))! 16 evaluate ($komega = min(100,$komega * 1.4)) ! 9 parameter bond (not resn tip3 and not name H*)(not resn tip3 and not name H*) $kbonds TOKEN angle (not resn tip3 and not name H*)(not resn tip3 and not name H*)(not resn tip3 and not name H*) $kangle TOKEN improper (all)(all)(all)(all) $kimpro TOKEN TOKEN !VAL: stereo CB improper (name HB and resn VAL)(name CA and resn VAL)(name CG1 and resn VAL)(name CG2 and resn VAL) $kchira TOKEN TOKEN !THR: stereo CB improper (name HB and resn THR)(name CA and resn THR)(name OG1 and resn THR)(name CG2 and resn THR) $kchira TOKEN TOKEN !LEU: stereo CG improper (name HG and resn LEU)(name CB and resn LEU)(name CD1 and resn LEU)(name CD2 and resn LEU) $kchira TOKEN TOKEN !ILE: chirality CB improper (name HB and resn ILE)(name CA and resn ILE)(name CG2 and resn ILE)(name CG1 and resn ILE) $kchira TOKEN TOKEN !chirality CA improper (name HA)(name N)(name C)(name CB) $kchira TOKEN TOKEN improper (name O) (name C) (name N) (name CA) $komega TOKEN TOKEN improper (name HN) (name N) (name C) (name CA) $komega TOKEN TOKEN improper (name CA) (name C) (name N) (name CA) $komega TOKEN TOKEN improper (name CD) (name N) (name C) (name CA) $komega TOKEN TOKEN end vector do (vx=maxwell($bath)) (all) vector do (vy=maxwell($bath)) (all) vector do (vz=maxwell($bath)) (all) dynamics verlet !was a time step of 0.004 nstep=$mdsteps.cool timest=0.004{ps} tbath=$bath tcoupling = true iasvel=current nprint=50 end evaluate ($bath = $bath - 25) end loop cool !final minimization: mini powell nstep 200 end {*==========================================================================*} {*======================= CHECK RESTRAINT VIOLATIONS =======================*} {*==========================================================================*} print threshold=$accept.noe noe evaluate ($v_noe = $violations) print threshold=0.5 noe evaluate ($v_noe_0.5 = $violations) print threshold=0.4 noe evaluate ($v_noe_0.4 = $violations) print threshold=0.3 noe evaluate ($v_noe_0.3 = $violations) print threshold=0.2 noe evaluate ($v_noe_0.2 = $violations) print threshold=0.1 noe evaluate ($v_noe_0.1 = $violations) evaluate ($rms_noe = $result) print threshold=$accept.cdih cdih evaluate ($rms_cdih=$result) evaluate ($v_cdih = $violations) coupl print thres = $accept.coup class c1 end evaluate ($rms_coup = $result) evaluate ($v_coup = $violations) sani print threshold = $accept.sani class rdc1 end evaluate ($rms_sani = $result) evaluate ($v_sani = $violations) vean print threshold = $accept.vean class vea1 end evaluate( $rms_vean = $result) evaluate( $v_vean = $violations) print thres=0.05 bonds evaluate ($rms_bonds=$result) evaluate ($v_bonds = $violations) print thres=5. angles evaluate ($rms_angles=$result) evaluate ($v_angles = $violations) print thres=5. impropers evaluate ($rms_impropers=$result) evaluate ($v_impropers = $violations) if ($v_noe > 0) then evaluate ( $accept=$accept + 1 ) end if if ($v_cdih > 0) then evaluate ( $accept=$accept + 1 ) end if if ($v_coup > 0) then evaluate ( $accept=$accept + 1 ) end if if ($v_sani > 0) then evaluate ( $accept=$accept + 1 ) end if if ($v_vean > 0) then evaluate ( $accept=$accept + 1 ) end if if ($accept = 0 ) then evaluate ( $label = "ACCEPTED" ) exit main else evaluate ( $label = "NOT ACCEPTED" ) evaluate ( $count = $count + 1 ) end if end loop main {*==========================================================================*} {*================ COMPUTE RMS DIFFERENCE BETWEEN OBSERVED =================*} {*========================= AND MODEL DISTANCES. ===========================*} {*==========================================================================*} constraints interaction (not resname TIP* and not resname ANI) (not resname TIP* and not resname ANI) end energy end remarks Structure $label remarks E-overall: $ener remarks E-NOE_restraints: $noe remarks E-CDIH_restraints: $cdih remarks E-COUP_restraints: $coup remarks E-SANI_restraints: $sani remarks E-VEAN_restraints: $vean remarks RMS-NOE_restraints: $rms_noe remarks RMS-CDIH_restraints: $rms_cdih remarks RMS-COUP_restraints: $rms_coup remarks RMS-SANI_restraints: $rms_sani remarks RMS-VEAN_restraints: $rms_vean remarks NOE Acceptance criterium: $accept.noe remarks NOE violations > $accept.noe: $v_noe remarks All NOE viol. (>0.5,0.4,0.3,0.2,0.1A): $v_noe_0.5 $v_noe_0.4 $v_noe_0.3 $v_noe_0.2 $v_noe_0.1 remarks CDIH Acceptance criterium: $accept.cdih remarks CDIH violations > $accept.cdih: $v_cdih remarks COUP Acceptance criterium: $accept.coup remarks COUP violations > $accept.coup: $v_coup remarks SANI Acceptance criterium: $accept.sani remarks SANI violations > $accept.sani: $v_sani remarks VEAN Acceptance criterium: $accept.vean remarks VEAN violations > $accept.vean: $v_vean write coordinates sele= (not resn TIP3) output =$pdb.out.file.1 end !write coordinates sele= (all) output =$pdb.out.file.1 end stop