J. Phys. Chem. 1995, 99, 1846-1855
1846
ARTICLES Molecular Dynamics Simulation Analysis of a Sodium Dodecyl Sulfate Micelle in Aqueous Solution: Decreased Fluidity of the Micelle Hydrocarbon Interior Alexander D. MacKerell Jr.* University of Maryland at Baltimore, School of Pharmacy, Department of Pharmaceutical Sciences, Baltimore, Maryland 21201, and Swarthmore College, Department of Chemistry, Swarthmore, Pennsylvania 19081 Received: August 8, 1994; In Final Form: October 25, 1994@
Structural and dynamic properties of a sodium dodecyl sulfate micelle were studied in aqueous solution via a molecular dynamics simulation using periodic boundary conditions. Results are presented for both the average structure and dynamic properties of the micelle. Over the course of the simulation the micelle remained spherical with a radius of gyration in agreement with experiment. Motions of individual lipid head groups were significant, with calculated changes of up to 7 8, occurring with respect to the micelle center of mass and 8 A parallel to the surface of the micelle. These motions were reminiscent of a piston in a cylinder or the movement of the head groups along the surface of the micelle. The micelle hydrocarbon interior is predicted to be less fluid than a pure alkane based on decreased dihedral transition rates and an increased free energy barrier to dihedral rotation of the aliphatic tails as compared to pure dodecane. This result contrasts calculations on a dipalmitoyl phosphatidylcholine lipid bilayer where the fluidity of the hydrocarbon interior was similar to that of pure hexadecane (Venable, R. M.; Zhang, Y.; Hardy, B. J.; Pastor, R. W. Science 1993, 262, 223-226). The predicted decrease in fluidity should be taken into account when micelles are used as model systems for lipid bilayers. The overall relative trans to gauche populations, however, are equivalent for the micelle and dodecane. Interactions between water and the micelle involve the sulfate head groups while the interior of the micelle is void of water. It is predicted, however, that the terminal methyl group of the hydrocarbon chain, in specific instances, may be located at the micelle surface and exposed to solvent. Interactions of the sodium ions with the micelle sulfate head groups occur primarily via the second hydration shell of the sulfate; no stable sodium to sulfate contact pairs were observed.
Introduction Micelles have been the subject of numerous experimental and theoretical studies in order determine the nature of their structure and dynamics. Interpretation of experimental results is often model dependent making theoretical studies of micelles relevant in order to obtain a molecular interpretation of experimental observables. Theoretical approaches initially applied to micelles employed lattice-based chain models. With increased computional resources, molecular dynamics (MD) simulations have recently been used, allowing for a more detailed microscopic understanding of the micellar properties to be ~ b t a i n e d . ~ MD -'~ studies have been performed on a variety of systems, including dodecyl surfactants, sodium octoanate, lysophosphatidylethanolamine, and n-decyltrimethylammonium chloride, in environments ranging from vacuum to fully solvated with counterions in the presence of periodic boundary conditions. Micelles have been exploited in experimental studies due to their potential as a model system for membrane bilayers. A wide variety of studies have used detergent micelles to explore the role of the lipid environment on peptide conformational properties, including structural determination via NMR approaches. Examples include glucagon,14the epidermal growth factor,15 and the bacteriophage Pfl coat protein.16 Several structural studies of peptides have been performed in SDS
* Address correspondence to the University of Maryland at Baltimore. @
Abstract published in Advance ACS Abstracts, January 1, 1995.
0022-365419512099-1846$09.00/0
micelle^.^^-^^ Additional incentive for the study of the structural properties of peptides in micelles is evidence for an essential role of lipids in the mediation of the receptor binding of peptide hormones.*l The relationship of peptide structure to environment and the essential role of lipids in physiological processes makes a molecular understanding of the structure and dynamics of micelles essential for the validation of their use as lipid model systems. A variety of physical measurements have been made to quantify the structure and dynamics of micelles in order to determine their validity as models of lipid bilayers. Experimental studies on dodecyl sulfate micelles include small angle neutron scattering experiments on lithium dodecyl sulfate micelles showing a radius of gyration of 15.4 A, corresponding to a sphere volume of 18.9 A, a maximum in the spatial distribution of terminal methyl groups at 10.5 A, water penetration to the first methylene group adjacent to the sulfate head group and the hydrocarbon interior to have a degree of order between that of a liquid and an all-trans state.** Solvent penetration to only the first methylene group has also been indicated by I9FNMR relaxation experiment^,^^ however, it has been suggested that the terminal methyl group has a high probability of being close to the micelle surface.24 Other NMR experiments yield correlation times of 34, 44, and 7 ps for the a, /3, and o C-H bond motions and predict an overall radius of the micelle of 20 8, based on the slow correlation time.25 X-ray scattering results show a paraffinic radius of 16.7 A and 0 1995 American Chemical Society
Molecular Dynamics of a Micelle
J. Phys. Chem., Vol. 99, No. 7, 1995 1847
a total radius of 22.3 Electron spin echo modulation To relax the numerous VDW contacts associated with the (ESEM) experiments on frozen sodium as well as lithium and all trans geometries of the alkyl chains in the generated structure trimethylammonium dodecyl sulfate micelles indicate no water a series of minimizations and MD simulations were performed. penetration into the aliphatic core, alkyl chains being disordered Vacuum calculations were performed using a truncation scheme to a lesser extent than in the liquid state and a broad distribution of 13.0, 12.0, and 10.0 A for the nonbond update list, the of the terminal methyl groups throughout the m i ~ e l l e In . ~ ~ ~ truncation ~~ of the electrostatic shift and VDW switch smoothing that study differences between the sodium, lithium, and tetrafunctions and the initiation of the VDW switch smoothing methylammonium dodecyl sulfate micelles are noted. Fourier function, re~pectively.~~ All lists were maintained on an atom transform infrared (FTIR) experiments have shown the number by atom basis. Following the fixing of the positions of the sulfur of defect states in the alkyl chains of SDS micelles and in atoms the micelle was minimized for 100 steepest-descent (SD) tridecane to be similar.29 Further clarification of the molecular steps with harmonic constraints of 5 kcal/(mol/8,) applied to interpretation of these experiments via theoretical studies is all atoms, followed by 100 adopted-basis Newton-Rapheson desireable. (ABNR) steps with harmonic constraints of 1 kcal/(mol/A) Presented below are results from a MD simulation of a applied to all atoms. The resulting structure was then subjected sodium dodecyl sulfate micelle in a periodic box of water. The to a MD simulation in the presence of harmonic constraints of micelle is comprised of 60 dodecyl sulfate monomers, in accord 100 kcal/(mol/8,) on the sulfur atoms. The initial portion of with the experimental aggregation number of 62,30.31and is the simulation involved gradual heating from 0 to 298 K over surrounded by 4398 water molecules and 60 sodium ions 2 ps using the Verlet algorithm with a 0.001 ps integration time rendering the system electrically neutral. Parameters used in step. SHAKE of all covelant bonds involving hydrogens was the c a l ~ u l a t i o nhave ~~,~ recently ~ been shown to yield excellent performed in this and all other MD simulation^.^^ This was agreement with experiment in MD simulations of a dipalmitoyl followed by a 5 ps equilibration of the system with scaling of phosphatidylcholine lipid bilayer,34a 1-palmitoyl-2-oleoyl-snvelocities every 100 steps if the temperature of the system had glycero-3-phosphatidylcholinelipid bilayer,35 and o-phosphodeviated more than f 5 K from 298 K. rylcholine and o-phosphosphorylethanolamine in aqueous soFollowing relaxation of the micelle in vacuum the system l ~ t i o n . For ~ ~ comparison of the properties of the micelle was immersed in a cube of water and 60 sodiums added to hydrocarbon interior results are presented from a simulation of render the system electrically neutral. Placing the micelle in a pure dodecane using periodic boundary conditions. Analysis waterbox involved overlaying the system with a 57.8 x 57.8 emphasizes comparison of the calculated results with available x 57.8 8, cube of water. The watercube was obtained from a experimental data, comparison of the micelle and dodecane alkyl previously optimized box of water38with image generation used chain structure and dynamics and details of the solvation of to obtain the desired size. Following the overlay any water the micelle and of the sodium ions. Emphasis will also be molecules whose oxygen was within 2.8 8, of a micelle heavy placed on comparison of the present results with those of atom were removed. With the micelle fixed the waters were Shelley, Watanabe and Kleing on a MD simulation of a sodium minimized for 100 SD steps followed by a 1 ps thermalization dodecyl sulfate micelle. This comparison allows for the at 298 K with velocity scaling every 100 steps if the temperature influence of the parameters used in the calculations on the of the system had deviated more than f 5 K from 298 K. For obtained results to be addressed. this and all remaining calculations a truncation scheme of 9.0, 8.5, and 7.5 8, for the nonbond update list, the truncation of Methods the electrostatic shift and VDW switch smoothing functions and All calculations were performed with version 22 of the the initiation of the VDW switch smoothing function, respecprogram CHARMM.37 Parameters were the CHARh4M22 allt i ~ e l y ,and ~ ~ a 0.002 ps integration time step were used. hydrogen lipid parameter^.^*,^^ The water model was the Following the water thermalization, water molecules within 12.0 CHARMM modified TIP3P m ~ d e l . ~ * - ~ ~ 8, of the geometric center of the micelle were deleted. This Micelle Simulation. Preparation of the micelle was initiated was done in agreement with experimental evidence indicating by building the micelle in vacuum, allowing the alkyl chains to the core of the micelle to be water free.22.23 equilibrate via vacuum MD simulations and overlaying the Placement of the sodiums was performed using a variation micelle with water and sodiums. The aggregation number of of a previously applied m e t h ~ d o l o g y . Following ~~ the therSDS micelles in 50 mM NaCl is approximately 62 and the malization of the waters, all water hydrogens were deleted. A critical micelle concentration is 0.008 M.30331To systematically single SDS molecule was then selected. In a sequential fashion, distribute the sulfate head groups over the surface of a sphere each water oxygen within 8.5 8, of the sulfate sulfur or oxygen the molecule C60 bu~kminsterfullerene~~ was selected as a atoms was individually replaced by a sodium ion and the total model. This offered a reproduceable lattice with approximately interaction energy of that sodium determined with the surroundthe correct number of sites corresponding to the SDS aggregaing environment, including the remaining water oxygens. Once tion number. The individual dodecyl sulfate molecules were the interaction energies for a sodium at each water oxygen built along vectors extending from the center of the fullerene position within 8.5 8, of the selected lipid head group was through the individual carbons of the C60 molecule. The determined the lowest interaction energy site was selected and dodecyl sulfate C12 methyl carbons were placed 3.0 8, from the original water oxygen was replaced by a sodium. With the the center of the fullerene and the remaining carbons were built first sodium in place a second SDS molecule was selected and in an outward fashion such that every other carbon, the thioester the procedure repeated. This cycle was iterated until 60 sodiums oxygen and one of the anionic oxygens were placed along the were positioned. Following placement of the sodiums the water vector in an all-trans configuration based on the geometries hydrogen atoms were replaced in their original positions. The included in the CHARMM22 all-hydrogen lipid parameter set. selected methodology ensures that a sodium is located within This procedpre led to the construction of an idealized micelle 8.5 8, of the sulfate head group of each SDS molecule. In with a 19.5 A radius based on the sulfur atoms. This radius is addition, as sodiums were added, previously placed sodium in good agreement with experiments which yield a paraffinic molecules were included in the interaction energy calculation. radius of approximately 16.7 8, and a total radius of apThe final system was comprised of 15 774 atoms, including 2520 proximately 22-23 A.26
MacKerell
1848 J. Phys. Chem., Vol. 99, No. 7, 1995 micelle atoms, 60 sodiums, and 4398 water molecules. Following placement of the sodiums the water molecule positions were allowed to relax by performing a minimization of 100 SD steps and 2 ps of MD thermalization in the presence of images with all micelle and sodium atoms fixed. This was followed by 2 ps of MD thermalization of the sodiums and waters in the presence of images with the micelle fixed. The final stage of the preparation was an "T simulation of the system in order to allow the total volume to relax to a value corresponding to approximately 1 ATM at 298 K. The NFT simulation43was performed for 6 ps using a temperature coupling constant of 0.1 and a pressure coupling constant of 0.1. During the NPT simulation the volume of the system contracted from an edgelength of 57.8 to a final value of 54.1 A. The production simulation was performed for 120 ps using constant volume and constant temperature via velocity scaling every 0.1 ps if the temperature was &5 K from 298 K. Analysis was performed on the final 100 ps of the production simulation using 0.1 ps time frames, unless noted. Statistical analysis of the results based on the 60 individual dodecyl sulfate molecules was performed as previously described.44 Calculations on the solvated micelle were performed on a CRAY-YMP at the National Cancer Institute Fredrick Biomedical Supercomputing Center. MD simulations of the full system, including periodic boundary conditions, required 4 CPU hours and 22MW of memory per 2 ps of simulation time. Dodecane Simulation. Preparation of dodecane for a periodic boundary simulation was performed as follows. A box of 125 dodecane molecules in the all trans configuration was built on a 5 x 5 x 5 lattice using an offset of 4 A in the y and z directions and 18 A in the x direction between the individual molecules. The system was minimized for 20 SD steps followed by 100 ABNR steps. For this and all other minimization and dynamics runs periodic boundary conditions, a 9.0,8.5, and 7.5 A for the nonbond update list, the truncation of the electrostatic shift and VDW switch smoothing functions, and the initiation of the VDW switch smoothing function, re~pectively,~~ and SHAKE4' of all covalent bonds involving hydrogens were used. To ensure converged populations of the trans versus gauche torsion angles a series of simulations were performed. A NVT simulation at 298 K using the temperature scaling procedure of B e r e n d ~ e nwas ~ ~ performed for 10 ps with an integration time step of 0.002 ps and a temperature coupling constant of 1.0. This was followed by a 50 ps "T simulation at 298 K and 1 ATM using a temperature coupling constant of 0.01 and a pressure coupling constant of 0.1. An NVT simulation was next performed by gradually heating the system to 1000 K over 5 ps with a temperature coupling constant of 1.0. At this stage of the preparation the density of the box of dodecane was lower than the experimental value of 0.749 g / m ~ l e .To ~ ~bring the system to the experimental density an additional NPT simulation was performed at 298 K with a reference pressure of 10,000 ATM, a temperature coupling constant of 0.01 and a pressure coupling constant of 10.0. This simulation was continued until the calculated volume of the system corresponded to the experimental density. The production simulation was performed for 200 ps using constant volume and constant temperature via velocity scaling every 0.1 ps if the temperature was &5 K from 298 K. Analysis was performed on the final 100 ps of the production simulation using 0.1 ps time frames. Values of dihedral transition rates have been shown to be dependent on the methodology used in their determination!6 In the present study a dihedral transition was defined as a change in a dihedral of between 60 and 150" within 0.3 ps. An additional check was made for recoils back to the original well
350
A I
100' , 750
3 700
d T
,
20
,
I
60
4b
*
,
, .
1
. . I
1
1
1
500',
d
.-I
80 100 120 140 160 180 200
B1
I
,
20
,
,
40
I
60
.
I
I
1
a j
60 100 120 140 160 180 200 Time, picoseconds
Figure 1. Dihedral energy versus time for the (A) micelle and (B)
dodecane simulations. and for overshoots beyond 150" within an additional 0.4 ps window. If such events did occur the previous transition was not included in the calculated transition rate. The presence of recoils and overshoots have previously been shown to not make significant contributions to correlation times for dihedral decays and, therefore, should not be included in the determination of transition rates.&
Results and Discussion Micelle structure may be defined based on the ensemble average of the conformations being sampled by the lipid molecules comprising the micelle. From MD simulations information on the ensemble average may be obtained from the time a ~ e r a g e . ~This ~ , ~equivalence * allows for the comparison of calculated and experimental data. Once agreement between the calculated and experimental macroscopic quantities is verified, further analysis of the calculated data allows for an understanding of the microscopic details of the structural and dynamic properties of the micelle. This information is invaluable as the majority of experimental observables are macroscopic in nature, such that microscopic interpretation of the data is generally model dependent. System Convergence. Due to the lack of a defined structure of a micelle efforts were made to build the initial structure in an unbiased fashion; C60 buckminsterfullerenea was selected as the model. The use of an all-trans conformations in both the micelle and dodecane initial structures required a series of NVT and NFT simulations prior to the production simulations to insure convergence in the trandgauche distributions of the aliphatic chain dihedral angles. Presented in Figure 1 are the dihedral angle energies as a function of time for the micelle and dodecane production simulations. In both cases there was an initial relaxation of the energy followed by stable values for the remainder of the simulations. The relaxation is complete by 20 ps in the micelle simulation while 100 ps were required for the dodecane simulation. Accordingly, the dodecane simulation was extended for 200 ps, and the final 100 ps used for analysis. The stability of the dihedral energy as a function of time indicates that the dihedral conformations had converged. Analysis of the total potential energies of both systems versus time (not shown) revealed a similar initial relaxation followed
Molecular Dynamics of a Micelle
J. Phys. Chem., Vol. 99,No. 7, 1995 1849
7
0.03-
6
0
u
-g
5
0.02-
4
0
3 0.01
2
0
0
1
0
5
10 15 20 Distance lo Micelle Center of Mass, A
25
30
Figure 2. Density distribution for the atom-to-micelle center-of-mass the sulfur atoms (O),the sodiums distance for the micelle carbons (O), (A), and waters (e).The sulfur and sodium densities have been multiplied by 5 for clarification of viewing.
0
5
10 15 20 Distance to Micelle Center 01 Mass, A
25
30
Figure 3. Probability distribution for the atom to micelle center-ofmass distance for the micelle C1 (O), C6 (O), and C12 (A)carbons, the sulfur atoms (*), and the sodiums (W). The probability is defined as the number of occurrences of the specified atom type in a 0.5 A window located the specified distance from the center of mass of the micelle per time frame.
by stable values through the remainder of the simulations. Average internal pressures from the production portions of the between 20 and 25 A. These density profiles are in good simulations were 669 f 196 and 704 f 273 ATM for the agreement with the experimental paraffinic radius of 16.7 A micelle and dodecane simulations, respectively, where the error and the total radius of 22.3 A from X-ray scattering.26 represents the root-mean-square (rms) fluctuations. In the Comparison of the carbon with the water and sodium density micelle simulation the radius of gyration of the micelle and the profiles show the interior of the micelle to be completely void average distance and standard deviation of the sulfur atoms to of solvent, again in good agreement with e ~ p e r i m e n t . More ~~,~~,~~ the micelle center of mass were monitored as a function of time detailed analysis of the micelle structure was obtained from the and observed to be stable following an initial relaxation during probability distribution of the C1, C6, C12, sulfur atoms, and the initial 20 ps of the simulation (not shown). One exception the sodium ions with respect to the micelle COM (Figure 3). was the rms distance of the sodium ions with respect to the The probability distributions of the carbon atoms in the micelle center of mass; this value was observed to increase hydrocarbon tails were observed to broaden as the distance from during the initial 50 ps of the production simulation, followed the head group increases. For the C1 atom the distribution is by a low frequency oscillation (not shown). This behavior similar to that of the head group sulfur atom. An increase in indicates that the sodium ions distribution has not totally the distribution occurs upon going to the C6 atom, with the equilibrated, possibly due to the initial placement of the sodium largest distribution occuring with the C12 methyl group. In ions within 8.5 8, of the sulfate head groups. addition, the C12 distribution has a tail in the region of the Micelle Structure. Generation of the starting structure of sulfur distribution which corresponds to a single lipid molecule. the micelle was based on the assumption of a spherical shape Due to the observation of this event in only a single lipid in accord with e~periment.4~ To check if the spherical geometry molecule statistical analysis is not possible, however, the was maintained during the simulation, the ratios of the principle occurrence of this event indicates that such a conformation is moments of inertia of the micelle heavy atoms were obtained feasible. These trends are similar to those calculated by Shelley every 0.2 ps. Average values and standard deviations of 1.02 et al., although the previous study showed the C6 probability f 0.02, 1.03 f 0.03, and 0.99 f 0.02 were obtained for the distribution to be more skewed and that of the C12 carbon did 11/13, 12/13, and 11/12 ratios, respectively, showing the geometry not extend out to a range similar to that of the sulfur atoms. of the micelle to remain spherical. From the micelle simulation Previous theoretical work based on a lattice model predicted a the radius of gyration was 16.02 iz 0.06 A, where the error maximum in the terminal methyl distribution near one-half of represents the root-mean-square (rms) fluctuation. The calcuthe radius of the micelle model and the possibility that chain lated value is slightly larger than the experimental value of 15.4 termini may lie near the micelle-solvent interface.' ExperiA for a lithium dodecyl sulfate micelle,?2 however, the agreemental work has shown the peak of the terminal methyl ment may be considered satisfactory. For the positions of the distribution to occur at 10.5 A and extend out to 25 A in lithium head groups, as designated by the sulfur atom, the average of dodecyl sulfate in good agreement with the present the rms distance from the center of mass was 19.70 f 0.08 A, calculated results. The present results are also consistant with which is in satisfactory agreement with an experimental value NMR relaxation studies;25order parameters associated with the of the 18.9 A.22 Care must be taken in comparing experimental C-H bonds were observed to decrease as the distance from results on micelles in the presence of different counterions; studies emphasizing both the differencesz7 and ~ i m i l a r i t i e s ~ ~ , ~the ~ polar head groups increased. This decrease may be associated with the increase in the width of the probability of various properties have been reported. Shelley et al? reported distributions upon going from the C1 to C6 to C12 atoms. The the average of the two largest ratios of the moments of inertia increased widths may indicate increased accessible conformafrom an SDS micelle MD simulation to be 1.13 and 1.07, tional space allowing for more extensive relaxation of the C-H indicating a slightly less spherical structure,13 and a diameter bonds. of 16.2 8, based on the sulfur atoms, which is smaller than the present results. Both differences are suggested to be attributable Structural aspects of the aliphatic chains were also analyzed to the smaller number of monomers, 42, in the previous study. in terms of the populations of the trans and gauche conformers. In the micelle there were minimal differences in the trans to The micelle structure was also analyzed in terms of the density of selected atoms as a function of distance from the gauche population as a function of dihedral; values of 0.58 f micelle center of mass (COM), as presented in Figure 2. From 0.03, 0.83 f 0.04, 0.69 f 0.05,0.73 f 0.05,0.82 f 0.04,0.79 the plot it may be seen that the carbon density rapidly decreases f 0.04, 0.80 f 0.04, 0.88 f 0.03, 0.81 f 0.04, and 0.91 & between 15 and 20 A, while the sulfur density rapidly decreases 0.03 for the trans population were obtained for the 0-C1-
1850 J. Phys. Chem., Vol. 99, No. 7, 1995
MacKerell
study a gauche-trans energy difference of 0.5 kcdmol was used in contrast to the present difference of 0.85 kcaymol, possibly contributing to the observed difference. Examination of the change in defects upon going from the micelle to dodecane reveals only minor differences. For the bend end conformation the calculations predict an increase in the neat alkane beyond the 2-fold increase expected due to the symmetry of the molecule, while differences in the double-gauche and kink states are negligable. Experiment shows all three states to be similar for the micelle and neat alkane. The similarities of the micellar and dodecane defect states are consistant with the results presented above showing the relative trans to gauche C2-C3, C1-C2-C3-C4...C9- 10-C11 -C 12 dihedrals, repopulations to be similar. spectively. Comparison of the dihedrals shows the differences Micelle Dynamics. Micellar dynamics may be analyzed with to be beyond the estimated error, however, no trends are respect to both the sulfate head groups and the aliphatic chains. observed. The largest difference was the larger gauche populaAnalysis of the motion of the sulfate head groups was performed tion for the dihedral (O-Cl-C2-C3) adjacent to the sulfate by calculating the distance of the individual lipid sulfur atoms head group. A similar result was obtained in a simulation of a to COM and the motion of the sulfur atoms parallel to the lysophosphatidylethanolaminemicelle.8 In the present calculasurface of the micelle as a function of time. Figure 4 contains tions the alteration of the relative trans-gauche population is plots for the sulfur-to-micelle COM distances and motion of due to a lowering of 0 1 -C1 -C2-C3 gauche energy relative the sulfur atoms parallel to the micelle surface for selected lipid to that of the trans conformation associated with the substitution molecules. For example, with lipid 18 (Figure 4a,D) the sulfur of an oxygen versus the normal alkane carbon atom. The overall to COM distance is initially 18 A at 20 ps. It then increased to trans and gauche population for the micelle, with the dihedral 24 8, at 80 ps and decreased back to 22 A at 120 ps. For the adjacent to the sulfate head group excluded, was 0.81 f 0.01 and 0.19 f 0.01. These values are both ~ i m i l d , and ~ , ~ ~ motion along the surface of the micelle lipid 18 moves almost 8 8, away from its initial position from 20 to 90 ps following dissimilar to previous micelle which is, in part, which it moves closer to its starting position. In several of the attributable to differences in the parametrization of the C-Clipids changes of up to 7 A in the sulfur to micelle COM C-C torsional surfaces. In the Shelley et al. study the gauche distance and 8 8, along the surface of the micelle relative to the ratio was in the range of 0.30, significantly larger than the starting position were observed. Presented in Figure 5 are time present value. This may be due to the smaller radius of the frames every 20 ps for the lipids included in Figure 4. micelle leading to an increased gauche population in order to Inspection of the figure shows two different types of motions. obtain the required packing in the hydrocarbon interior as well Lipids 8 and 60 are predominantly in all-trans conformations as differences in the parametrization. Analysis of the pure with a motion reminiscent of a piston in a cylinder. Such a dodecane simulation yielded an identical trans and gauche motion is dominated by translations along a vector between the population as in the micelle. Thus, the environment on the head group and the COM, with minimal motion along the interior of the micelle does not significantly alter relative gauche surface of the micelle. Lipids 5 and 18 have large motions along to trans populations as compared to pure dodecane. the surface of the micelle. Following the initial time frame, Additional analysis of the aliphatic chain conformations was the motion of lipid 5 is dominated by the head group, while performed via determination of defects in the chains. Such relatively small movements occur in the aliphatic tail. This is defects, or bent conformations, are experimentally accessible contrasted by lipid 18, where the motion of the head group along via IR by correlating different carbon-hydrogen bending modes the surface of the micelle is accompanied by movement of the with c ~ n f o r m a t i o n .In ~ ~such studies three conformations are aliphatic tail. In lipids 6 and 39 the movement contains a larger accessible based on characteristic carbon-hydrogen bending contribution from motions of the head group along the surface frequencies, including the bent-end methyl, double gauche and of the micelle, although motion with respect to the COM is kink (gauche-trans-gauche) conformations. Resented in Table significant. Inspection of lipids 8 and 60 show the aliphatic 1 are results from the micelle and dodecane simulations and chains to be in an extended conformation which may be related from e ~ p e r i m e n t . The ~ ~ experimental neat alkane work was to the minimal motion along the surface of the micelle, performed on tridecane rather than dodecane as used in the suggesting a possible relationship between the nature of the lipid present study, however, the additional methylene group will be motion and the hydrocarbon tail conformation. Further study assumed to make a negligible contribution. Comparison of the is required to clarify such a relationship. These results suggests experimental and calculated results reveals good agreement for that the motion of lipid molecules may be placed into two classes the bent-end methyl state, while the double-gauche and kink associated with motion along a vector from the micelle COM state calculated populations are significantly smaller than the to the sulfate head group and motion of the head group along experimental values. The magnitude of the difference and the the surface of the micelle. standard error reported for the calculated results indicates these Examination of the flexibility of the aliphatic chains was differences to be significant. Considering the good agreement performed via calculation of the energy surface from the transof the bent-end methyl populations, the discrepancies are gauche population analysis and direct determination of the difficult to explain. Shelley et al. reported bent-end methyl, aliphatic dihedral transition rates. Analysis of the free energy double-gauche, and kink populations of 0.34, 0.21, and 0.84. The increased kink population as compared to the present study surfaces in Figure 6, calculated from histograms of the dihedral is consistent with the increased population of gauche dihedrals populations and assuming a Boltzmann distrib~tion:~shows and is in better agreement with experiment, however, the differences between the micelle and dodecane simulations. In population of double-gauche state is similar to the present study. accord with the identical trans and gauche populations (see Interpretation of the experimental data was based on the above) the energy differences between the two conformers are rotational isomeric state model of Fl01-y.~~ In the experimental identical for the two systems. In contrast, the banier to rotation
TABLE 1: Populations of Aliphatic Chain Defect States per Molecule for the Micelle and Dodecane Simulations structure bent-end methyl double-gauche kink micelle calculated 0.27 f 0.04 0.23 f 0.06 0.23 i0.04 experimental" 0.36 0.77 0.68 dodecane 0.18 f 0.02 0.24 f 0.02 calculated 0.71 f 0.03 experimentalb 0.68 0.64 0.77 Experimenal results for tridecane. a Experimental data from ref29.
Molecular Dynamics of a Micelle
-22
J. Phys. Chem., Vol. 99, No. 7, 1995 1851
-
A
181
4
/
,
,
,
,
.
.
,
,
.
I
42
0 8
C 4 1
4
,
14,
,
,
,
,
,
,
,
,
,
.
,
,
,
,
,
C
2
D
H
0 ~~
22
1
-
8
F
F 4-
4
20
1
,
,
,
,
40
60
.
.
, 80
, 100
.
1
0
120
Time, picoseconds
,
'
l
'
I
.
,
'
Time, picoseconds
Figure 4. (a, left) Sulfur to micelle center-of-mass distance versus time and (b, right) sulfur motion parallel to the surface of the micelle surface for lipids (A) 5 , (B) 6, (C) 8, (D)18, (E) 39, and (F) 60. Sulfur motion parallel to the surface of the micelle is defined as the motion of the sulfur atom perpendicular to the sulfur to the micelle center-of-mass vector obtained at the 20 ps time frame.
5
5
60
5
60
Figure 5. Stereodiagram of times frames at 20 ps intervals for lipids (a, top) 5 , 8, and 39 and (b, bottom) 6, 18, and 60. Oxygen atoms are represented by dashes and hydrogens have been removed for clarity. Line thickness increases as a function of time with the thin line representing the 20 ps time frame. The left and center diagrams represent the cross-eyed stereopair and the center and right diagrams represent the wall-eyed stereopair. is significantly larger in the micelle. Shown in Figure 7 are the transition rates as a function of the carbon-carbon bonds in the SDS hydrocarbon chain and in dodecane. In agreement with the energy surface, the transition rates are significantly larger in dodecane than in the micelle; error bars calculated for the transition rates show the differences to be statistically significant. The increase in the transition rates occur throughout the aliphatic chain, being approximately 3-fold larger in dodecane at the termini and approximately 2-fold larger in the center of the chain. Calculation of the transition rate from the energy surfaces presented in Figure 6, assuming an Arrhenius relationship with the preexponential factor identical for the two
systems, yields a transition rate 2.4 times slower in the micelle as compared to dodecane. Although the trans-gauche populations are identical the hydrocarbon interior of the micelle is predicted to be less fluid than that of pure dodecane. The decreased fluidity is due to an increased free energy barrier to rotation of the aliphatic chain dihedral angles. The decreased fluidity of the micelle aliphatic chains contrasts calculations on a lipid bilayer using the same parameter set.34 In the dipalmitoyl phosphatidylcholine bilayer the dihedral transition rates for the hydrocarbon tails were similar to those in neat hexadecane. The transition rates for dodecane in the present study are similar to those reported for hexadecane by
MacKerell
1852 J. Phys. Chem., Vol. 99, No. 7, 1995
6 1 5
E $ 4
P $3
t
i2 1
0 0
60
120
180 240 Dihedral Angie. degrees
300
Carbon Numbei
360
Figure 6. Free energy versus the hydrocarbon C-C-C-C dihedral angle from the micelle (U) and dodecane ( 0 )simulations. The free energy surfaces were determined from histograms of the dihedral populations assuming a Boltzmann distribution.
Figure 8. Hydration number of the carbon atoms in the hydrocarbon tails. Hydration numbers were determined based on all carbon to water oxygen interactions less than 3.5 A.
above in Figures 2 and 3, respectively, and several radial distribution functions presented below. The water density in 30 Figure 2 is zero within 12 8, of the micelle COM and gradually increases to 0.034 at 23 8,. This is followed by a broad plateau out to 30 A, during which a gradual decrease in the density to 0.033 k3occurred. The present density is in good agreement with the expected water density of 0.0334 k3. Although by no means conclusive, the gradual decrease in the water density from 23 to 30 8, suggests an influence of the solute micelle on the adjacent bulk solvent properties. Similar observations have been previously reported in macromolecular s i m ~ l a t i o n s .The ~ ~ head . ~ ~ groups, as defined by the sulfur atoms (Figures 2 and 3), are distributed between 14 and 26 8,. 01 Comparison of the head group and water density distributions 1 2 3 4 5 6 7 8 9 10 11 Dihedral Number show both to go to zero approximately 14 8, from the micelle COM, suggesting no significant penetration of water into the Figure 7. Dihedral transition rates for the C-C bonds from the micelle (U) and dodecane (0)simulations. Error bars represent the standard hydrocarbon interior. From the densities in Figure 2 the length error.44On the x axis 2 represents the Cl-C2-C3-C4 dihedral and of the interface region may be calculated based on the location so on. where the hydrocarbon and water densities pass through 10% of their bulk densities.13 This approach was used to determine Venable et al.,34 while the transition rates in the micelle are an interface length of 4.0 8,13 in the micelle simulation b significantly lower. Previous theoretical studies' have also Shelley et aL9 In the present study an interface length of 4.5 calculated a decrease in the mobility of the hydrocarbon interior of the micelle (see Introduction). Experimental data from N M R between 15 and 19.5 8, was calculated. This suggests that the interaction of water molecules with the micelle is dominated relaxation experiments lends further support to the decreased by interactions with the head group. In the Shelley et al. study fluidity of the micelle hydrocarbon interior. For a SDS micelle additional water penetration below 10% of the bulk density correlation times of 34, 44, and 7 ps for the a,/3, and w C-H occurred, moving approximately 6 8, further into the hydrocarbond motions have been measured.25 These values may be compared to correlation times of 41, 32, 20, 10, and 7 ps for bon core. This difference may be attributed to the less spherical the C2, C3, C4-13, C14, and C15 methylene groups, respeccharacter of the micelle in that study.13 Analysis of the sodium tively, of a dipalmitoyl phosphatidylcholine lipid b i l a ~ e ras~ ~ . ~ to ~ micelle COM distribution (Figure 2) shows the inner penetration of the ions into the micelle to be less than that of reported by Venable et al.34 Comparison shows the ,b C-H water, emphasizing that interactions with water are required for correlation times in the micelle to be longer than the C3 and sodium ions to interact with the micelle surface (see below). C4-13 times in the lipid bilayer. The terminal methyl group correlation times are identical for the two systems while that Experimental studies using fluorinated hydrocarbons suggest that interactions between water and the aliphatic tails are of the methylene group adjacent to the head group (C2) is shorter observed only at the C1 although it has been in the lipid bilayer. Influence of the head group on the motions of the C2(a) methylene groups limits interpretation of those suggested that the terminal methyl group may also be exposed.24 correlation times, however, the increased /3 methylene correlaFrom Figures 2 and 3 it may be seen that overlap between the water distributions and the carbons 1, 6, and 12 does occur. tion time supports the decreased fluidity of the micelle interior The overlap is greatest for the C1 carbon, followed by the C12, as compared to the dipalmitoyl phosphatidylcholinelipid bilayer with the C6 carbon having the least overlap. To closer analyze as well as the pure alkane. Preliminary analysis of the micelle trajectory in order to obtain C-H correlation times was hindered interactions between the aliphatic chains and water, the hydration by difficulties in obtaining satisfactory fits of the relaxation data. numbers of the individual carbons were determined, as shown Work is in progress on the development of a model to use the in Figure 8. The hydration number decreases rapidly upon going micelle trajectories (MacKerell, Jr., A. D.; Pastor, R. W., work from C1 to C3 with values close to zero for C4 to C11. A small increase occurs with the C12 methyl group, which is in progress) for the interpretation of the available N M R relaxation data.25 consistant with the extension of the C12 probability distribution out to 25 A in Figure 3. These results are similar to those Micelle-Solvent Interactions. Analysis of the interaction reported by Shelley et al.; however, the increased hydration of of the micelle with the surrounding solvent was performed using the density profiles and the probability distribution presented the terminal methyl group was not observed. Thus, the present
4
I
1
I
!
1
Molecular Dynamics of a Micelle
J. Phys. Chem., Vol. 99, No. 7, 1995 1853
4.5
1
n
4
3.5 3
-m
2.5
I
1.5 1
0.5
0
2
5 Distance,
6
i
i
A
Figure 9. Radial distribution functions between the sulfate head group sulfur and water oxygen (W), sodium ions (O),and sulfate sulfur atoms (A).Rdfs were calculated from time frames obtained from the micelle simulation every 0.5 ps for the sulfur-water rdfs and 0.2 ps for the sodium and sulfur to sulfur rdfs. Rdf s are normalized based on the density of 0.0334 A-3 for water and the average density of the sodium and sulfur atoms in the simulation system. TABLE 2: Hydration Numbers Determined from the Presented Radial Distribution Function integrated hydration fnst peak to no. position interacting pair 14.82 3.85 sulfur to water 5.25 0.71 5.35 sulfur to sodium 6.65 0.65 5.55-5.65 sulfur to sulfur 6.85 6.10 2.25 sodium to water 3.05 0.10 3.45-3.55 4.65 sodium to sodium 0.32 4.45 sodium to anionic oxygen 5.35 2.78 2.75 anionic oxygen to water 3.25 17.21 4.85 ester oxygen to water 5.85 Distances in angstroms. results predict the interactions between the hydrocarbon tail and water to be primarily associated with the C1, C2, and C3 atoms, but that some exposure of the C12 methyl group may occur. These results support previous theoretical descriptions of micellar s t r u c t ~ r ewhere ~ ~ ~ ~penetration of water into the hydrophobic core does not occur, in agreement with experimental r e s ~ l t s . ~Results ~ ~ * ~from a lattice-based model suggest that interactions of the terminal methyl groups with water may 0ccur.l It should be noted that the removal of waters from within 12 8, of the micelle geometric center during preparation of the micelle and the limited duration of the simulation make these conclusion preliminary in nature. Solute-Solvent Interactions. Study of the solvation of the sulfate head groups was performed via sulfur to water, sulfur to sodium, and sulfur to sulfur radial distribution functions (rdf) as shown in Figure 9 and the associated hydration numbers and first peak positions presented in Table 2. The sulfur to water oxygen rdf shows a peak at 3.85 8, and a coordination number of 14.8. The sulfur to sodium and sulfur to sulfur rdfs have overlapping peaks in the region 5-6 8, with coordination numbers of 0.7 and 0.7, respectively. In the startin micelle structure the first sulfur to sulfur peak occurred at 8.65 During the simulation some of the head groups moved closer to each other. Visual inspection of the final time frame of the simulation indicated the 5.75 8, peak in the rdf to be associated with solvent-separated head groups. In the solvent-separated peaks intervening waters would typically assume an orientation where one of the water protons would hydrogen bond to both of the adjacent sulfate groups or where the two water protons individually hydrogen bonded to adjacent sulfate head groups. These results portray a model where the first hydration shell of the sulfate head groups is almost exclusively comprised of water
1.
,
0 2
3
4
5 Distance,
I
6
,
7
8
A
Figure 10. Radial distribution functions between water and the anionic (0)and ester (W) sulfate oxygens. Rdf's were calculated from time frames obtained every 0.5 ps from the micelle simulation and are normalized based on a water density of 0.0334 A-3, molecules. Close inspection of the sulfur to sodium rdf indicates occasional penetration of sodiums into the first hydration shell; however, this contribution is negligable. The second hydration shell is again dominated by water molecules with contributions from both the sodium ions and the adjacent sulfate head groups. Of note is the lack of a significant number of sodiums in the first hydration shell of the head groups. This result contrasts previous empirical force-field studies on sodium octoanate4,7~10~22,23 and SDS9 micelles where such interactions were observed. In a study on a sodium octoanate micelle1° it was reported that the inclusion of induced polarization in the force field eliminated the presence of contact ion pairs. In the present study there is no polarization term present, however, no stable sulfur to sodium contact pairs were observed (Figure 9), although low probablilty occurrences of contact pairs are evidenced by the nonzero behavior in the region 3.4-4.2 8, in the sulfur to sodium rdf. These low-probability events included an interaction less than 4.2 8, that lasted for 3 ps at approximately 30 ps of the production simulation and again at 108 ps; the second close interaction lasted for approximate1 12 ps; however, it included three excursions beyond 4.2 during that period. These results suggest that induced polarization in the potential energy function is not required to avoid stable ion contact pairs. In the optimization of the CHARMM22 all-hydrogen p a r a m e t e r ~ ~emphasis ~ 2 ~ was placed on the proper balance between solute-solvent, solute-solute and solventsolvent portions of the intermolecular force field. The present results suggest that proper balance of the intermolecular parameters yields an adequate treatment of the ion-solvent interactions, although additional improvements upon inclusion of induced polarization are expected. More detailed analysis of the sulfate head group to water interactions was performed via the sulfate anionic oxygen to water and sulfate ester oxygen to water interactions (Figure 10). For the anionic oxygens the first hydration peak occurs at 2.75 8, with a coodination number of 2.8. This peak correspond to the sulfur-water peak at 3.85 8, in Figure 9. In contrast, no distinct peak was obtained between the ester oxygen and water, indicating the lack of a distinct first hydration shell at that site. The large hydration number associated with the ester oxygen is due to the position of the first peak such that all waters solvating the individual head groups are included. These results suggest that the hydration of the sulfate head groups is dominated by interactions with the anionic oxygens while the exclusion of water from the micelle interior is already occuring at the head group ester oxygen. Presented in Figure 11 are rdfs between sodium and water, the sulfate anionic oxygens, and sodium. The first solvation
f
1854 J. Phys. Chem., Vol. 99, No. 7, 1995
MacKerell
:iIt
"il I 2
3
4
5
6
i
8
Distance, A
Figure 11. Radial distribution functions between the sodium ions and water (W), sodiums (0).and the anionic sulfate oxygens (A).Rdf s were calculated from time frames obtained every 0.5 ps from the micelle simulation for the sodium-water rdf and every 0.2 ps for the sodiumsodium and sodium-anionic sulfate rdfs and are normalized based on the density of 0.0334 A-3 for water and the average density of the sodium and sulfate oxygen atoms in the simulation system.
shell of sodium is dominated by water, in accord with the sodium-water peak at 2.25 8, and a coordination number of 6.1 (Table 2). The sodium-water rdf as well as the calculated coordination number are in good agreement with previous ~ o r k . ~ * -In ~ Oa few instances sodium-sodium contact pairs occurred, as evidenced by the peak at 3.6 8, in the sodium to sodium rdf. From the sodium to sodium coordination number of 0.1 (Table 2) such occurences were infrequent. Sodium to sodium contact peaks have been previously observed in MD simulations of micelles and DNA?,60The second hydration shell of sodium is comprised of either water or anionic sulfate oxygens as indicated by the presence of sodium-water and sodium-anionic oxygen peaks in the vicinity of 4.3 8,. The coordination number of sodium by the anionic oxygens is 0.3 when integrated to 5.35 8, (Table 2), indicating that even secondshell interactions between sodiums and the sulfates occur for only a minority of the sodium ions. This is supported by the probablity distributions of the sodiums from the micelle COM (Figure 3) showing some sodiums to be located 10 or more A from the head group sulfur atoms. A second sodium to sodium peak in the rdf occurs at 6.6 A, corresponding to the solventseparated pair. Thus, in the present calculations the sodium ions are coordinated primarily by water molecules, while interactions with the micelle via the sulfate head group are solvent mediated. The present analysis of the solvent-solute interactions must be tempered by the truncation scheme used in the calculations. Ideally, MD simulations should be performed under conditions similar used in the parametrization. The use of a shift truncation scheme37 for the electrostatic interactions alters both the energetics and shape of the interaction energy surface as well as excluding interactions beyond the truncation distance. While such a truncation approach has been used extensively in MD simulations and shown to yield reasonable results,6l its influence on the present results cannot be excluded. In addition and as discussed above, omission of induced polarizability in the potential energy function may lead to additional complications. Further studies on model systems will aid in the clarification of these points.
Conclusion The presented calculations portray a picture of a sodium dodecyl sulfate in agreement with a variety of experimental results. On a microscopic level the micelle is spherical. However, significant mobility is observed in the individual lipids
with changes in the head group to micelle COM distances of up to 7 8, and movements of up to 8 A of the head groups along the surface of the micelle occumng during the simulation. The interior of the micelle is void of solvent molecules, however, the present results predict that the terminal methyl groups may, in certain instances, occupy a location near the surface of the micelle where they can be exposed to solvent. Conformations of the aliphatic chains are similar to that of neat dodecane, however, a decrease in the fluidity of the chains in the micelle is predicted. This results contrasts previous calculations using the same parameter set on a lipid bilayer34 where the chain fluidity was similar to that of a neat alkane. Considering the role of SDS and other micellar systems as model systems for the experimental study of membranes and for peptidemembrane interactions this finding warrants further study. Results from the simulation were analyzed with respect to the interactions of water with the micelle and the sodiums ions. The sulfate head groups are almost exclusively solvated by water with contributions from sodium ions and adjacent sulfate groups in the second hydration layer. First shell hydration of the sodiums is dominated by water while the second hydration shell includes interaction with water and micelle sulfate head groups. Previous work on a sodium octoanate micelle reported that induced polarization in the potential energy function eliminated the presence of ion contact pairs.1° Stable ion contact pairs were not observed in the present study, in which induced polarization was not included in the potential energy function. This is in contrast to the SDS micelle simulation by Shelley et aL9 but in agreement with the intepretation of results from N M R experiment^.^^,^^ While the inclusion of polarization will undoubtly lead to a more exact treatment of the intermolecular energetics of the system, the present results emphasize the role of the parameters in empirical force-field studies of condensedphase systems. In particular, it is essential that the delicate balance between the solvent-solvent, solvent-solute, and solute-solute interactions be properly maintained. The current CHARMM22 all-hydrogen parameter set was designed to be used with the modified TIP3P water potential. As such, the partial charges and Lennard-Jones parameters in the force field were optimized with respect to that solvent model. Substitution of alternate water models with the present parameters should be attempted with caution.
Acknowledgment. Financial support from the University of Maryland at Baltimore, School of Pharmacy, in the form of a Biomedical Research Support Grant and computer time provided by the NCI's Fredrick Biomedical Supercomputing Center are acknowledged. Thanks to Drs. Ronald D. Guiles, Richard W. Pastor, and Richard M. Venable for helpful discussions. References and Notes (1) Dill, K. A.; Flory, P . J . Proc. NatZ. Acad. Sci. U.S.A. 1981,78, 616. (2) Dill, K. A. J. Phys. Chem. 1982,86, 1498. (3) Dill, K. A.; Koppel, D. E.; Cantor, R. S.;Dill, J. D.; Bendedouch, D.; Chen, S.-H. Nature 1984,309,42. (4) Watanabe, K.;Ferrario, M.; Klein, M. L. J. Phys. Chem. 1988,92, 819. (5) Haile, J. M.; O'Connell, J. P. J. Phys. Chem. 1984,88, 6363. (6) Woods, M. C.; Haile, J. M.; O'Connell, J. P. J. Phys. Chem. 1986, 90, 1875. (7)Jonsson, B.;Edholm , 0.;Teleman, 0. J. Chem. Phys. 1986,85, 2259. ( 8 ) Wendoloski, J. J.; Kimatian, S. J.; Schutt, C. E.; Salemme, F. R. Science 1989,243,636. (9) Shelley, J.; Watanabe, K.; Klein, M. L. Int. J. Quantum Chem., Quantum Bid. Symp. 1990,17, 103. (10) Shelley, J. C.; Sprik, M.; Klein, M. L. Lmgmuir 1993,9, 916.
J. Phys. Chem., Vol. 99, No. 7, 1995 1855
Molecular Dynamics of a Micelle (11) Bocker, J.; Brickmann, J.; Bopp, P. J. J. Phys. Chem. 1994, 98, 712. (12) Laaksonen, L.; Rosenholm, J. B. Chem. Phys. Let. 1993,216,429. (13) Shelley, J. C. Ph.D. Thesis, University of Pennsylvania, 1992. (14) Braun, W.; Bosch, C.; Brown, L.; Go, N.; Wutrich, K. Biochim. Biophys. Acta 1981, 667, 377. (15) Kohda, D.; Inagaki, F., 1992,31,677Biochemistry 1992,31, 677. (16) Shon, K.; Kim, Y.; Colnago, L.; Opella, S. J. S. Science 1991, 252, 1305. (17) Yang, J. T.; Bewley, T. A,; Chen, G. C.; Li, C. H. Proc. Natl. Acad. Sci. U.S.A. 1977, 74, 3235. (18) Arseniev, A. S.; Barsukov, I. L.; Bystrov, V. F.; Lomize, A. L.; Ovchinnikov, Y. A. FEBSLett. 1985, 186, 168. (19) Zetta, L.; DeMarco, A.; Becchio, G.; Mazzola, G.; Longhi, R.; Motta, A.; Pastore, A,; Goud, N. A,; Castiglione Morelli, M. A. Biochemistry 1991, 30, 10444. (20) Yan, C.; Sarma, S.; MacKerell, A. D., Jr.; Guiles, R. G., Manuscript in preparation. (21) Moroder, L.; Romano, R.; Guba, W.; Mierke, D. F.; Kessler, H.; Delporte, C.; Winand, J.; Christophe, J. Biochemistry 1993, 32, 13551. (22) Bendedouch, D.; Chen, S.-H.; Koehler, W. C. J. Phys. Chem. 1983, 87, 153. (23) Ulmius, J.; Lindman, B. J. Phys. Chem. 1981, 85, 4131. (24) Mukejee, P.; Yang, A. Y. S. J. Phys. Chem. 1976, 80, 1388. (25) Soderman, 0.; Carlstrom, G.; Olsson, U.; Wong, T. C. J. Chem. SOC.,Faraday Trans. 1 1988, 84, 4475. (26) Itri, R.; Amaral, L. Q. J. Phys. Chem. 1991, 95, 423. (27) Jones, R. R. M.; Maldanado, R.; Szajdzinska-Pietek, E.; Kevan, L. J. Phys. Chem. 1986, 90, 1126. (28) McManus, H. J. D.; Kang, Y. S.; Kevan, L. J. Phys. Chem. 1992, 96, 5622. (29) Holler, F.; Callis, J. B. J. Phys. Chem. 1989, 93, 2053. (30) Croonen, Y.; Geladt, E.; Van der Ziegel, M.; Van der Auweraer, H.; Vandendriessche, F. C.; DeSchryver, F. C. J. Phys. Chem. 1983, 87, 1426. (31) Attwood, D.; Florence, A. T. In Surfactant Systems; Chapman and Hall: New York, 1983; p 71. (32) Schlenkrich, M.; Brickmann, J.; MacKerell Jr., A. D.; Karplus, M., manuscript in preparation. (33) Roux, B. Ph.D. Thesis, Harvard University, 1990. (34) Venable. R. M.: Zhang, Y.: Hardy. B. J.; Pastor, R. W. Science 1993,262, 223. (35) Heller, H.; Schaefer, M.; Schulten, K. J. Phys. Chem. 1993, 97, 8343. (36) Woolf, T. B.; Roux, B. J. Am. Chem. SOC. 1994,116, 5916. (37) Brooks, B. R.; Bruccoleri, R. E.; Olafson, B. D.; States, D. J.; Swaminathan, S.; Karplus, M. J. Comput. Chem. 1983, 4, 187. I
(38) Jorgensen, W. L.; Chandrasekhar, J.; Madura, I. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926. (39) Reiher ID,W. E. Ph.D. Thesis, Harvard University, 1985. (40) Kratschmer, W.; Lamb, L. D.; Fostiropoulos, K.; Huffman, D. R. Nature 1990, 347, 354. (41) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. J. Comput. Phys. 1977, 23, 327. (42) Jung, S.-H. Ph.D. Thesis, Harvard University, 1989. (43) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; D ~ o l a , A.; Haak, J. R. J. Chem. Phys. 1984, 81, 3684. (44) Loncharich, R. J.; Brooks, B. R.; Pastor, R. W. Biopolymers 1992, 32, 523. (45) Handbook of Chemistry and Physics; CRC Press: Boca Raton, FL, 1988. (46) Zhang, Y.; Pastor, R. W. Mol. Sim. 1994, 13, 25. (47) McQuarrie, D. A. Statistical Mechanics; Harper and Row: New York, 1976. (48) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford University Press: New York, 1989. (49) Hayashi, S . ; Ikeda, S. J. Phys. Chem. 1980, 84, 744. (50) Stilbs, P.; SMerman, 0.;Walderhaug, H. J. Magn. Reson. 1986, 69, 411. (51) Marconcelli, M.; Qi, S. P.; Straws, H. L.; Snyder, R. G. J. Am. Chem. SOC.1982, 104, 6237. (52) Flow.P. J. Statistical Mechanics of Chain Molecules;. Wilev: - New York, '1969.. (53) Brown. M. F.: Ribeiro. A. A.: Williams. G. D. Proc. Natl. Acad. Sci.~U.S.A. 1983, 80, 4325. (54) Brown, M. F. J. Chem. Phys. 1984, 80, 2808. (55) Wong, C. F.; McCammon, J. A. J. Am. Chem. SOC.1986, 108, 3830. (56) Damodaran, K. V.; Merz, K. M., Jr.; Gaber, B. P. Biochemistry 1992, 31, 7656. (57) Menger, F. M. N. Nature 1985, 313, 603. (58) Impy, R.; Madden, P. A,; McDonald, I. R. J. Phys.Chem. 1983, 87, 5071. (59) Smith,P. E.; Pettitt, B. M. J. Am. Chem. SOC.1991, 95, 8430. (60) Mohan, V.; Smith, P. E.; Pettitt, B. M. J. Phys. Chem. 1993, 97, 12984. (61) Loncharich, R. J.; Brooks, B. R. PROTEI": Struct. Funct. Genet. 1989, 6, 32. (62) Wennerstrom, H.; Lindman, B. Rev. Sect. Phys. Lett. 1979, 52, 1. (63) Lindman, B.; Wennerstrom, H. Topics in Current Chemistry; Springer-Verlag: New York, 1980; Vol. 87. JP9420887
