Acyl Chain Conformation and Packing in ... - American Chemical Society

J. Phys. Chem. B 2002, 106, 6273-6288

6273

Acyl Chain Conformation and Packing in Dipalmitoylphosphatidylcholine Bilayers from MD Simulation and IR Spectroscopy Robert G. Snyder,*,† Kechuan Tu,‡ Michael L. Klein,† Richard Mendelssohn,§ Herbert L. Strauss,† and Wenjun Sun† Department of Chemistry, UniVersity of California, Berkeley, California 94720-1460, Department of Chemistry, UniVersity of PennsylVania, Philadelphia, PennsylVania 19104-6323, and Department of Chemistry, Newark College of Arts and Science, Rutgers UniVersity, Newark, New Jersey 07102 ReceiVed: June 5, 2001; In Final Form: October 24, 2001

The lipid chain conformation predicted by a recent molecular dynamics simulation (MD) of the fluid and gel phases of a dipalmitoylphosphatidylcholine (DPPC) phospholipid bilayer is described in detail and, where possible, is compared to the results of IR measurements. Conformational statistics for the MD simulated fluid-phase chains are compared with those for an unconstrained liquidlike n-C16H34 system estimated using a rotational isomeric state model. The constraints imposed on the bilayer chains are found to reduce the concentration of gauche bonds, suppress sterically unfavorable conformational sequences, and increase trans/ gauche bond segregation. The concentrations of short (1-3 bond) conformational sequences for the MDsimulated bilayer are found to be nearly the same as those for the unconstrained RISM n-C16 system with the same gauche bond concentration as the constrained chains. For long sequences, however, especially long all-trans sequences, there are large differences between the systems. The bilayer-normal distribution of CC bonds across the gel phase bilayer is found to be largely determined by longitudinal chain displacement, whereas for the fluid phase, the distribution is determined by conformational disorder. Order-disorder heterogeneity is observed in the MD-simulated fluid-phase bilayer in the form of a low concentration of more or less well-defined, highly transient microdomains consisting of the order of 20 all-trans or nearly all-trans chains. The packing of these chains is similar to that found for the simulated gel phase. Both chaintilt and chain order/disorder correlations between opposing monolayers were found. The two most reliable IR methods for determining chain conformation are critically reviewed. Previously reported IR measurements on the gel- and fluid-phase DPPC bilayer using these methods are reinterpreted. The distribution of gauche bonds in the gel is found to be determined primarily by longitudinal displacements of the chain. By combination of the measured concentrations of bond-pair conformations tt, gt, and gg in the fluid phase, the overall concentration of gauche bonds is estimated to be 0.14 ( 0.04. This value is significantly lower than most previously reported values and much lower than the 0.28 derived from the MD simulation.

I. Introduction Nothing approaching an accurate molecular-level description of the conformation and packing of the hydrocarbon chains in a fully hydrated phospholipid bilayers in the fluid phase is currently available. Such information would be useful for a number of reasons, among them our ability to identify and measure local changes in the biomembrane bilayer when entities such as proteins are introduced. Our ignorance is due in part to the immense complexity of the fluid-phase bilayer resulting from the high degree of conformational disorder and in part from the paucity of effective methods for measuring transient chain conformation and packing. Among the experimental techniques currently available for determining chain conformation, vibrational spectroscopy, particularly infrared spectroscopy (IR), remains perhaps the most incisive. A number of informative IR studies on the conformation of the acyl chains in model bilayers have been reported. The results of these independent studies are clearly inter-related but have not been correlated. Similarly, the conformational measurements have not been related in a quantitative way to * To whom correspondence should be addressed. † University of California, Berkeley. ‡ University of Pennsylvania. § Rutgers University.

the increasingly reliable bilayer structures generated by molecular dynamics (MD) methods. To some extent, progress in this direction has been impeded by the absence of detailed quantitative descriptions of the chain conformations predicted by the MD simulations and by the lack of a common basis or format for a direct comparison of the spectroscopic and MD derived conformational statistics. This paper has two aims. One is to describe in relevant detail the conformation and packing of the chains in a fluid-phase bilayer generated from an MD simulation. The modeled bilayer we employ is the MD-simulated fully hydrated dipalmitoylphosphatidylcholine (DPPC) bilayer reported by Klein and coworkers for both the gel1 and fluid2 phases. Their simulation successfully reproduce a number of experimentally determined structural parameters observed for the DPPC bilayer, including, for the gel phase, the distribution of carbon atoms across the bilayer3 and the electron density profile,4,5 and, for the fluid phase, the surface area per lipid.6 We have characterized the conformation and packing of the chains in the simulated fluidphase bilayer and, to a lesser extent the gel phase, in a variety of ways: chain-length distributions of conformational sequences involving trans and gauche bonds, distributions CC bonds across the bilayer, and distributions of disorder in both the lateral and

10.1021/jp012145d CCC: $22.00 © 2002 American Chemical Society Published on Web 05/29/2002

6274 J. Phys. Chem. B, Vol. 106, No. 24, 2002 transverse directions. This subject constitutes roughly the first half of our paper and, in a sense, is complete in itself. Our second aim concerns the IR determination of bilayer conformation, which is expressed largely in terms of the concentrations of bonds pairs having tt, gt, and gg conformations. The methods used in these determinations are critically reviewed, and the interpretation of the spectral quantities measured are in some ways revised. The bilayer measured statistics are compared with both those derived experimentally for the n-alkanes and those predicted from the MD bilayer simulation. This paper is necessarily rather detailed. To help the reader, an outline follows. We begin with a brief account of the MD calculation. The concentrations and concentration ratios of both short conformational sequences and long trans-bond sequences found for the MD-simulated fluid-phase bilayer are compared with those for the liquidlike n-C16 system whose chains are unconstrained. The differences between the two systems are noted and discussed in terms of bilayer packing constraints. The distribution of CC bonds in the direction perpendicular to the bilayer plane is generated for both the MD-simulated gel and fluid phases. The distributions are interpreted in terms of longitudinal chain displacement and conformational disorder. We then discuss the unexpected presence of lateral order/ disorder heterogeneity in the MD-simulated fluid phase. Finally, localized structural correlations between the monolayers are characterized. The remaining part of this paper concerns the IR measurement of lipid chain conformation in the DPPC bilayer for the gel and fluid phases. Two IR methods, originally developed for n-alkanes and related systems and more recently applied to the DPPC bilayer, are described and critically reviewed. For the gel phase, the distribution of gauche bonds along the acyl chains, derived from earlier IR measurements, is interpreted as resulting from longitudinal displacements of the chains. For the fluid phase, the concentrations of the bond-pair conformations tt, gt, and gg were determined by combining results from the two IR methods. These were used to estimate the concentration of gauche bonds and the degree of segregation of the trans and gauche bonds along the chains. The large discrepancy between the experimentally determined value the gauche concentration and the MD-derived value is discussed. II. Organization of the Acyl Chains in a Fluid-Phase MD-Simulated DPPC Bilayer A. MD Simulation Method. The full account of the MD modeling of the DPPC bilayer presented by Tu et al.1,2 will be summarized here. The simulation box (repeat unit) consists of two lipid monolayers, each monolayer having 32 (4 × 8) DPPC molecules, giving a total of 64 lipids and 128 acyl chains. The initial configuration of the molecules in the bilayer was constructed in accordance with the X-ray-determined structure of related crystals. After a series of setup and equilibration stages, constant-pressure and -temperature MD simulations were performed on the lipid bilayer system at given external temperatures (19 °C for the gel phase and 50 °C for the fluid phase) and at zero pressure in a fully flexible simulation box using the hybrid algorithm developed by Martyna et al.7 Periodic boundary conditions were applied in three dimensions so as to generate an infinite, multilamellar system. The calculations were executed in parallel/vector mode on the Cray C90 supercomputer at the Pittsburgh Supercomputing Center using the CHARM program.8 The simulation time was 1070 ps for the gel phase and in excess of 1500 ps for the fluid. An MD simulation of

Snyder et al. liquid n-C16 was carried out under the same conditions. The distribution of dihedral angles associated with the CH2CH2 bonds show sharp maxima at values corresponding to trans, gauche, and gauche′ conformations. B. Unconstrained n-C16 Reference System. To determine the effects of packing constraints on the conformation of the chains in the fluid-phase bilayer, we have used an unconstrained reference system closely related to the bilayer, namely, liquid n-C16. The effect of the constraints were identified by comparing the conformational statistics for the MD-simulated bilayer with those of the unconstrained system. The conformer populations of n-C16 can be estimated using the rotational isomeric state model (RISM).9 A few general comments concerning the relative degree of conformational disorder associated with the gel and fluid phase bilayer and the liquid and crystalline n-alkanes are in order. The acyl chains in the gel-phase bilayer are highly ordered. Nearly all the CH2CH2 bonds are trans so that the chains are predominantly planar zigzag, except near the methyl chain-ends, where there is a significant increase in the concentration of gauche bonds. The gel phase is a little less ordered than crystalline n-alkanes, whose chains are essentially entirely alltrans. Conformational disorder in the fluid-phase bilayer is greater than that in the gel but less than that in a liquid-phase n-alkane. The use of the RIS model greatly simplifies the statistical mechanical calculations. The model assumes the CH2CH2 bonds of the PM chain have three stable conformations (t, g, and g′) that correspond to the minima in the CC bond rotation potential. It is also assumed that the conformations of the chains are not measurably affected by interchain interaction. Evidence supporting the latter assumption follows from a comparison of the room-temperature IR spectra of neat liquid n-C16sin a frequency region (1300-700 cm-1) that is highly sensitive to chain conformationswith the spectrum of a 20% room-temperature solution of n-C16 in CCl4. The two spectra are nearly identical. The small differences may be due to solvent effects on the IR band intensities. Our results suggest that conformer concentrations differ by less than 5% (unpublished measurements carried out at UC Berkeley) .Similar results derived from Raman spectra are reported in ref 10. The RISM derived concentrations of the interconverting conformers were obtained using a statistical weight for each chain based on the sum of the energies of the individual gauche bonds and gg′ pairs multiplied by the number of occurrences. All possible conformers of n-C16 were included without regard for chain self-overlap. A temperature of 50 °C is assumed here and, we note, elsewhere unless otherwise indicated. The RIS model we used involves two energy parameters (Eg and Egg′) whose values are given relative to that of the trans bond. The exact value of Eg, the gauche bond energy, is not known precisely. However, nearly all experimental determinations reported for n-alkanes fall in the range 0.50-0.65 kcal/ mol. The high-level ab initio calculation carried out by Smith and Jaffe11 indicates Eg for n-butane is 0.59 kcal/mol, which is in the middle range of the experimental values. Therefore, we have used 0.575 kcal/mol for Eg, with 0.50 and 0.65 kcal/mol taken as lower and upper limits. The second parameter, Egg′, is associated with the stearically unfavorable bond pair combination gg′. Its value is assumed to be 2.5 ( 0.3 kcal/mol.9 The RISM-calculated gauche-bond concentration Xg for an unconstrained n-C16 system is highly dependent on Eg. When comparing the conformational statistics of the MD-simulated bilayer to those for the unconstrained n-C16 system, we have

Dipalmitoylphosphatidylcholine Bilayers TABLE 1: List of RISM Calculations RISM calculation RM1 [0.575]a RM2 [0.81] RM3 [1.1]

Xg

Eg set to reproduce value of:

0.35 0.29 0.22

measured value of Xg for liquid n-C16 Xg for MD-simulated liquid n-C16 Xg for the MD-simulated fluid bilayer

a The bracketed quantity is the value of Eg in kcal/mol used in the RISM calculation; Egg is assumed to be 2.5 kcal/mol; T ) 50 °C.

set Xg for the n-C16 system to the bilayer value by adjusting Eg in the RISM calculation. The resulting n-C16 system is then referred to as “equivalent” to the MD system. Different situations require different n-C16 RISM equivalents. The three most frequently used equivalents are designated RM1[0.575 ( 0.075], RM2[0.81], and RM3[1.1], where the value of Eg (in kcal/mol) employed in the RISM calculation is bracketed. These systems are briefly described in Table 1. The following conventions are used. “Conformationally unconstrained” means that the conformations of the chains are not measurably affected by the conformation of the other chains. The chains in the fluid-phase bilayer are “constrained” in that their conformation is in some degree affected by other chains. An “ordered” chain consists of all-trans chains. An “alkyl chain” refers to the alkyl part of the acyl chain. The concentrations of trans and gauche bonds, or of specific conformational sequences of bonds, are expressed as mole fractions. It is important to note that the gauche bond concentration for liquid n-alkanes modeled by the MD technique is significantly smaller than the measured value. For MD-simulated liquid n-C16 (Table 1), for example, the MD generated value of Xg is 0.29, whereas the experimentally based value is in the range 0.33-0.38. This difference must be taken into account when MD derived conformational populations are compared to values obtained in other ways. The low MD-derived value of Eg results from the use of an all-atom potential. Such a force field is not well suited to reproduce the CC bond rotational potential, since the latter has its origin in the electronic structure of the carbon sp3 bonds. On the other hand, the concentrations of the shorter conformational sequences in the n-C16 system calculated using the MD and RISM methods are remarkably close, provided that the value of Xg for the n-C16 system is set (by adjusting Eg) equal to the MD value. The similarity in the conformational statistics of the two systems is evident in Table 2. However, the differences between the systems, which are small for short sequences, tend to increase dramatically with sequence length. C. Effect of Bilayer Constraints on Conformation and Packing. In going from the unconstrained n-C16 to the fluidphase bilayer, the conformational statistics change. The value of Xg decreases. For example, for unconstrained n-C16 (RM2[0. 81]) Xg is 0.29, whereas for the MD-simulated DPPC bilayer it is 0.225, about 25% lower. The concentrations of all the conformational sequences change, mostly the way expected from a decrease in the gauche bond concentration (Table 2). However, there are exceptions. These are readily identified by comparing the MD calculated conformational statistics for the bilayer to the RISM calculated (RM3[1.1]) values for the equivalent unconstrained n-C16 system (Tables 2). The anomalies occur because the geometry of certain conformational sequences is sterically unfavorable to parallel chain packing. For example, the concentration ratio Xgtg′/Xgtg for unconstrained n-C16 is unity because the conformational energies for gtg′ and gtg are the same. For the bilayer, this is not the case. The gtg′ sequence is much preferred because its presence in an otherwise all-trans chain does not change the direction of the chain, whereas this

J. Phys. Chem. B, Vol. 106, No. 24, 2002 6275 is not true for gtg. As a result, the value of Xgtg′/Xgtg in the bilayer is not 1.0, but rather 1.6. D. Distribution of Specific Conformational Sequences along the Chain. In this section, various chain-length distributions of trans and gauche bonds are used to characterize the MD-simulated fluid-phase bilayer. The differences between the bilayer distributions and those for the unconstrained, equivalent n-C16 system are discussed. 1. Single-Bond Distributions. We first consider the distribution of gauche bonds for the unconstrained n-C16 system that is the equivalent (calculation RM2[0.81]) of the MD-simulated bilayer. This distribution, shown in Figure 1, is essentially flat except at the chain ends, where there is a small increase. The increase occurs because the nearest-neighbor interaction energy for the penultimate gauche bonds is about one-half that for an interior gauche bond. This comes about because for a given interior gauche bond, which is designated as g, there are two bond pairs involving it, g′g and gg′, that are in essence energy-forbidden. As result, the overall gauche bond concentration of the system is reduced. However, for the terminal CH2CH2 bond (penultimate CC bond), there is only one energetically forbidden pair, egg′. (The CH2CH3 bond, which is designated e, acts like a trans bond.) Therefore, the terminal CH2CH2 gauche bonds are suppressed half as much as the interior ones. As a result, the gauche concentrations for the terminal bonds are somewhat higher than those for the interior bonds. The gauche bond distribution for the MD modeled fluid-phase bilayer, unlike that for unconstrained n-C16, is not flat, although it tends to be so in the central region, from bonds 5-8 (Figure 1). There is a shallow minimum at bond 7, after which there is a monotonic increase to the chain end. The gauche concentration for the last CH2CH2 bond (bond 13) is slightly enhanced by the end effect discussed above. At the acyl-group end of the chain, the gauche concentrations for bonds 1 and 2 are anomalously low (less than half the average value for the other bonds), probably due to steric effects associated with the ester group. The distributions for the SN1 and SN2 chains are nearly alike. (SN1 and SN2 refer to glycerol positions 1 and 2.) 2. Bond-Pair Distributions. The chain-length distributions of the bond pair sequences, tt, gt, and gg for the bilayer are of special interest because their concentrations are measurable by IR methods, as will be discussed below. Figure 2 shows the chain-length distributions of tt, gt, and gg bond pairs for the fluid-phase MD-simulated bilayer along with the corresponding distributions for the unconstrained equivalent n-C16 system. The distributions for the n-C16 system are flat, except for a small change at the chain ends due to end effects analogous to the one discussed above for single gauche bonds. The values of the pair concentrations in the flat region are marked on the farright side of Figure 2. For the fluid-phase bilayer, the MD calculated bond-pair distributions are not flat, although the deviations from the n-C16 values are relatively small. The shapes of the gt and gg distributions are similar to that for single gauche-bonds. The distribution of tt pairs complements that for the gt pair since, to a good approximation, Xtt + Xgt ) 1. 3. Trans-Bond Sequence Distributions. By far the most outstanding differences in sequence concentrations between the MD-simulated bilayer and the n-C16 system occur for the alltrans sequences. The differences between the two systems become much greater as the length of the trans sequence increases. This trend is clearly evident in Figure 3 even for the shortest sequences. In this figure, the ratio (R(mt), bilayer to n-C16) of the concentrations of the sequences made up of mt

6276 J. Phys. Chem. B, Vol. 106, No. 24, 2002

Snyder et al.

TABLE 2: MD and RISM Calculated Average Concentration of Sequences for Unconstrained n-C16 and the Fluid-Phase DPPC Bilayer liquid n-C16 seq

RISMa lower

RISMa upper

t g tt gt gg gg′ ttt tgt ttg tgg gtg gtg′ tttt tttg ttgt ttgg tgtg tg′tg gttg gttg′

0.624 0.376 0.360 0.536 0.100 0.004 0.208 0.191 0.310 0.142 0.058 0.058 0.120 0.179 0.220 0.082 0.082 0.082 0.033 0.033

0.669 0.331 0.418 0.506 0.077 0.003 0.261 0.190 0.317 0.115 0.048 0.048 0.164 0.199 0.238 0.072 0.072 0.072 0.030 0.030

DPPC bilayer (fluid phase)

diff

MD simul

RISMb (equiv)

+0.045 -0.045 +0.058 -0.030 -0.023 -0.001 +0.053 -0.001 +0.007 -0.026 -0.010 -0.010 +0.044 +0.020 +0.018 -0.010 -0.010 -0.010 -0.003 -0.003

0.708 0.292 0.486 0.454 0.053 0.008 0.339 0.182 0.297 0.078 0.037 0.038 0.247 0.212 0.228 0.047 0.050 0.069 0.014 0.027

0.708 0.292 0.478 0.464 0.056 0.002 0.324 0.183 0.314 0.089 0.038 0.038 0.219 0.213 0.248 0.060 0.060 0.060 0.026 0.026

diff

MD simul

RISMc equiv

diff

0.000 0.000 +0.008 -0.010 -0.003 +0.006 +0.015 -0.001 -0.017 -0.011 -0.001 0.000 +0.028 -0.001 -0.020 -0.013 -0.010 +0.009 -0.012 +0.001

0.775 0.225 0.595 0.368 0.032 0.005 0.473 0.153 0.238 0.054 0.026 0.038 0.383 0.134 0.169 0.041 0.025 0.066 0.013 0.011

0.776 0.224 0.589 0.379 0.030 0.001 0.447 0.161 0.288 0.052 0.023 0.023 0.339 0.218 0.244 0.039 0.039 0.039 0.018 0.018

-0.001 +0.001 +0.006 -0.011 +0.002 +0.004 +0.026 -0.008 -0.050 +0.002 +0.003 +0.013 +0.044 -0.084 -0.074 +0.002 -0.014 +0.073 -0.005 -0.007

a Sequence concentrations based on rotational isometric state model (RISM). Lower and upper limits cover the acceptable range of experimentally determined values of Eg: lower RM1[0.50]; upper RM1[0.65] (see text). b RISM calculation RM2[0.81], which yields the MD value for Xg. c Same as footnote b, except this uses calculation RM3[1.1].

Figure 1. Distribution of single gauche bonds along the alkyl chains in the MD-simulated fluid-phase DPPC bilayer and in the equivalent unconstrained n-C16 system, RM3[1.1]. The concentration of gauche bonds for CC bond number i is plotted against i (i ) 1, 2, ..., 13).

Figure 2. Distribution of bond pairs tt, gt, and gg along the alkyl chains of the MD-simulated fluid-phase DPPC bilayer and of the equivalent, unconstrained n-C16 system, RM3[1.1]. Bond pair number i goes from 1 to 12.

trans bonds (mt ) 1, 2, 3, and 4) are plotted as a function of the position of the sequence along the chain. The sequence position refers to the position of the center of the sequence. The CC bonds are numbered beginning at the ester-group end. The plots in Figure 3, which have been smoothed, show that R(mt) exceeds unity at nearly all positions along the chain. Thus, the populations of trans-bond sequences in the bilayer generally exceed those for the equivalent unconstrained n-C16 system. There are two maxima. The one near the ester end of the chain is probably associated with steric effects involving the head and ester groups. A second maximum, located very near the center of the alkyl chain between bonds 6 and 7, marks the region of highest conformational order. The region of least order is at the interfacial region between monolayers, that is, in the region

where the terminal methyl groups of the chains are most concentrated. Here, the conformational order is about the same as that for the unconstrained n-C16 system. E. Other Distributions. 1. Gauche Bond Distribution among Chains. Conformational differences between the fluid-phase bilayer and unconstrained n-C16 systems are also manifest in the manner in which the gauche bonds are distributed among (not along) the chains constituting the system. This can be seen in Figure 4, in which the fraction X(mg) of alkyl chains that contain exactly mg gauche bonds is plotted against mg, the number of gauche bonds. The four distributions in the upper half of the figure are for unconstrained n-C16 systems. Distribution a is derived from conformational statistics generated using the RISM method with Eg equal to 0.50 kcal/mol, the lowest

Dipalmitoylphosphatidylcholine Bilayers

J. Phys. Chem. B, Vol. 106, No. 24, 2002 6277

Figure 5. Side view (yz projection) of the hydrocarbon chains in the MD-modeled gel-phase DPPC bilayer. The centers of CH2-CH2 bonds are designated ] and 4, depending on whether the bonds are trans or gauche.

Figure 3. Concentrations of trans-bond sequences as a function of the position i of the sequence along the chains compared for the MDsimulated fluid-phase DPPC bilayer and the unconstrained equivalent n-C16 system (RM3[1.1]). The sequence concentration ratios R(mt) (bilayer to n-C16) are plotted against the sequence position (sequence center) i along the chain. mt is the number of trans bonds in the sequence. Distributions are shown for mt equal to 1, 2, 3, and 4.

Figure 4. Distribution of gauche bonds among the alkyl chains for the MD-simulated fluid-phase DPPC bilayer and for the unconstrained equivalent n-C16 system. The distributions are obtained by plotting X(mg), the concentration of chains having exactly mg gauche bonds, against mg.

acceptable value. This calculation gives Xg ) 0.38. Distribution b is like a, except Eg ) 0.65 kcal/mol, the highest acceptable value of Eg, in which case Xg ) 0.33. Distributions c and d also represent unconstrained n-C16 systems, both with Xg equal to 0.29. Distribution c was obtained using the RISM model (calculation RM2[0.8]). Distribution d was generated using the MD technique.

The dependence of the distribution on Xg for unconstrained systems is indicated in a comparison of distributions a and b, for which Xg equals 0.38 and 0.33, respectively. The maximum of distribution b occurs at a lower value of mg, as would be expected. The shapes of the two distributions are, however, essentially the same. A comparison of the shapes of the two distributions c and d for which Xg is the same (0.29) is of interest because distribution c was RISM generated whereas d was MD generated. The two distributions nearly coincide, affirming our earlier observation in section II.B that the MD and RISM simulations yield essentially the same conformational statistics, if Xg is the same in both cases. Finally, we found the shape of the distribution for the fluidphase bilayer to be substantially different from the equivalent unconstrained C16 system, as is apparent in the lower half of Figure 4. The bilayer distribution is flatter; that is, the populations of bilayer chains with the fewest gauche bonds and those with the most gauche bonds have increased at the expense of chains with an intermediate number. This change may be thought of as a manifestation of trans/bond segregation in the bilayer, a topic we will return to later. 2. Distribution of Acyl Chain CC Bonds across the Bilayer. This distribution is reflects bilayer randomness associated with longitudinal chain displacement in the gel phase and conformational disorder in the fluid phase. After we describe how the distributions were obtained, structural implications derived from the CC bond peak positions are considered, first for the gel phase whose MD-simulated structure is shown in Figure 5 and then for the fluid phase. Then the peak-width distributions are similarly considered. Distributions for individual CC bonds consist of plots of Xik versus k, where Xik is the fraction of the CC bonds with chainlength position number i that are located in the kth lateral slice of the bilayer. The CC bonds are numbered i ) 1, 2, 3... 15, beginning with the ester (O)C-CH2 bond. The lateral slices are cut parallel to the bilayer plane and are therefore perpendicular to the z axis. They are normally 1.0 Å thick and are numbered beginning at the headgroup side of the upper monolayer. (The two monolayers are distinguished by tagging them individually as “upper” and “lower”.) The CC bond distributions for the gel and fluid phases are shown in Figures6 and 7, where the axis system in Figure 5 is used. a. Peak Positions in the CC Bond Distribution. First the gel phase. Figure 5 shows a side view (yz) projection of the MD modeled DPPC gel reported by Tu et al.1 The position of the center of each CH2CH2 bond is indicated by ) or 4, depending on whether the CC bond is trans or gauche. We note that the sign of the chain tilt-angle reverses between adjacent mono-


Figure 6. Bilayer-normal (z axis direction) distribution of acyl chain C-C bonds across the MD-simulated gel-phase DPPC bilayer. The SN1 and SN2 chains are considered separately.

layers, whereas it is known from X-ray diffraction studies that the sign is the same. This difference would not be expected to affect chain conformation. The CC bond population profiles across the gel-phase bilayer are shown in Figure 6 for each bond of the SN1 and SN2 chains. As noted earlier by Tu et al.,1 the longitudinal positioning of the SN1 and SN2 chains differs, the ester-group end of the SN2 chain being somewhat closer to the bilayer surface than that of the SN1 chain. As a result, the distribution peaks for the CC bonds of the SN1 chains are slightly displaced from those of the SN2 chains. In going down the chain, the shapes of the distribution peaks do not change much because the gel phase

Snyder et al. is uniformly ordered, except that, as the methyl end of the chain is approached, there is a small increase in distribution peak width due to an increase in the number of gauche bonds. The corresponding distributions for the fluid phase are shown in Figure 7. Relative to the gel, the peaks for the fluid are broader, and consequently, the distinction between the SN1 and SN2 chains is blurred. Also, in contrast to the gel phase, the fluid-phase peaks tend to broaden in going down the chain. Figure 8A displays, for each CC bond for both the gel and fluid-phase bilayer, the z axis peak positions. In the gel-phase, the peak position selected is that of the dominant peak for the SN2 chain because these peaks are better defined than those for the SN1 chains. The use of SN1 peaks would result in a uniform shift of the data points toward the mid-plane, leaving the slope unchanged. The fluid-phase peak positions were obtained by weight-averaging the subpeaks that make up the main peak. For the fluid phase, all bond positions are plotted, but for the gel only the even numbered bonds. The thickness of the gel-phase hydrocarbon layersdefined as the distance between the center of the first bond in the upper monolayer and the corresponding bond in the lower layersis about 38 Å. If the chains were not tilted, the value would be about 46 Å. For the fluid, this distance is 28 Å. A uniform CC bond-density for both phases is indicated, since the plots are essentially linear for both the gel and fluid. The gel-phase plot shows the expected gap at the bilayer midpoint (monolayer interface). No corresponding gap appears for the fluid because the distribution peaks associated with terminal or near terminal CC bonds of the opposing monolayers are overlapped. This suggests a degree of interlayer interdigitation, which, if present, would account for some of the decrease in bilayer thickness observed in going through the gel-to-fluid transition. b. Peak Widths in the CC Bond Distribution. Figure 8B displays for both phases the widths of the CC bond distribution peaks for each monolayer. Base widths are used as a relative measure of peak widths, because the peaks appear to have a triangular shape, so base widths can be measured more accurately than halfwidths. In Figure 8B, the base widths equal the lengths of the vertical line-segments. The peak widths for the gel phase remain nearly constant in going along the chain. This consistency implies the widths are determined by variations in the longitudinal positions of the chains, since such variations affect the positions of all the CC bonds of a given chain in the same way. It is relevant that longitudinal chain displacement plays an important role in determining the distribution of gauche bonds in the gel-phase bilayer, a subject discussed in section IV. The situation for the fluid phase is different. The distribution peak widths for the fluid are considerably larger than those for the gel and tend to increase in going down the chain. These differences from the gel can be attributed to a much greater degree conformational disorder in the fluid phase. Unlike longitudinal displacement, conformational disorder tends to have a cumulative effects so that widths increase in going down the chain. As an aside, we note that account must be taken of the bilayer normal distribution of CC bonds when interpreting spectral data obtained from localized probes placed at a specific sites along the chain. An example of a localized probe is the CD2 methylene group which is commonly used in both NMR and IR spectroscopy to determine the dynamics and conformation of polymethylene chains at specific methylene positions. Because the position of a given methylene is spread along the bilayer normal,



Figure 7. Bilayer-normal (z axis direction) distribution of acyl chain C-C bonds across the MD-simulated fluid-phase DPPC bilayer.

the probe position is not well-defined, especially in the fluid phase. We estimate that a CD2 group samples a region of the fluid-phase bilayer that extends over about 5 Å along z, a distance that spans at least 3 static CC bonds. F. Fluid-Phase Order-Disorder Heterogeneity. 1. Snapshots of the Simulated Bilayer. The presence of heterogeneity in the conformation and packing of the hydrocarbon chains in the fluid phase is apparent in Figure 9. This unexpected phenomenon is manifest in the form of small, more or less welldefined domains of high conformational order. Figure 9 displays a molecular scale view (xy projection) of the bilayer that shows the alkyl chains in the upper and lower monolayers of the MDsimulated fluid-phase DPPC at a simulation time of 1500 ps. The positions of the centers of the CC bonds are indicated. To

give an overview, most of the nine repeat units (simulation boxes) are included. The central repeat unit is outlined. Figure 10 is an expanded scale view of Figure 9 that covers just one repeat unit. The trans and gauche bonds are indicated by the symbols ] and 4 respectively. The most highly ordered chains are located in the upper monolayer near the center of the repeat unit. In the ordered domains, a few of the chains are all-trans, and the rest are nearly so. They pack with an average chain-tilt angle near that found for the MD-simulated gel phase. Chain ordering in the lower monolayer is less pronounced as may be seen in the side-view (xz) projection shown in Figure 11. To simplify this projection, the repeat unit has been divided in two by an xy plane. The “front” half (Figure 11A) shows highly ordered regions in the upper layer that correspond to


Snyder et al.

Figure 8. Part A: Locations along the bilayer normal (z axis) of the concentration maxima for each CC bond (bonds are numbered i ) 1, 2, ..., 15) along the acyl chains for the MD-simulated DPPC bilayer in its gel and fluid phases. Part B: Locations along the z axis of the ends of the line segment whose length equals the basewidth of the distribution peak associated with bond i, plotted against i.

the ordered regions in the xy projection in Figure 10. The back half (Figure 11B) displays no regions of exceptionally high order in either the upper or lower monolayer. To distinguish the ordered and disordered regions in a more quantitative way, we have constructed order/disorder (O/D) maps that indicate the approximate location of each chain and provide some information about its conformation. These maps, for example the one shown in Figure 12, are xy projections of the monolayer repeat unit, upon which the approximate position of each ester group is marked. The marker consists of a pair of integers, rt/ng. The first (rt) is the number of trans bonds making up the longest trans-bond sequence in this chain. The second integer (ng) is the total number of gauche bonds in the chain. The more ordered regions are enclosed by contour lines drawn between those chains with rt g 8 and those with rt < 8. The boundaries of the ordered regions are imprecise because the cutoff value of 8 trans bonds is (necessarily) arbitrary and because the trajectories of the chains are not known. Nevertheless, the ordered regions designated in the O/D maps (Figure 12) have shapes that correspond quite well to those found in xy projections (Figures 10 and 11). The boundaries between the most highly ordered regions and the surrounding disordered regions tend to be sharp, and in that way the ordered regions are domain-like. A good example is the previously discussed highly ordered microdomain in Figure 10, which also appears in the O/D map in Figure 12. In crossing the boundary between the highly ordered region and the adjoining region, the average number (rt) of trans bonds in the longest trans sequences drops from about 12 to 5. The relation between simulation time and chain ordering was explored in a limited way. Monolayer O/D maps at simulation

Figure 9. Part A: Top view (xy projection) showing the acyl chain CH2-CH2 bond centers for the upper monolayer of the MD-simulated fluid-phase DPPC after a simulation time of 1500 ps. The rectangular box in the center represents a single repeat unit. Part B: Same as Part A, but for the lower monolayer.



Figure 11. Side view (xz projection) of the acyl chain CC bond centers for the repeat unit shown in Figure 10. The repeat unit is divided into halves A and B, which extend from y ) 0 to +20 Å and from y ) 0 to about -20 Å, respectively. Figure 10. Projection of the acyl chain CC bond centers onto the xy plane of the fluid-phase bilayer. Same as Figure 9, except only one repeat unit is shown. The trans and gauche bonds are marked ] and 4, respectively.

times of 0, 500, 1000, and 1500 ps are shown in Figure 13. (As noted earlier, “zero simulation time” occurs after a series of setup and equilibrium steps have been performed.2) After 500 ps, the shapes of highly ordered regions bear little resemblance to the initial shapes. It appears that the total fraction of the bilayer identified as ordered in the O/D maps does not change much with time, although it is possible that these regions become slightly more compact. Changes in the shape of the ordered regions monitored at shorter time intervals of 10 ps, between 1460 and 1500 ps, were also examined. During this interval, the shapes change significantly but remain recognizable. 2. EVidence for Heterogeneity from Anomalies in Trans Sequence Populations. Evidence for order/disorder heterogeneity in the fluid phase, in addition to that found in the O/D maps, is found in anomalies in the distribution of trans bond sequences obtained by plotting the relative frequency of occurrence of alltrans sequences against sequence length. Such distributions are displayed in Figure 14 for n-C16 and bilayer systems. In this figure, X(rt) is the fraction of chains whose longest all-trans sequence consists of exactly rt trans bonds. For reference, it is important to know the dependence of the shape of these distributions on the gauche concentration Xg. To this end, time-averaged distributions a and b are shown in Figure 14A for two n-C16 systems with gauche bond concentrations of 0.35 and 0.22, respectively. Due to the time averaging, anomalies associated with heterogeneity would not show in these distributions. As expected, the peak of the distribution associated with the system with the higher value of Xg occurs at a shorter sequence

length, and the population of long all-trans sequences is lower. The distribution for the MD-simulated fluid bilayer c is also displayed. The difference between this distribution and that of the equivalent n-C16 system b (also with Xg ) 0.22) is accountable to the bilayer packing constraints. The main difference is that the concentration of long trans segments is higher for the bilayer. The high concentrations of long transsequences in the chains making up the ordered domains do not stand out in the bilayer distributions in Figure 14A because these chains represent a relatively small fraction of the bilayer chains. (We note there is a minimum in the populations of long trans chains in the distributions in Figure 14A for low values of Xg. The minimum occurs because, at low gauche concentrations, the formation of the longest trans-sequences requires the gauche bonds to be located at or near the ends of the chains. This requirement together with the substantially higher energy associated with gauche bonds relative to trans bonds can, for low Xg, reduce the populations of the longest trans-sequences to values below that of the all-trans chain.) The four bilayer distributions in Figure 14B were measured at simulation times 0, 500, 1000, and 1500 ps with no time averaging. Heterogeneity is evident in the form of anomalous spikes in the long trans-sequence region. Most of the spikes are associated with fluctuations attributable to the small sample size of 64 chains per repeat unit. However, some are of a magnitude far outside the fluctuation range. For example, the anomalously large spike at rt ) 13 for the upper monolayer in Figure 14B that is present at a simulation time of 1500 ps is associated with the extended chains located in the previously identified highly ordered region shown in the bilayer projections in Figures 9-11. G. Correlations between Monolayers. Structural correlations between the upper and lower monolayers were found for the


Figure 12. Order/disorder maps (xy projections) for the simulated fluidphase DPPC bilayer at 1500 ps. The integer combinations rt/ng, which are placed at the position of the attached ester group, pertain to the chain associated with the ester group: the integer rt is the number of trans bonds in the longest trans run in the chain; ng is the total number of gauche bonds in the chain. The contour lines enclose chains with rt g 8.

fluid-phase MD-simulated bilayer. In the most ordered regions, the tilt angle tends to change sign in going from one monolayer to the other (Figure 11) and thus behaves in the same way as the tilt angle in the MD-simulated gel phase (Figure 5). A second kind of correlation is found in the tendency for an exceptionally ordered region of one monolayer to face an exceptionally disordered region in the opposing monolayer (Figure 11A). III. IR Methods for the Determination of Alkyl Chain Conformation Many techniques, among them NMR, ESR, Raman, and IR spectroscopies, as well as calorimetry and dilatometry, have been used to detect and estimate conformational disorder in model biomembranes. However, quantitative determinations have proved extremely difficult. Of the techniques, vibrational spectroscopy is uniquely well suited to measure conformation because the time scale associated with the method is sufficiently fast to distinguish conformers and because certain bands, mostly IR bands, represent well-defined conformational sequences. For the quantitative measurement of conformation, infrared spectroscopy is generally better suited than Raman spectroscopy because infrared bands tend to be associated with vibrations that are spatially localized. Spatial localization is essential to ensure that the intensity of a conformationally sensitive band

Snyder et al.

Figure 13. Order/disorder maps for the MD-simulated fluid-phase bilayer at simulation times 0, 500, 1000, and 1500 ps.

is linearly related to the concentration of the conformational sequence it represents.13-15 The earliest quantitative estimates of conformational disorder in lipid bilayers were measured using Raman spectroscopy. As a result, the earliest estimates of disorder tended to be on the high side because bands representing highly delocalized CC stretching modes were used.16 There exist two IR methods suitable for the quantitative measurement of specific conformational sequences. These methods, originally developed to study disordered n-alkane systems,17-20 are described below. A. Wag-Band Method. This method employs certain IR methylene wagging bands that are associated with bond pairs eg and gg.17 The concentrations of eg and gg sequences can be estimated from the intensities of the wagging bands at 1341 and 1354 cm-1, respectively, since the vibrations represented by these bands are localized within the sequences. For calibration, the spectra of liquid n-alkanes are used since estimates of the concentrations of eg and gg for these systems are available from the RISM model. The methyl umbrella (U) band near 1375 cm-1 is used to normalize the measured intensities to a perchain quantity. This band is well suited as an internal standard because its intensity is derived almost entirely from the methyl group at the end of the alkyl chain and is therefore essentially independent of chain conformation.



Figure 14. Fraction X(rt) of alkyl chains whose longest all-trans sequence consists of exactly rt trans bonds, plotted against rt. Part A: (a) Liquid n-C16 with Xg equal to 0.35, RM1[0.575]; (b) unconstrained n-C16 with Xg equal to 0.22, RM3[1.1]; (c) the time-averaged MD-simulated bilayer with Xg equal to 0.22. Part B: Same as Part A but for the simulated fluid-phase bilayer at simulation times 0, 500, 1000, and 1500 ps.

Values for Xeg and Xgg can be obtained from the equation

XζB ) XζA(IζB/IζA)

(1)

where ζ represents eg or gg and where A and B refer to the n-alkane and the fluid bilayer, respectively. The quantities IζA and IζB are measured intensities. The data used and the results obtained are summarized in Table 3. Two values are listed there for both Xeg,A and Xgg,A. These RISM derived concentrations were calculated using the lower-limit value of Eg (0.50 kcal/ mol) and once again using the upper-limit value (0.65 kcal/ mol), so as to cover the acceptable range of Eg. B. CD2-Band Method. This method exploits the conformational sensitivity of the frequency of the IR active methylene

rocking vibration of an isotopically isolated CD2 group. The frequency of the CD2 rocking vibration depends on the conformation of the pair of CC bonds adjoining the CD2 group and, to a lesser but still measurable extent, on the conformation of the next adjoining pair of CC bonds. The CD2 rocking vibrations are highly localized, and therefore, their intensities are good measures of the concentrations of the pair of CC bonds flanking the CD2 group. The method is site specific since the CD2 can be placed at any specific position along the chain. We have previously used this method to determine the distribution of gauche bonds along the chains in crystalline n-alkanes.18-20 If the two CC bonds adjoining the CD2 group are both trans (tt), the frequency of the CD2 rocking band is about 620 cm-1 (Figure 15). If one bond is gauche and the other trans (gt), the


Snyder et al.

TABLE 3: Data and Results Associated with the IR Determination of Xgg and Xeg for the Fluid-Phase DPPC Bilayer frequencies intensities

concentrations

observed (cm-1) Iζ (normalized measureda n-C16 (liquid) DPPC bilayer Iζ,DPPC/Iζ,C16 (measured) Iζ,C16 (calculated)b Eg ) 0.50 kcal/mol ) 0.65 kcal/mol Xζ,DPPC (determined)c Eg ) 0.50 kcal/mol ) 0.65 kcal/mol av

gg

eg

1354

1341

0.771 0.377 0.49

0.356 0.271 0.76

0.100 0.077

0.424 0.376

0.049 0.037 0.043

0.325 0.287 0.31

a ζ represents gg or eg. b Two values of Eg, representing upper and lower acceptable limits, were used to calculate (using the RIS model) the n-C16 reference system values. c Based on the measured value of Iζ,DPPC/Iζ,C16. See eq 1 in text.

Figure 15. Infrared spectra showing the CD2 rocking bands associated with an isotopically isolated CD2 group in a polymethylene chain. Part A: Fluid-phase DPPC bilayer with a CD2 group at carbon number 6; liquid n-C21H42D2-4,4-d2 (redrawn from refs 12 and 20, respectively). Part B: Liquid n-C7H14D2-4,4-d2 at 30 and -70 °C.

frequency is about 650 cm-1. There is also a band at a higher frequency, probably between 655 and 670 cm-1, that is associated with the gauche-gauche (gg) pair. We have observed this band at 665 cm-1 in the spectrum of crystalline cyclohexadecane-1,1-d2.21 However, this band is seldom in evidence in the spectra of disordered PM chains because it is relatively weak (the concentration of gg pairs is low) and tends to be broad due to the high sensitivity of its frequency to changes in the CC

torsional angles associated with the two adjoining gauche CC bonds.18 The 650 cm-1 band can therefore be assumed to consist of both the gg and gt bands. In addition, there is an extraneous band near 660 cm-1 associated with a CCD bending mode of the CHD methylene groups that are present in CD2 substituted samples. Usually this band is weak and can be neglected. The relation between the concentrations of tt, gt, and gg pairs and the intensities of the bands that represent them is

(Igt + Igg)/Itt ) (kXgt + k′/Xgg)/Xtt

(2)

where k and k′ are constants whose values have been shown to be near one.20,21 Since both the gt and gg pairs are represented by the 650 cm-1 band, (Igt + Igg)/Itt ) I650/I620. We will designate the intensity ratio I650/I620 as R. In summary then

R ) I650/I620 ) (Igt + Igg)/Itt ) (Xgt + Xgg)/Xtt

(3)

In view of the difficulties that plague quantitative determinations of chain conformation, it is appropriate to comment on the reliability of the CD2-band method. Perhaps the most incisive evidence for its accuracy rests on our analysis of the IR spectrum of a sample of crystalline cyclohexadecane with each ring having one CD2 group.21 The conformation of cyclohexadecane in the crystalline phase is ggttg′g′ttggttg′g′tt, which defines a squareshaped ring with alternating gg and g′g′ pairs at the corners. The numbers of tt, gt, and gg pairs per ring are therefore 4, 8, and 4, or, expressed as ratios, 1:2:1. The IR spectrum of this sample exhibits three well-separated CD2 rocking bands with frequencies near 618, 645, and 665 cm-1 that represent tt, gt, and gg pairs. The observed intensity ratios of these bands are 1:2:1 to within about 10%. Therefore, the intrinsic intensities associated with the 3 different bond pairs are nearly equal, that is, k ) k′ ) 1.0 ( 0.1. A bond moment calculation yields similar values for k and k′.20 There is one more check available. To account for the value of I650/I620 observed for CD2 substituted n-alkanes, it is necessary to assume that the intrinsic intensities of the tt, gt, and gg bands are nearly equal.18,20 The equivalent of this check is found in the remarkably close agreement between the value of R measured for liquid n-C21H42-4,4-d2 and the RISM value calculated of R for liquid n-C16 obtained using the experimentally measured value of 0.575 kcal/mol for Eg. Both values of R are 1.6. The accuracy of the CD2-band method must be due in large part due to the similarity of the normal coordinates of the CD2 rocking modes associated with the tt, gt, and gg pairs and to the common environment these bond pairs share. As mentioned earlier, the CD2 rocking band intensities can be used to estimate the concentration ratio, Xtgt/Xggt.18,20,22 The bands at 652 and 646 cm-1 representing the tgt and ggt sequences make up the 650 cm-1 band that represents the gt pair. If this gt pair is preceded by a trans bond, the CD2 rocking band associated with the resulting t(gt) sequence appears at 652 cm-1. If preceded by a gauche bond, the CD2 band associated with g(gt) appears at 646 cm-1. The intensity ratio I652/I646 is then proportional to Xtgt/Xggt. The tgt and ggt component bands appear in the room-temperature spectra of both liquid n-C21H42D24,4-d2 and the fluid-phase CD2 substituted bilayer, represented here by DPPC-6,6-d2. Both spectra are shown in Figure 15A. The assignments of 652 and 646 cm-1 bands to tgt and ggt sequences have been established from normal coordinate calculations.18,20 This assignment is consistent with the temperature behavior of the I652/I646 ratio for n-C7H14-4,4-d2, as shown in Figure 15B.23 The ratio increases as the sample is


J. Phys. Chem. B, Vol. 106, No. 24, 2002 6285 TABLE 4: Measured Values of the Intensity Ratio I650/I620 for Fluid-Phase DPPC Bilayer at 50 °C at 4 Different (CD2) Positions along the Chains carbon numbera 4 6 10 13

bond pair numberb 2 4 8 11 ave

I650/I620

Xgc

d

0.125 0.183 0.120 0.106 0.14 ( 0.04

0.261 0.477d 0.245d 0.205e 0.30 ( 0.08

a The carbon numbered one is the ester carbon. b The bond pair numbering is that used in Figure 2. c Calculated from eq 8, with Xgg equal to 0.043, which was determined by the wag-band method. See text. d From ref 12. e From ref 25.

disorder,19 we conclude that the bilayer distribution is also determined by longitudinal displacement. This explanation is supported by the fact that the MD-simulated gel-phase bilayer structure exhibits longitudinal displacement disorder. V. Fluid-Phase Bilayer Chain Conformation Figure 16. Distributions of gauche bonds along the chains of DPPC gel-phase bilayer and hexagonal n-C21H44. Data points 0 are MD calculated values for the gel-phase bilayer. Points 9 are the experimental values for the gel-phase bilayer at 35 °C.12 Filled circles (b) are experimental values for crystalline hexagonal n-C21 at 39°.19

cooled, in keeping with the fact that the tgt band at 652 cm-1 is the more stable, since it has only one gauche bond. IV. Distribution of Gauche Bonds in the Gel-Phase Bilayer Comparison of the measured chain-length distribution of gauche bonds for the gel-phase bilayer with the distributions for the MD-simulated gel-phase bilayer and for crystalline n-C21H44 in its hexagonal phase19 indicates that gel distribution is determined by random longitudinal displacements of the chains. Figure 16 shows gauche bond distributions: (i) for the MD modeled gel-phase bilayer, (ii) for the bilayer determined from CD2 rocking measurements, and (iii) for crystalline n-C21H44 in its hexagonal phase, likewise determined by the CD2 method. (The hexagonal phase of an n-alkane is a hightemperature phase that is significantly more conformationally disordered than the low-temperature phase. For n-C21, the transition from the low-temperature phase to the hexagonal phase occurs about 10 °C below the melting point.) To take into account the differences in the chain length between the systems, we have employed a “reduced length” defined as (m - 1)/(M - 1), where M is the number of CH2CH2 bonds and m is the CC bond position along the chain. The first and last bonds in the chain are then fixed at 0 and 1, independent of chain length. We have used the MD calculated gauche bond distribution for the gel phase reported by Tu et al.1 The three experimentally determined values of Xg provided by the CD2 measurements of Mendelsohn et al.12 were used for the gel. The values of Xg for hexagonal n-C21H44 are those of Maroncelli et al.19 found using the CD2 method. The gauche-bond distribution for the MD-simulated gel-phase bilayer has a shape very similar to the experimentally determined distributions for both the bilayer and the hexagonal phase of n-C21. (The differences between the bilayer and n-alkane are mostly due to headgroup differences.) It seems likely therefore that the distributions have a common origin. Since the observed gauche bond distribution for hexagonal n-C21 has been quantitatively accounted for in terms of longitudinal-displacement

A. IR Measurements in General. Casal and McElhaney24 and Senak et al.22 used the wag method to estimate the concentrations of eg and gg bond pairs for the DPPC fluidphase bilayer. We have remeasured the wag band intensities of a series of liquid n-alkanes used to calibrate the method. Our values for Xeg and Xgg derived using the new n-alkane data are very near those reported earlier. The parameters used in these determinations are listed in Table 3. The CD2 method was employed by Mendelsohn and coworkers to measure the intensity ratio I650/I620 for four DPPC fluid-phase bilayer samples, each with a CD2 group at a specified position along the alkyl chains.12,25,26 These data, listed in Table 4, were used by the authors to estimate the gauche bond concentrations Xg at the various position. Unfortunately, the equation used, which was meant to be the one we derived in ref 18, was mistranscribed with the consequence that the values reported for Xg in refs 12, 25, and 26 are too high by roughly a factor of 2. The CD2 measurements they reported are reinterpreted another way in the next section. We have determined the concentration ratio Xtgt/Xggt for the bilayer using the measured value of the intensity ratio (I652/ I646)B and the equation

(Xtgt/Xggt)B ) (I652/I646)B(Xtgt/Xggt)A(I646/I652)A

(4)

where B refers to the bilayer and A to liquid n-C7H14-4,4-d2. The value of (Xtgt/Xggt)A was estimated from RISM calculation RM1[0.575] and that for I652/I646)A from IR measurements on n-C7H14-4,4-d2. The accuracy of the intensity ratio for the bilayer may be compromised in some degree by the uncertainty in the separation of the overlapping 652 and 646 cm-1 bands. B. Determination of the Gauche Bond Concentrations Xg and Xeg. 1. Interior Bonds (Xg). At present, there is no experimental technique capable of determining Xg in a direct way. However, Xg can be evaluated from the measured values of the bond-pair concentrations Xtt, Xgt, and Xgg. An equation relating Xg to these experimentally determined quantities can be derived by combining the following three equations. The first corresponds to eq 3:

R ) (Xgt + Xgg)/Xtt

(5)

The second equation is

Xtt + Xgt + Xgg ) 1

(6)


Snyder et al.

TABLE 5: Experimental and Calculated Concentrations and Concentration Ratios of Relevant Conformational Sequences for the MD-Simulated Fluid-Phase Bilayer and Unconstrained n-C16 DPPC fluid phase bilayer (50 °C) MD-simulated measured directly Xg Xeg Xgg Xtt Xgt Rk Xtgt/Xggt

0.31 ( 0.02a 0.043 ( 0.006a 0.30 ( 0.08b 2.1 ( 0.6b

n-C16 equiv of the MD simulation

liquid n-C16 (50 °C)

derived from meas

corrf

(uncorr)

corrg

(uncorr)h

RISM calci

MD calcj

0.14 ( 0.04c 0.28 ( 0.02d

0.27 0.36 0.054 0.46 0.43 1.10 1.8

(0.22) (0.30) (0.032) (0.59) (0.37) (0.063) (2.8)

0.27

(0.22)

0.29

0.047 0.51 0.44 0.86

(0.030) (0.59) (0.38) (0.69)

0.35 0.40 0.090 0.38 0.52 1.61l 1.50

0.77 ( 0.04e 0.19 ( 0.04e

0.053 0.49 0.45 0.92 2.3

a Measured using the wag-band method. b Measured using the CD -band method. c From eq 8. d Corrected for “end effect” (see text). e From the 2 equation Xtt ) (1 + R-1)-1 - Xgg. f MD calculated values corrected for the underestimation of Xg. Corrections made according to eq 9 (see text). g The n-C equivalent of the corrected MD bilayer (RISM with E ) 0.90 kcal/mol). h The n-C equivalent of the uncorrected MD bilayer (RM3[1.1]). 16 g 16 i Calculation RM1[0.575], in which the value of Eg is determined experimentally. j MD-simulated under the same conditions used for the fluidphase bilayer. k R ) (Xgt + Xgg)/Xtt (see text). l Measured for n-C21H42-4,4-d2.

which would be exact if Xgg′, whose value is less than 0.01, was included. The last is

Xt ) Xtt + Xgt/2

(7)

Xg ) Xgg + Xgt/2

(7a)

or its equivalent

Equations 7 and 7a express the contributions to Xt and Xg from the tt, gt, and gt pairs that adjoin the CD2 group. Each tt pair contributes one trans bond to Xt and each gt one-half a trans bond, and similarly for Xg. Combining eqs 5-7 and using the definition of R given in eq 3 leads to the needed relation

Xg ) 0.5[(1 + R-1)-1 + Xgg]

(8)

This equation replaces our earlier one, Xg ) (1 + 2R-1)-1, in which the gg pair contribution to Xg was ignored.19 To verify eq 8, we can use the conformational statistics generated for n-C16 using the RISM model that are listed in the first column in Table 2 (calculation RM1[0.50]). The listed values of Xg, Xtt, Xgt, and Xgg are, respectively, 0.376, 0.360, 0.536, and 0.100. Using these values and eq 5 we find R to be 1.77. Now, we will reverse the calculation by assuming R and Xgg are equal to 1.77 and 0.100 and then use these values in eq 8 to calculate Xg. This gives a value of 0.369 for Xg, which is very near the correct value of 0.376, the small difference being due to the neglect of Xgg′. For comparison, the value of Xg obtained using our earlier eq 18 with R ) 1.77 is 0.469, which is substantially higher than the correct value, 0.376. To evaluate Xg for the fluid phase, DPPC bilayer, R is assumed to be 0.30 ( 0.08, the average of the four values of R measured at different sites along alkyl chains (12,26) (Table 4). The value used for Xgg, is the 0.043 ( 0.006 obtained from the wag-band measurements (Table 3). Equation 8 then gives Xg ) 0.14 ( 0.03. This value is significantly lower than nearly all those previously reported12,27,28 and is well below the respected estimate of 0.24 for the fluid-phase DPPC bilayer determined by Nagle and Wilkinson,29 who used a combination of calorimetry and dilatometry. Our bilayer value of 0.14 is less than one-half of that (0.35) for liquid n-C16 (Table 5). The large difference between the bilayer and n-alkane values is reflected in the large difference in the relative intensities of the 650 and 620 cm-1 bands in the CD2 spectra of these systems, as may be seen in Figure 15A.

Equation 8 indicates that the intensity ratio R, (I650/I620), essentially determines the value of Xg, since Xgg is much smaller than R. For the bilayer, the 620 cm-1 band is dominant, giving R a value of 0.30. For the n-alkane, however, the 650 cm-1 is dominant, giving R a value of 1.6, which is about 5 times the bilayer value. As a result, Xg for liquid n-C16 is found to be about twice that of the bilayer. Evidence supporting the reliability of the CD2 method was discussed above. In addition we have found no evidence indicating environmental factors or band shape distortion due to high band absorptivities have significant effects on the measured intensities. The frequencies, widths, and shapes of the CD2 rocking bands observed in the bilayer spectra are essentially the same as those previously observed in the spectra of CD2 substituted n-alkanes, in which case background absorption is much lower.19-21 It is also noteworthy that the CD2-band method has been successfully used to determine conformational disorder in the fluid-phase DPPE bilayer.30 In this study, the CD2 was positioned at carbon 4. The measured value of R was about 0.20, which is near the 0.26 measured for DPPC with the CD2 group also at the carbon 4 position. The method has also been used to measure the effect of cholesterol on the conformational ordering in the fluid-phase DPPC bilayer26 and to determine the changes in conformational disorder for the gel and fluid phases resulting from the addition of gramicidin.31 2. Terminal Bond (Xeg). The value of the gauche concentration Xeg associated with the CH2CH2 bond nearest the methyl endgroup is 0.31, a value about twice that for the internal CC bonds. There are at least two factors that tend to increase the value of Xeg relative to Xg. One is the chain-end effect discussed earlier. On the basis of RISM calculations on liquid n-C16, we estimate this effect increases Xeg by about 10%. A second factor, which is responsible for a much larger increase, stems from the fact that the methyl end of the chain resides in the most highly conformationally disordered region of the bilayer, that is, in the interface between the two monlayers. We have, in fact, encountered this situation in our analysis, above, of the distribution of conformational disorder in the gel-phase DPPC bilayer, in which case we observed a spectacular increase in the gauche concentraiton in the interfacial region (Figure 16). C. Trans/Gauche Bond Segregation. Trans/gauche (t/g) bond segregation in the fluid-phase bilayer is revealed by comparing the measured values of Xgg for the bilayer and its unconstrained n-C16 system equivalent. For the bilayer Xgg is 0.043, whereas for the n-C16 system it is 0.010. Thus, on the



average, each chain in the bilayer has the equivalent of onehalf of a gg pair, whereas for the unconstrained chain, the corresponding value is one-tenth. In general, an increase in t/g segregation results in a decrease in R and an increase in Xgg. The sensitivity of R and Xgg to segregation is illustrated by comparing the values of these parameters for two conformers with disparate degrees of t/g segregation. The chain we consider here has m CH2CH2 bonds, half of them trans and the other half gauche. The values of R and Xgg for the conformer (gt)m, whose t/g segregation is the least possible, are 8 and 1, respectively. In contrast, the values of R and Xgg for the most segregated conformer possible, tmgm, are 0 and 1/2. Chain-end gauche aggregation in the fluid-phase bilayer may be responsible for the fact that the values of R,12,26 measured at 4 positions along the alkyl chain, do not increase as the chain end is approached (Table 4), as would be expected since disorder increases. However, chain-end gauche aggregation can negate the expected increase if the gauche bond concentration is low, as is the case for the bilayer. To illustrate this, we note that since Xg for the bilayer is 0.14, the average number of gauche bonds per chain is about 2.0. If each gauche bond is flanked by trans bonds, as is the case, for example, for the conformer t2gt6gt3, then ntt, ngt, and ngg equal 8, 4, and 0, respectively. Then R ) 0.5. However, if there is aggregation and the two gauche bonds are adjacent and positioned at the end of the chain to form the t/g aggregated conformer g2t10, then ntt, ngt, and ngg equal 10, 1, and 1. Then R ) 0.20, which is less than half the value estimated for the less aggregated conformer. D. Comparison of Conformational Statistics for the FluidPhase Bilayer. For the comparison between the experimentally determined conformational statistics and the statistics derived from the MD simulation to be meaningful, we must account for the fact that the gauche concentration for the MD-simulated liquid n-C16 is significantly lower than the experimental value and therefore the MD bilayer value must also be too low. As a result, all bilayer statistics will be adversely affected. To allay this problem, first-order corrections to the MD-derived conformational statistics for the bilayer were made using scaling factors Tη based on liquid n-C16 statistics obtained from the RIS and MD models. The scaling factors are defined

T η ) Xη

RISM

/Xη

MD

(9)

where η is a specific conformational entity of interest and where XηRISM and XηMD are the calculated concentrations of the entity in liquid n-C16. The corrected and uncorrected values are listed in Table 5. The MD-derived statistics for the fluid-phase bilayer are in poor agreement with those derived from experiment (Table 5). The paramount difference is in the concentration of gauche bonds. The uncorrected MD value of Xg is 0.22, while the measured value is 0.14. If the “corrected” MD value (0.27) is used in the comparison, the agreement is even worse. This anomaly underlines the need to avoid corrections entirely by modifying the force field so as to bring the MD calculated value of Xg for liquid n-C16 into agreement with the experimental value. It is of interest to compare the MD-derived value of Xg for the fluid bilayer with the MD value for liquid n-C16. Comparable degrees of conformational disorder are indicated. The value of Xg for the fluid bilayer is about 20% less than that for liquid n-C16, irrespective of whether corrected or uncorrected MD values are used. The corrected values for the fluid bilayer and liquid n-C16 are 0.22 and 0.27, respectively; the uncorrected

values are 0.29 and 0.35. The similarity, if real, would appear to belie the fact that the conformational structures of the chains in the two systems are substantially different. The difference is in the distribution of trans and gauche bonds along the length of the chains. This difference is lost in the averaging. As previously discussed, the concentrations of short conformationally-specific sequences are, with few exceptions, nearly the same for fluid bilayer and liquid n-C16 systems having the same value of Xg. The near equal concentrations are a consequence of the near equal conformational energies determined by short-range intrachain forces, which are essentially the same for the fluidphase bilayer and liquid n-C16 systems. Thus, for example, the concentration (Xtttt) of sequences consisting of runs of 4 trans bonds is given closely, for both the bilayer and n-C16, by the term (Xt)4 derived on the (apparently excellent) assumption that the probability for the occurrence of a trans bond is independent of the conformation of adjoining bonds. The value of Xtttt given by (Xt)m for liquid n-C16, for which Xt is 0.708 (Table 1) is 0.251, a value that is essentially identical to the value of 0.247 obtained from the MD simulation. For the fluid bilayer, for which Xt ) 0.775, the (Xt)m value of Xtttt is 0.36, in surprisingly close agreement with the calculated MD value of 0.38. For long conformational sequences, however, the situation changes. Differences between the bilayer and n-alkane concentrations of conformational sequences increase significantly with increasing sequence length due to a tendency, driven by packing constraints, for the bilayer chains to favor trans bond sequences. The degree of this type of ordering might be a highly sensitive function of the interchain potential. The preference for trans sequences is evidenced by the appearance of highly ordered microcrystalline like regions in the MD-simulated fluid bilayer. We note that the value of Xg in the ordered regions is 0.12 ( 0.04, a value comparable to the 0.14 determined experimentally for the fluid bilayer. The ordered regions in the simulated structure could be easily overlooked since they involve only a small fraction, perhaps 15%, of the chains, and consequently have only a small effect on the conformational statistics. Further studies are clearly needed to explore the possibility of such highly ordered regions. Acknowledgment. We gratefully acknowledge support of this research from the National Institutes of Health (GM27690 to the University of California at Berkeley, GM40712 and GM14463 to the University of Pennsylvania, and GM29864 to Rutgers University at Newark). We thank Kay Newhouse for help in drawing the O/D map figures. References and Notes (1) Tu, K.; Tobias, D. J.; Blaise, J. K.; Klein, M. L. Biophys. J. 1996, 70, 595. (2) Tu, K.; Tobias, D. J.; Klein, M. L. Biophys. J. 1995, 69, 2558. (3) Zaccai, G.; Buldt, G.; Seelig, A.; Seelig, J. J. Mol. Biol. 1979, 134, 693. (4) Wiener, M. C.; Suter, R. M.; Nagle, J. F. Biophys. J. 1989, 55, 315. (5) Sun, W.-J.; Suter, R. M.; Knewtson, M. A.; Worthington, C. R.; Tristram- Nagle, S.; Zhang, R.; Nagle, J. F. Phys. ReV. E. 1994, 49, 4665. (6) Nagle, J. F.; Zhang, R.; Tristram-Nagle, S.; Sun, W.; Petrache, H. I.; Suter, R. M. Biophys. J. 1996, 70, 1419. (7) Martyna, G. J.; Tobias, D. J.; Klein, M. L. J. Chem. Phys. 1994, 101, 4177. (8) Brooks, B. R.; Bruccoleri, R. E.; Olafson, B. D.; States, D. J.; Swaminathan, S.; Karplus, M. J. Comput. Chem. 1983, 4, 187. (9) Flory, P. J. Statistical Mechanics of Chain Molecules; WileyInterscience: New York, 1969. (10) Snyder, R. G. J. Phys. Chem. 1982, 76, 3342. (11) Smith, G. D.; Jaffe, R. L. J. Phys. Chem. 1996, 100, 18718.

6288 J. Phys. Chem. B, Vol. 106, No. 24, 2002 (12) Mendelsohn, R.; Davis, A. M.; Brauner, J. W.; Schuster, H. F.; Dluhy, R. A. Biochemistry 1989, 28, 8934. (13) Pink, D. A.; Green, T. J.; Chapman, D. Biochemistry 1980, 19, 349. (14) Snyder, R. G. Macromolecules 1990, 23, 2081. (15) Snyder, R. G.; Cameron, D. G.; Casal, H. L.; Compton, D. A. C.; Mantsch, H. H. Biochim. Biophys. Acta. 1982, 684, 111. (16) Gabor, B. P.; Peticolas, W. L. Biochim. Biophys. Acta 1977, 465, 260. (17) Snyder, R. G. J. Chem. Phys. 1967, 47, 1316. (18) Snyder, R. G.; Poore, M. W. J. Chem. Phys. 1973, 6, 708. (19) Maroncelli, M.; Strauss, H. L.; Snyder, R. G. J. Chem. Phys. 1985, 82, 2811. (20) Maroncelli, M.; Strauss, H. L.; Snyder, R. G. J. Phys. Chem. 1985, 89, 4390. (21) Shannon, V. L.; Strauss, H. L.; Snyder, R. G.; Elliger, C. A.; Mattice, W. L. J. Am. Chem. Soc. 1989, 111, 1947.

Snyder et al. (22) Senak, L.; Davis, M. A.; Mendelsohn, R. J. Phys. Chem. 1991, 95, 2565. (23) Mendelsohn, R.; Snyder, R. G. In Biological Membranes; Mertz, K. M., Rox, B., Eds.; Birkhouser: Boston 1996; pp 145-174. (24) Casal, H. L.; McElhaney, R. N. Biochemistry 1990, 29, 5423. (25) Mendelsohn, R.; Davies, M. A.; Schuster, H. F.; Xu, Z.; Bittman, R. Biochemistry 1991, 30, 8558. (26) Davies, M. A.; Schuster, H. F.; Brauner, J. W.; Mendelsohn, R. Biochemistry 1990, 29, 4368. (27) Schindler, H.; Seelig, J. Biochemistry 1975, 14, 2283. (28) Petersen, N. O.; Chan, S. I. Biochemistry 1977, 16, 2657. (29) Nagle, J.: Wilkinson, D. A. Biiophys. J. (1978, 23, 159. (30) Davies, M. A.; Hubner, W.; Blume, A.; Mendelsohn, R. Biophys. J. 1992, 63, 1059. (31) Davies, M. A.; Brauner, J. W.; Schuster, H. F.; Mendelsohn, R. Biochem. Biophys. Res. Commun. 1990, 168, 85.

Acyl Chain Conformation and Packing in ... - American Chemical Society

Recommend Documents