Molecular Dynamics Simulations of Cardiolipin Bilayers - American

Aug 20, 2008 - chrome c oxidase with two different P-P distances and conformations are shown in Figure 7. The choice of angles and distances above ena...
0 downloads 0 Views 609KB Size
Download PDF
J. Phys. Chem. B 2008, 112, 11655–11663

11655

Molecular Dynamics Simulations of Cardiolipin Bilayers Martin Dahlberg* and Arnold Maliniak DiVision of Physical Chemistry, Arrhenius Laboratory, Stockholm UniVersity, S-106 91 Stockholm, Sweden ReceiVed: April 19, 2008; ReVised Manuscript ReceiVed: June 19, 2008

Cardiolipin is a key lipid component in the inner mitochondrial membrane, where the lipid is involved in energy production, cristae structure, and mechanisms in the apoptotic pathway. In this article we used molecular dynamics computer simulations to investigate cardiolipin and its effect on the structure of lipid bilayers. Three cardiolipin/POPC bilayers with different lipid compositions were simulated: 100, 9.2, and 0% cardiolipin. We found strong association of sodium counterions to the carbonyl groups of both lipid types, leaving in the case of 9.2% cardiolipin virtually no ions in the aqueous compartment. Although binding occurred primarily at the carbonyl position, there was a preference to bind to the carbonyl groups of cardiolipin. Ion binding and the small headgroup of cardiolipin gave a strong ordering of the hydrocarbon chains. We found significant effects in the water dipole orientation and water dipole potential which can compensate for the electrostatic repulsion that otherwise should force charged lipids apart. Several parameters relevant for the molecular structure of cardiolipin were calculated and compared with results from analyses of coarse-grained simulations and available X-ray structural data. Introduction The lipid bilayer is the basic structural component of biological membranes. Membrane protein function is often strongly affected by the chemical and physical properties of the bilayer, which in turn depend on the composition of lipids and the complex molecular interactions in lipid systems.1 Simplified models of lipid bilayers have become useful tools for understanding the fundamental properties of biological membranes.2 The inner mitochondrial membrane (IMM) is an example of a highly folded and dynamic biological structure.3 The topology and membrane properties are strongly linked to mitochondrial function in the energy metabolism of the cell,4,5 and even though proteins have been found to control the folding morphology of the IMM cristae,6 the lipid components are central for energy homeostasis, apoptosis, and transport processes across the IMM.7,8 The main lipid components of the IMM in humans are cardiolipin (CL), phosphatidylcholine (PC) (Figure 1), and phosphatidylethanolamine (PE) at the approximate ratio of 1:2:2 by weight.9 Cardiolipin is a phospholipid with two phosphate groups and four acyl chains of varying length and saturation. The CL content outside of the IMM is very low. Due to its role in the mitochondrial metabolism and structure, CL has become the focus of previous research efforts.7,10-13 The lipid cocrystallizes with several IMM proteins, and the electrostatic interactions between proteins and the charged CL headgroup are thought to be important, for example, for cytochrome c localization at the outer face of the IMM.14,15 Several experimental studies of IMM models with unilamellar vesicles or supported planar bilayers have been reported.16-22 Also, an attempt to describe the IMM with theoretical models of PC and PE has been reported recently.23 Obtaining a detailed molecular understanding of this system has thus become of considerable interest. Other lipid mixtures have been examined with computer simulations on a variety of levels of approxima* To whom correspondence [email protected].

should

be

addressed.

E-mail:

tion,24-31 but the number of reports on mixtures of charged and zwitterionic lipids relevant for the IMM is still limited, as are simulations of CL.32,33 In this article we investigate the molecular details of CL and its effect on bilayer properties using molecular dynamics (MD) computer simulations. Bilayers with three compositions, pure CL, pure PC, and a mixture of CL and PC that approximates the composition found in vivo (∼20 and 9% CL by weight and mole, respectively) in human mitochondria, were modeled at the atomistic level. Many species of CL with respect to the acyl chain composition have been observed in vivo, ranging from short saturated chains to long polyunsaturated chains.13,34 For the present work, tetraoleoyl-CL was chosen because of its low melting transition temperature and its relatively high degree of unsaturation, which is consistent with the human IMM CL species. In line with previous coarse-grained (CG) simulations,32 we chose a doubly deprotonated CL headgroup. These lipids, with two negative charges on the phosphate groups, were found to form bilayers rather than the inverted hexagonal phase (HII). The actual charge state in vivo is not well established, although there is evidence that the second pKa value is unusually high (between 7.5 and 9.5) in CL.35 The connection between protonation and intramolecular geometry in CL is also not determined but is most likely dependent on headgroup conformation. Thus, as a basis for further studies, we compared recently available structures of CL from structures deposited in the RSCB Protein Data Bank (PDB) with those obtained from the atomistic simulations. In addition we compared the results from our simulations with experimental and theoretical data obtained in studies of similar model systems, such as the previous CG simulations of CL. Methods The lipid interaction model was based on the force field latest described by Berger et al.36-39 The POPC parameters and coordinates were obtained from previous work by Zhao et al.40 The R/R/R stereoisomer of tetraoleoylcardiolipin (CL) (Figure 1) was constructed from a phosphatidylglycerol (POPG) model41

10.1021/jp803414g CCC: $40.75  2008 American Chemical Society Published on Web 08/20/2008

11656 J. Phys. Chem. B, Vol. 112, No. 37, 2008

Dahlberg and Maliniak

Figure 1. 1,1′,2,2′-Tetraacylcardiolipin (CL) and 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC). The asterisk (*) marks the central glycerol carbon atom.

TABLE 1: Details of the Simulated Bilayersa simulation

composition

water per lipid

simulated time [ns]

0-CL 9-CL 100-CL

118 POPC 118 POPC + 12 CLb 64 CL

45.7 45.7 45.7

145 300 (500) 300

area per lipid [nm2] 0.64 ( 0.01 0.62 ( 0.01 0.99 ( 0.01

KA [mN/m]

thickness [nm]

350 340 1100

3.76 ( 0.05 4.1 ( 0.1 4.72 ( 0.04

a Compositions and simulation time of the POPC (0-CL), POPC+CL (9-CL), and CL (100-CL) bilayers. Area per lipid, compressibility modulus (KA), and bilayer thickness calculated as the distance between phosphate groups in adjacent leaflets of the bilayer. b 9.2% CL by mole.

with the following modifications: (i) the palmitoyl chain was exchanged for an oleoyl chain with the same parameters as the oleoyl in POPC, (ii) a small repulsive potential, 10-7 kJ mol-1 nm12 in the C12 term for the r-12 interaction, was added to the hydroxyl group hydrogen to avoid excessive intramolecular hydrogen bonding to the phosphoester oxygen atoms, in analogy with a previous modification of the DOPE model,42 (iii) standard 1,4-Lennard-Jones parameters based on the ffgmx force field were used for 1,4 interactions with the hydroxyl oxygen, and (iv) the connecting glycerol moiety was made symmetric with respect to charges and atom types. After energy minimization of a single CL lipid, 16 CL were assembled into a small bilayer patch. Subsequently, about 46 water molecules and two sodium ions were added per lipid. The amount of water was then comparable to the previously observed hydration levels in experiments on tetramyristoyl-CL.43 The system was energy minimized and equilibrated in a short MD run. A larger bilayer with 64 CL (32 per leaflet) was constructed by copying and translating the smaller 16 CL patch. This system was then equilibrated for 100 ns before the production simulation, which lasted 300 ns. Thereafter, 12 CL molecules (6 in each leaflet) were extracted from the pure CL bilayer system and superimposed on the equilibrated pure POPC bilayer. With a mix of 12 CL and 118 POPC, the lipid composition was 9.2% CL by mole. POPC molecules overlapping CL were removed, and the system was energy minimized before a 60 ns equilibration simulation, which allowed the bilayer area to relax. The production simulation of this system was also 300 ns. An additional 200 ns was then simulated for investigations of ion binding. The pure POPC bilayer was simulated for 145 ns, with a 15 ns equilibration where the area per lipid stabilized. Compositions and simulation times are summarized in Table 1. In the following, the three investigated systems will be denoted as 0-CL, 9-CL, and 100-CL for pure POPC, 9.2% CL in POPC, and pure CL, respectively. The Berendsen thermostat and semianisotropic barostat44 were used with coupling times of 0.1 and 1 ps, respectively, with the temperature set to 310 K and the pressure to 1 bar. A time step of 5 fs was used throughout, and bonds in lipids and water

were constrained with the LINCS and SETTLE algorithms, respectively.45,46 Results Molecular dynamics simulations were carried out on the three bilayers with explicit water and counterions. The resulting trajectories were analyzed with respect to the area per lipid, width, and compressibility, which give succinct information on the compactness of the bilayer. Binding and association of ions and water to the bilayer and the individual lipid components were investigated with geometric definitions of bonds. Furthermore, the density and dipole potential along the bilayer normal for various key atom groups and the systems as a whole were calculated. The polarizing effects of the bilayer and ions on water were investigated by monitoring the angle between the dipole and the bilayer normal. As a measure of the changes in lipid packing, the order of the acyl chains was compared in the different bilayers. Finally, parameters related to the lipid geometry, such as the intramolecular P-P distance and the P-C*-P angles in the headgroup, were compared to previous CG simulations and the PDB coordinates of CL. Area per Lipid and Compressibility. The area per lipid and compressibility modulus are collected in Table 1. The pure POPC reference bilayer showed an area per lipid of 0.64 nm2, which is in good agreement with previous simulations of POPC.41 In the mixed bilayer system, the area per lipid species is not available directly from the simulation. The canonical method described by Edholm and Nagle could be used if additional lipid compositions were studied.47 Instead we directly calculated the average area per lipid, which in the 9-CL system was 0.62 nm2. Thus, even though CL has four chains and two phosphate groups, the total area per lipid in the 9-CL system decreased relative to pure POPC. The excess area per lipid, calculated as the difference between the ideal combination of the pure phases and the area measured in the mixed system, was -0.05 nm2, which is in reasonable agreement with experiments on mixed POPC/CL monolayers.20

Cardiolipin Bilayers

J. Phys. Chem. B, Vol. 112, No. 37, 2008 11657 TABLE 2: Sodium Ion to Lipid Carbonyl Oxygen Coordination Numbera

Figure 2. Distribution of distances between phosphorus atoms in the direction normal to the bilayer for 0-CL, 9-CL, and 100-CL. Snapshots from the respective systems are also shown.

The area per lipid for the 100-CL bilayer was 0.99 nm2, which is remarkably low when compared to uncharged phospholipids. The area per acyl chain was smaller, though comparable to the area of the closely related anionic lipid POPG (0.530 nm2 per lipid).41 Charged lipids have indeed been reported to have significantly reduced areas per lipid when compared to similar zwitterionic lipids.41,48,49 Tetramyristoyl-CL, with all saturated chains, was previously estimated to have an area per lipid of 1.0 nm2 in the lamellar (LR) phase,43 which is in good agreement with the current result for tetraoleoyl-CL. From the total area and its fluctuations, the area compressibility modulus, KA, can be defined: KA ) kT〈A〉(N〈(A - A0)2〉)-1, where A is the area per lipid, N the number of lipids per leaflet, and A0 is the equilibrium area per lipid.50 For 0-CL, 9-CL, and 100-CL, KA was found to be 350, 340, and 1100 mN/m, respectively. Whereas the 0-CL and 9-CL systems were very similar, the 100CL system showed a significant increase in compressibility modulus consistent with high rigidity. Interestingly, there is experimental evidence that small additions of CL reduce the cohesive strength of PC vesicles, from 153 mN/m at 0% CL to 115 mN/m at 9.2% CL together with a decreased lysis tension.18 However, in the same work the elastic modulus changed very little upon adding the charged lipids to POPC vesicles, which is in agreement with the 9-CL and 0-CL simulations having similar compressibility modulus. The systematic deviation between the MD simulation and experiments on CL bilayers can be explained by the difference in length scales. For our pure POPC bilayer, the compressibility modulus can be compared with 250 mN/m,2 and 178 mN/m from vesicle experiments.51 Long undulation and buckling modes, which tend to decrease the compressibility modulus in experimental settings, are not seen in the simulation because of the relatively small dimensions of the box.52 Because of the inverse dependence on the variance of the area, KA will also tend to be systematically overestimated for short simulations. Thickness variations in the bilayers also give some information about bilayer flexibility. The z-component of the interphosphate distances is presented in Figure 2. The 100-CL bilayer was highly planar and in plane movement was limited on the 300 ns time scale. The width of the distribution is connected to all deviations from the planar average structure, and from the changes upon adding CL we can conclude that CL and counterions induce bilayer rigidity. Representative snapshots from the trajectory, also shown in Figure 2, show the thickening of the bilayers with increasing CL contents. The mean distance between the phosphate groups is given in Table 1. The significant increase in bilayer thickness is a consequence both of the increased fraction of the 18 carbon oleoyl chains, as compared to the 16 carbon palmitoyl chains on the sn-1 position of POPC, and of the increased packing in the x-y plane. Note the significant increase in thickness (more than 0.3 nm) between 0-CL and 9-CL. Ion Binding. In Table 2, the number of lipid carbonyl oxygens in the first coordination shell (within 0.33 nm) of the

N

9-CL

100-CL

0 1 2 3 4

2(3 4(4 15 ( 7 42 ( 9 37 ( 7

35 ( 2 14 ( 3 22 ( 3 23 ( 3 6(1

a

Coordination number, N, from the 9-CL and 100-CL systems (in %). Uncertainties are standard deviations.

sodium ions is shown. The carbonyl oxygen was the main site for sodium binding. There is experimental evidence that the primary binding site for metal ions is the phosphate group.53-55 However, in simulations with the current choice of force field,41,56 sodium ions interact with the carbonyl oxygen, with a preference for the sn-1 chain over the sn-2 chain. Recently, MD simulations of ion binding to DPPC bilayers with a different choice of ion parameters was reported, and significant binding of ions to both phosphate and carbonyl oxygens was found.57 With the observed degree of binding, at least 65% of the positive ions are on average located below the negatively charged phosphate groups in the 100-CL system. This strong association between ions and lipids implies that the anionic surface is significantly screened. In the 9-CL bilayer, the number of ions with bonds to lipid oxygens (either from CL or POPC) was even larger, 98%. Notable differences between the 100CL and 9-CL systems are the very low amount of free sodium ions in the aqueous subphase in 9-CL, and the corresponding increase in sodium ions with 3 or 4 coordinating carbonyls. The charge from the CL anionic headgroups and ions would now be just 2% of the nominal surface charge of deprotonated CL, which is considerably smaller than the 16% ionization arrived at by Nichols-Smith et al.17 Still, the experimentally observed surface charge in CL mixtures was significantly smaller than the expected ionization given the nominal surface charge of two negative charges per CL. Because our simulation models have a fixed protonation state, and the size of the system prohibits reliable modeling of protons at the relevant concentrations, we are unable to determine the role of protonation to the electrostatic potential of the bilayer. When we compare 9-CL with 100-CL, however, we see that charge density is not the sole determinant of ion binding to the bilayer. Rather, for low CL content (low lipid charge density) the sodium ion binding was almost complete (98%), whereas the pure CL bilayer, with a much higher lipid charge density, bound 65% of the ions. The high degree of ions within the headgroup region suggests that the observed low long-range electrostatics is not necessarily because of protonation of the phosphates, a conclusion also reached in experimental work.17 Because CL has two pairs of carbonyl groups, an ion can bind to several carbonyls from the same CL. The occurrence of such ions was monitored in the 9-CL and 100-CL simulations, and we found that there was indeed a significant population of ions with at least two bonds to the same CL. In the 9-CL system, 3.6 ( 0.9 ions out of 24 were bound to more than one carbonyl of the same CL. This means that 30% of the CL molecules were engaged in double coordination of ions. The corresponding number in the 100-CL simulation was 12.8 ( 2.6 ions out of 128 and engages 20% of the CL molecules. Double coordination of ions has been suggested for calcium ion binding to CL but was then proposed to occur between the phosphate groups.58 Multiple coordination at the level of the carbonyl groups adds another possible mechanism for the reduction of lipid area upon adding salt to CL bilayers.

11658 J. Phys. Chem. B, Vol. 112, No. 37, 2008

Dahlberg and Maliniak

Figure 3. Time dependence of number of bonds with CL carbonyl oxygens in the 9-CL system: hydrogen bonds from water (top), and ions with an indicated number of carbonyl oxygens within 0.33 nm (bottom).

In Figure 3 the time dependence of the binding preference to CL in the 9-CL bilayer can be seen. Interestingly, there was a substantial variation in the binding preference between POPC and CL on the time scale of the simulation (500 ns). In contrast, the total number of ions bound to either POPC or CL was very high (98%, Table 2) and varied little during the simulation (data not shown). On average, the binding of water and ions to CL was complementary, which can be seen in the approximately mirrored changes in ion binding and hydrogen bonds between water and the carbonyl group. There was always an excess of binding to CL, counted per carbonyl group, with an average ratio of 1.5 ( 0.3 between CL and POPC (with variance given as standard error over the 500 ns). The slow lateral movement of lipids and strong association of ions to the bilayer make changes in coordination intrinsically slow. The ions were found to go from the water phase into the lipid headgroup region immediately after the start of the equilibration simulation. This behavior is similar to the results seen previously in simulations.59 In the production simulation there was some exchange between ions bound to lipids and ions in the water phase, although most ions remained bound to the bilayer for the duration of the simulation. Density and Electrostatic Potential. The electrostatic dipole potential in the direction perpendicular to the bilayer was calculated by integrating the negative charge density twice with respect to the z coordinate. There were two main changes in the dipole potential, shown in Figure 4, as we increased CL concentration. First, the total dipole potential, calculated as the difference between the middle of the bilayer and the middle of the aqueous phase, was increased in magnitude by 22 mV between 0-CL and 9-CL. The shape of the potential was overall the same but shifted along the bilayer normal. Even though the change in the total dipole potential was small, there were substantial changes in the contributions from the different components of the bilayer. The increase in the choline contribution to the total dipole potential was 0.48 V, whereas the glycerol and carbonyl contribution changed 0.12 V. The pure POPC bilayer, 0-CL, induced a water dipole potential of -5.4 V, which can be compared to -4.6 V in the 9-CL system. A similar change in the water dipole potential has been seen in previous work when adding salt to POPC membranes.56 The total dipole potential and the increased height of the potential barrier make crossing the interface more difficult for positive ions. The increase in the total dipole potential was in part due to a decreased tilt away from the bilayer normal by the large P-N dipole. Second, we observed a lower potential barrier for passing the interface in the 100-CL than in the 9-CL and 0-CL systems, with 170 and 120 mV, respectively. A significant surface

Figure 4. Electrostatic potential (A) and mass density distributions (B-E) are shown as a function of distance from the middle of the bilayer. The -40 mV surface potential of 100-CL is indicated in A.

potential appeared at 2.9 nm from the bilayer center, just outside the phosphate groups. The surface potential was -40 mV, calculated as the difference between the water midplane potential and the potential at 2.9 nm. Experimental studies of the electrostatics of bilayers with CL report surface potentials comparable in magnitude to the potential found in the 100-CL simulation. For a 1:9 CL:POPC mixture, a surface potential of -55 mV was calculated from Ca2+ binding studies,60 and -81 ( 28 mV was reported for a 9.2% CL bilayer with a surface force apparatus technique.17 The 9-CL simulation, however, showed no significant surface potential. The first local potential minimum from the bilayer center, seen in the 0-CL and 9-CL potentials, disappeared in 100-CL. The total dipole potential of the system did not change significantly (about -2 mV) between 9-CL and 100-CL. Interestingly, there was a large change in the dipole potential generated by the water molecules, which decreased from -4.6 to -6.5 V. Thus, this polarization screened the other components in the system. From the mass density along the bilayer normal, displayed in Figure 4, we see that the highest ion density was found closer to the bilayer center that the phosphate groups in both the 100CL and 9-CL systems, which is in line with a substantial association between ions and carbonyl groups. The bimodal ion density distribution is similar to that found previously for POPG,41 where the outer local maximum was found to be ions not bonded to the lipid oxygens. Whereas a significant number of ions stayed in the water phase in the 100-CL simulation, we noticed that no ions in the 9-CL system stayed in that region. Hydration and Water Dipole Distribution. Water molecules were found deep in the lipid headgroup region, which resulted in an effective shielding of the charges from the lipids and counterions. All oxygens in the lipids showed some degree of

Cardiolipin Bilayers

J. Phys. Chem. B, Vol. 112, No. 37, 2008 11659

Figure 5. The averages of the cosine of the angle between the dipole vector of H2O and the outward bilayer normal for the 0-CL, 9-CL, and 100-CL systems as a function of distance from the bilayer center.

Figure 6. The local order parameter, SCD, for 100-CL (black solid), 9-CL POPC oleoyl (blue solid), 9-CL CL oleoyl (black dash), and 0-CL (green solid). POPC palmitoyl from 9-CL (blue dash), and 0-CL (green dash).

solvation by water. Although the interface between the lipids and water became narrower with increasing CL concentration (see water density distributions in Figure 4), the amount of hydrogen bonds from water actually increased. As a measure of the hydration level, we introduce the number of hydrogen bonds from water to lipids per bilayer area. This reduced number is introduced to make comparisons between lipids of different sizes easier. For the 0-CL, 9-CL, and 100-CL systems, the hydration level was 10.1, 10.3, and 13.8 bonds/nm2, respectively. The hydration level of the membrane thus increased with CL concentration, although the difference was quite small between the 0-CL and 9-CL. Earlier results, based on fluorescence lifetime measurements,61 indicated that the amount of bound water at the lipid-water interface increases considerably in the presence of CL, which is in agreement with the hydration level reached here. The most drastic change in hydration of POPC occurred at the sn-2 carbonyl oxygen, which lost 35% of its hydrogen bonds to water when comparing 0-CL (1.65 ( 0.01 bonds) and 9-CL (1.08 ( 0.02). Conversely, the sn-2 carbonyl oxygen of CL increased its number of hydrogen bonds to water by 59% going from 100-CL (0.51 ( 0.01) to 9-CL (0.8 ( 0.2). The water bound to the phosphate group oxygens in CL decreased from 1.9 to 1.7 hydrogen bonds per phosphate oxygen in 100-CL and 9-CL bilayers, respectively. This change in hydration of CL is consistent with the free choline groups of POPC, which are located in the outermost regions of the bilayer, and interact electrostatically with CL phosphate charge. Counted per acyl chain, CL has a substantially smaller headgroup than phosphatidylcholines and gives less sterical hindrance for water binding. This also affects the headgroup hydroxyl group, which was dehydrated from 1.48 ( 0.01 to 1.39 ( 0.01 hydrogen bonds to water when going from 100-CL to 9-CL. Also, water molecules clearly experience a very different environment close to a highly charged surface when compared to a zwitterionic surface. The effect is clear in the average orientation of the water dipoles relative to the outward bilayer normal. The average cosine of the dipole angles is shown in Figure 5. There was a minimum in the distribution at 1.8 nm for 0-CL was shifted outward with increasing amounts of CL. The distributions were overall similar to those obtained previously for POPG.41 We also saw the effect of the changes in charge density at the interface, causing oscillations in the water dipole orientation distribution at around 2 nm from the bilayer center. Water molecules close to the midplane between the interfaces retained a nonisotropic orientation in 100-CL. This implies that no real bulk water was present in the simulation, even though 45.7 waters per lipid were used, which is 50% more water than the saturation limit for zwitterionic lipids.62,63 Finally, the small amount of water passing through or into the bilayer gave poor statistics within about 1 nm of the center of the bilayer. Order Parameters. The order parameter, SCD, for the acyl chains was calculated using SCD ) |1/2〈3 cos2 θ - 1〉|, where θ

is the angle between the CD vector and the bilayer director, taken to be parallel to the z-axis of the simulation reference frame. This vector was chosen because it corresponds to the principal component of the electric field gradient (EFG) at the deuterium nucleus and reflects the motion that determines the 2H NMR spectrum.64 Tetrahedral symmetry was assumed when calculating CD for the methylene groups. For unsaturated carbons, the hydrogen atoms were assumed to be in plane with the closest connected carbons, and in the direction of the bisector to the angle formed by Cn-1 to Cn, and Cn+1 to Cn. The effect of CL and ions on the area per lipid and the increased rigidity was apparent in the order parameter of the hydrocarbon chains, shown in Figure 6. Compared to POPG,41 the overall behavior for the oleoyl chain was similar, although the order was slightly larger in the CL chains. This implies that although the amount of monounsaturated lipid increased from 50% in POPG to 100% in CL, the disorder due to kinks from the cis bonds was compensated by other factors, one of which is likely the mobility constraints due to the unique headgroup in CL. Even though the fraction of CL was small, the 9-CL simulation showed a significantly increased order for POPC, relative to the pure POPC bilayer. The increased order of POPC with 9-CL can be explained by the high degree of counterion binding, not only to CL but also to the carbonyl groups of POPC. The 0-CL POPC order parameters are in good agreement with the experimentally derived counterparts,64 also for the two C-H vectors connected to the unsaturated carbons. The effect of increased order when salt is added to POPC bilayers has been observed previously in simulations.56 In our simulations, the ordering effect was larger, even though the amount of sodium ions was comparable. Thus, counterion binding alone was not responsible for the increased order. The relatively small headgroup of CL, when compared to its chain volume, likely adds to the reduction in area per lipid and thus increases the order parameter. The oleoyl chains in CL and POPC had very similar order parameters, implying that there was a similar environment for the chains on either lipid. This conclusion seems reasonable, given that no lipid domains could be found when visually inspecting the trajectory. Molecular Geometries. We begin the description of the geometry of CL by investigating the distance between phosphate groups in the headgroup. The distributions from the present united-atom (UA) model, the previously studied CG model, and the PDB structures are shown in Figure 7. In our united-atom simulations the distributions showed two main states, with P-P distances of 0.55 and 0.63 nm. There was also significant population of intermediate states. Compared to the 100-CL system, the average intramolecular P-P distance increased slightly as POPC was introduced, which can be seen as a shift in population from the site at 0.55 nm to the site at 0.63 nm in 9-CL. This expansion is connected to the larger area per lipid seen in the 9-CL bilayer compared to the 100-CL bilayer, but

11660 J. Phys. Chem. B, Vol. 112, No. 37, 2008

Figure 7. Distributions of CL geometries: intramolecular P-P distances and molecular structures of two representative members from the PDB (A), P-C*-P angle (B), and acyl-C*-acyl angle (C); 100CL (black), 9-CL (blue), CG (gray), and PDB (dash). Order parameters for the vector from particle n to n + 1 for coarse-grained (CG) and united-atom (UA) models of 100% CL bilayers are shown in the inset (C).

it is also interesting to note that whereas the general driving force for a decreased area per lipid is based on tighter binding through ions between lipids, the relative compactness seen in 100-CL also has a slight contribution from changes in the intramolecular geometry. The increased sodium ion concentration close to, but not bound to, the phosphate groups can also screen the Coulombic interactions of the latter, allowing a more compact headgroup. Because detailed information of CL structure is scarce, but several known protein structures contain these lipids,65 we compared our model geometries to those of CL found in the PDB (the names of the structures can be found in the Supporting Information). The total number of structures was 94, although several of these represent different reported structures of the same protein. An indication of the accuracy in these structures is given by the average resolution 0.24 nm (min 0.16 nm, max 0.31 nm) and the average temperature (Debye-Waller) factor, which is 0.82 nm2 (min 0.112 nm2, max 1.39 nm2), the latter corresponding to displacements of about 0.1 nm around the equilibrium position. The distribution of P-P distances is clearly bimodal, with a maximum at about 0.55 nm. At 0.7 nm there is another cluster of structures. The interval between 0.55 and 0.7 nm is most populated at the shorter distances. All PDB data were derived from crystal structures solved by X-ray diffraction. Because the lipid was cocrystallized with protein, these lipid conformations do not necessarily represent bilayer conformations. However, the information does give some insight into the molecular geometry. First, the distribution showed a large spread, where the entire range sampled in the simulations is covered. Second, the large collection of structures around 0.7 nm was not represented in either the 9-CL or 100-CL bilayers

Dahlberg and Maliniak from simulations, although such geometries are not forbidden; the simulations sample as far as 0.78 nm with low population. There are at least two possible reasons for this discrepancy. The force field has not been optimized for the unique CL headgroup structure, which means that charges and dihedral potentials could be poorly chosen. We find it unlikely that standard bonded parameters, such as bond lengths and angles, gave a contribution large enough to shift the maximum of the P-P distance from 0.63 to 0.7 nm. Additionally, specific binding to proteins can induce changes in the lipid conformation. The largest contribution to lipids with 0.7 nm P-P distance are from a group of structures with a P-C*-P angle close to 170°. These are found wedged at the interface between two protein subunits in the cytochrome c oxidase models.66,67 Such conformations were not sampled in the UA simulations and do not represent viable bilayer conformations. Examples of CL structures from cytochrome c oxidase with two different P-P distances and conformations are shown in Figure 7. The choice of angles and distances above enables comparisons with the previous study of CL with a CG model,32 and the distributions of these parameters are also included in Figure 7. There are several significant differences between the present united-atom (UA) MD simulations and the CG study. With the CG model, the P-P distance was overestimated relative to the UA and PDB models, although the width of the distributions were similar. The headgroup angle was also overestimated with CG, when comparing it with the UA model. Although the phase properties of the CG model were shown to agree with experimental findings,32 we find it likely that the UA approach provides a better model of the geometry in the CL headgroup given the explicit description of bond distances and valence angles. The splayed chain conformations found when nearby CG lipid aggregates, bilayers, or micelles came into contact with each other were not found in the UA simulations. The high hydration level used here prevents close contact between bilayers, and the splayed chain conformation then has a highenergy penalty. As a measure of the overall lipid geometry, the angle distribution between the headgroup glycerol group and the outer terminal tail beads was calculated and shown in Figure 7. Groups of three to four atoms in the UA trajectories were mapped onto CG particles by taking the group center of mass as the position of the CG particle. Different mappings were tried but did not change the qualitative outcome of the angle distributions. The UA model had a significantly smaller average spread angle than the CG model. Important differences between the two models are the intramolecular distances and the site of counterion binding to the bilayer. The P-P distance, for instance, was exaggerated in the CG model, and the current UA results are more in line with the geometries found in the PDB. The CG model also lacks the detail to allow specific binding at the carbonyl groups. The combination of a smaller headgroup and deep counterion binding result in higher order and tighter chain packing. A significant difference in the area per lipid should thus be noted, 0.99 nm2 for the UA model compared to 1.29 nm2 in the CG model. To further explore the differences in overall lipid structure, the mapped trajectories were also used to calculate the order parameters for the vectors between particles n and n+1 in the acyl chains. A comparison of the CG order parameters and the current UA derived order parameters is shown in Figure 7. The CG model of CL was previously found to have a higher order than zwitterionic lipids within the same CG model set. Here, the UA model showed a substantially stronger ordering than

Cardiolipin Bilayers the CG model, and the effect was especially clear near the center of the bilayer. With order parameters close to unity, the chains were essentially parallel with the bilayer normal. Discussion From the simulations in this work we conclude that CL tends to increase the order of POPC bilayers, and that the effect is significant even with 9.2% CL by mole, the concentration found in the inner membrane of eukaryotic mitochondria. Recent work on monolayers, using grazing incidence X-Ray diffraction also suggests that the order of PC hydrocarbon chains is increased when CL is introduced.68 Previous fluorescence anisotropy experiments indicate that inclusion of up to 20 mol% CL in PC vesicles does not, however, change bilayer fluidity.61 Shibata et al.69 reported similarly, that up to 20 mol% CL does not cause a phase change in the acyl chain region in phosphatidylcholine bilayers, as judged by FTIR measurements. We agree with the conclusion that a phase transition is not reached, even for the pure CL system, but that there is a tendency to form structurally more compact bilayers with CL. The negative excess area per lipid found in our simulations was comparable to previous work.70 However, in that work a continuous increase in the area per lipid for increasing CL content was observed, which is at odds with our slight decrease of the total area per lipid at 9.2% CL. Nichols-Smith et al.18 found small positive excess areas at 10% CL in egg-PC/CL monolayers but also reported that the excess area turned negative as salt was added. From other experimental work on a similar mixture of lipids, DPPC and CL, small deviations from ideality in the area per molecule has been reported.71 These observations highlight the sensitive balance of forces responsible for the area variations with composition, and that ion interactions with the lipids are at least in part responsible for the condensing effect CL shows on bilayers. Because lipid mobility is low in bilayers with bound ions, the time needed for equilibration and sampling is long. We observed long time fluctuations in the number of ions bound to CL when it is mixed with POPC, even though the overall preference for binding to carbonyl oxygens of CL was higher than for POPC. This is what is expected for an anionic lipid, although the site for interaction was not the charged phosphate but the carbonyl groups. It is likely that the choice of ion parameters affects the depth of penetration into the headgroup region. Recent simulations of POPC bilayers with NaCl solution, using CHARMM parameters for added ions, indeed show a shift from binding primarily at the carbonyl groups toward binding to the phosphate group.72 Similarly, new ab initio-derived charges for the lipid headgroup have been used in simulations of phosphatidic acid (PA), which resulted in ion binding primarily to the phosphate group.73 In both cases the overall bilayer properties were not drastically changed by the change in ion position. Preliminary results from simulations of the 9-CL system with OPLS parameters for the sodium counterions also indicate that the position of the ions shifts toward the phosphate groups, but that a significant number of bound ions still coordinate the carbonyl groups (approximately 60%). The reason for preferring one parameter set over another is, however, currently weak.57,74 It is notable that practically all ions bound to the bilayer in the 9-CL system, making the bilayer effectively neutral. For the 100-CL system, ion binding was not complete and there was evidence of a negative surface potential, which is expected of a bilayer with anionic headgroups. The 0.3 nm increase in bilayer width upon introducing 9.2% CL into POPC bilayers can in principle be verified with X-ray or neutron scattering experiments.

J. Phys. Chem. B, Vol. 112, No. 37, 2008 11661 Mixtures of the main IMM components, PE/CL and PC/PE/ CL, have previously been found to form lipid domains.18,19,21,22,71,75 Unlike these mixed systems, PC/CL has been found not to associate into domains, although the existence of superlattice structures has been hypothesized.76 We did not see domain formation in the simulations, but this should not be regarded as strong evidence against domain or superlattice formation, because of the limited lateral diffusion of lipids on the 100-ns time scale. Also, it is uncertain if the hypothesized superlattice structures would be thermodynamically stable on the length scale of the simulated systems. Unlike zwitterionic lipids like PC and PE, which are able to shield their phosphate groups by aligning the positively charge amine group, CL cannot shield its charges by rearranging itself. This leaves the CL headgroup sensitive to the composition of the aqueous subphase. Indeed, apart from the ions discussed previously, there were an increased number of water molecules forming hydrogen bonds with the 100-CL bilayer, and a large change in the water dipole potential between 9-CL and 100-CL simulations. The conformations obtained from simulations showed partial overlap with the structures of CL obtained from PDB structures. The differences can be due to (i) conformational changes in the lipid upon binding to protein, (ii) poor representation of small molecules in the PDB, (iii) force field parameters which do not allow the fully stretched conformation (180° between P-C*-P). We find it likely that conformational changes are induced by protein binding, although further study of the nature of this interaction is needed. The poor representation of small molecules in the protein structures deposited in the PDB77,78 certainly also limits direct comparisons between simulations and X-ray structures of CL. We find it likely that the clearly bimodal P-P distance distribution from the PDB structures really does represent two different structural states, given the clear difference in the P-C*-P angle, and the position of the lipid between two protein subunits. Further investigations into the structure and charge state of the CL headgroup are necessary to establish the connection between CL in the bilayer and CL bound to proteins. However, the present model provides some insight into the association of ions to bilayers with CL and the strong order it induces. The present model also provides a straightforward way of improving coarse-grained models, which hopefully will result in more accurate mesoscale models of the IMM. Acknowledgment. We thank Samuli Ollila, Siewert Jan Marrink, and Alex H. de Vries for helpful discussions. This work was supported by the Swedish Research Council. Supporting Information Available: PDB structures used for structure analysis. This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Evans, D. F.; Wennerstro¨m, H. In The Colloidal Domain: Where Physics, Chemistry, Biology, and Technology Meet; Wiley - VCH: New York, 1999. (2) Nagle, J. F.; Tristram-Nagle, S. Structure of Lipid Bilayers. Biochim. Biophys. Acta 2000, 1469, 159–195. (3) Mannella, C. A. Structure and Dynamics of the Mitochondrial Inner Membrane Cristae. Biochim. Biophys. Acta 2006, 1763, 542–548. (4) Mannella, C. A. The Relevance of Mitochondrial Membrane Topology to Mitochondrial Function. Biochim. Biophys. Acta 2006, 1762, 140–147. (5) Chvanov, M. Metabolic Control of Elastic Properties of the Inner Mitochondrial Membrane. J. Phys. Chem. B 2006, 110, 22903–22909. (6) John, G. B.; Shang, Y.; Li, L.; Renken, C.; Mannella, C. A.; Selker, J. M.; Rangell, L.; Bennett, M. J.; Zha, J. The Mitochondrial Inner Membrane Protein Mitofilin Controls Cristae Morphology. Mol. Biol. Cell 2005, 16, 1543–1554.

11662 J. Phys. Chem. B, Vol. 112, No. 37, 2008 (7) Schlame, M.; Rua, D.; Greenberg, M. L. The Biosynthesis and Functional Role of Cardiolipin. Prog. Lipid Res. 2000, 39, 257–288. (8) Gohil, V. M.; Hayes, P.; Matsuyama, S.; Schagger, H.; Schlame, M.; Greenberg, M. L. Cardiolipin Biosynthesis and Mitochondrial Respiratory Chain Function are Interdependent. J. Biol. Chem. 2004, 279, 42612– 42618. (9) Daum, G. Lipids of Mitochondria. Biochim. Biophys. Acta 1985, 822, 1–42. (10) Hoch, F. L. Cardiolipins and Biomembrane Function. Biochim. Biophys. Acta 1992, 1113, 71–133. (11) McAuley, K. E.; Fyfe, P. K.; Ridge, J. P.; Isaacs, N. W.; Cogdell, R. J.; Jones, M. R. Structural Details of an Interaction between Cardiolipin and an Integral Membrane Protein. Proc. Natl. Acad. Sci. U.S.A. 1999, 96, 14706–14711. (12) Haines, T. H.; Dencher, N. A. Cardiolipin: A Proton Trap for Oxidative Phosphorylation. FEBS Lett. 2002, 528, 35–39. (13) Schlame, M.; Ren, M.; Xu, Y.; Greenberg, M. L.; Haller, I. Molecular Symmetry in Mitochondrial Cardiolipins. Chem. Phys. Lipids 2005, 138, 38–49. (14) Rytomaa, M.; Mustonen, P.; Kinnunen, P. K. Reversible, Nonionic, and pH-Dependent Association of Cytochrome c with Cardiolipin-Phosphatidylcholine Liposomes. J. Biol. Chem. 1992, 267, 22243–22248. (15) Gorbenko, G. P.; Molotkovsky, J. G.; Kinnunen, P. K. Cytochrome C Interaction with cardiolipin/phosphatidylcholine Model Membranes: Effect of Cardiolipin Protonation. Biophys. J. 2006, 90, 4093–4103. (16) Nichols-Smith, S.; Bao, T.; Kuhl, T. Electrostatic Interactions in Model Mitochondrial Membranes: Cardiolipin and Egg Pc Mixtures. Biophys. J. 2002, 82, 156A–157A. (17) Nichols-Smith, S.; Kuhl, T. Electrostatic Interactions between Model Mitochondrial Membranes. Colloids Surf., B 2005, 41, 121–127. (18) Nichols-Smith, S.; Teh, S. Y.; Kuhl, T. L. Thermodynamic and Mechanical Properties of Model Mitochondrial Membranes. Biochim. Biophys. Acta 2004, 1663, 82–88. (19) Lupi, S.; Perla, A.; Maselli, P.; Bordi, F.; Sennato, S. Infrared Spectra of Phosphatidylethanolamine-Cardiolipin Binary System. Colloids Surf., B 2008, 64, 56–64. (20) Domenech, O.; Sanz, F.; Montero, M. T.; Hernandez-Borrell, J. Thermodynamic and Structural Study of the Main Phospholipid Components Comprising the Mitochondrial Inner Membrane. Biochim. Biophys. Acta 2006, 1758, 213–221. (21) Domenech, O.; Morros, A.; Cabanas, M. E.; Montero, M. T.; Hernandez-Borrell, J. Thermal Response of Domains in Cardiolipin Content Bilayers. Ultramicroscopy 2007, 107, 943–947. (22) Domenech, O.; Redondo, L.; Picas, L.; Morros, A.; Montero, M. T.; Hernandez-Borrell, J. Atomic Force Microscopy Characterization of Supported Planar Bilayers that Mimic the Mitochondrial Inner Membrane. J. Mol. Recognit. 2007, 20, 546–553. (23) Ponnuswamy, A.; Nulton, J.; Mahaffy, J. M.; Salamon, P.; Frey, T. G.; Baljon, A. R. Modeling Tubular Shapes in the Inner Mitochondrial Membrane. Phys. Biol. 2005, 2, 73–79. (24) Pandit, S. A.; Bostick, D.; Berkowitz, M. L. Mixed Bilayer Containing Dipalmitoylphosphatidylcholine and Dipalmitoylphosphatidylserine: Lipid Complexation, Ion Binding, and Electrostatics. Biophys. J. 2003, 85, 3120–3131. (25) de Joannis, J.; Jiang, Y.; Yin, F.; Kindt, J. T. Equilibrium Distributions of Dipalmitoyl Phosphatidylcholine and Dilauroyl Phosphatidylcholine in a Mixed Lipid Bilayer: Atomistic Semigrand Canonical Ensemble Simulations. J. Phys. Chem. B 2006, 110, 25875–25882. (26) de Joannis, J.; Jiang, F. Y.; Kindt, J. T. Coarse-Grained Model Simulations of Mixed-Lipid Systems: Composition and Line Tension of a Stabilized Bilayer Edge. Langmuir 2006, 22, 998–1005. (27) Sugar, I. P.; Biltonen, R. L. Lateral Diffusion of Molecules in TwoComponent Lipid Bilayer: A Monte Carlo Simulation Study. J. Phys. Chem. B 2005, 109, 7373–7386. (28) Shi, Q.; Voth, G. A. Multi-Scale Modeling of Phase Separation in Mixed Lipid Bilayers. Biophys. J. 2005, 89, 2385–2394. (29) Niemela, P. S.; Ollila, S.; Hyvonen, M. T.; Karttunen, M.; Vattulainen, I. Assessing the Nature of Lipid Raft Membranes. PLoS Comput. Biol. 2007, 3, e34. (30) Baoukina, S.; Monticelli, L.; Amrein, M.; Tieleman, D. P. The Molecular Mechanism of Monolayer-Bilayer Transformations of Lung Surfactant from Molecular Dynamics Simulations. Biophys. J. 2007, 93, 3775–3782. (31) Ollila, S.; Niemela, P.; Rog, T.; Hyvonen, M. T.; Karttunen, M.; Vattulainen, I. Differences between Lateral Pressure Profiles of OneComponent and Many-Component Lipid Membranes. Chem. Phys. Lipids 2007, 149, S15-S16. (32) Dahlberg, M. Polymorphic Phase Behavior of Cardiolipin Derivatives Studied by Coarse-Grained Molecular Dynamics. J. Phys. Chem. B 2007, 111, 7194–7200. (33) Mukhopadhyay, R.; Huang, K. C.; Wingreen, N. S. Lipid Localization in Bacterial Cells through Curvature-Mediated Microphase Separation. Biophys. J. 2008,

Dahlberg and Maliniak (34) Krebs, J. J. R.; Hauser, H.; Carafoli, E. Asymmetric Distribution of Phospholipids in the Inner Membrane of Beef-Heart Mitochondria. J. Biol. Chem. 1979, 254, 5308–5316. (35) Kates, M.; Syz, J. Y.; Gosser, D.; Haines, T. H. PH-Dissociation Characteristics of Cardiolipin and its 2′-Deoxy Analogue. Lipids 1993, 28, 877–882. (36) Berger, O.; Edholm, O.; Jahnig, F. Molecular Dynamics Simulations of a Fluid Bilayer of Dipalmitoylphosphatidylcholine at Full Hydration, Constant Pressure, and Constant Temperature. Biophys. J. 1997, 72, 2002– 2013. (37) Chiu, S. W.; Clark, M.; Balaji, V.; Subramaniam, S.; Scott, H. L.; Jakobsson, E. Incorporation of Surface Tension into Molecular Dynamics Simulation of an Interface: A Fluid Phase Lipid Bilayer Membrane. Biophys. J. 1995, 69, 1230–1245. (38) Egberts, E.; Marrink, S. J.; Berendsen, H. J. Molecular Dynamics Simulation of a Phospholipid Membrane. Eur. Biophys. J. 1994, 22, 423– 436. (39) Chiu, S. W.; Clark, M. M.; Jakobsson, E.; Subramaniam, S.; Scott, H. L. Optimization of Hydrocarbon Chain Interaction Parameters: Application to the Simulation of Fluid Phase Lipid Bilayers. J. Phys. Chem. B 1999, 103, 6323–6327. (40) Zhao, W.; Rog, T.; Gurtovenko, A. A.; Vattulainen, I.; Karttunen, M. Atomic-Scale Structure and Electrostatics of Anionic Palmitoyloleoylphosphatidylglycerol Lipid Bilayers with Na+ Counterions. Biophys. J. 2007, 92, 1114–1124. (41) Zhao, W.; Rog, T.; Gurtovenko, A. A.; Vattulainen, I.; Karttunen, M. Atomic-Scale Structure and Electrostatics of Anionic Palmitoyloleoylphosphatidylglycerol Lipid Bilayers with Na+ Counterions. Biophys. J. 2007, 92, 1114–1124. (42) de Vries, A. H.; Mark, A. E.; Marrink, S. J. The Binary Mixing Behavior of Phospholipids in a Bilayer: A Molecular Dynamics Study. J. Phys. Chem. B 2004, 108, 2454–2463. (43) Lewis, R. N.; Zweytick, D.; Pabst, G.; Lohner, K.; McElhaney, R. N. Calorimetric, X-Ray Diffraction and Spectroscopic Studies of the Thermotropic Phase Behavior and Organization of Tetramyristoyl Cardiolipin Membranes. Biophys. J. 2007, 92, 3166–77. (44) Berendsen, H. J. C.; Postma, J. P. M.; Van Gunsteren, W. F.; Dinola, A.; Haak, J. R. Molecular-Dynamics with Coupling to an External Bath. J. Chem. Phys. 1984, 81, 3684–3690. (45) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M. LINCS: A Linear Constraint Solver for Molecular Simulations. J. Comput. Chem. 1997, 18, 1463–1472. (46) Miyamoto, S.; Kollman, P. A. Settle - an Analytical Version of the Shake and Rattle Algorithm for Rigid Water Models. J. Comput. Chem. 1992, 13, 952–962. (47) Edholm, O.; Nagle, J. F. Areas of Molecules in Membranes Consisting of Mixtures. Biophys. J. 2005, 89, 1827–1832. (48) Elmore, D. E. Molecular Dynamics Simulation of a Phosphatidylglycerol Membrane. FEBS Lett. 2006, 580, 144–148. (49) Petrache, H. I.; Tristram-Nagle, S.; Gawrisch, K.; Harries, D.; Parsegian, V. A.; Nagle, J. F. Structure and Fluctuations of Charged Phosphatidylserine Bilayers in the Absence of Salt. Biophys. J. 2004, 86, 1574–1586. (50) Marrink, S. J.; deVries, A. H.; Mark, A. E. Coarse Grained Model for Semiquantitative Lipid Simulations. J. Phys. Chem. B 2004, 108, 750– 760. (51) Shoemaker, S. D.; Vanderlick, T. K. Intramembrane Electrostatic Interactions Destabilize Lipid Vesicles. Biophys. J. 2002, 83, 2007–2014. (52) Duncan, S. L.; Larson, R. G. Comparing Experimental and Simulated Pressure-Area Isotherms for DPPC. Biophys. J. 2008, 94, 2965– 86. (53) Tatulian, S. A. Binding of Alkaline-Earth Metal-Cations and some Anions to Phosphatidylcholine Liposomes. Eur. J. Biochem. 1987, 170, 413– 420. (54) Huster, D.; Arnold, K.; Gawrisch, K. Strength of Ca2+ Binding to Retinal Lipid Membranes: Consequences for Lipid Organization. Biophys. J. 2000, 78, 3011–3018. (55) Garidel, P.; Blume, A.; Hubner, W. A Fourier Transform Infrared Spectroscopic Study of the Interaction of Alkaline Earth Cations with the Negatively Charged Phospholipid 1,2-Dimyristoyl-sn-Glycero-3-Phosphoglycerol. Biochim. Biophys. Acta 2000, 1466, 245–259. (56) Bo¨ckmann, R. A.; Hac, A.; Heimburg, T.; Grubmu¨ller, H. Effect of Sodium Chloride on a Lipid Bilayer. Biophys. J. 2003, 85, 1647–1655. (57) Cordomi, A.; Edholm, O.; Perez, J. J. Effect of Ions on a Dipalmitoyl Phosphatidylcholine Bilayer. A Molecular Dynamics Simulation Study. J. Phys. Chem. B 2008, 112, 1397–1408. (58) Shah, D. O.; Schulman, J. H. Binding of Metal Ions to Monolayers of Lecithins, Plasmalogen, Cardiolipin, and Dicetyl Phosphate. J. Lipid Res. 1965, 6, 341–349. (59) Bo¨ckmann, R. A.; Grubmu¨ller, H. Multistep Binding of Divalent Cations to Phospholipid Bilayers: A Molecular Dynamics Study. Angew. Chem., Int. Ed. 2004, 43, 1021–1024.

Cardiolipin Bilayers (60) Macdonald, P. M.; Seelig, J. Calcium-Binding to Mixed Cardiolipin Phosphatidylcholine Bilayers as Studied by Deuterium Nuclear-MagneticResonance. Biochemistry 1987, 26, 6292–6298. (61) Chen, Q. P.; Li, Q. T. Effect of Cardiolipin on Proton Permeability of Phospholipid Liposomes: The Role of Hydration at the Lipid-Water Interface. Arch. Biochem. Biophys. 2001, 389, 201–206. (62) Ho¨gberg, C. J.; Lyubartsev, A. P. A Molecular Dynamics Investigation of the Influence of Hydration and Temperature on Structural and Dynamical Properties of a Dimyristoylphosphatidylcholine Bilayer. J. Phys. Chem. B 2006, 110, 14326–14336. (63) Marrink, S. J.; Berkowitz, M.; Berendsen, H. J. C. MolecularDynamics Simulation of a Membrane Water Interface - the Ordering of Water and its Relation to the Hydration Force. Langmuir 1993, 9, 3122– 3131. (64) Seelig, J.; Waespesarcevic, N. Molecular Order in Cis and Trans Unsaturated Phospholipid Bilayers. Biochemistry 1978, 17, 3310–3315. (65) Palsdottir, H.; Hunte, C. Lipids in Membrane Protein Structures. Biochim. Biophys. Acta 2004, 1666, 2–18. (66) Muramoto, K.; Hirata, K.; Shinzawa-Itoh, K.; Yoko-o, S.; Yamashita, E.; Aoyama, H.; Tsukihara, T.; Yoshikawa, S. A Histidine Residue Acting as a Controlling Site for Dioxygen Reduction and Proton Pumping by Cytochrome c Oxidase. Proc. Natl. Acad. Sci. U.S.A. 2007, 104, 7881–7886. (67) Tsukihara, T.; Shimokata, K.; Katayama, Y.; Shimada, H.; Muramoto, K.; Aoyama, H.; Mochizuki, M.; Shinzawa-Itoh, K.; Yamashita, E.; Yao, M.; Ishimura, Y.; Yoshikawa, S. The Low-Spin Heme of Cytochrome c Oxidase as the Driving Element of the Proton-Pumping Process. Proc. Natl. Acad. Sci. U.S.A. 2003, 100, 15304–15309. (68) Etienne, F.; Roche, Y.; Peretti, P.; Bernard, S. Cardiolipin Packing Ability Studied by Grazing Incidence X-Ray Diffraction. Chem. Phys. Lipids 2008, 152, 13–23. (69) Shibata, A.; Ikawa, K.; Shimooka, T.; Terada, H. Significant Stabilization of the Phosphatidylcholine Bilayer Structure by Incorporation of Small Amounts of Cardiolipin. Biochim. Biophys. Acta 1994, 1192, 71– 78.

J. Phys. Chem. B, Vol. 112, No. 37, 2008 11663 (70) Domenech, O.; Sanz, F.; Montero, M. T.; Hernandez-Borrell, J. Thermodynamic and Structural Study of the Main Phospholipid Components Comprising the Mitochondrial Inner Membrane. Biochim. Biophys. Acta 2006, 1758, 213–221. (71) Sennato, S.; Bordi, F.; Cametti, C.; Coluzza, C.; Desideri, A.; Rufini, S. Evidence of Domain Formation in Cardiolipin-Glycerophospholipid Mixed Monolayers. A Thermodynamic and AFM Study. J. Phys. Chem. B 2005, 109, 15950–15957. (72) Gurtovenko, A. A.; Vattulainen, I. Effect of NaCl and KCl on Phosphatidylcholine and Phosphatidylethanolamine Lipid Membranes: Insight from Atomic-Scale Simulations for Understanding Salt-Induced Effects in the Plasma Membrane. J. Phys. Chem. B 2008, 112, 1953–1962. (73) Dickey, A. N.; Faller, R. Examining the Contributions of Lipid Shape and Headgroup Charge on Bilayer Behavior. Biophys. J. 2008, doi: 10.1529/biophysj.107.128074. (74) Patra, M.; Karttunen, M. Systematic Comparison of Force Fields for Microscopic Simulations of NaCl in Aqueous Solutions: Diffusion, Free Energy of Hydration, and Structural Properties. J. Comput. Chem. 2004, 25, 678–689. (75) Pinheiro, T. J. T.; Duralski, A. A.; Watts, A. Phospholipid Headgroup-Headgroup Electrostatic Interactions in Mixed Bilayers of Cardiolipin with Phosphatidylcholines Studied by H-2 NMR. Biochemistry 1994, 33, 4896–4902. (76) Somerharju, P.; Virtanen, J. A.; Cheng, K. H. Lateral Organisation of Membrane Lipids. the Superlattice View. Biochim. Biophys. Acta 1999, 1440, 32–48. (77) Marsh, D.; Pali, T. Lipid Conformation in Crystalline Bilayers and in Crystals of Transmembrane Proteins. Chem. Phys. Lipids 2006, 141, 48– 65. (78) Kleywegt, G. J.; Henrick, K.; Dodson, E. J.; van Aalten, D. M. F. Pound-Wise but Penny-Foolish: How Well do Micromolecules Fare in Macromolecular Refinement. Structure 2003, 11, 1051–1059.

JP803414G