J . Phys. Chem. 1985,89, 3655-3663
3655
Hydration Analysis of the Intercalated Complex of Deoxydlnucleoside Phosphate and Proflavin: Computer Simulations Kwang S. Kim* and E. Clementi IBM Corporation, Department 48B, Mail Station 428, Neighborhood Road, Kingston, New York I2401 (Received: September 24, 1984)
Using the metropolis Monte Carlo method and previously reported ab initio atom-atom pair interaction potentials, we have studied the networks of water molecules hydrating a 2:2 complex of proflavin cations and deoxycytidylyl-3',5'-guanosine. Simulations were performed with different numbers of water molecules in the crystal, and the energetics of these water molecules were then analyzed on the molecular level. We compare the positions of the water molecules with the X-ray data, and the hydration structures in the crystal are discussed in detail. We have investigated both the degree of hydration for each hydration site and the hydrogen bnding of each water molecule. This analysis indicates 60 to 66 water molecules per half-unit cell. The latter case is analyzed in detail.
Introduction The role of water in stabilizing and mediating the structure and conformation of nucleic acids has been studied and speculated upon for some time. This feature has been of particular interest among theoreticians in order to elucidatt the conformation of water molecules hydrating ions and biomolecules.14 A number of X-ray diffraction studies have recently focused on drug complexes with dinucleoside monophosphates1*20 in order to improve our understanding of biomolecular mechanism of drug action on nucleic acids. In particular, the network of water molecules in a crystal of a deoxydinucleosideproflavin complex was studied with X-ray diffraction by Neidle, Berman, and Shieh.loJ1 Their study identifiad 50 water molecules per half-unit cell of the crystal. However, our recent theoretical studies2' by computer simulation
(1) E. Clementi, 'Determination of Liquid Water Structure, Coordination Numbers for Ions and Solvation for Biological Molecules", in "Lecture Notes in Chemistry", Vol. 2, Springer-Verlag, Berlin, 1976. (2) 8. Clementi, "Computational Aspects for Large Chemical Systems", in "Lecture Notes in Chemistry", Vol. 19, Springer-Verlag. Berlin, 1980. (3) E. Clementi and G. Corongiu, In?.J. Quuntum Chem., 22, 595 (1982); Biopolymers, 21, 193 (1982). (4) E. Clementi in "Structure and Dynamics: Nucleic Acids and Proteins", E. Clementi and R. H. Sarma, Eds., Adenine Press: New York, 1983, pp 321-364. (5) P. J. Rossky and M. Karplus, J . Am. Chem. SOC.,101, 1913, (1979). (6) A. T. Hagler and J. Moult, Nature (London), 272, 222 (1978). (7) J. Hermans and M. Vacatello in "Water in Polymers", S.P. Rowland, Ed., American Chemical Society, Washington, DC, 1980, ACS Symp. Ser. No. 127. (8) H. J. C. Berendsen, J. P. M. Postma, W. F. van Gunsterer, and J. Hermans, 'Intermolecular Forces", B. Pullman, Ed., Reidel, Dodrecht, 1981, (9) J. M. Goodfellow, J. L. Finney, and P. Barnes, Proc. R. SOC.London, Ser. E. 214. 213 11982). (10) H. S . Shieh, H. M. Berman, M. Dabrow, an4 S. Neidle, Nucleic Acids Res., 8, 85 (1980). (11) S. Neidle, H. M. Berman, and H. S. Shieh, Nature (London).288. 129 (1981). (12) H. M. Sobell, C. C. Tasai, S. C. Jain, and J. Gilbert, J . Mol. Eiol., 114, 333 (1977). (13) H. M. Berman, S.Neidle, and R K. Stodola, Proc. Natl. Acad. Sci. U.S.A.,75, 828 (1978). (14) S. Neidle, Prog. Med. Chem., 16, 151 (1979). ( 1 5 ) S.Neidle, G. Taylor, M. Sanderson, H. S.Sbieh, and H. M. Berman, Nucleic Acids Res., 5, 4417 (1978). (16) E. Westhof, S.T. Rao, and M. Sundarlingam, J. Mol. Eiol., 142,331 (1980). (17) E. Westhof and M. Sundarlingham, Proc. Natl. Acad. Sci. U.S.A., 77, 1852 (1980). (18) H. S. Shieh, H.M. Berman, M. Dabrow, and S. Neidle, Nucleic Acads Res., 8, 85 (1980). (19) A. Aggarwal, S.A. Islam, R. Kuroda, and S. Neidle, Biopolymers, in press. (20) P. Swaminathan, E. Westhof, and M. Sundarlingam, Acta Crystallogr., Sect. B, 38,515 (1982).
favored more water molecules in the crystal. Beveridge et aLzz also reported a related simulation study using 50 water molecules; they found that the computed oxygen atom positions agreed with those from the X-ray data within 4040%. Neidle and BermanZ3 now speculate that there could be as many as 54 to 58 water molecules, because of mobile or disordered water molecules in the crystal. For this reason, it has been of interest to investigate theoretically the most energetically favorable number of water molecules in the crystal. These theoretical analyses of the water networks in the crystals are especially interesting because difficulties in X-ray studies can be complemented by theoretical studies using computer simulations. X-ray experiments can easily determine crystal structures, while computer simulations encounter difficulty because of the need for complicated interaction potentials. On the other hand, X-ray experiments can barely supply information about the orientations of water molecules or about hydrogen bonds, even in well-ordered crystal systems. Moreover, the X-ray studies cannot detect mobile or disordered water molecules, while computer simulations are well suited to obtaining this type of information. Our recent study based on (1) the total internal energy change of the water-crystal system due to hydration and (2) the pattern comparison of the water networks between the X-raylo*" and Monte Carlo (MC) simulations has shown that the most favorable number of water molecules per half-unit cell is between (1) 70 and (2) 65, respectively. The total internal energy change of the system is very sensitive to interaction potentials. Because our study shows more water molecules compared with the X-ray study, it is necessary to seek reasonable interpretations for an apparent discrepancy between the theoretical and experimeptal conclusions. Therefore, in this paper, we study the crystal system via hydration analysis. Because this hydration analysis is less dependent on the theoretical model (such as interaction pair potentials) than analyses of total energies or water network patterns, it should help in determining the correct number of water molecules in the crystal and in obtaining a reasonable interpretation for the discrepancies between computer simulations and X-ray experiments. From MC simulations, the effects of varying the number of water molecules (NW) from 50 to 82 per half-wit cell will be investigated in two aspects: (i) the degree of hydration of each hydration site of the solute and (ii) the hydrogen bonding and the energetics of each water molecule. In this way, the mobility of each water molecule and the role of each hydration site in the formation of the water (21) K. S. Kim, G. Corongiu, and E. Clementi, J . Eiomol. Struct. Dyn., 1, 263 (1983).
(22) M. Mezei, D. L. Beveridge, H. M. Berman, J. M. Goodfellow, J. L. Finney, and S. Neidle, J . Eiomol. Struct. Dyn., 1, 286 (1983). (23) Private communication from S. Neidle and H. M. Berman to K.S.K.
0022-3654/85/2089-3655$01.50/0 0 1985 American Chemical Society
3656
The Journal of Physical Chemistry, Vol. 89, No. 17, 1985
Kim and Clementi
Z(K)
I
1
I
I
I
I
I
I
I
I
I
I
I
1
I
b X(A7
1
bXW)
22.0A
11.0 .
0 .
-11.0 .
-22.0
-
1
I
-33
-24.75 -16.5 -8.25
I
I
0
8.25
16.5 24.75
33
Figure 1. Projection of 16 asymmetric unit cells with the P2,2,2 space group symmetry onto the x-z (top) and x-y (bottom) planes. Our MC simulation considers three such sets, stacked in three layers along the z axis. (*) For the analysis, the crystal is partitioned into four portions, denoted as 1, 2,3, and 4. Note the periodicity along the z axis with the Pf stocks in subvolumes 2 and 4 and the C-G stacks in subvolumes 1 and 3.
network can be determined. Thus, this analysis should detect those water molecules that cannot easily be determined from X-ray experiments. The crystal investigated is a self-complementary duplex of two dinucleoside phosphate strands with Watson-Crick hydrogen bonds in the 2:2 complex of proflavin and deoxcytidylyl-3',5'guanosine, Le., P f d (CpG), as shown in Figure 1. The group symmetry is P2,2,2. A unit cell consists of four asymmetrical units arranged in such a way as to define major and minor grooves. The length, width, and height of a unit cell are x = 32.991 A, y = 21.995 A, and z = 13.509 A, respectively. In order to facilitate the analyses, the crystal is partitioned along the z axis (height of a unit cell: 13.509 A) into four subvolumes indicated by 1, 2, 3, and 4 seen in the top insets of Figure 1. Subvolumes 2 and 4 show proflavin (Pf) stacks; subvolumes 1 and 3 show cytosine-guanine (C-G) stacks.
MC Simulations In the X-ray data, the water molecules are connected continuously along the z axis and are also connected between M (a major groove) and m (a minor groove) along the x axis (see Figure 1). However, they are not connected between M and m along the y axis (Le., the nearest ones are more than 5.5 A apart according to the X-ray data). Therefore, for the water-water interactions in our M C simulation, we can consider a half-unit cell as a primitive cell for the water molecules, while utilizing suitable symmetry operations with these water molecules for a periodic boundary condition. This one-half unit cell corresponds to the portion inside the rectangular box in Figure 1. If a primitive cell is reduced to a quarter unit cell as in the solute of the crystal, some water molecules are near the symmetry axis, and thus, the thermal statistical effect for these molecules cannot be properly
considered. Also, such symmetry constraint on water would be, in general, too strong and artificial for real complex-crystal systems. In our M C simulations, the d(CpG) and Pf species are assumed rigid at the positions given from the X-ray data; the attached hydrogen atoms for these species have been added with the standard bond lengths and angles. The boundary condition is not limited to a volume which would include only a single unit cell, as is often done for neutral solutes. The presence of the negatively charged groups (PO,-) and of the positively charged proflavins ( P P ) compells one to consider a much larger volume in order to take into account the long-range ionic field which extends further than the dimension of a unit cell. If either d(CpG) or Pf would be partitioned into more than two fragments, then the water molecules would experience a strong monopole electric field (an artifact of the partitioning) rather than a dipole electric field. Therefore, in our M C simulation we considered three sets of the 16 asymmetric unit cells of Figure 1, stacked in three layers along the z axis. For the water molecules in the major groove, the left 12 X 3 asymmetric unit cells were considered as the periodic boundary conditions, while for the water molecules in the minor groove, the right 12 X 3 asymmetric unit cells were used. In this way, the boundary condition was imposed in a very symmetrical manner as it should be, and the discontinuity of the potential near the border of major and minor grooves was really negligible in our M C study. In the M C simulations the previously reported a b initio atomatom pair interaction potential^*^*^^-^^ have been used. The sim(24) K. S. Kim and E. Clementi, IBM Technical Report, IS & TG, POK-36. (25) K. S. Kim and E. Clementi, J . Am. Chem. SOC.,in press.
The Journal of Physical Chemistry, Vol. 89, No. 17, 1985 3657
Proflavin-Guanosine Complexes
c
Figure 2. Hydrogen-bonding patterns for a simulation with NW = 50X,in which the oxygen positions are taken from X-rays, while the hydrogen atoms are free to rotate. The insets show x-y projections of subvolumes 1 through 4. The darkened circles are either oxygen or nitrogen atoms, showing hydrophilicity. The water molecules are numbered from 1 to 50 (cf. Table 11). Note the symmetric conjugation of each odd-numbered water molecules with the next-nearest even-numbered one, except for 49 and 50 which are self-conjugated. Note that P = phosphate, F = proflavin, C = cytosine, and G = guanine.
ulations were performed a t 0 OC for N W = 50, 57, 61, 69, 72, 77, and 82, following the method of Metropolis et aL31s3* The details of the simulation methods are the same as those presented in ref 25. For these simulations the water molecules were allowed to have all translational and rotational degrees of freedom; the solute molecules were assumed to be rigid. For the comparison of the M C simulations with the X-ray data, another simulation performed with only the orientational degree of freedom (namely, fixing the oxygen atoms to the X-ray positionsI0*") is labeled in this paper with an X, Le., N W = 50 X. For this simulation, 200000 M C steps were performed until a reasonable equilibration was obtained, followed by 200 000 M C steps for an analysis. Figure 2 shows the simulated result which is denoted as N W = 50X to distinguish it from the simulation denoted as N W = 50, where the oxygen positions are not frozen. The initial conformation of water molecules for the N W = 50 simulation was chosen with the equilibrium conformation of N W = 50X. However, since the water molecules for our system are crystallized and caged, multiminima problems should be considered. Therefore, we generated several different initial conformations by selecting a conformation every 20 000 M C steps in an M C simulation performed at 3000 K so as to erase the memory of the initial conformation. These conformations were very different from each other as well as from the X-ray data, with their total energies positive. After 200000 M C steps for each different simulation with a different initial conformation were performed at 0 OC, each total energy was almost converged. These total energies for each case were slightly different, while their conformations were rather different. By choosing the conformation with the minimum energy, an M C simulation was performed up (26) 0. Matsuoka, E. Clementi, and M. Yoshimine, J . Chem. Phys., 64,
1351 (1976). (27) G.Corongiu and E. Clementi, Gam. Chim. Ital., 108, (506), 273 (1978). (28) E. Clementi, G.Corongiu, and F. Lejl, J . Chem. Phys., 70, 3726 (1979). (29) R.Scordamaglia, F.Cavallone, and E. Clementi, J . Am. Chem. SOC., 99,5545 (1977). (30) G.Corongiu and E. Clementi, Int. J Quantum Chem., 9,213 (1982). (31) N. Metropolis, A. W. Rosenbluth, M. N. Rosenbluth, A. H. Teller, and E. Teller,J . Chem. Phys., 21, 1078 (1953). (32) J. A. Barker and R. 0.Watts, Chem. Phys. Lett., 3,144-145 (1968).
to 800000 MC steps; the last 200000 M C steps were chosen for the analysis. After equilibrium, seven water molecules migrated to a vacant region in the major groove; thus, there were at least seven vacant places in which there were no water oxygens and solute atoms within the radii of about -3.5 A. Therefore, seven water molecules were added in each vacant position. Then, this initial conformation of water molecules for N W = 57 was used to obtain several different samples, just as for N W = 50. After equilibration, this still gave four vacant places, in which no water oxygens and solute atoms were within the radius of 2.5 A. With the same procedure as in N W = 57, we generated several different initial conformations of N W = 61 and 66. When N W = 66, the most energetically favorable final equilibrium conformation did not have any reasonable vacant places. However, several water molecules were again added in the largest remaining vacant places. In this way, the initial conformations of water molecules were chosen for N W = 69, 72, 77, and 82, then several different initial confomrations for each case were obtained from a 3000 K randomization. All these MC simulations were performed in the same way as in N W = 50. In this way, we located the global minimum for each simulation by carefully choosing many different initial samples. This procedure gives a reasonable assurance that we have located the absolute minimum (or at least one of the most significant conformations). The final total equilibrium energies of the watersolute system were -60.2, -71.4, -69.0, -67.4, -66.3, -65.3, -64.7, -61.8, and -59.3 kJ/mol of water for N W = 50, 57, 61, 66, 69, 72, 77, and 82, respectively. The standard deviations for each case were less than 0.1 kJ/mol (see ref 25 for the total solvation energies of the water-solute crystal system). Hydration Sites of the Solute As N W changes, the M C simulations indicate that some sites of the solute molecules are always well hydrated, while the hydrations of some other sites are strongly dependent upon the variation in NW. In the following discussion, a phosphate group is denoted by P,a proflavin cation by F, cytosine by C, guanine by G, a sugar residue by S,and the terminal atoms (end site) of the sugar by E (see Figure 2 ) . Figures 2-6 show the hydration structures of the M C simulations for N W = 50X, 50, 57, 61, and 66, respectively. In the figures water molecules have been num-
3658 The Journal of Physical Chemistry, VoE. 89, No. 17, 1985
Kim and Clementi
A
Figure 3. Hydrogen-bondingpatterns for the sublimation with NW = 50 (cf. Figure 2).
P
4
Figure 4. Hydrogen-bonding patterns for the simulation with NW = 57 (cf. Figure 2).
bered: numbers 1 to 50 correspond to the 50 water molecules identified in the case N W = 50X. Numbers 51 to 66 are reserved for crystals with more than 50 water molecules; N W = 66 is used as our reference system. Each odd-numbered water molecule is symmetrically conjugated with the next highest even-numbered one, except for the water molecules numbered as 49, 50,5 1, and 52. Water molecules 49 and 52 are self-conjugate molecules which do not have any symmetrically corresponding partner. Water molecule 50 has 51 as partner; except for the simulations with N W = SOX and N W = 50, where 50 is self-conjugate. Table I shows that some hydration sites are always hydrated by the same water molecules in the same configuration regardless of NW. For example, the phosphate group, P 1 / 0 (see notation of Figure l ) , is hydrated by 39 P 2 / 0 by 33, and P 4 / 0 by 13, 15, and 21. The proflavin group, F2/NH,/’ is hydrated by 33, F2’/NH2/’ by 39, F 3 / N H by 11, F3/NH; by 17, and F3’/NH; by 13. The cytosine group C1/NH2 is hydrated by 7; for the
guanine groups, G l / O is hydrated by 4, G3-0 by 19, and G3/N7 by 23. The site where the sugar is connected to guanine (SG2/0) is hydrated by 47; and the end site which is connected to S G 2 / 0 (EG2/0) is hydrated by 41. On the other hand, Table IB shows that some hydration sites are strongly dependent on the value of NW. The site G1/N7 is hydrated by 5 for N W = SOX and N W 1 61 but is not hydrated for N W I57. The F l / N H and F l / N H sites are hydrated by 25 and 31 for N W = SOX and N W 1 61, respectively, while they are barely hydrated in the simulation with N W I57. In general, as long as N W is larger than 6 1, the hydration appears to be in accordance with that reported in the X-ray analysis (cf. N W = SOX). However, for N W 2 50 the Fl/NH;, F1/NH2”, and C3/NH2 sites are hydrated by 54, 56, and 57, respectively; while for N W = 50X they are not hydrated. These observations tend to support the conjecture that the X-ray experiment has likely missed six water molecules: the three discussed above and the
The Journal of Physical Chemistry, Vol. 89, No. 17, 1985 3659
Proflavin-Guanosine Complexes
4
Figure 5. Hydrogen-bonding patterns for the simulation with NW = 61 (cf. Figure 2).
P
3
4
I
Figure 6. Hydrogen-bondingpatterns for the simulation with N W = 66 (cf. Figure 2).
three related ones, 53, 5 5 , and 58. Reference 25 has reported that the X-ray position of 50 would be the ensemble-averaged values of doubly occupied water molecule positions. Including molecule 5 1 (Le., the conjugate hydration site of 50) which hydrates both G2’/N3 and G2/NH2, seven water molecules seem to be missing from the X-ray data. Although the hydration energies of the sites discussed above are weak compared with those of the phosphate groups, the hydrophilicity of these sites is generally believed to be strong enough to be hydrated. A similar situation is observed for the sites S C l / O , S C 2 / 0 , and EG4/0H, which are hydrated by 39, 33, and 29, respectively. These water molecules are reported in the X-ray data, but they do not appear to hydrate these sites (cf. Figure 2). The resolution limit of the X-ray data is approximately 0.8 A. If the experimental error in the positions of the oxygen atoms is considered, the water molecules near the sites S C 1 / 0 , S C 2 / 0 , and E G 4 / 0 H in the
simulation N W = 50X can hydrate the solute, in agreement with other simulations. Another discrepancy between the experiment and our simulations exists for the sites P 3 / 0 and E C l / O H . P 3 / 0 is hydrated by 9 and 17 for N W 5 51, and begins to be weakly hydrated by 61 for 61 I N W I 66, and finally is hydrated by all three for NW 2 72. ECl/OH is not hydrated for NW < 66, but begins to become hydrated by 61 for N W = 66, and finally is fully hydrated by 61 for N W = 72. That is, molecule 61 hydrates (i) neither sites of P 3 / 0 and E C l / O H for N W < 61, (ii) weakly hydrates only the P 3 / 0 site for N W = 61, (iii) weakly hydrates both sites for N W = 66, and (iv) fully hydrates both sites for N W 2 72. This result may be attributed to the hydrogen orientation of the terminal site of the sugar, which was not optimized to be at the most energetically favorable position. In a more stable orientation, molecule 61 might not even be present. However,
3660 The Journal of Physical Chemistry, Vol. 89, No. 17, 1985
Kim and Clementi
TABLE I: Hydration Analysis for Hydration Sitesa A. No Hydration Change in spite of Change of NW P1/0 F2/NH2/’ CI/NH2 SG2/0
39 33 7 47
P2/0 F2’/NH2/’ G1/0 EG2/0
33 39 4 41
P4/0 F3/NH G3/0
13, 15, 21 11 19
F3/NHi G3/N7
17 23
F3’/NH2’
13
B. Change of Hydration due to Change of NW 50X G1 / N 7
F ~ ~ N H F1 /NH2 Fl/NHi Fl/NH2/’ C3/NH2 sc 1/o sc2/0 EG4/0H P3/0 ECl/OH G2/N3 GZ’/NH, c4/0 G4/N3 G4/NH2 SG4/0
50
5 (25) (31) 54 30b 57 39 33 29
* *
5OC 50‘ (43h (45)
*
57
* * * 56
*
*
* * * * 50 50 43 45 43.49
61
5 25 31
* * * * * *
66
*
72
*
* * $
* * * *
*
*
9, 17, 61 61
* * * *
* *
9, 17, (61)
*
*
$
1
*
“See Figure 2 for the notation. An asterisk denotes that this value is the same as the value in the column immediately to the left. Parentheses denotes very weak hydration. b30 actually corresponds to the position of 56. c N o symmetric conjugate for the number, i.e., self-conjugate.
Figure 7. Probability distribution (the top four insets) and isoenergy density maps (the bottom four insets) projected onto the x-y and x-z planes for the simulated water molecules of N W = 66. The 0 and H denote the probability distribution maps of oxygen and hydrogen, respectively. The standard deviations of the oxygens and hydrogens are 0.36 and 0.43 A, respectively. The St and W denote the isoenergy density maps for the interaction energy of each water molecule at a given position with the solute (St) and with all the other water molecules (W), respectively.
The Journal of Physical Chemistry, Vol. 89, No. 17, 1985 3661
Proflavin-Guanosine Complexes
TABLE 11: Average Total Energies and DR of Each Water Molecule for Each NW Corresponding to That of NW = 50X’ av total energies, kJ/mol of water X-ray coord of 0 W no. X Y 2 50X 50 57 61 66 69 72 1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49 50
-4.5 -7.4 -6.0 -3.3 -0.6 -6.1 -5.5 -3.8 -1.7 -4.9 -3.7 -2.3 -5.9 -3.2 -2.0 -1.5 -1.1 0.3 1.4 2.5 3.5 5.9 6.0 7.0 8.2 8.2
-1.0 1.7 -2.6 1.1 1.3 2.4 -1.7 1.o 2.7 2.1 -0.9 3.5 2.3 1.6 -1.1 2.9 -3.9 0.2 2.4 4.8 0.7 0.6 0.6 0.1 0.2 0.0
1.2 13.3 13.3 13.3 13.1 10.1 10.3 10.4 9.5 6.7 6.7 6.9 3.6 3.1 3.3 4.2 4.0 1.9 3.3 2.5 3.5 5.0 8.1 10.5 4.0 1.1
-3 1 -46 -42 -70 -86 -64 -82 -80 -8 5 -46 -89 -43 -6 3 -29 -52 -3 5 -100 -42 -44 -97 -46 -65 -49 -68 -49 -62
-51* -67 -56* -6 1 -98 -67 -97 -83* -98 -65 -94 -77* -54* -57* -53 -45 -1 25 -53* -46 -121 -56 -6 5 -6 8 -6 9 -6 5 -69
-46* -60 -57* -7 2 -90 -65 -93 -89 -102 -73 -95 -64 -53 -57’ -43 -48 -1 25 -71. -53 -120 -58 -6 1 -64 -70 -67 -54
-4 1 -58 -43 -76 -89 -69 -9 1 -92 -9 3 -58 -103 -55 -6 2 -44 -55 -54 -125 -60* -4 5 -118 -56 -74 -4 5 -66 -69 -5 5
-40 -6 1 -40 -7 1 -89 -70 -9 7 -86 -90 -64 -103 -57 -60 -48 -56 -5 1 -122 -40 -48 -123 -58 -57 -6 8 -7 3 -67 -5 2
-4 3 -57 -48 -73 (-74 -68 -94 -8 2 -86 -69 (-82 -6 1 -56 -47 -44 -5 1 -123 -4 1 -59 (-108 -59 -56 -69 -7 1 -64 -5 3
-56* -60 -40 -57 -77 -7 3 -99 (-58 -95 -66 -9 7 -78 -59 -50 -54 -4 5 (-97 -41* -48 -116 -6 1 -6 1 -60 -68 -67 -52
77
82
-51 -61 -39 -64 -74 -62 -91 -59 -93 -57 -83 -77 -64 -44 -44 -42 -94 -36 -49 -109 -48 -58 -61 -67 -73 -55
-35 -57 -44 -61 -78 -62 -95 -56 -86 -45 -82 -44 -50 -38 -48 -46 -94 -38 -50 -110 -44 -53 -60 -51 -62 -53 ~
W no. 1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49 50
50 2.9* 0.9 2.7* 1.3 1.6 0.4 0.3 1.3* 0.7 1.4 0.8 1.5* 1.2* 2.3* 1.8) 1.2 0.7 9.2* 0.9 0.8 0.7 0.6 1.1 0.2 0.4 1.3
57 2.5*) 0.7 2.9*) 1.2 1.4 0.5 0.5 1.5 1.1 1.1 1.o 1.5 0.9 2.6.) 1.1 1.1 0.8 5.1* 0.7 0.8 0.5 0.6 0.8 0.4 0.3 1.5
61 0.4 0.7 0.6 0.2 0.8 0.4 0.4 0.5 0.3 0.4 0.8 1.o 0.3 0.8 0.8 1.7 0.8 2.0*) 1.4 0.8 0.5 0.6 1.4 0.2 0.3 1.5
66 0.4 0.7 0.7 0.3 0.7 0.4 0.4 0.5 0.3 0.5 0.7 1.o 0.3 0.6 0.7 0.7 0.7 0.5 0.5 0.8 0.4 0.6 0.9 0.4 0.5 1.5
69 0.8 (1.6 0.9 0.4 1.2 0.5 1.2 0.9 0.4 1.2 1.3 1.3 1.3 0.9 1.3 1.4 0.8 1.3 0.9 0.6 0.6 0.6 0.9 0.5 0.6 1.5
72 (2.0* 1.1 0.6 1.2 0.9 0.8 0.3 1.4 1.o 1.3 1.1 1.9 1.1 1.1 1.3 1.2 1.1 (1.8* 0.6 0.7 0.2 0.1 1.o 0.5 0.3 1.5
77 0.9 1.4 0.9 1.1 0.9 1 .o 0.8 1.4 0.9 1.2 1.4 1.8 1.1 1.1 1.o 0.8 0.9 1.3 1.o 0.6 0.9 0.6 1.1 0.8 0.4 1.5
82 0.9 1.1 0.8 0.9 1.2 1.2 0.5 1.5 1 .o 1.5 1.4 1.1 1.2 1.2 0.6 1.3 1.o 0.8 0.6 0.7 0.9 1.1 1.1 0.9 0.2 1.5
hydration sites X X
0 0 0
0
G1‘/0 G1/N7 Cl/NH2 P3/0 F3/NH P4/0 P4/0 P3/0 G3/0 P4/0 G3/N7 Fl/NH,
X
0
(EG4/0H) F1 / N H Pl/O,FZ/NH,”,SCZ/O
X
0 0 0 0
P2/0,F2’/NH,”,SCl/O EG2/0 C4/O,G4/N3,G4/NH, C4/O,G4/N3,G4/NH, SG2/0
0
G2/N3,G2’/NH2
Units of coordinates and DR are A. The values corresponding to the symmetrically conjugate water molecules are averaged. No symmetric conjugate for 49. No symmetric conjugate for N W = 50X and 50; otherwise the conjugate is 51. An asterisk indicates no corresponding position to the X-ray data. (, sudden jump from left to right. ), sudden jump from right to left. x, in general, no agreement with X-ray data except 60 < N W < 70. 0,in general, very good agreement with X-ray data for any NW. this extra water molecule hardly affects the energetics and pattern analysis since this position does not correspond to one of the X-ray positions. Consideration should also be given to water molecule 65, which is near 45. The average positions of these two molecules (45 and 65) in N W = 66 correspond to position 45 of the X-ray data. For N W = 66, the S G 4 / 0 (and partially G4/N3) is hydrated by 65, while G4/N3 and G 4 / N H 2 are hydrated by 45. Since this phenomena did not occur for N W < 66, water molecule 65 may not be a reliable one, though it is energetically preferred. For the N W = 66 case there are 16 extra water molecules compared to the X-ray data. Seven water molecules (51 and 53 to 58) hydrate the solute; two molecules (61 and 62) may not be favorable ones; and another two molecules (65 and 66) are still
questionable. The remaining five water molecules (52,59, 60, 63, and 64) do not hydrate, but make bridges between their neighboring water molecules. For example, in subvolumes 1 and 3, some sites of the bases are hydrated by the hydrogens of water molecules, while in subvolumes 2 and 4, proflavins are hydrated by the oxygens of water molecules. Since each subvolume generally has small space, the water molecules hydrating the solute sites at the positive y values in each subvolume are mostly oppositively oriented to those hydrating the solute sites at the negative y values (e.g., 3 and 4, 25 and 26, 11 and 12, 13 and 14, etc.). Also, most of these water molecules are oppositely oriented to the water molecules in the neighboring parallel subvolumes. That is, most of the water molecules in each subvolume are oppositely oriented to the neighboring water molecules along y and z axes.
Kim and Clementi
3662 The Journal of Physical Chemistry, Vol, 89, No. 17, 1985
NW=50X
a
NW=
a
NW=66
b
b
Figure 8. Hydrogen-bonding patterns (x-z projections) for the simulationsof NW = SOX, 57, 61, and 66. The b insets represent the networks of the water molecules. Note that the symmetrically conjugate networks are also possible. The a insets represent the 50 oxygen positions obtained by averaging the symmetrically conjugate ones that best meet one-to-one correspondence to the X-ray data. If the possibility for hydrogen-bonding exists, the oxygen atoms in the a insets are connected. The filled circles indicate the average positions of double occupancies.
This makes the water-water interaction repulsive. In addition, the distance between the centers of two parallel subvolumes is about 3.4 A, which is somewhat longer than the oxygen-oxygen distance for pure water. Since the water molecules are not near the center of the subvolumes, the insertion of the extra water molecules between the subvolumes tends to be energetically favorable by hydrogen bonding of these molecules with the neighboring water molecules. This situation is shown in Figure 4 with molecule 52. Molecule 52 is located between subvolumes 1 and 2. Each hydrogen of 52 is aligned toward the oxygens 3 and 4 in subvolume 1, while the hydrogens 53 and 54 of subvolume 2 are oriented to the oxygen of 52. Therefore, 52 plays the role of bridging the adjacent water molecules by hydrogen bonding, making the resulting configuration energetically favorable. Note that even the
shortest distance from 52 to the neighboring water molecules is more than 2.8 A, which implies the available space is large enough to accommodate one water molecule. Likewise, the molecules 59, 60, 63, and 64 play similar roles as 52. Therefore, although those extra nonhydrating water molecules seem to be energetically favorable, we may not be certain that all of them are really necessary. In conclusion, the hydration analysis favors 60 to 66 water molecules. Bonding Structure and Energetics of Water Molecules In this section the interactions between the water molecules and the solute, and among the water molecules themselves are studied. The relative stability of this system is investigated by analyzing the energy and positions of the water molecules as a function of NW. Table I1 shows the simulation results for each
Proflavin-Guanosine Complexes water molecule for different N W values (nine cases, including N W = SOX). The table also reports the X-ray oxygen coordinates, total interaction energy, distance difference (DR) between the X-ray oxygen positions and the corresponding simulated ones,and the hydration sites for each water molecule. These values are given as the average value of the two symmetrically conjugate molecules with the exception of 49. For 50, the exception is applied only for the N W = 50X and 50 cases. First, consider the water molecules that hydrate the phosphates, Le., 9, 13, 15, 17, 21, 33, and 37 (7 molecules per asymmetric unit cell, or 14 molecules per half-unit cell). These are the most characteristic molecules for water network formation and also for the energetics. The interaction energies are very strong (about 100 kJ/mol) and their distance deviations from the X-ray values (DR) are small. In particular, 33 and 39 partially hydrate even the proflavin and the sugar group, and for this reason their distance deviations are constant for any value of N W . These two water molecules behave almost like rigid molecules, with their hydrogen orientations nearly fixed. In addition, the interaction energies are very large (about 120 kJ/mol). As N W increases, most of the water molecules hydrating the phosphates show an energy jump near N W = 70 (marked by open parentheses in Table 11). This implies that additional water molecules tend to share the hydration of the phosphates, and thus, an excess of water molecules near N W = 70. Second, water molecules 1, 27, 35, 37, and 49 do not hydrate any sites of the solutes. The oxygen positions for 1, 27, and 35 are often very different from the X-ray positions, except for the case of N W = 66. The positions of 37 and 49 are somewhat fixed; though these molecules do not hydrate the solute, they are hydrogen bonded to other water molecules that are strongly bound to the solute, in,effect making 37 and 49 also rather nonmobile. In particular, the available space for 49 is small just so as to accommodate only one water molecule. Third, 5 does not hydrate G1/N7 unless N W I61, which directly contradicts the X-ray data. This result again implies that there may be too little water when N W < 61. Fourth, for N W = 50, seven water molecules (marked with asterisks) are deviated considerably from the X-ray positions and, consequently, show no correspondence between the simulated and the X-ray data. Sometimes both positions of the asymmetrically conjugate molecules are far from the X-ray positions, as in 35 and 36. However, we often find that one of the two symmetrically conjugate molecules corresponding to the aforementioned seven positions is in rather good agreement with the X-ray experiment, while the other is far from agreement. For N W = 57 there are four such water molecules deviated from the X-ray positions; for N W = 61 there is one such molecule, for N W = 66 there are no such molecules, and for N W = 72 there are two such molecules. In Figure 7 we report the probability and isoenergy maps for the N W = 66, because N W = 66 is one of the most favorable
The Journal of Physical Chemistry, Vol. 89, No. 17, 1985 3663 values among our simulated systems. The upper four insets show the probability distribution of oxygen and hydrogen, and the lower four insets show the isoenergy density maps of the interaction energy of each water molecule (at the given position) with the solute (St) and with all the other water molecules (W). The strongly hydrating water molecules have large solute-water interaction energies E(St-W) and small mobility. From the isoenergy density maps, the well-hydrated water molecules are darker in the St(x-y) and St(x-z) maps, while the water-bound molecules are darker in the W(x-y) and W(x-z) maps. Discussion and Conclusion The hydration analyses of our M C simulations indicate that there is notable agreement between computer simulations and X-ray experiments, especially for the cases of N W = 61 and 66. Figure 8 brings about this conclusion very clearly. For these cases the deviations from the experiments are minor, and those deviations which exist have been reasonably explained. We recall from ref 25 that the total energy of the system favored 70 water molecules and the pattern comparison of the simulated water networks with the X-ray data favored about 65 water molecules. Also recall that the pair interaction energy was somewhat overestimated. The pair potentials with the counterpoise c o r r e c t i ~ n ~were ~ - ~estimated ~ to reduce the most energetically favorable N W by about 5. Therefore, all analyses are now in good agreement with a value of N W between 60 to 66. The uncertainty in the number of water molecules in the crystal from our M C simulations arises from water molecules 61 and 65, and possibly one or two of the nonbound mobile water molecules. The uncertainty of 61 is attributed to the orientation of the terminal hydroxyl connected to the sugar group. Since we did not optimize those hydrogen positions of the solute, we cannot be certain whether 6 1 is really energetically favorable. However, the X-ray suggested that the disordering of the terminal groups of the solutes is natural and does not change the pattern of water molecules. This may be inferred from the argument that water molecules can hydrate either oxygens or hydrogens of the terminal hydroxyl. Therefore, the orientation of the terminal hydroxyl changes only the orientation of attached water molecules, yet generally would not change the number of water molecules.
Acknowledgment. K.S.K. thanks the National Foundation for Cancer Research (NFCR) for financial support. We also acknowledge the efforts of Drs. M. Probst and S. Chin for their critical review of the manuscript. ~
~
~~
(33) S.F. Boys and F. Bernardi, Mol. Phys., 19, 553 (1970). (34) W.Kolos, Theor. Chim. Acta, 51, 219 (1979). (35) A. Johansson, P. Kollman, and S. Rothenberg, Theor. Chim. Acta, 29, 167 (1973).
