13898
2009, 113, 13898–13900 Published on Web 09/30/2009
NMR Scalar Coupling Constant Reveals That Intraprotein Hydrogen Bonds Are Dynamically Stabilized by Electronic Polarization Chang G. Ji†,‡ and John Z. H. Zhang*,‡,§ Institute of Theoretical and Computational Chemistry, School of Chemistry and Chemical Engineering, Nanjing UniVersity, Nanjing 210093, China, State Key Laboratory of Precision Spectroscopy, Department of Physics, East China Normal UniVersity, Shanghai 200062, China, and Department of Chemistry, New York UniVersity, New York, New York 10003 ReceiVed: August 18, 2009; ReVised Manuscript ReceiVed: September 16, 2009
Molecular dynamics simulations based on the standard nonpolarizable AMBER force field and on quantumderived polarized protein-specific charge (PPC) are performed to compute NMR scalar coupling constants across hydrogen bonds for three benchmark protein systems: ubiquitin, the GB1 domain of protein G, and the SMN Tudor domain. Direct comparison of the simulation result with experimental data gives strong evidence that intraprotein hydrogen bonds are significantly stabilized by electronic polarization, both in terms of NMR scalar coupling constants and X-ray determined geometries of hydrogen bonds. Without the polarization effect in the force field, hydrogen bonds are found to be “too loose”, which leads to less stable or even unstable local structures of proteins. Hydrogen bonds are of fundamental importance in stabilizing biomolecular structure and play a key role in many of protein’s functions.1 Thus, accurate description of hydrogen bonds has been a long-standing goal in both experimental and theoretical studies in structural biology. Recent observation of J couplings across hydrogen bonds by NMR experiments2 provided direct observation of static and dynamic characters of hydrogen bonding. These couplings are sensitive to geometries of hydrogen bonds and therefore provide considerable insight in the determination of secondary and tertiary structure in biological systems. Validation of molecular force fields is often done through direct comparison of experimental and calculated scalar coupling constants.3-6 Previous MD simulation of NMR scalar coupling constants shows that hydrogen bond length is overestimated from explicit water MD simulation based on the nonpolarizable CHARMM force field, despite significant improvement over that from implicit water model simulation.5 In view of the difference between simulation and experimental data, one should realize that the current force field lacks electronic polarization of proteins. Specifically, the atomic charges in the standard force field are based on fitting of the electrostatic potential (ESP) of individual amino acids, without the polarization effect of the protein. Since electrostatic interaction plays a major role in protein interaction, the lack of protein polarization may be responsible for observed deficiencies in proteins’ structural and dynamical properties extracted from MD simulations. Recently, polarized protein-specific charge (PPC) has been developed to provide new atomic charges of proteins for better * Author to whom correspondence should be addressed. E-mail:
[email protected]. † Nanjing University. ‡ East China Normal University. § New York University.
10.1021/jp908002n CCC: $40.75
description of protein dynamics near native structure.7 Fully quantum mechanical calculation of electronic structure of proteins in solution can be achieved through the MFCC8-10 (molecular fragmentation with conjugated caps) method, in which the protein is partitioned into fragments and the Poisson-Boltzmann equation for protein in solvent is numerically solved. PPC is then derived by fitting atomic charges of protein fragments to their electrostatic potentials calculated for protein native structure in solution.7 Thus, effective atomic charges of PPC correctly represent the polarized electronic state of the particular protein at a given structure (native structure). As is understood, the ability to represent the inhomogeneous electrostatic environment of protein by protein’s atomic charge is essential for accurate description of protein dynamics through MD simulation. Previous studies demonstrated that PPC gives better conformation sampling for accurate calculation of free energy change in the protonation process in proteins7a and also prevents the protein system from being driven away to unphysical structures in protein-ligand binding of PPARγ.7b Here, we perform quantum fragment calculations to derive PPCs for three benchmark proteins and then employ both PPC and Amber99SB11 to simulate NMR scalar coupling constants across hydrogen bonds in these benchmark proteins for explicit comparison with experimental measurement. Such direct comparison should provide explicit information on the effect of electronic polarization on the structure and dynamics of hydrogen bonding in proteins. In our theoretical approach, the solute (protein) is represented by a charge distribution F(r) embedded in a cavity surrounded by a polarizable medium with dielectric constant ε. The solute charge distribution F(r) polarizes the dielectric medium and creates a reaction field which reacts back to polarize the solute untilequilibriumisreached.BycombiningthePB(Poisson-Boltzmann) solution for the electrostatic field and the MFCC calculation 2009 American Chemical Society
Letters
J. Phys. Chem. B, Vol. 113, No. 42, 2009 13899
Figure 1. Comparison of calculated and experimental h3JNC′ values for all three proteins (1mhn, 1pgb, 1ubq). The computed results are from MD simulations using, respectively, AMBER and PPC.
for protein, we obtain the MFCC-PB version for quantum calculation of protein in solution as described in ref 7. The calculated electron densities of protein fragments are fitted to ESP (electrostatic potential). The quantum chemistry calculation of protein is performed at the level of DFT/6-31G*. Three proteins used in our study are ubiquitin (1ubq: 1.8 Å), the GB1 domain of protein G (1pgb: 2.1 Å), and the SMN Tudor domain (1mhn: 1.8 Å). For MD simulation in explicit water, the protein was solvated in an octahedron-like box and the system was neutralized by adding counterions. After heating and equilibration, another 10 ns production run was performed at 300 K (NPT). For each protein, separate simulations were performed using, respectively, Amber99SB and PPC. In PPC simulation, the corresponding AMBER force field parameters11 were retained except for the atomic charges that are replaced by PPC. The charge fitting philosophy used in PPC was the same as that used in the Amber force field, and this guarantees that PPC charge is consistent with other parameters of the Amber force field. Protein conformations were saved at every 0.1 ps for the calculation of the J scalar couplings. On the basis of DFT and finite perturbation theory calculation, the geometric dependencies for J couplings on H-bonds of proteins can be parametrized by the formula given by Barfield:3
JNC′ ) 〈(-357 Hz) exp(-3.2rHO /Å) cos2 θ〉
h3
(1)
where θ is the H · · · OdC angle and rHO is the distance between the hydrogen and oxygen atoms. For the purpose of direct comparison with experimental results, only backbone hydrogen bonds were studied here, with a total of 76 NsH · · · OdC H-bond pairs in all three proteins combined. Figure 1 shows the comparison between the calculated and experimentally measured coupling constant h3JNC′ for all three proteins (pdb id: 1mhn, 1pgb, and 1ubq) combined. As is clearly shown in Figure 1, the scalar coupling constants h3 JNC′ computed from MD sampling of conformations under PPC agree much better with experimental values12 than those under AMBER. Most of the h3JNC′ computed under Amber are larger than the experimental values with an average rmsd of 0.246, vs 0.125 under PPC simulation, for all three proteins computed together. It is understood from eq 1 that the value of h3JNC′ strongly correlates with the geometry and dynamics of the hydrogen bond. Thus, systematic deviation of h3JNC′ calculated under the Amber force field indicates that hydrogen bonds in MD simulation under Amber are “too loose”. This conclusion can be easily verified from the distribution of the average H · · · O bond distance shown in Figure 2 for the three proteins studied. As shown in Figure 2, the hydrogen bond length rHO under PPC
Figure 2. Distribution of H · · · O bond distance in MD simulations under PPC and AMBER force fields for three proteins. The experimental values are indicated by XRD.12
Figure 3. (a) Distribution of H-bond geometry from MD simulation (black circle: experimental value). (b) Distribution of H · · · O bond length and H · · · OdC angle of the hydrogen bond of Ala28 · · · Asp32 in 1ubq.
is populated with the maximum at XRD (experimental) value for all three proteins. In contrast, shifts of 0.08, 0.12, and 0.18 Å toward large rHO distance are detected for 1mhn, 1pgb, and 1ubq, respectively, in Amber, indicating weaker hydrogen bonding strength in the force field due to the lack of polarization. The result is highly consistent with our previous finding for protein-ligand binding where the lack of electronic polarization caused breaking of critical hydrogen bonds and is responsible for partial collapse of the stable protein-ligand binding complex.7b To further investigate dynamical properties of individual H-bonds, we plot the contours of statistical distribution of rHO and bond angle θ for a particular backbone H-bond (Ala28 · · · Asp32 in 1ubq) under both AMBER and PPC simulations, as shown in Figure 3a. As we can see from the figure, heavily sampled geometries under Amber simulation are not distributed around experimental values but off from them. This can cause systematic shifts of the equilibrium structures of proteins under the standard Amber force field. The direction of the shift indicates that AMBER underestimates H-bond strength. Figure 3b also shows that larger flexibility was observed under Amber simulation for the H-bond of Ala28 · · · Asp32 in 1ubq. Previous work by Case13 concluded that this “over-flexibility” under the Amber force field is also the origin of systematic error in predicted NMR order parameters from simulation. In contrast,
13900
J. Phys. Chem. B, Vol. 113, No. 42, 2009
Figure 4. Comparison of dipole moments of individual residues calculated from PPC and those from AMBER99 charges for three proteins.
simulation under PPC gives excellent agreement with NMR scalar coupling constants and hydrogen bond geometrics, as shown in Figure 3b, indicating proper description of hydrogen bonding strength. Without proper inclusion of electronic polarization in the force field, protein’s conformational distribution in vivo could not be accurately represented from MD simulation. This defect of the classical force field has been recognized for years.13 For this reason, direct artificial constraints are sometimes applied on proteins during MD simulation in order to obtain better agreement with experimental data. For example, systematic effort through extra constraints based on experimental NOE distances and NMR order parameters have been carried out in MD simulation for years.4 However, these empirical approaches do not solve the fundamental problem. Here, we demonstrated that PPC provides an excellent choice to solve this problem without introducing any artificial constraint in MD simulation. A significant challenge in computational biology is to correctly predict free energy change in important biochemical processes (such as protein-ligand and protein-protein bindings, etc.). Although the free energy landscape can be constructed directly from the distribution of conformational ensemble from MD simulation, computed free energy changes critically depend on the accuracy of the underlying molecular force field. Thus, the accuracy of the interaction force field is of fundamental importance for correct conformational sampling and for accurate prediction of protein structure and dynamics. Our previous work on MD simulation demonstrated that PPC performs very well in predicting free energy change during the protonation process in proteins.7a To explore further differences between PPC and Amber charge, dipole moments of individual residues in proteins were analyzed. Figure 4 shows that PPC predicts generally larger dipole moments for individual amino acids than the AMBER force field. This can be interpreted as the fact that PPC are more polarized than Amber charges and electrostatic interaction represented by PPC is relatively stronger than the one represented by unpolarized Amber charge. Since PPC is derived from quantum electron structure calculation of protein in solution, it is fundamentally different from AMBER or other standard force fields in which atomic charges are derived for individual amino acids in an inhomogeneous electrostatic environment. High resolution crystal structure experiments14 verified that amino acids have different polarities under different environments. It should be mentioned that a number of similar quantum electron
Letters structure methods for protein calculation have also been developed by several groups.15-17 These new generation quantum methods for protein structure calculation should enable us to provide new insight as well as more accurate description of protein structure and dynamics in computational biology. In summary, we presented molecular dynamics simulation study for NMR J coupling across hydrogen bonds for three benchmark proteins. Much improved agreement between experimental structural data and simulated result using PPC is obtained. The present result demonstrates the important effect of electronic polarization on stabilizing intraprotein hydrogen bonds and protein structure. A note is in order here for possible further improvement of PPC. Since proteins are dynamic and exist as an ensemble of structures, it may be desirable to fit PPC charges from multiple configurations such as the NMR-derived structural ensemble. Another possible approach is to include solvent effects using the explicit water model instead of the implicit PB method and thus can include the solvent fluctuation effect in multiconfiguration fitting of PPC.18,19 These ideas need to be explored in future studies. Acknowledgment. We thank the National Basic Research Program of China (Grant No. 2004CB719901), the National Natural Science Foundation of China (Grant No. 20773060 and 20933002), and Shanghai PuJiang program (09PJ1404000) for financial support. References and Notes (1) (a) Pauling, L.; Delbruck, M. A. X. Science 1940, 92, 77. (b) Jeffrey, G. A.; Saenger, W. Hydrogen Bonding in Biological Structures; Springer: New York, 1991. (2) Cordier, F.; Nisius, L.; Dingley, A. J.; Grzesiek, S. Nat. Prot. 2008, 3, 235. (3) Barfield, M. J. Am. Chem. Soc. 2002, 124, 4158. (4) Gsponer, J.; Hopearuoho, H.; Cavalli, A.; Dobson, C. M.; Vendruscolo, M. J. Am. Chem. Soc. 2006, 128, 15127. (5) Sass, H. J.; Schmid, F. F. F.; Grzesiek, S. J. Am. Chem. Soc. 2007, 129, 5898. (6) Markwick, P. R. L.; Sprangers, R.; Sattler, M. J. Am. Chem. Soc. 2003, 125, 644. (7) (a) Ji, C. G.; Mei, Y.; Zhang, J. Z. H. Biophys. J. 2008, 95, 1080. (b) Ji, C. G.; Zhang, J. Z. H. J. Am. Chem. Soc. 2008, 130, 17129. (8) Zhang, D. W.; Zhang, J. Z. H. J. Chem. Phys. 2003, 119, 3599. (9) Mei, Y.; Zhang, D. W.; Zhang, J. Z. H. J. Phys. Chem. A 2005, 109, 2. (10) Mei, Y.; Ji, C. G.; Zhang, J. Z. H. J. Chem. Phys. 2006, 125, 7. (11) Hornak, V.; Abel, R.; Okur, A.; Strockbine, B.; Roitberg, A.; Simmerling, C. Proteins: Struct., Funct., Bioinf. 2006, 65, 712. (12) (a) Markwick, P. R. L.; Sprangers, R.; Sattler, M. J. Am. Chem. Soc. 2003, 125, 644. (b) Cornilescu, G.; Ramirez, B. E.; Frank, M. K.; Clore, G. M.; Gronenborn, A. M.; Bax, A. J. Am. Chem. Soc. 1999, 121, 6275. (c) Cordier, F.; Grzesiek, S. J. Mol. Biol. 2002, 317, 739. (13) Case, D. A. Acc. Chem. Res. 2002, 35, 325. (14) Lario, P. I.; Vrielink, A. J. Am. Chem. Soc. 2003, 125, 12787. (15) (a) Xie, W. S.; Gao, J. L. J. Chem. Theory Comput. 2007, 3, 1890. (b) Xie, W. S.; Song, L. C.; Truhlar, D. G.; Gao, J. L. J. Phys. Chem. B 2008, 112, 14124. (16) Fedorov, D. G.; Kitaura, K.; Li, H.; Jensen, J. H.; Gordon, M. S. J. Comput. Chem. 2006, 27, 976. (17) Wang, B.; Merz, K. M. J. Chem. Theory Comput. 2006, 2, 209. (18) Gao, J. L.; Luque, F. J.; Orozco, M. J. Chem. Phys. 1993, 98, 2975– 2982. (19) Sanchez, M. L.; Aguilar, M. A.; delValle, F. J. O. J. Comput. Chem. 1997, 18, 313–322.
JP908002N
