`
`2773
`
`Molecular Dynamics of Proteins with the OPLS Potential
`Functions. Simulation of the Third Domain of Silver Pheasant
`Ovomucoid in Water
`
`Julian Tirado·Rives and \Villiam L. Jorgensen•
`
`Contribution from the Department of Chemistry, Purdue Unieersity,
`West Lafayerte, Indiana 47907. Received }.fay 15. 1989
`
`Abstract: A molecular dynamics simulation using the OPLS nonbonded potential functions has been carried out for the third
`domain of silver phe.asant ovomucoid in aqueous solution. Insights have been obtained on the quality of the force field, the
`convergence of such calculations, differentts in the protein's structure in the crystal and in aqueous solution, protein hydration,
`and the dynamics of water molecules near a protein. The simulation covered JOO ps at 25 °C, which allowed complete equilibration
`prior to averaging and analysis of the results. Continuous monitoring of the potential energy, root·mean·square deviations
`from the crystal structure, and other properties indicated that convergence to a stable structure was achieved after 30-40 ps.
`The RMS deviation of the instantaneous structure from the crystal structure after 100 psis 1.43 A for the backbone atoms
`of residues 8-56 and 1.61 A for all residues. There is substantial reorgani1.ation of hydrogen bonds that do not invoke secondary
`structure in comparing the crystal and solution structure>, though in the simulation Ala·44 is displaced from the a·helix and
`L)S·29, Thr·30, Tyr·3 l, and Gly·32 arc moved out of hydrogen·bonding distance in the triple·stranded antiparallel i'.i·shcct.
`Analyse> of the protein-water hydrogen bonding were also carried out and are compared with results from previous simulations
`and :'\MR experiments.
`
`The rational design or modification of biomolecule£, including
`the development of selective inhibitors for enzymes, requires
`detailed knowledge of the structure, dynamics, and corresponding
`energetics, Importantly, the continuous improvement of crys·
`tal\ographic techniques has made possible the precise determi·
`nation of the structures of many proteins, as reflected in the more
`than JOO entries now deposited in the Brookhaven Protein Data
`Bank.I However. the structures obtained in this fashion represent
`an ordered crystalline state, while biological processes normally
`occur in solution. Furthermore, the data obtained from crystal
`structures arc static in nature, although some dynamic information
`can be obtained from the temperature (8) factors. 2 Recent
`advances in nuclear magnetic resonance spectroscopy (NMR),
`particularly the use of 2-D nuclear Overhauser effects, have been
`very valuable, since sets of distance constraints are obtained that
`can be transformed into three-dimensional structures of proteins
`in solution.l Although this methodology allows the direct ob·
`servation of proteins in their native solution state, the structures
`obtained reflect conformational averaging and are not unique
`solutions for the data sets. Nevertheless, the combination of the
`experimental structural results with molecular dynamics (1'-.·tD)
`calculations is proving to be a powerful approach to the detailed
`characterization of the structure, dynamics, and energetics of
`proteins.'
`Since the pioneering MD simulation of bovine pancreatic tryp:;in
`inhibitor (BPTI) in vacuo,5 there have been numerous molecular
`dynamics calculations of proteins.~ However, the aqueous medium
`has rarely been represented in molecular detail. Some exceptions
`are for BPTl,6-1 avian pancreatic polypeptide hormone (APP),9
`
`(I) Bernncin, F. C.; Koetzle, T. F.; Williams, G. J.B.; ?.feyer, E. f., Jr.;
`Brice, M. D.; Rodgen, J. R.; Kennard, O.; Shimanouehi, T.; Tarnmi. M. J.
`Mo/. Bio/, 1977, 111, 5J5.
`(2) Peuko, G. A.; Ringe, D. A1111u. Rev. Blophys. Bl~ng. 1984, /J, 331.
`Fraunfc!der, H.; Hartmann, ll; Karp!ui, hf.; Kuntz, I. D., Jr.: Kuriyan, J,;
`Parak, f.; Petsko, G. A.; Ringe, D.; Tilton, R. A., Jr.; Connolly, M. L.; Max.
`N. Biod1fmiJ1ry 1987, 16, 254.
`(3) \\'iithrich, K. NMR of fro1ti111 and Nucleic AddJ; John \\'iky &
`Sons: New York, 1986.
`(4) For reo::ent reviews, see: (a) Karplus, M.; McCammon, J. A. CRCCri1.
`Rev. Bixhem. 19&6. 9, 293. {b) McCammon, J. A.; Har.,.ey, S. C. Dynomii:s
`of ProliilU and Nudeic Acids; Cambridge University Pre>s: Cambridge,
`England. 1987.
`(5) McCammon, J. A.; Gelin, B. R.: Karplu,, M. Nature 1977, 167, SSS.
`(6) 1·an Gun1teren, W. F.; Bcrend~n. H.J. C. J. Mo/. Biol. 1984, 176, 559.
`(7) Ghosh, I.; McCammon, J. A. J, Phys. Chem. 1987, 91, 4878.
`(8) Le1itt, M.; Sharon, R. Proc. Not/. Acod. Sci. US.A. 1988, 8S. 1551.
`
`and parvalbuminlO in aqueous solution and BPTI in its full
`crystalline environment. 11 Other recent calculations have been
`more focused toward the mcdeling of enzyme-inhibitor complexes
`in water, including trypsin-beniamidine12 and thermolysin(cid:173)
`phosphonamidate.ll
`\Vith the exceptions of the parvalbumin10 and one of the BPTI
`calculations,8 all of these simulations were run for very short times.
`Total times including the equilibrium periods were 15-30 ps for
`some of the BPTI and the APP simulations and 45 ps for the
`trypsin-benzamidine complex, as compared to 106 and 210 ps
`for the parvalbumin and the longest of the BPTI calculations,
`respectively. It is not clear that the shorter simulation times allow
`the systems to achieve equilibrium and to remove the biases from
`the starting conformation, typically obtained from a crystal
`structure. In this setting, the present study was undertaken to
`follow a molecular dynamics simulation for a protein in water for
`a Jong enough time to assess the COn\'ergcnce issue, to further test
`the performance of the OPLS force field, 14 and to obtain insights
`on protein hydration and possible differences between the solution
`and crystal structures. The protein selected for the study was the
`third domain of silver pheasant ovomucoid (OMSVP3). B
`Ovomucoids make up about 10% of the protein in avian egg
`whites, in which they arc the dominant inhibitors of serine pro(cid:173)
`teases. They are members of the Kaza! family of pancreatic
`secretory inhibitors. which are generally important in controlling
`the premature activation of pancreatic zymogcns. The complete
`ovomucoid consists of three homologous, tandem domains, each
`
`(9) K1!iger, P.; Stra.Bburger, \\'.; \\'ollmer, A.: van Gunsteren, \V, F. Eur.
`Biophys. J, 1985, JJ. 77.
`(10) (a) AhhtrOm, P.; Teleman, 0.: JOns.wn, B.: ForsCn,S. /.Am. Chtm.
`Soc. 1987, 109, 1541. (b) Ahlstrom, P.; Teleman, O.; JOn~son, B. J, Am,
`Chem. Soc. 1988, 110, 4198.
`(11) (a) ~·an Gun1tercn, W. F.; Karplus, M. BlrxhonfSfry 1981. 21, 2259.
`(b) S'>l-aminathan, S.; khi>·e, T.; van Gumteren, W. F.; Kup!us, M. Bio·
`chemistry 1982, 11, 5230. (c) van Gunstcrcn, W, F.; Berendsen, H.J. C.;
`Hermans, J.; Hol, \V. G. J.; Postma, J.P. M. Proc. Noll. A rod. Sd. U.S.A.
`1983, 80, 4315.
`(12) Wong, C. F.; McCammon, J. A. !Jr. J, Chtm. 1986, 17, 211; J. Am.
`Chtm. Soc. 1986, 103, 3g30,
`(13) Bash, P.A.; Singh, U. C.: Brown, F. K.: Langridge, R.: Kollman, P.
`A. Sdena 1987, 215, SN.
`(14) Jorgensen, \\'. L.; TiradirRlvc1, J. J. Am. Chtm. Soc. 1988, 110,
`1657.
`(!5) A preliminary report on this work was provided in the Proceedings
`or the !988 Nobd S)'m)X!Sium: Jorgcn.s<:n, \\', L; TiradirRive>, J. Chrm. Sa.
`1989, 19A, 191.
`
`0002· 7863/90/ 1512-2773502.50/0
`
`© 1990 American Chemical Society
`
`Breckenridge Exhibit 1026
`Breckenridge v. Novartis AG
`
`
`
`2774 J. Am. Chem. Soc., Vol. 112, /'r'o. 7, 1990
`
`Tirado· Rives and Jorgensen
`
`of which may inhibit a protease. The third domain contains 56
`amino acid residues that can be detache<i by controlled proteolysis.
`These fragments are typically at least as active as the full ovo·
`mucoid.
`Sequences have now been determined for ovomucoid third
`domains from o~·er I 00 avian species by Laskowski and co·
`workers. 16 They have also determined association e-0nstants for
`many of these inhibitors with a·chymotrypsin (AC), elastase
`(HLE), subtilisin, and Streplon1yces griseus proteases A and B
`(SGPA and SGPB). 11•18 Many of the sequences differ by only
`one or a few residue changes, so an unusually complete struc·
`lure/activity data base is being constructed.
`Besides the association and hydrolytic constants, some structural
`data are also available. Crystal structures have been determined
`for complexes of turkey ovomucoid third domain with SGPB, 19
`HLE,"° and Ac.ii for two isolated third domains, silver pheasant22
`and Japanese quail, 23 and for their corresponding hydrolyzed
`forms. 24 However, structural data in solution are limited to two
`studies by 2·D N~fR, spedfically, for turkey ovomueoid third
`domain in both its nativeH and hydrolyzed forms. 26
`The typical structure of an ovomucoid third domain contains
`three disulfide bridges, a 10-11 residue long a-helix, and a tri·
`pie-stranded antiparalle\ P'·sheet. The combination of small size,
`a large body of experimental data, and interesting structure makes
`the ovomucoids unusually attractive for a series of molecular
`dynamics investigations. Such simulations could provide detailed
`structural and thermodynamic information, which would be of
`great assistance in the interpretation of the biophysical data and
`the development of selective inhibitors. The third domain of silver
`pheasant ovomucoid was chosen for this initial study owing to its
`greater sequence homology with other ovomueoids and the better
`resolution of its crystal structure- than that of Japanese quail.
`Computational Procedure
`The entire simulation was conducted on Sun-4 romputers using the
`AMBER 10 program,11 with minor local modifications to irnproYe its
`use of the UNIX environment. The OPLS nonbonde..:I parameters" were
`used for the prlltein atoms in eonjuction with the TIPJP model for
`water. 21 As spt-eified in the OPLS model, 1he dielectric Cllnstant was
`kept fhed at 1.0, and the scaling factors for the l,4-nonbonded interac·
`tions were 8.0 for the Lennard.Jones and 2.0 for the el~trootatic inter·
`actions.14 The energetics for angle bending and torsional motion were
`describ(d with the AMBER united-atom force fic!d. 21 During the sim·
`ulation, all bond lengths and the H-H distances in water were kept
`constant by using the SHAKE algorithm a with a tolerance of0.0004 A,
`
`(16) La skol' ski, M., Jr.: Kato, I.; Ardelt, W.: Cook, J.: Denton, A.; Empie,
`M. W.: Kohr, \V. J.; Park, S. J.; Parks, K.: Schatiley, B. L.; Seh(}(nberger,
`0. l.; Tashiro, M.: Vichot, G.; Whatley, H. E.; Wieczorek, A.: Wieciorek,
`M. BtochtmfJllY 1987, 16, 102. Kato, I.; Kohr, W. J.: Lasko111ski, M., Jr.
`Biochtmillry 1987, 16, 193.
`{17) Empie, M. \V.; Laskowski, M., Jr. Bfocl1tmfrrry 1987, 11. 2214.
`(!8) Luko11o1ki, M., Jr.: T1.1hiro, M.; Empie, M. W.: Puk, S.J.: Kato, I.:
`Ardell, W.: Wiccrorek, M. In Prottlnau l11hiblt(Jls: Medic-al tmd Biological
`Aspecu; Katunuma, N., Umeuwa, H., Hozer, H., Eds.; Japan Scientific
`Societies Press: Tokyo, 1983: p SS.
`(19) Fujinaga, M.: Read, R.: Sie!ccki, A. R-; Ardell, \V.: Laskowski, M.,
`Jr.; James, M. N. G. Pux. Nail. Acad. Sd. U.S.A. 1982, 79, 4868. Read,
`R,; Fujinaga, M.; Sidccld, A. R.; James, M. N. G. 8iochtml11ry 1983, }2,
`4420.
`(20) BOOe. W,; Wd, A.·Z.: Huber, R.: Meyer, E.: Travis. J.; Neumann,
`S. EMBO J. 1986, 10. 24S3.
`(21) Read, R.; Fujinaga, M.; Sielecki, A. R.: Arde!t, W.: Lasko,,.~ki, M.,
`Jr. Aaa CryJtaf/agr., Sect. A, Suppl, 1984, 40, C-SO.
`(22) BOOe, W.; Epp, 0.: Huber, 0.; Laskowski, M., Jr.; Ardell, \V. Eur.
`J. Bfrxhtrn. 1985, 147, 387.
`(23) \Vcber. E.; Papamokos. E.: Bode, W.: Huber, R.: Kato, I.; lad:o-•nki,
`M., Jr. J. /.fol. Biol. 1981, 149, 109. Papamokcii, E.; \Vebcr, E.; Bode, W.:
`Huber, R.: Empie, M. W.; Kato, J.: Laskowski, M., Jr. J, /.fol. Biol, 1982,
`158, sis.
`(24) ~hail, D.; Bode, \V.: Ma)'f, I.; Huber, R.; Laskowski, M., Jr,; Lin.
`T.·Y.; Ardell, W. Unpublished re~ults.
`(25) Rob-ert~on, A.: \Veitkr, W. M.; Markley, J. L. Biochemistry 1988.
`17, 2519.
`(26) Rhyu, G. I.: Markley, J. L. Bioclumfs/Ty 1988, 27, 2529.
`(27) \Veiner, S. J.: Kollman, P.A.: Case, D. A.; Singh, U. C.: Ghio, C.:
`Alagona. G.: Profeta, S.: Weiner, P. J. J. Am. Chtm. Soc. 1984, 106, 165.
`{28) Jorgensen, W. L.; Chandrasekhar. J.: Madura, J. D.: lmpey, R. W.:
`Kldn, M. L J. Clitm. Phys. 1983. 70, 926.
`(29) nn Gunstcren, W. F.: Bcrcndscn, H.J. C. Mof. Phy1. 1917, J4, 131 !.
`
`·1~(.l))·r
`
`-191(1") ~
`
`~:: .
`
`I
`
`I ::~1, ~//ij!~iJ/o 1.
`
`11~
`
`l 1~1~ I
`
`O
`
`N
`
`M
`40
`THNiplJ
`Figure I. Potential energy variation during the course of the MD sim·
`ulation.
`
`W
`
`100
`
`.1.ll~
`
`~.I!-
`
`:.<)-V
`
`[Ill ~~~~~~~~-~~
`>'Oii
`111.•0
`IH'
`f.'111
`
`'1111
`
`Figure 2. Rt.fS de.,.iations OCtv..oen the instantaneous computed structure
`and the crystal structure for all residue; as a fonction of simulation time.
`
`which allowed the use of a time step of 2 fs. The temperature and
`pressure were kept constant at 298 K and l bar (0.987 atm). A non·
`banded pair list \\'as used to accelerate the cakulations and was updated
`eYery 10 steps. This list was generated by using a residue-based cutoff
`(9 A) to avoid spliuing dipoles. All the calculations utilized periodic
`boundary conditions to a\·oid edge effects.
`The initial coordinate• for the third domain of sih·er pheasant ovo·
`mucoid were obtained from the crystal structure 21 depo>ited in the
`BrooKhaven PrQtein Data Bank.uo The protein molecule, \1.1\hout any
`of the erystallographically located water molecules, was centered in a
`rectangubr bo.'t of water obtained by p.::ricdk translations in the :c, y, and
`z directions of a cube of water previously equilibrated Yia Monte Carlo
`c-alculations. Any water molecule el00;er than 1.5 A to any protein atom
`or farther away than 6 A from the closest protein atom in any one
`Cartesian direction was then deleted to give an initial system containing
`the solute plus 1721 water molecules (5676 total atoms) in a rectangular
`box of dimensions 43.8 x 42.0 x J4.J A.
`The initial preparation of the system consisted of JOO steps of steepe-st
`descent energy minimization, followed by a short ( l pi) constant volume
`molecular dynamics run. The resulting structure was then used as the
`starting point for the MD simulation at constant temperature and pres·
`sure. Initial atomic velocities were assigned from a Maxwellian distri·
`bution corresponding to a temperature of 298 K. The values of the
`potential energy, R~iS deYiation from the crystal structure, and volume
`were monitored eontinuou5Jy in order to follow the equilibration of the
`system. A total time of 100 p> was C{IYercd in the simulation, during
`which the coordinates, velocities, and energies were saved every 50 time
`sleps (0.J ps) for further analy>is.
`Results
`Cornergence Ueha1lor. The potential energy and Rl\·1S devi·
`at ions of the main-chain atoms (N. C", C, and 0) of residues 8-56
`for each instantaneous structure arc plotted in Figures I and 2,
`respectively, as a function of simulation time. The first seven
`
`(30) Entry 20VO, \"enion of No.,. 8, 198S.
`
`
`
`hfo/ecu/ar Dynamics Simu/alion for Pheasam Ovonwcoid
`
`J. Am. Che111. Soc., Vol. !12, No. 7. 1990 2775
`
`residues were ex.eluded from the RMS deviations, since they form
`a pendant tail that is less well resolved in the crystal structure.22
`The R~1S deviations are typicaUy 0.2 A higher when these residues
`are included. The system appears to have reached an equilibrium
`state after 30-40 ps, though the ultimate longevity of this state
`could only be determined in a much longer run. After this initial
`period, the simulation was continued to 100 ps in order to acquire
`enough data for averaging and analysis. In the 210-ps simulation
`for the other protease inhibitor, BPTJ, in water, an equilibrium
`state was achieved after ca. 50 ps.i
`Protein Structure and Dynamics. (A) Comparison "Ith the
`Crystal Structure. The calculated mean structure of the protein
`in solution, obtained by directly averaging the Cartesian coor(cid:173)
`dinates of all the saved structures from 30 to 100 ps, shows an
`RMS deviation from Lhe crystal structure of 1.28 A for the
`backbone atoms of residues 8-56 and 1.49 A for all residues.
`Results from previous simulations in water include 1.05 A for the
`C" of APP, averaged from 5 to 15 ps,9 and 0.77 A for all the
`backbone atoms of residues 1-56 of BPTI, averaged from 105
`to 210 ps.8 Oe\·iations from simulations in vacuo are typically
`larger by at least a factor of 2.s.11
`The value of the R~fS deviation of the instantaneous structure
`at 100 ps from the crystal structure is 1.43 A for the backbone
`atoms of residues 8-56 and 1.61 A for all residues. These values
`are in the same range as those obtained in previous, shorter sim·
`u!ations in water, e.g., 1.5 A for only the Ct atoms of BPTI after
`20 ps6 and 1.72 A for the C"" atoms of the trypsin-benzamidine
`complex at 45 ps. 11 Clearly, the AMBER/OPLS force field is
`providing a comparatively reasonable representation of the protein.
`An approximate base line for the RMS deviations can be deduced
`from the comparison of coordinates from proteins whose crystal
`structures have been solved in different crystalline forms or that
`have different molecules in the asymmetric unit. For five such
`cases. Chothia and Lesk obtained RMS deviations ranging from
`0.25 to 0.40 A, with a mean of 0.33 A.ll
`A comparative plot of the backbone atoms of the crystal
`structure with the instantaneous structure at the end of the MD
`It can be seen that the tertiary structure
`run is given in Figure 3.
`is well preserved in the simulation. However, several striking
`differences with the crystal structure are revealed by more detailed
`graphical analyses:
`( l) In the crystal structure, the amino acid side chains are
`mostly folded onto the surface of the protein, while in the solution
`simulation they ex.tend more into the solvent. This effect likely
`renccts the lower water content in the crystal and the removal
`of the interprotein interactions. It also allows polar side chains
`to be better stabilized by hydration.
`(2) The last residue of the o-helix in the crystal, Ser-44, does
`not form the required hydrogen bond with Ala·40 to be part of
`the helix in solution. This terminal hydrogen bond is also absent
`in the crystal structures of isolated Japanese quai113 and hydrolyzed
`silver pheasantH ovomucoids and in the complex.es of turkey
`ovomucoid with a-chymotr)'psin11 and S, gr/se11s protease B. 19
`(3) As discussed in the next section, some of the residues of
`the outermost strand of the /j-sheet (Lys-29, Thr-30, Tyr-31, and
`Gly-32) are tilted out of the plane of the other two strands in the
`simulation. This displacement inhibits the formation of the in·
`terstrand hydrogen bonds and yields greater hydrogen bonding
`between these residues and water mo!«:ules.
`(4) The hydrogen bond in the crystal structure between the side
`chains of Glu-19 and Thr-17, which surround the scissile peptide
`bond between Met-18 and Glu-19, is not found in the simulation
`results. Instead, it is replaced by a hydrogen bond between the
`guanidinium fragment of Arg-21 and the carbox.ylate in Glu-19,
`while the hydroxyl group of Thr-17 is h)'drogen-bonded to solvent
`molecules. This effect appears to be related to the loss of in(cid:173)
`terprotein interactions. In the crystal lattice, the charged side
`charged side chain of Arg-21 is close to the carbonyl oxygens of
`Glu-10 and Pro-12 from a protein molecule in a neighboring unit
`cell. Interestingly, the same interprotein contacts were found in
`
`:\10, IOOps
`
`Crystal Structure
`Figure 3. Compariwn of backbone atoms in the inst.antanwus structure
`at the end of the simulation (I = 100 ps) and in the crystal structure
`
`the crystal structure of Japanese ~uail ovomucoid, which has a
`very different crystalline form. 23· 2
`In addition, the Glu-19
`earboxylate group in the crystal is hydrogen-bonded to the
`side-chain ammonium group of Lys-13 in the same neighboring
`cell. Figure 4 compares the different interactions in this region
`for the crystal structure and the instantaneous structure at 80 ps.
`The overall conformation of a protein can be expressed in terms
`of the backbone torsional angles ¢; (C1-1-N1-C/-C1) and ~1
`(N1-Ci-CrN1+1)· The average values for these angles during the
`30--100-ps period are compared with the corresponding values for
`the crystal structure in Figure 5. Consistent with the statements
`above, the biggest differences are found in the region around the
`C-terminus of the a-helix (Glu-43, Ser-44, and Asn-45), where
`the middle residue is twisted away from the helix, and in the region
`near the turn connecting the central and outer strands of the
`.8-sheet (Ser-26, Asp-27, and Asn-28). Some of the other dif·
`fercnces are of little consequence to the overall conformation, since
`compensation occurs when there are differences. of opposing signs
`in an angle "11 and the subsequent oti1+1• This is the case for J..1et-18
`and Glu-19 and for the Thr-47, Leu-48, and Thr-49 region. The
`
`(31) Chothia, C.; led;, A. M. EMBOJ. 1986, .5, 823.
`
`(32) SilHr phe.uant oYomucoid crystallizes in the C2 5pace gioup, while
`for Japanese quail, the cryital~ belong in the tetragonal P42 12 ip.ace group.
`
`
`
`2776 J. Am. Chem. Soc., Vol. 112, iVo. 7, 1990
`
`Tirado- Rives and Jorgensen
`
`.__,,,
`
`0
`
`~II
`
`I
`''
`
`' ,,
`
`Figure 4. Hydrogen bonds near the scissi\e bond in the instantaneou5
`structure at the end of the simulation (I = 100 ps) and in the crystal
`structure.
`
`R;\fS differences from the crystal structure computed for each
`backbone dihedral angle over the entire protein are 39°, 48°, and
`11° for q,, '},and w, respectively. These differences are in the
`same range as the results from previous simulations, e.g., 31°,
`37°, and 8° for !Ji, ~. and w in trypsin 12 and 26° and 33° for 4>
`and Y, for BPTI in Mvan dcr \Vaals water". 1h
`A more global impression of the overall conformation of the
`protein can be obtained from the plots of¢ vs Vt (Ramachandran
`maps) in Figure 6. The general trend observed in the m~ps is
`that residues in the a-helix stay close to the crystallograph1ca\ly
`observed angles, while those outside the helix are shifted toward
`values more typical of the C5 and ~ conformations. This dis·
`placement is consistent with the results.of our previous ene~gy
`minimiiation studies on the conformations of N-acetylglycme·
`N-methytamide and N-acetylalanine·N·methy\amide, which
`showed that these conformations arc lower in energy than the
`corresponding a. "R· or al alternatives. 14
`(8) Comparison wilh Nl\1R Data. The NMR data in solution
`on native and hydrolyzed turkey ovomucoid should be relevant
`to the present case, since turkey and silver pheasant ovomucoi~s
`differ by only one residue (18: ~fet/Lcu). In general, the main
`effort in the NMR investigations was devoted to the complete
`assignment of resonances, and only qualitative information re·
`garding secondary structure was obtaincd.H,16
`Although the principal conclusion of these studies was that the
`solution structure had to be very similar to that in the crystal,
`some details were explored further. In particular, the NOESY
`conn~tivities that establish the antiparallel triple-stranded ,B-sheet
`were given. Figure 7 provides a comparison of the segments
`comprising the .B·sheet in the crystal and In the mean calculated
`solution structure, obtained by direct averaging of the Cartesian
`
`Figure 5. Plots of the b3ckbone angles¢ and~ for the crystal structure
`and the averages from the 1fD simulation.
`
`coordinates from 30 to 100 ps in the MD sintulation. The in·
`terproton distances gh·en with the crystal fragments are measured
`in the crystalll but correspond to observed NOE pairs in the
`solution Nl\1R. Only the distances that show significant deviations
`from the l\fD averages are given. As alluded to above, the largest
`differences are found in the separation between the central and
`the outermost strands of the .8-sheet. The two distances of 6-7
`A calculated in the simulation are beyond the normally accepted
`range for NOE detection. A problem with the force field or the
`preparation of the system could be indicated. However, it is also
`possible that the protein is exploring the expanded phase of a
`low-frequency Mbrcathing mode", and the separations may docrease
`at later times in the simulation. 34
`(C) F1uctuations. The overall mobility of the different atoms
`during an MD simulation can be expressed as their R~fS fluc·
`tuations, (w-2) 1/l = ((rf - (r1))2}112, computed o;·er the averaging
`period. An examination of the accumulated average atomic
`fluctuations as the time span of their evaluation increases shows,
`as expected, a monotonic increase for about 40 ps that levels
`toward a plateau after ca. 50 ps. These fluctuations can be
`compared to the average movement observed in the crystallo·
`graphic determination as expressed in the B thermal factors.
`Crystallographic RMS fluctuations can then be derived via the
`relation (ar2) 112 = (3B/8r2) 112 from the B-values.lS Figure 8
`compares the RMS t1uctuations calculated during the last 25 ps
`of the simulation with the fluctuations derived from the thermal
`factors. 22
`
`(33) Any mining h)drogen atoms in the X·ray and the MD structures
`were added by using the SYBYL program from TRIPOS Asroo:ia1es.
`(34) Suczaki, Y.; Go, N. l/11. J. Ptpl. P10tti11 Res. 1975, 7, 333.
`(35) \Vi!lis, B. T. M.; Pryor, A. 'ti.'. ThtrntfJl Vibrations in Crptalfogra·
`phy. Cambridge Univenity Piess: Cambridge, F.ngland, 1915.
`
`
`
`,\fo/ec11lar Dynamics Simulation for Pheasant Ovon111coid
`
`J. Am. Chem. Soc .. Vol. 112, A'o. 7. 1990 1777
`
`\~lj-,-----·---------
`
`. . <+:++-:..
`+-
`
`+.
`
`1
`
`I
`1:0-
`
`'
`-WI~
`i
`
`-1!0 -,
`
`'" i__::_,_. -~----~--~-~
`
`-1,\:1)
`
`·120
`
`-(ll
`
`1,£1
`
`1,0
`
`'"'
`
`lW
`
`x
`'X« ~x xf
`x x
`'
`x1.:x)(
`x x
`
`x
`
`'\.
`xft
`
`x
`
`r
`
`ro 1
`'
`O·
`i
`
`·t0--1
`
`-120
`
`·l!:O
`-1&0
`
`.1io
`
`.w
`
`rh·:··,
`
`x
`
`x
`
`Cry:ml
`
`,;,
`
`ISO
`
`Figure 7. Backbone atoms of the /3·shcet segments for the average MD
`structure and in the cnytal structure. Some interproton distances are
`sho'.\'n.
`
`~:-U('J
`Flgurt 6. Ramachandran maps for the MD results and the crystal
`structure.
`
`23
`
`Table I. Average RMS Fluctuations (A) for Different Atom Types
`averaging time, ps
`30----70
`l.16
`1.24
`I .JO
`I.JO
`1.27
`J.96
`
`atoms
`c·
`c•
`('
`C'
`all protein
`0"'1"
`
`30----50
`0.73
`0.82
`0.90
`0.92
`0.86
`2.76
`
`30--100
`1.32
`1.40
`1.50
`1.45
`1.44
`5.17
`
`-MD
`- - X·R•y
`
`~
`'
`~
`" ~
`"
`
`,,
`
`u
`
`05
`
`~ ~ I! I I
`LO \ J~Ut, ,~,~~,1 ~~wl,~~\ll~ ~
`
`1 1 ·
`
`1111
`
`\: I ,
`
`~
`
`It can be seen from the plot that qualitative agreement in the
`width and location of the largest peaks is obtained. The most
`noticeable difference is the comparativ~ly high fluctuations for
`the N-terminus in the simulation, which decrease until Cys-8 is
`reached. The first seven residues form a pendant chain whose
`increased mobility in solution can be attributed to the formation
`of a four-stranded /3-channel with the N-tcrminal segments of
`neighboring molecules in the crystal.ll The remaining qualitati\•e
`differences reflect greater motion in the simulation for the side
`chains of Glu-10, ~tet-18, Asp-27, Lys-29, GJu-43, and Asn-45,
`which are more exposed to the solvent in the simulation than in
`the crystal structure.
`The Ri\1S fluctuations for the different atoms are correlated
`to their distances from the backbone of the protein. Table I lists
`the average R~·IS fluctuation~ for different carbon atoms for three
`averaging periods. The values obtained for the ci atoms have
`larger statistical uncertainties, since there are only 30 of these
`atoms in the protein. The fluctuations calculated in the initial
`20 ps of the averaging period are comparable to those obtained
`in the simulations of aqueous trypsin-bcnzamidine during 19 ps
`(C" = 0.52, C" = 0.59, Cl= 0.69, C6 = 0.70, all= 0.68).12 BPTI
`in "van der \Vaats water" for 25 ps (Ca = 0.54, Cll = 0.69, Cl
`
`oO
`
`0
`
`''° "'
`'°
`A1om r-.·orr.OC1
`Fii1Jrt 8. Comf'lri>0n of RMS nuetuations from the MD ealculltion and
`from the crsytallographie B·factors.
`
`360
`
`450
`
`540
`
`= 0.89, ct= 1.15, all = 0.78), lb and aqueous APP for 10 ps (0"
`=- 0.53, cs= 0.62, Cl= 0.73, Cl= 0.79, all= 0.75).9 Howe\•er,
`the fluctuations calculated over the full averaging period of 70
`ps are considerably larger than those reported by Levitt and
`Sharon for the final 105 ps of residues 2-56 of BPTI (C" = 0.42,
`cs= 0.50, D = 0.54, c~ = 0.67). 8 Greater experience is needed
`to ascertain if these differences are associated more with the
`proteins, the forox fields, or the details of the simulation procedures,
`though further comment on this issue is made below.
`(D) Hydrogen 8-0ndlng. Analyses of the intra protein hydrogen
`bonding were carried out on the 700 coordinate sets saved during
`the last 70 ps of the MD simulation. The criteria used to define
`a hydrogen bond were purely geometric. A list of all potential
`donors and acceptors (hydrogens attached to heteroatoms and the
`heteroatoms themselves) was generated at the beginning of the
`analysis. For each coordinate set, every potential donor-acceptor
`
`
`
`2778
`
`J. Am. Chem. Soc., Vol. I J 2, No. 7, 1990
`
`Tirado- Rives and Jorgensen
`
`Table II. Protein-Protein Hydrogen Bonds and Their Frequencies
`donor
`acceptor
`crystal
`re>idue distance, A
`atom
`residue
`atom
`Backbone to Backbone
`3.4
`C)'S·8
`2.1
`Gly·32
`Tyr-31
`2.9
`2.8
`His·52
`Lys-29
`2.8
`Thr-49
`2.8
`Gly-25
`2.9
`Leu-23
`2.8
`Asn-33
`3.0
`2.8
`Lys-34
`Cys-35
`3.0
`Asn-36
`3.1
`Phe-37
`2.9
`Cys-38
`3.1
`Asn-39
`2.9
`Ala-40
`2.9
`Val-42
`2.9
`Val-41
`2.8
`Scr-44
`3.2
`Cys-24
`2.9
`Cys-24
`3.2
`Pro-22
`3.1
`
`MD
`frequency, %
`
`22.6
`99.6
`
`99.4
`
`96.6
`
`94.\
`100.0
`98.9
`52.3
`70.6
`98.6
`85.3
`0.3
`
`8.0
`
`55.J
`77.7
`53.4
`
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`N
`
`T~·r·l l
`0
`Arg-21
`0
`0
`Leu-23
`Cys-24
`0
`G!y-25 D
`0
`Scr-26
`Asn·28
`0
`T)'r-31
`0
`0
`Phe-37
`Cys-38
`0
`Asn-39
`0
`Ala-40
`0
`Val-41
`0
`Va!-42
`0
`Glu-43
`0
`Ser-44
`0
`Asn-45
`0
`0
`Gl)'·46
`Thr-47
`0
`Scr-51
`0
`His-52
`0
`C!y-54
`0
`
`RMS f\c,!ntiun !A}
`Fig~re 9. Distribution of RMS fluctuations for heavy atoms of the
`protein and water oxygens from the MD simulation.
`
`conserved in the t\iD simulation. This is also reasonable, since
`the backbone atoms show smaller RMS deviations from the crystal
`structure than do the side chains. A total of 27 backbone hydrogen
`bonds were ohserved during the period from 30 10 100 ps for more
`than 10% of the structures, with a mean occurrence of 64o/o.
`Among these, 14 of the 22 hydrogen bonds detected in the crystal
`structure are conserved at the 10% or greater level. Some of the
`previous MD simulations have reported conservation of 13 out
`of l6 cr)'stallographic backbone hydrogen bonds in the simulation
`of BPTI in a truncated octahedron of water,6 16 out of 20 in the
`APP simulation,9 and 20 out of 21 in Levitt and Sharon's simu·
`\ation of BPTl in water with similar hydrogen· bond criteria. 8•36
`The latter authors also reported a mean occurrence of 68% for
`the backbone hydrogen bonds during the simulation. Since BPTI
`and Ot\1SVP3 ha,·e similar size and secondary structural clements,
`the smaller fluctuations and greater hydrogen·bond conservation
`in the study of Levitt and Sharon suggest that the Lifson force
`field is stiffer than the OPLS/Al\fBER model and keeps the
`protein closer to the starting crystal structure. Another indication
`of the relative rigidity of their force field may be found in the
`fluctuations calculated during their MD simulation in vacuo,
`which, although larger than their corresponding values in solution,
`are still smaller than the nuctuations observed in the crystal. This
`is not surprising, since the Lifson force field was specifically
`designed to reproduce crystal structures.11
`The backbone hydrogen bonds that are absent from the present
`solution calculation fall in two major groups: those for residues
`44-40, 45-42, and 47-44, formed in the crystal by the last segment
`of the a· helix. and the turn needed to connect back to the .B·sheet,
`and the hydrogen bonds formed in the crystal by the third strand
`of the .B·sheet, namely, 23-31, 25-29, 28-25, and 31-23. These
`structural features, as mentioned before, arc not present in the
`structure calculated in solution.
`Sohent Structure and Dynamics. In the course of the molecular
`dynamics simulation, the solvent exhibits, as expected, greater
`mobility than the protein. The average Rl\IS nuctuation of the
`pro