Molecular DynamicsIntegrating equations of motionVerlet AlgorithmLeapfrog, Velocity VerletCHARMMLoop flexing in triose phosphate isomeraseUnfolding of ß2-microglobulinCHARMM resourcesScriptsSlide 10Non-bonded cutoffsSimulating Water EnvironmentExplicit watersWater MoleculesSlide 15Slide 16Thermodynamic ensemblesMolecular Dynamics•Verlet - integrating equations of motion•Langevin dynamics – simulating water•Monte Carlo•Brownian•Newton’s 2nd law•coupled first-order ODEs (IVP)–many-body problem; simulation•Runge-Kutta (4th-order) –4th-order: error decays with step size as O(h5)Integrating equations of motioncooling of a sphere:Verlet Algorithm•write forward and backward Taylor expansions to 3rd-order:•calculate a from forces (deriv. of poten. energy):•estimate velocities from positions:•conservation of energy: •error is O(t2)•nice property: symplectic integrator–conserves time-evolution of (perturbed) Hamiltonian–conserves volume in phase spaceLeapfrog, Velocity Verlet•Velocity Verlet: calculate v from a (mid-step):•try to reduce errors to O(t4)•LeapfrogCHARMM•procedure–build (H’s), solvate?–minimize–heat–equilibrate–production run, 100ns / ~1fs = 100,000 iterations–(quenching)–analyze trajectory•minima (MC)•events and energy changes•time-correlations, spectra•Monte Carlo simulations–statistical mechanics, thermodynamic averages–kinetic energyLoop flexing in triose phosphate isomerase•Derreumaux and Schlick (1998)•changes in energy correlated with 2 H-bonds...•importance of solvation (waters) –1200K LD in vac., 1 opening event, 25 kcal/mol; interactions with water reduces the energy barrierUnfolding of ß2-microglobulin•Ma and Nussinov (2003)•amyloid formation•determine which sub-structures (beta-strands) dissociate first•compare stability of N-terminal truncation and proteolytic cleavagenative, 400Knative, 500Ktrunc, 400KCHARMM resources•tutorials–http://www.ch.embnet.org/MD_tutorial/–http://www.esi.umontreal.ca/accelrys/life/insight2000.1/charmm_principles/Ch04_mol_dyn.FM5.html–http://www.mdy.univie.ac.at/lehre/charmm/course1/course1-1.html–http://www.charmmtutorial.org/index.php/CHARMM:The_Basics•documentation (current version c35b1) –http://www.charmm.org/html/documentation/c34b1/index.html•180-page online book by Robert Schleif–http://gene.bio.jhu.edu/book.pdfScripts• charmm < minimize.inp > minimize.out•PAR/TOP param files; PSF files (PDBs have coords, not connectivity)•read write files: units (fortran)•read the log files! (fix any errors)•selection syntax• copy coor comp•trajectory files•grep for “DYNA>”; monitor energy, temp in log file output•time-correlation statistics (CORRel) – dihedrals, rms, energy...Non-bonded cutoffs•ignore VDW, ES interactions beyond ~10Å–most forces are small at that range–greatly reduces run-time–may causes artifacts (non-smooth); may leak energy•Particle-Mesh Ewald summation (PME)–allows calculation of total electrostatic energy (including long-range interactions) assuming periodic boundary conditions–model as lattice of protein-water boxes (like crystal)–use Fast Fourier-transform (in reciprocal space) to calculate product of potential at different points in unit cell with partial charges at all other pointsNBOND VSWITCH...CUTNb 16.0 CTOFnb 12. CDIE EPS 1. CTONnb 8.Simulating Water Environment•implicit: distance-dependent dielectric (rdie=1)–extra 1/r in E:–dampens strength•Langevin dynamics–simulate friction of water via dissipative + random forces on surface atoms– /m= collision frequency, relaxation time–fluctuation-dispersion theorem:–in the limit of =0, get Brownian dynamicsExplicit waters•explicit waters–each water adds 6-9 degrees of freedom–water molecules strike side-chains on surface, impart random forces, exchange energy–dipoles interaction with surface atoms (H-bond)–add counter-ions to neutralize system•water box–periodic boundary conditions–pressure/volume changesWater Molecules•Solvation–polarity (H-bonds)–density? velocity? residence times•Important for surface side-chains to not be in vacuum; helps them flex properly and reduce self-interactions •SPC model–3-points, H’s rigid (no bond stretch or bend)–6-dof rather than 9–VDW interaction site centered on O•TIP3P model (Jorgensen et al., 1983)–partial charges: q(O)=-0.83, q(H)=+0.415•TIP4P–separate center of O charge from VDW centerTIP4PiceSPC/E(Caleman)diffusion, delectric permittivity, vibration spectra...•quenching: lower temp or minimize at end of trajectory•SHAKE: constrain bonds to H’s, decreases dof•temperature drift (losses due to NB approx.) –temperature scaling (artificial)–couple to heat bath –Nose-Hoover (Evans and Holian, 1985)•re-write equations of motion so T is a constant rather than E (scale momenta and time)Thermodynamic ensembles•must generate proper ensemble for calculating thermodynamic averages (Monte Carlo)–standard simulations conserve energy (NVE), not temp–but proper (canonical) ensembles (NVT) have constant temp and variable energy; must couple to a heat bath•NVT: canonical ensemble –constant temperature, Helmholz free energy•NVE: micro-canonical ensemble –constant energy (microstate)–omega is count of microstates (related to entropy)•NPT (cpt) simulations –constant pressure (isobaric) –allow size of box to expand/contract (additional parameter)–Gibb’s free energy,
View Full Document
Unlocking...