J. Phys. Chem. B 2007, 111, 11237-11243
11237
Molecular Modeling of Oligopeptide Adsorption onto Functionalized Quartz Surfaces Giuseppe Forte,*,† Antonio Grassi,† and Giovanni Marletta Laboratory for Molecular Surfaces and Nanotechnology (LAMSUN), Dipartimento di Scienze Chimiche, UniVersity of Catania and CSGI, Viale A. Doria 6, Catania, Italy ReceiVed: December 21, 2006; In Final Form: May 25, 2007
The adsorption of an EAK 16-II oligopeptide sequence in aqueous medium onto functionalized quartz surfaces has been studied by using force field calculations and molecular dynamics methods. Two different surfaces have been simulated respectively involving fully methylated and fully silanolic quartz surfaces. Geometry optimization and molecular dynamics simulations showed that the adsorption process is mainly governed by the electrostatic interactions between SiO- surface groups and the charged residues of the oligopeptide sequence. In particular, it was found that strong electrostatic interactions (a) prompt the parallel orientation of the oligopeptide with respect to the hydrophilic charged surface, resulting in an effective physisorption process and (b) stabilize the β-sheet configuration of the physisorbed molecules. In particular, the end-on oligopeptide orientations are demonstrated to progressively lie back onto the hydrophilic surface, but this does not happen onto the hydrophobic surface. In any case, no physisorption process was observed for the fully methylated surface, where the molecule is seen to move away from the surface during the simulation time.
Introduction Geometry optimization and dynamical properties of single amino acids and oligopetides have been extensively studied in both the gas and solution phases by using ab initio and molecular dynamics techniques.1-7 Moreover, due to the wide application in technology, the silicon-based surfaces have been the object of many works and, in the last years, a growing interest has been directed to the study of the biological molecules adsorbed onto inorganic surfaces, due to their applications in several areas of outstanding technological interest such as biosensors, biomaterials, etc.8 In previous works9,10 we have analyzed in detail the nature and the strength of the interaction between the l-lisyne and R-quartz surface in an aqueous medium with (1) a fully hydrophilic surface (SiO-/SiOH groups), (2) a partially hydrophilic surface (SiO-/SiCH3 groups in 1:5 ratio), and (3) a fully hydrophobic surface (only -SiCH3 groups). In these previous papers, both geometry optimization and dynamical behavior of the interaction were studied, and the results showed that the nature of the interactions with fully or partially hydrophilic surfaces is essentially electrostatic, with a relevant contribution due the H-bond, whereas for fully methylated surfaces no noticeable interaction between the amino acid and the surface was observed. Furthermore, it was also found that the termination of the various surfaces strongly affects the relative ordering of the water shells.9 In recent years, a great attention has been devoted to study the structure-properties relationships in solution of a newly synthesized oligopeptide sequence, indicated as EAK 16-II,11,12 consisting of the following series of 16 amino acid residues:
AEAEAKAKAEAEAKAK where A ) alanine, E ) glutamic acid, and K ) lysine. * Corresponding author. Tel.: +390957385066. Fax: +39095580138. E-mail:
[email protected]. † Faculty of Pharmacy at the Univeristy of Catania and CSGI.
EAK 16-II sequences exhibit very interesting biological properties, involving the capability to stimulate specific cell adhesion and proliferation processes, particularly relevant to advanced tissue engineering applications.13,14 Such properties have been basically attributed to the occurrence, in saline solutions, of complex extensive supramolecular assembly processes able to promote the formation of biocompatible membranes, basically related to the occurrence of ordered β-sheet configuration as the basis for the more complex organization processes.11 Recent experimental results show that the supramolecular organization also occurs not only in saline11 but also in nonsaline solution.15 One of the most stimulating challenges now involves the attempt to transfer the self-assembling process of oligopeptide sequences from solution to surfaces, to obtain selective cell interactions on covered synthetic prosthesis and tissue engineering scaffolds.14 It should be noted that preliminary experimental results from our laboratory revealed the occurrence of massive formation of fibrils onto various inorganic surfaces.15 To analyze complex processes involving two or more molecules of oligopeptide, one needs first to understand both the basic interactions and forces that occur between the molecule and the surface. In the present paper both geometry optimization and molecular dynamics simulations have been performed in the framework of the molecular mechanics methods to investigate the nature and the strength of the primary interactions in an aqueous medium between the EAK 16-II sequence and two model functionalized quartz surfaces, respectively, showing fully hydrophilic and fully hydrophobic character. In particular, the calculations aimed to establish the dependence of the physisorption processes upon the surface features, as polarity, and the relative weight of the different terms of interaction forces. Furthermore, the simulations are aimed to clarify the nature of the anchoring sites for the different surfaces and their relationship with the possible occurrence of the peculiar β-sheet configuration, needed for the EAK 16-II self-assembling processes11
10.1021/jp068803h CCC: $37.00 © 2007 American Chemical Society Published on Web 09/06/2007
11238 J. Phys. Chem. B, Vol. 111, No. 38, 2007
Forte et al.
Figure 1. EAK 16-II.
Figure 2. Initial configurations and relative orientation end-on (Ia), end-on (Ib), side-on (Ic).
Methods All calculations were performed in the framework of the molecular mechanics approximations, using the Material Studio (MS) package.16 In both geometry optimization and molecular dynamics simulations, the Discover package included in MS software was adopted and the Consistent Valence Force Field (CVFF)17-19 parametrization was used in all simulations for all the components of the system (oligopeptide, surface and water molecules). Periodic boundary conditions (PBC), using the Ewald summation20 method with a dielectric constant equal to 1, was applied and the electric neutrality condition was maintained using Li+ as a counterion. This choice is due to the experimental observation that lithium salts are the most effective, inducing the self-assembling processes,11 so that it is important to check its influence on the oligopeptide-surface interaction. In the geometry optimization, a minimization energy procedure was performed by using the Smart Minimizer algorithm, starting from the selected initial structure. This method combines the features of other available methods in sequence, starting with the steepest descent method (SD), followed by the conjugated gradient method (CG), ending with a NewtonRaphson gradient technique (NR).16 The convergence criterion adopted for the value of maximum force, were 1000.0 (kcal/ mol)/Å for SD, 10 (kcal/mol)/Å for CG, and 0.001 (kcal/mol)/Å for NR. TABLE 1: Energy Difference between the Two Structures at the Lowest Energy (∆Emin, kcal/mol), Average of the Ten Sampled Structures (Eav, kcal/mol), and Related Standard Deviation (σ) geometry
∆Emin
Ia Ib Ic
24.2 28.3 26.5
Ia Ib Ic
28.2 36.7 25.1
S1
S2
Molecular dynamics simulations were run under NVT conditions at 298 K (Berendsen thermostat21 was used with a decay constant of 1 ps) with a 1 fs time step. The initial structure of silica cluster was built from an R-quartz unit cell by replicating it 3 × 3 × 3 along the a b c directions, the supercell obtained was cleaved along the 001 plane between two oxygen layers, and the oxygen atoms in the bottom layer were saturated with hydrogen. In the methylated surface the oxygen atoms were replaced by carbon atoms and saturated with hydrogen atoms, so that the dimensions of the silica slab were 6.15 × 3.72 × 0.96 nm. The simulations have been performed at pH 7.16, which is the calculated isoelectric point for EAK 16-II,22 i.e., the pH value when all the sites on the silanolic surface are deprotonated, due to the fact that the isoelectric point of SiO2 is found at pH ) 2.23 Such a deprotonated surface will be henceforth be noted as the S1 surface, and the fully methylated surface will be indicated as S2. The EAK 16-II (C68H118N20O25, MW ) 1615.8) oligopeptide, is a sequence of amino acids AEAEAKAKAEAEAKAK (- + + - - + +), where A ) alanine, E ) glutamic acid, and K ) lysine (see Figure 1). A has a neutral hydrophobic group, E and K have negatively and positively charged hydrophilic groups, respectively, and in parentheses is reported the charge distribution summarized by the Roman number II, which indicates the number of the same TABLE 2: Electrostatic and Van der Waals Contributions (kcal/mol) of Binding Energy for Ia, Ib and Ic Forms onto S1 and S2 Surfaces
Eav
σ
geometry
energy
electrostatic contribution
dispersive (VdW) contribution
-50042.1 -50604.7 -52634.4
12.9 13.4 13.2
Ia Ib Ic
-50016.3 -50573.7 -52610.2
S1 -50921.5 -51861.4 -54456.7
-1014.6 -1010.5 -1020.8
-44536.1 -44399.5 -45719.9
12.0 13.2 12.9
Ia Ib Ic
-44507.4 -44365.3 -45689.1
S2 -46614.3 -46343.8 -47312.6
-399.8 -401.5 -390.4
Molecular Modeling of Oligopeptide Adsorption
J. Phys. Chem. B, Vol. 111, No. 38, 2007 11239
Figure 3. Dihedral angles Ψ and Φ between LYS-8 and ALA-9.
Figure 4. Time evolution of dihedral angles Ψ and Φ of EAK 16-II onto S1 at 298 K.
kind of charges grouped together. To obtain an initial β-sheet configuration, the Ψ and Φ angles were set to 180°. In all simulations the oligopeptide was always considered as a zwitterion so that the four glutamic acids and four lysines were negatively and positively charged, respectively.11,24,25 The following three initial oligopetide orientations with the respect to the surface (see Figure 2) have been chosen: (a) an orthogonal orientation with the -NH3+ terminal group pointing toward the surface (end-on, Ia); (b) an orthogonal orientation with the -COO- terminal group toward the surface (end-on, Ib); (c) a TABLE 3: ∆E (Electrostatic and Van der Waals, kcal/mol) between the System Oligopeptide Surface, without Water Molecules, and Isolate Oligopeptide and Surface configuration
∆E electrostatic elec elec (E Ic+surf - E elec pep - E surf )
∆E VdW VdW VdW (E Ic+surf - E VdW pep - E surf )
Ia-S1 Ia-S2
-219.9 -91.3
-47.3 -26.2
configuration
∆E electrostatic elec elec (E Ia+surf - E elec pep - E surf )
∆E VdW VdW VdW (E Ia+surf - E VdW pep - E surf )
Ib-S1 Ib-S2
-215.6 -88.4
-46.2 -24.8
configuration
∆E electrostatic elec elec (E Ib+surf - E elec pep - E surf )
∆E VdW VdW VdW (E Ib+surf - E VdW pep - E surf )
Ic-S1 Ic-S2
-842.2 -396.8
-275.7 -98.8
parallel orientation of the backbone with respect to the surface (side-on, Ic). Structures Ia, Ib, and Ic were placed in a box of 6.5 × 4.0 × 12.0 nm, together with 1024 water molecules; hence we have analyzed the following six configurations: Ia-S1, Ib-S1, IcS1, Ia-S2, Ib-S2, and Ic-S2. For each configuration an equilibration protocol was performed, consisting of 1500 iterations of Smart Minimizer optimization applied to the whole system, followed by a 90 ps MD simulation at 298 K. The RDF and the potential energy fluctuation analysis (not reported here) provided the evidence that the system was equilibrated after the above procedure. After this simulation period, the system was run for other 10 ps, during which 10 structures were randomly sampled for each configuration, then an energy minimization procedure was applied at each of these structures and the structure corresponding to the lowest energy value was used for the successive molecular dynamics simulations. It should be noted that, for each of the six configurations, we have found a structure having a much lower energy than the nine other structures. In Table 1 are reported the energy differences between the two structures at lowest energy (∆Emin), the average energy of the ten sampled structures (Eav), and the related standard deviation (σ) for each configuration. As one can see from the reported data, the lowest energy structure is meaningfully more stable than the nine others, for all the six configurations here investigated. Accordingly, we
11240 J. Phys. Chem. B, Vol. 111, No. 38, 2007
Forte et al.
Figure 5. Time evolution of dihedral angles Ψ and Φ of EAK 16-II onto S2 at 298 K.
have chosen these minimum energy structures as the starting point for molecular dynamics simulations that were run for 2 ns while the coordinates were sampled every 0.05 ps for successive analysis. Results and Discussion A. Analysis of the Optimized Configurations. The structures identified, for each configuration, according to the procedure described in the previous section have been analyzed in term of the various energy contributions. Table 2 reports the energy values (kcal/mol) as well as the corresponding energy partitions in terms of electrostatic and Van der Waals contributions. The other energy contributions (torsional, stretching, bending, etc., not reported in Table 2) are positive, as one can determine from the difference between the total energy and the sum of the Van der Waals and Coulombic terms. This fact is possibly due to the expected strong intramolecular or molecule-surface electrostatic interactions, among the charged groups respectively on the oligopeptide sequence and on the surface, which may well determine conformational tensions. From Table 2, the Ic-S1 configuration, corresponding to the oligopeptide backbone parallel to the S1 surface, appears to be the more stable among all the configurations here considered, with an energy difference of 2593.9 kcal/mol for Ic-Ia and 2036.5 kcal/mol for Ic-Ib, respectively. Moreover, from the geometry optimization point of view, for the Ic-S1 configuration the mean distance of the oligopeptide chain from the surface is about 0.3 nm, suggesting that the oligopeptide is adsorbed on the surface. This distance was calculated as the average among the minimal distances of each atom in the oligopeptide and the closer surface atoms. At variance with this, the calculated energy differences between the Ic-Ia and Ic-Ib configurations with respect to the S2 surface are 1181.7 and 1323.8 kcal/mol, respectively. In this case, the mean distance of the most stable Ic structure from S2 is 0.5 nm, suggesting a weaker interaction. The above hypothesis is confirmed by the energy differences calculated according to elect elect elect elect E int pep-surf ) E pep+surf - E pep - E surf
(1)
elect is the electrostatic energy of the system where E pep+surf elect oligopeptide + surface and E elect pep and E surf are the electrostatic energies of the isolate oligopeptide and surface, respectively (the same approach was followed for the calculation of the Van der Waals interaction term). Equation 1 is the Coulombic contribution due only to the oligopeptide-surface interaction and has been calculated by using the geometry corresponding to the minimum energy and subtracting the terms due to the water molecules. The results, reported in Table 3, show that the electrostatic as well as Van der Waals energy interactions the between the oligopeptide and the surface is higher for the Ic-S1 configuration with respect to all the other configurations. Indeed, this term does not change appreciably for both the Ia-S1 and Ib-S1 configurations and the trend of the values relative to the two surfaces remains essentially the same. In conclusion, the analysis of Table 3 shows that, from a minimum energy configuration point of view, the Ic form may easily be adsorbed on the S1 surfaces, because the high-energy stabilization due to the electrostatic and Van der Waals interaction energy. B. Molecular Dynamics Results. To get a deeper understanding of the behavior of the EAK 16-II sequence on the model surfaces, we have performed molecular dynamics simulations, starting from the optimized geometries described above. At first, we have studied the relative rigidity of the oligopeptide backbone with respect to both S1 and S2 surfaces monitoring the dynamical behavior of the 14 couples of Ψ and Φ angles for the nine different configurations discussed above (see Figure 3). In particular, it turned out that the Ψ and Φ couple for the eighth amino acid from the N end in the EAK 16-II sequence, i.e., the central LYS unit, shows the lowest mobility with respect to the other angle couples, allowing us to choose this couple as the reference to evaluate the global mobility of the whole oligopeptide sequence. The dynamical analysis of the various angle values has been performed by averaging every 1000 steps. Such analysis for the three forms Ia-Ic onto S1 and S2 surfaces, respectively, reported in Figures 4 and 5, showed that the Ic form undergoes the smallest fluctuations around the characteristic β-sheet value
Molecular Modeling of Oligopeptide Adsorption
J. Phys. Chem. B, Vol. 111, No. 38, 2007 11241
Figure 6. Ramachandran Φ, Ψ distribution of EAK 16-II onto S1 obtained after 2 ns at 298 K.
Figure 7. Ramachandran Φ, Ψ distribution of EAK 16-II onto S2 obtained after 2 ns at 298 K.
Figure 8. Time fluctuation of the distance of Ia, Ib, and Ic forms from (a) S1 and (b) S2 surfaces.
for Ψ and Φ angles onto the S1 surface. In other words, the Ic form substantially maintains the β-sheet configuration on this surface. At variance of this, the Ic form is more flexible on the S2 surface, as is the case for the Ia and Ib forms for both S1 and S2 surfaces. The observed trend is confirmed by the analysis of the Ramachandran plots, reported in Figures 6 and 7. These plots show indeed the predominance of the β-sheet configuration for the Ic forms onto S1, as well as the apparent loss of such a β-sheet configuration for the S2 surface. The stabilization of the β-sheet configuration onto the S1 surface should be attributed to the strong electrostatic forces characteristic of the oligopeptide interactions with hydrophilic sites on this surface. Indeed, the simulation of the behavior of a single oligopeptide molecule in solution (results not shown) demonstrates that no β-sheet configuration is found in that case. In other words, the strong intermolecular interactions with the hydrophilic S1 surface are the stabilizing factor for the β-sheet configuration.
Once we determine the relative configurational stability of the various forms on the S1 and S2 surfaces, we study in detail the time evolution of the distance of the oligopeptide sequence from the S1 and S2 surfaces (see Figure 8a,b). For the side-on orientation this distance is an averaged minimal distance that is obtained by calculating the shortest distance between the surface and each atom of the oligopeptide sequence, in a given configuration at time [t], and then averaging over 0.4 ps. This value was adopted to obtain the fluctuation, over 2 ns, of the distance between the surface and the oligopeptide. The same procedure was followed for the end-on orientations but in this case the minimal distance from the surface of the terminal atom of the oligopeptide, on one hand N for end-on Ia and, on the other hand, C for end-on Ib, was considered. Figure 8a shows that for the S1 surface the molecule keeps a constant average distance of about 0.3 nm from the surface, indicating that, in agreement with the above-described results of geometry optimization (section A above), the oligopeptide can be considered just physisorbed on the surface. At variance with these results,
11242 J. Phys. Chem. B, Vol. 111, No. 38, 2007
Forte et al.
Figure 11. Variation of the H-bond distance among the nearest O- of the S1 surface and the amidic hydrogen of the residue 9-ALA. Figure 9. Time evolution of Ia onto S1 surfaces (a) at the beginning of the simulation, (b) after 1 ns and (c) after 2 ns. Solvent molecules are omitted for clarity.
respect to the one expected for a true H-bond formation, i.e., 0.13-0.27 nm,26 and it remains constant during the system evolution. This fact seems to rule out any significant H-bond formation, supporting our previous considerations about the dominance of electrostatic interactions with respect to the β-sheet configuration stabilization and oligopeptide physisorption on the various surfaces. Conclusions
Figure 10. Time evolution of Ib onto S1 surfaces (a) at the beginning of the simulation, (b) after 1 ns and (c) after 2 ns. Solvent molecules are omitted for clarity.
the continuously increasing distance with simulation time found for the oligopeptide from the S2 surface (Figure 8b) clearly indicates that the molecule freely diffuses away from the fully methylated surface, so that no effective physisorption can be considered for Ia, Ib, and Ic. The simulation performed by starting with a different initial conformation of the Ic form onto S2, involving all the methyl groups oriented toward the surface, showed the same behavior (results not shown). It should be noted that molecular dynamics simulations demonstrate that for Ia-S1 and Ib-S1 the initial end-on orientation tends to be converted to a side-on one by means of a peculiar mechanism involving the progressive lie back of the molecule onto the surfaces, as shown in Figures 9 and 10. Accordingly, the molecular dynamics simulations also have shown that at the interface between the oligopeptide and S1 surface there is no statistically meaningful number of water molecules, whereas a continuous water layer has been found at the S2-oligopeptide interface, this last effect being obviously related to the mentioned outdiffusion of EAK 16-II sequence. In the force field used in the present work no special hydrogen bond function is used18,19 and hydrogen bonds are a natural consequence of the standard Van der Waals and electrostatic parameters; therefore nonbond terms account for the possible interactions among the amidic hydrogens in the oligopeptide and the Si-O- groups on the S1 and S2 surfaces. To analyze the relative importance of this contribution, we have monitored the distances between the SiO- groups on the S1 surface and each amidic hydrogen in the oligopeptide, reporting for the sake of clarity only the behavior of the amidic hydrogen in 9-ALA residue in the sequence as representative of the behavior of the others. The results, reported in Figure 11, clearly show that the calculated average distance 0.29-0.30 nm is very high with
In the present paper we have reported the results of geometry optimization and molecular dynamics simulations concerning the interaction between EAK 16-II oligopeptide sequence and two kinds of functionalized R-quartz surfaces, respectively consisting in a hydrophilic surface, containing silanolic groups, indicated as S1, and a hydrophobic one, terminated with methyl groups, indicated as S2. Two end-on and one side-on orientation of the oligopeptide sequence were taken into account with respect to each surface. The energy optimized values for the various molecule-surface configurations indicated that the side-on orientation is by far more stable than the end-on ones onto hydrophilic S1 surface, prompting an effective physisorption in that case. This finding is also supported from the fact that the end-on orientations are demonstrated to progressively lie back on such surface. On the contrary, in the case of hydrophobic S2 surface, also if the sideon orientation is again the most stable, no effective physisorption could be found for the end-on as well as the side-on orientations, every form moving away from the S2 surfaces in the employed simulation time. Furthermore, we have also found that, within the molecular dynamics simulations scheme, the average distance between the oligopeptide and the hydrophilic S1 surface is relatively shorter (0.30 nm) than the one calculated for the highly hydrophobic S2 surface. These results suggest that pure electrostatic interactions are the driving factor for effective oligopeptide physisorption. Indeed, these interactions promote a larger contact area between the surface and the oligopeptide, strongly favoring the side-on orientation adsorption onto S1 surface. Finally, the electrostatic interactions for the S1 surface are also demonstrated to be responsible for the observed high rigidity of the backbone oligopeptide, resulting in a very stable β-sheet configuration on this surface, as is confirmed by the time evolution analysis of Φ and Ψ dihedral angles. It is to stress that the β-sheet configuration is completely absent for the oligopeptide interacting with hydrophobic S2 surface. The results obtained in the present paper are a firm basis to further extend the study to the important case of interacting sequences onto surfaces, which is mandatory for understanding in detail the primary steps of the experimentally observed fibril
Molecular Modeling of Oligopeptide Adsorption and membrane assembly process. The attention should be focused in particular on the way the orienting and ordering effects, mainly due to the molecule-surface electrostatic interactions, may affect the growth mechanism of supramolecular structures. Acknowledgment. We acknowledge the financial support of the National Program on “Self-assembling films of synthetic oligopeptides for biomimetics surfaces” (COFIN 2005-MIURRome) as well as the grant for “Ricerca di Ateneo 2006” (University of Catania). References and Notes (1) Csa`sza`r, A. G. J Am. Chem. Soc. 1992, 114, 9659. (2) Tun˜o´n, I.; Silla, E.; Ruiz Lo´pez, M. F. Chem. Phys. Lett. 2000, 433. (3) Csa`sza`r, A. G. J Phys. Chem. 1996, 100, 3541. (4) Shirazian, S.; Gronert, S. J. Mol. Struct. (THEOCHEM) 1997, 397, 107. (5) Tortonda, F. R.; Tun˜o´n, I.; Silla, E.; Ruiz Lo´pez, M. F.; Rinaldi, D. Theor. Chem. Acc. 2000, 104, 89. (6) Pertsemlidis, A.; Safena, A. M.; Soper, A. K.; Head-Gordon, T.; Glaeser, R. M. Proc. Natl. Acad. Sci. U.S.A. 1996, 93, 10769. (7) Carravetta, V.; Monti, S. J. Phys. Chem. B 2006, 110, 6160. (8) Horbett, T. A. Protein Adsorption on Biomaterials. In Biomaterials: Interfacial Phenomena and Applications; Copper, S. L., Peppas, N. A., Eds.; American Chemical Society: Washington, DC, 1982; p 233. (9) Gambino, G. L.; Marletta, G.; Grassi, A.; Lombardo, G. M. J. Phys. Chem. B 2004, 108, 2600.
J. Phys. Chem. B, Vol. 111, No. 38, 2007 11243 (10) Gambino, G. L.; Marletta, G.; Grassi, A. J. Phys. Chem. B 2006, 110, 4836. (11) Zhang, S.; Holmes, T.; Lockshin, C.; Rich, A. Proc. Natl. Acad. Sci. U.S.A. 1993, 90, 3334. (12) Zhang, S.; Cook, R.; Lockshin, C.; Rich, A. Biopolymers 1994, 34, 663. (13) Holmes, T. C.; de Lacalle, S.; Su, X.; Liu G.; Rich, A.; Zhang, S. Proc. Natl. Acad. Sci. U.S.A. 2000, 97, 6728 (14) Holmes, T. C. Trends Biotechnol. 2002, 20, 16 (15) (a) Hong, Y.; Protzker, M. D.; Legge, R. L; Chen, P. Colloids Surf. B: Biointerfaces 2005, 46, 152. (b) Cataldo, S.; Satriano, C. and Marletta, G. Private communication. (16) Accelrys Inc., MS Modeling Getting Started, San Diego: Accelrys Inc., 2003. (17) Dauber-Osguthorpe, P.; Roberts, V. A.; Osguthorpe, D. J.; Wolff, J.; Genest, M.; Hagler, A. T. Proteins: Struct., Funct. Genet. 1988, 4, 3147. (18) Hagler, A. T.; Dauber, P.; Lifson, S. J. Am. Chem. Soc. 1979a, 101, 5131 (19) Hagler, A. T.; Dauber, P.; Lifson, S. J. Am. Chem. Soc. 1979b, 101, 5122. (20) Ewald, P. Ann. Phys. 1921, 64, 253. (21) Berendsen, H. J. C.; Postma, J. P. M.; Van, Gunsteren, W. F.; Di Nola, A.; Haak, J. R. J. Chem. Phys. 1984, 81, 3684. (22) Theorically calculated on www.embl-heidelberg.de. (23) Parks, G. A. Chem. ReV. 1965, 65 (2), 177. (24) Fung, S. Y.; Keyes, C.; Duhamel, J.; Chen, P. Biophys. J. 2003, 85, 537. (25) Jun, S.; Hong, Y.; Imamura, H.; Ha, B. Y.; Bechhoefer, J.; Chen, P. Biophys. J. 2004, 87, 1249. (26) Buemi, G. In Hydrogen BondingsNew Insights; Grabowski, S. J., Ed.; Springer: Berlin, 2006; Chapter 2 (in press).
