What are the true values of the bending modulus of simple lipid bilayers?
Values of the bending modulus KC are reviewed, and possible causes for the considerable differences are discussed. One possible cause is the use of glucose and sucrose in the classical micromechanical manipulation and shape analysis methods. New data, using the more recent low angle x-ray method, are presented that do not support an effect of glucose or sucrose on KC. Another possible cause is using an incomplete theory to interpret the data. Adding a tilt term to the theory clearly does not affect the value obtained from the shape analysis method. It is shown that a tilt term, using a value of the modulus K\xce\xb8 indicated by simulations, theory, and estimated from order parameters obtained from NMR and from the wide angle x-ray method, should also not affect the value obtained using the micromechanical manipulation method, although it does require a small correction when determining the value of the area compressibility modulus KA. It is still being studied whether including a tilt term will significantly affect the values of KC obtained using low angle x-ray data. It remains unclear what causes the differences in the experimental values of KC for simple lipid bilayers.
1. Introduction<!>2. What is the effect of sugar concentration cs on the bending modulus?<!>3. New x-ray results for possible sugar effect<!>4.1. Standard Theory<!>4.2. Refined theory<!>4.3. Refined theory - MM method<!>4.4. Refined theory applied to SA and x-ray methods<!>5. On obtaining values for the tilt modulus K\xce\xb8<!>6. Discussion and Conclusions<!>
<p>The bending modulus KC is a most important membrane mechanical property. Accordingly, it has been measured many times for many different lipid bilayers. Although uncertainties are typically reported at the 10% level, values obtained in different labs and with different measuring techniques typically differ by as much as a factor of two for the same lipid at the same temperature (Nagle, 2013, Marsh, 2006). As biological processes often involve transition states with curved membranes, that part of the activation energy would differ by a factor of two. For thermally activated processes, the predicted rate constant, depending as it does on the exponential of the activation energy, could easily be kinetically incompetent for the larger KC value, while being quite feasible for the smaller value (Nagle, 2013). It is therefore of some biophysical importance to obtain more accurate values of KC, as well as to alleviate the embarrassment that membrane researchers lack accepted values for a quantity that is recognized to be central.</p><p>The two most common methods for measuring KC are micromechanical manipulation (MM) of giant unilamellar vesicles (GUV) (Rawicz et al., 2000, Henriksen and Ipsen, 2004, Vitkova et al., 2006, Shchelokovskyy et al., 2011, Evans and Rawicz, 1990), often called the pipette aspiration method, and fluctuating shape analysis (SA) of GUV (Meleard et al., 1998, Henriksen and Ipsen, 2002, Meleard et al., 1997, Pecreaux et al., 2004, Gracia et al., 2010) and many earlier references (Bouvrais, 2012, Vitkova and Petrov, 2013). A few results have been obtained by pulling cylindrical tethers from GUV (Heinrich and Waugh, 1996, Sorre et al., 2009, Tian et al., 2009) and a variety of other techniques are reviewed by (Dimova, 2014). Here we will focus on some MM and SA results as well as results from an x-ray method. Analysis of low angle diffuse x-ray scattering from oriented stacks of membranes has also been employed more recently (Lyatskaya et al., 2001, Liu and Nagle, 2004, Salditt et al., 2003, Li et al., 2006, Pan et al., 2008, Pan et al., 2009), and the values so obtained agree well with those reported in a classic MM paper (Rawicz et al., 2000). However, ignoring for the moment differences between different labs using the same method, the values obtained using SA are generally larger than those obtained using MM or x-ray methods (Nagle, 2013, Marsh, 2006). One possible reason may be related to different length scales of the measurements and the theory involved, as we review in Section IV. A more mundane possibility regards an experimental aspect of the two GUV methods to which we first turn in the next section.</p><!><p>The MM method typically uses a sugar solution to conserve the volume of the GUV (Vitkova et al., 2006, Evans and Rawicz, 1990). Often a sucrose solution inside the GUV and a glucose solution outside the GUV is used to enhance visual contrast. A classic study used cs = 200mM sugar with the same molarity on both sides to ensure flaccid GUV with zero imposed surface tension or pressure; Fig. 1 shows the reported value (red right pointing triangle) of KC/kT for DOPC, that has two double bonds in the two oleoyl (DO) hydrocarbon chains, and the value (blue right pointing triangle) for SOPC, that has a saturated stearoyl (S) chain that makes the SOPC bilayer a bit thicker and stiffer than DOPC (Rawicz et al., 2000). Also shown in Fig. 1 are the results from the first study that reported an effect of sugar concentration on SOPC (Vitkova et al., 2006). This study utilized the SA method for small sucrose concentrations (two open blue down triangles in Fig. 1) and the MM method (four solid down blue triangles in Fig. 1) for larger sucrose concentrations (fluorescent dye was used for contrast instead of glucose outside). These results, when interpolated at 200mM, are about a factor of two smaller than the earlier results (Rawicz et al., 2000), perhaps attributable to differences in the way the two labs interpreted MM data. More importantly, an exponential decay with sugar concentration was indicated, as shown by the lower dashed blue line in Fig. 1. The same group, using the SA method exclusively, also reported a decreasing KC with increasing sucrose concentration, although the decrease was only about half as large as in Fig. 1 (Genova et al., 2006). Subsequently, decreasing KC was reported for other sugars (Genova et al., 2007), although, contrarily, maltose was reported not to decrease KC (Genova et al., 2010). Figure 1 also shows the most recent SOPC value (open square) with no sugar obtained after further development of the SA method (Genova et al., 2013).</p><p>An MM study of DOPC found that KC at small cs = 8mM sucrose/8mM glucose was twice as large as with 100mM sucrose/110 mM glucose (Shchelokovskyy et al., 2011). As shown by the red dashed line in Fig. 1, that is also consistent with a strong exponential decay with sugar concentration; extrapolation to 200mM gives a value three times smaller than the earlier DOPC result of (Rawicz et al., 2000). Interestingly, the same exponential dependence connects the earlier MM value for SOPC of (Rawicz et al., 2000) with the values of KC reported by (Henriksen and Ipsen, 2004) using the MM method and by (Henriksen and Ipsen, 2002) using the SA method (half closed squares in Fig. 1). Finally, Fig. 1 shows that our x-ray results at zero sugar (upright triangles) agree well with the earlier MM results (right pointing triangles).</p><p>Not surprisingly, considering the current state of KC results, there is the feeling that more experiments should be done (Genova et al., 2013) (Dimova, 2014). However, instead of only exhorting others to do more experiments, we have been addressing this issue using the x-ray method.</p><!><p>Briefly, four samples (i–iv) were made by mixing DOPC and sugar (first solubilized in heated trifluoroethanol or methanol) in excess 1:1 chloroform/(trifluoroethanol or methanol) organic solvent with mole ratios of (i) 0.12 glucose/DOPC, (ii) 0.12 sucrose/DOPC, (iii) and (iv) 0.22 glucose/DOPC. The mixtures were deposited on Si wafers using the rock and roll technique to achieve superior alignment in the stack of about 2000 bilayers (Tristram-Nagle, 2007). The dry sample was then hydrated in a humidity chamber in situ on the x-ray beamline. Hydration was conveniently even more rapid with sugar than for pure DOPC. Compared to the repeat spacing of fully hydrated DOPC, D=63.1 Å, the repeat spacing increased to D ≈ 68 Å for the lower concentrations in samples (i) and (ii), and for the higher concentrations in samples (iii) and (iv) D increased to 74 Å. Fig. 2 shows a grayscale image of the x-ray scattering. KC/kT was obtained using our analysis procedure (Liu, 2003).</p><p>New x-ray results are shown in Fig. 3 along with some of the literature results already shown in Fig. 1. In order to facilitate comparison between various results, KC has been normalized to 1 for zero sugar, except for the MM result for DOPC at 200mM for which the x-ray result was used for the normalization. The x-ray KC result with no sugar has been shifted to a slightly negative concentration to allow it to be seen distinctly from the SOPC result; the small x-ray uncertainty is the estimated uncertainty in the mean obtained from 27 DOPC samples studied over a period of ten years. The standard deviation for the distribution of KC for these 27 samples is 5.2 times larger, and that is the assumed uncertainty shown by the error bars for the new data. However, we emphasize that for samples (iii) and (iv), several D spacings were achieved for the same sample by systematically varying the relative humidity. This gave a range of sugar concentrations for the same sample, so the same uncertainty applies to all the KC values for that sample. That is to say, KC for all the concentrations has an uncertainty given by the error bar, but the concentration dependence shown by the slope of the line is not affected by this uncertainty. As that slope is essentially zero, samples (iii) and (iv) show negligible dependence of KC on glucose concentration. More importantly, the values for each of those samples are remarkably close to the average value for the 27 samples with no sugar.</p><p>The amount of sugar added in samples (i) and (ii) was designed to correspond to 200mM sugar in fully hydrated DOPC, calculated from the result that DOPC has nW=30 water molecules/lipid between the bilayers with lamellar repeat spacing D=63.1 Å (Nagle and Tristram-Nagle, 2000, Kucerka et al., 2008) and area A = 67.4 Å2 (Kucerka et al., 2008). There are two reasons that the concentrations displayed in Fig. 3 for samples (i) and (ii) are smaller. The first reason is that the repeat spacing was larger, D = 68 Å, for both glucose (i) and sucrose (ii). As the bilayer dimensions changed very little, nw should be increased to 35.5, decreasing these nominal sugar concentrations to 164 mM.</p><p>The effective sugar concentration can be further reduced if there is binding of sugar molecules to the lipid, over and above the number of sugars per lipid that would dissolve in the nw intercalated waters. If such binding occurs, then the number ns of sugars/lipid, which is ns=0.11 from the weights of sugar and lipid in samples (i) and (ii), is the sum of bound, nsb, and unbound, nsub, sugars, and the concentration that must be compared to MM and SA experiments is the nominal concentration of 164mM times nsub/ns. Whether there is net binding or the opposite, net exclusion, has been controversial, with recent work (Andersen et al., 2011) concluding that there is net binding for small concentrations and net exclusion for larger concentrations. Using the data in Fig. 3 of (Andersen et al., 2011), we estimate that nsub/ns ~ 0.64 for samples (i) and (ii), which multiplied by 164mM gives 105mM; this is the concentration shown in Fig. 3 for those samples. However, for our higher concentrations, the same figure indicates no net binding or exclusion, so no nsub/ns factor was applied to the concentrations of samples (iii) and (iv).</p><p>There is yet another consideration that could reduce the effective sugar concentration. During the evaporation of the organic solvent used to mix lipid and sugar, a fraction f of the sugar could have accumulated into defect regions free of lipid instead of being intercalated between the dry oriented bilayers. Assuming that this sugar did not migrate when the sample was hydrated would decrease the number of sugars between adjacent bilayers to (1-f)ns. Nevertheless, it is clear that some sugar resided between the bilayers because the fully hydrated D spacing increased when sugar was added. (In passing, note that the value of the D spacing reflects a balance mostly between the attractive van der Waals interaction and the repulsive fluctuation pressure, the hydration force being small near full hydration. As the fluctuation pressure depends on KC, and that value did not change with sugar, the total sugar, intercalated plus bound, likely reduces the van der Waals attraction.) In the extreme case, if sugar had been totally excluded from between the bilayers, then hydration would have formed pools of sugar/water in contact with the stack of bilayers. Such pools would have exerted an osmotic stress on the water between the bilayers, sucking it out and thereby reducing D compared to that of fully hydrated bilayers with no sugar, rather than increasing D as observed. For the relevant case of a fraction (1-f) of the sugar present between the bilayers and f excluded, then upon hydration, as water is more permeable by four orders of magnitude than sugar, water would enter both between the bilayers and into the pools, resulting in equal sugar concentrations in both locations. Then the total volume of the pools would be f/(1-f) times the volume of the water between the bilayers; the latter volume is 45% of the volume of the DOPC bilayer stacks when D=68Å and A=67.4Å2 (Kucerka et al., 2008). We estimate that f=1/2 is an upper bound because it is likely that even this much volume of defect pools would disorient the sample, but no difference in mosaic spread was observed when sugar was added. (It might also be noted that the presence of defect pools impacted the Luzzati gravimetric method for measuring area/lipid, (Nagle and Tristram-Nagle, 2000, Koenig et al., 1997) with about 20% of the water in such pools for DOPC at full hydration. A difference is that the relative size of the pools to the intercalated water could be shrunk by applying osmotic pressure with a polymer, but here that ratio remains constant.)</p><p>Although we can only say that the sugar concentrations cs for the x-ray results shown in Fig. 3 are upper bounds, we estimate that an additional reduction in the effective concentration cs by even a factor of 2 is unlikely. Therefore, our new results do not support a dependence of KC on sugar concentration. However, after initial submission of this paper, we noticed a recent result, included only in a book chapter, consistent with our result in Fig. 3 that glucose has a negligible effect on KC. Contrary to (Genova et al., 2007), (Vitkova and Petrov, 2013) reported that glucose negligibly decreases KC in eggPC bilayers, whereas sucrose continued to have an effect similar in magnitude to that reported for SOPC bilayers (Genova et al., 2006). Our result for sucrose at a single concentration in Fig. 3, while not indicating an effect on KC, could be consistent with one. We plan to extend our sucrose studies to higher concentrations.</p><!><p>The SA, MM and x-ray methods all measure time averaged displacements of the bilayer from the flat condition. If the flat bilayer is positioned in the r = (x,y) plane with z=0, then thermal fluctuations produce out-of-plane deviations z(x,y) that time average to zero, <z>=0, but that have non-zero mean square values <z2>≠0. It is assumed that the local energy of bending symmetric bilayers to positions z(x,y)=z(r) is given by
(1)Ebend(r)=(KC/2)(∂2z/∂x2+∂2z/∂y2)2≡(KC/2)C(r)2,
where C(r) is the local curvature. (Helfrich, 1973) Straightforward statistical mechanics enables calculation of the quantities of interest by Fourier transforming into reciprocal space with qr = (qx,qy). This gives the well-known undulation spectrum
(2)S(qr)=kTKCqr4+Î³qr2
for the mean square amplitudes of the modes with wavelength 2π/qr, where qr=|qr| and γ is the surface tension. The surface tension is identically zero for the x-ray method where no difference in osmotic pressure can be applied to different sides of the bilayersm whereas it is systematically varied in the MM method as aspiration pressure is changed. The SA method directly measures the amplitudes of the smallest qr modes with the longest wavelengths (λ=2π/qr, typically of order 10µm), so surface tension is kept small to be able to obtain KC.</p><p>The MM method measures changes in the projected area Δα of a GUV as γ is varied over a wide range. As is well known (Rawicz et al., 2000, Evans and Rawicz, 1990), there are two contributions to Δα. The obvious contribution is the expansion of a flat membrane given by ΔαA=γ/KA, where KA is the area compressibility modulus which is typically 250–300mN/m (Evans et al., 2013). However, the largest initial increase in projected area when small γ is applied comes from reducing the area αU that was hidden in the undulations; αU is obtained by integrating S(qr)qr3dqr from the smallest qm limited by the radius R≈10µm of the GUV to the largest qM ~ π/a where a is an intermolecular in-plane distance of order 8Å, so the integration range spans a little more than 4 orders of magnitude. When γ=0, the integral is logarithmic, so each decade of length scales contributes equally to the area. The area pulled out from undulations is then
(3)Δαu=αu(0)−αu(γ)=(kT/8πKC)ln(1+γ/KCqm2)
in the γ range of interest.</p><p>The x-ray method requires pair correlation functions <zn(r)zm(r')>, not only between points in the same bilayer, but also points in different n and m bilayers, which takes some additional calculation (Lyatskaya et al., 2001); the major physical difference is that fluctuations in the spacing between neighboring bilayers means that fluctuations in the z direction and the corresponding qz modes must also be considered. Also, γ=0, so Eq. (2) becomes
(4)S(qr,qz)=kTKCqr4+Bqz2,
where B is the interbilayer compressibility modulus. The interaction term Bqz2 reduces the amplitude of the small qr modes within the same bilayer, thereby reducing <zn(r)zn(r')> for the larger r−r' length scales.</p><!><p>All three methods in the preceding paragraphs should provide the same value for KC for the theory in Eq. (1). However, there is a more refined theory (Hamm and Kozlov, 2000) that adds to the bending energy in Eq. (1) a local free energy that depends on molecular tilt t(r) with respect to the local bilayer normal,
(5)Etilt(r)=Kθ|t(r)|2/2,
and it also redefines the curvature in Eq. (1) by the addition of gradr(t(r)). This theory has a spectrum (May et al., 2007)
(6)S(qr)=kTKCqr4+kTKθqr2.
Estimates of the tilt modulus Kθ of order 50mN/m will be discussed in Section V. Assuming this value here and a typical value of KC = 10−19 J = 24kT, the qr−4 term dominates S(qr) at the large length scale, so the true KC is obtained by the SA method. However, the second term in (6) becomes larger than the first term when qr exceeds qθ = (Kθ/KC)1/2 ≈ 0.7 nm−1, corresponding to a tilt length scale ξθ ≡ 1/qθ ≈ 1.4 nm. From there down to the molecular length scale of 0.8 nm comprises part of the decade of largest qr.</p><!><p>In order to quantify what the preceding theory implies for the MM method, it is necessary to consider the theory with surface tension γ. Just as for the derivation with no surface tension, one adds a term consisting of γ times the excess surface area. Following (May et al., 2007), this already results in the following rather more complicated expression
(7)S(qr)=kT(1+KCqr2/Kθ)KC(1+γ/Kθ)qr4+γqr2.
In addition, it can be argued that applying surface tension reduces the energy for tilting, in which case Kθ in Eq. (7) would be replaced by Kθ−γ to first order in γ. Recently a formula has been published for a model with even more additional parameters than the tilt modulus (Watson et al., 2013); when the extra terms are eliminated and the surface tension is retained only to leading order, which make only small differences for the experimentally accessible range, that formula reduces to Eq. (7) with Kθ replaced by Kθ−γ.</p><p>Starting from the modified Eq. (7), the same procedure as for the conventional theory in Eq. (2) then predicts the following result for the area pulled out from undulations in an MM measurement
(8)Δαu(γ)=kT8πKC[(1−γKθ)2ln(1+γ(1−γ/Kθ)KCqm2)+2γKθ(2−γKθ)ln(qMqm)]
Plotting this in Fig. 4 shows that there is little difference in the predicted MM measurements for the theory with tilt compared to the conventional theory in the small tension regime where KC can be obtained from the intial slope. For larger γ the theory with tilt has a larger Δα because there is an additional degree of freedom that gives rise to additional area αu that is then pulled out with increasing γ. However, for small γ it is the pulling out of large length scale (small qr) undulations that dominates Δα, giving the same form as the conventional theory in Eq. (3).</p><p>Fig. 5 emphasizes the large γ regime where Δα(γ) is dominated by the lateral area expansion ΔαA(γ) = γ/KA. However, even in the conventional case, the slope of the Δα(γ) curve is greater than ΔαA(γ) because ΔαU(γ) continues to increase as the undulations continue to be pulled out. Although a naïve approach to obtaining KA would just fit the slope of Δα(γ) in the large tension regime in Fig. 5, the appropriate way is to subtract ΔαU(γ) first and then one re-obtains the exact KA which was set to 250mN/m. When the tilt term is added, Fig. 5 shows that Δα(γ) is greater for all γ because of the extra tilt degree of freedom. Importantly, the slope is also greater at large γ.</p><p>The first interesting question now is, how much would the apparent KA value differ from the true value if one fits data, affected by a tilt degree of freedom, with the conventional model that has no tilt degree of freedom. For the example in Fig. 5, the fitted KA is 233mN/m, 7% smaller than the exact 250mN/m. The second interesting question is, are the residual errors in this fit large enough to diagnose that the conventional theory is incorrect? The answer is no; the largest residual to the fit is only 0.0002 for the example in Fig. 5, not even large enough to see in these plots and far smaller than experimental uncertainty. In other words, the conventional theory is adequate to fit MM data, but MM data are not sufficient to detect the presence or absence of a tilt degree of freedom, at least for typical values of the moduli.</p><p>The preceding analysis shows that not including a tilt degree of freedom in the analysis of MM data still obtains the true value of KC. Therefore, this does not account for the experimental differences compared to results obtained from the SA method. Importantly, the MM method can not provide experimental evidence for or against including the tilt degree of freedom in modern models of lipid bilayers. Furthermore, as the preceding analysis shows that the MM method underestimates KA if there is a tilt degree of freedom, determination of KA by the MM method would require a value of Kθ from other sources. However, the error of not including this term appears to be less than 10%. This point should nevertheless be kept in mind when evaluating KA from simulations; a value smaller than the true value should be expected when evaluating KA using A(∂γ/∂A)T. Interestingly, a recent simulation obtained a smaller value of KA (277±10 mN/m) using this method compared to the KA (321±37 mN/m) obtained from the area fluctuations σA using the formula 2AkT/NσA2 (Braun et al., 2013).</p><!><p>In the SA method, individual GUV have non-zero surface tension. The method observes fluctuations for q less than about 1 µm−1, so the surface tension has to be small, of order 10−4 mN/m in order for the spectrum in either Eq. (2) or (7) to depend measurably upon KC ≈ 10−19J. For such small γ the ratio of S(q) in Eq. (7) to S(q) in Eq. (2) differs from 1 only by about 10−6, so the SA analysis would give the same value for KC regardless of tilt, and the method would also have no possibility of determining the tilt modulus Kθ.</p><p>The x-ray method could be different. To model the samples which consist of bilayer stacks, the interbilayer compressibility energy must be included. This gives a height-height spectrum similar to Eq. (7) with each factor there of γqr2 replaced by Bqz2,
(9)S(qr,qz)=kT(1+KCqr2/Kθ)KCqr4+Bqz2(1+KCqr2/Kθ).
This is considerably more complicated than Eq. (4), and even for that conventional model the computation of x-ray intensities is non-trivial, because there are both in-plane qr and out-of-plane qz modes, and pair correlation functions, not just S(qr,qz), are required. With the extra complexity introduced by tilt in Eq. (9), data analysis becomes even more computationally difficult. Would such a development be useful, or would it, like the MM method, not make any significant difference? A difference compared to the MM method is that the undulations of bilayers in a stack are constrained by neighboring bilayers, and this reduces the amplitudes of the small qr undulations preferentially as the wavelength exceeds the lateral correlation length ξ ≡ (KC/B)1/4 ≈ 5nm. This is comparable to the tilt length scale ξθ = (Kc/Kθ)1/2 ≈ 1.4 nm below which tilt becomes dominant. Therefore, in contrast to MM data which are affected by four decades of length scales, the x-ray data are determined predominantly by the smallest decade, so it is more likely that tilt might affect the x-ray values of KC. Another reason that incorporation of a tilt modulus into the x-ray data analysis could be valuable is that, with signal/noise becoming higher, small but systematic differences have been appearing between the fits and the data using the conventional theory, as shown in Fig. 2.</p><!><p>Although the primary focus of this paper is on the bending modulus, the possible involvement of a tilt degree of freedom may make it of some interest to consider the value of its modulus in this section. A theoretical estimate of the tilt modulus gave Kθ ~ 100mN/m (Hamm and Kozlov, 2000) and analysis of inverted hexagonal phase data led to Kθ ~ 80mN/m (Hamm and Kozlov, 1998). (Note the factor of 2 is required to convert monolayer values to bilayers.) It is interesting that a model different from that of Hamm and Kozlov also provides a similar theoretical estimate for Kθ. This alternative model supposes that a tilted chain continues to have the same overall length and the same cross-sectional area perpendicular to its tilt axis, thereby requiring the interfacial area to increase by ΔA~t2/2 to conserve volume. The increase in hydrophobic free energy would then scale as γowt2/2 for each monolayer, where the oil-water surface tension γow≈50mN/m would then be half the value of the tilt modulus for a bilayer. We note that this model motivated the replacement after Eq. (7) of Kθ by Kθ−γ because an applied surface tension γ decreases the interfacial free energy of tilting. This model differs from the original theory of Hamm and Kozlov which does not rescale Kθ upon tilting, instead maintaining the interfacial area and requiring the tilted chains to stretch in order to maintain constant volume. We have been informed by Evan Evans that a polymer brush model (Rawicz et al., 2000) also predicts the same value of Kθ and predicts that tilted chains stretch, unlike the alternative model above, but only 1/3 as much as the Hamm and Kozlov model. It would be of interest to examine simulations to evaluate these different models by correlating chain stretching with chain tilting, both of which are easily defined in simulations. It is interesting that all three models give the same value of Kθ.</p><p>Values ranging from 50–110mN/m have been reported from spectral analysis of coarse grained simulations (May et al., 2007, Watson et al., 2012). The most realistic simulation that had a large enough length scale for spectral analysis has been an atomistic simulation for DMPC (Brandt et al., 2011). When analyzed using a direct Fourier method, there appeared to be no qr−2 term in the S(qr) spectrum, corresponding to infinite Kθ. However, analysis of the same simulation by imposing a real space envelope showed a qr−2 term, consistent with Kθ = 56mN/m (Watson et al., 2012); the difference between these two methods of analysis is the subject of ongoing study (Albert et al., in preparation). A different simulation method that fitted the angular distribution to a potential of mean force obtained a monolayer Kθ ≈ 6.7kT/A for DOPC (Khelashvili and Harries, 2013); for bilayers this is doubled to give 80 mN/m when A is approximated as 70 Å2.</p><p>We next consider the use of chain order parameter data to estimate Kθ. The generic order parameter is defined as
(10)S=(3<cos2β>−1)/2,
where β is the angle away from the bilayer normal. Most commonly, deuterium NMR measures averages over the orientational distribution of individual CD bonds to obtain the SCD order parameters as a function of carbon number. The value of SCD in the plateau region can then be converted into the molecular order parameter SNMR=2|SCD| (Seelig and Seelig, 1980). An order parameter Sx-ray has also been obtained from x-ray scattering in the wide angle range (Mills et al., 2008), distinct from the low angle x-ray scattering in Fig. 2 that has been considered so far. Sx-ray averages over the angular distribution of local bundles of tilted chains and is larger than SNMR by about a factor of 1.4. In principle, Sx-ray might be more appropriate in that it addresses structure at the same length scale as the tilt parameter instead of at the level of individual methylene groups, but the differences do not substantially affect the following estimates of Kθ.</p><p>In terms of β, the absolute value of the tilt parameter t is tanβ which for small β can be approximated as sinβ. Therefore the total average chain tilt energy from Eq. (5) is
(11)Etilt≈(Kθ/2)<sin2β>(ACNC/2)=(Kθ/2)ACNC(1−S)/3,
where AC is the mean area/chain and NC is the total number of chains in both monolayers. Assuming independent tilt modes, as was assumed by (Khelashvili and Harries, 2013) and (May et al., 2004), equipartition of NC chain tilt degrees of freedom then requires Etilt = NCkT/2. Combining with Eq. (10) then gives
(12)Kθ≈3kT/AC(1−S).
For DOPC at T=30° C, AC=0.34nm2, Sx-ray = 0.27 (Mills et al., 2008), yielding Kθ = 48mN/m. For DMPC at 30° C, AC=0.30 nm2, Sx-ray=0.40 (Pan et al., 2009), yielding Kθ = 70 mN/m. The simulated value of Kθ = 56 mN/m for united atom DMPC (Watson et al., 2012) is in reasonable agreement, but in the simulation the tilt parameter was defined for each lipid molecule rather than each chain; that replaces AC in Eq. (12) with AL=2AC which would reduce the order parameter estimate of Kθ to 35 mN/m. This raises the issue of what the basic unit, chains vs. molecules, should be for the tilt parameter.</p><p>Another issue is whether the basic unit, either chains or molecules, should be considered to be statistically independent. If tilt is correlated between chains, then the value of NC in the assumption Etilt = NCkT/2 would be reduced and that would reduce Kθ. For the tilted DPPC gel phase, the correlation between chains persists to 2900 Å (Sun et al., 1994) indicating large numbers of chains in each independent unit. Of course, such long range persistence does not occur in the disordered fluid phase where the width of the wide angle scattering indicates a lateral correlation length L that is only about 6 Å (Mills et al., 2008), thereby supporting nearly independent chains, although that L is the overall positional correlation length whereas the tilt-tilt correlation length could be larger. The tilt-tilt correlation length is an interesting quantity, one that could be determined by mining simulation data. Simulated tilt-tilt correlation lengths of cholesterol have been reported (Khelashvili and Harries, 2013) to be less than 1 nm and splay-splay correlation lengths of lipids have been reported to be less than 0.5 nm (Khelashvili et al., 2014 in press). Also, (May et al., 2004) suggest tilt-tilt correlations are negligible. Whatever reduction in NC would be indicated, it is important to note that the assumption of independent modes gives an upper bound to Kθ and that the smaller the estimate for Kθ, the more important it is to include a tilt mode in the analysis of data and simulations.</p><p>It should also be noted that the experimental orientational order parameter Sexp is a product of the local tilt order parameter Stilt, which is what should be used in Eq. (12), and an undulational order parameter Sund that accounts for undulations tilting the local bilayer normal from the laboratory z axis (Petersen and Chan, 1977). For small ϕ deviations, Sund ≈ 1 – 3<ϕ2>/4. For the classical theory with no tilt, <ϕ2> ≈ 0.028 giving Sund ≈ 0.98 (Nagle and Tristram-Nagle, 2000). With tilt,
(13)<ϕ2>≈kT2πKC[ln(πξa)+12(π2KCa2Kθ)]
where the x-ray analysis gives ξ ≈ 40Å. The extra term on the right hand side would reduce Sund close to 0.9 using a=8 Å but a larger value of a closer to the bilayer thickness is probably more appropriate, returning Sund closer to 0.98. In either case, this appears to introduce only a secondary numerical correction to Eq. (12).</p><p>It is interesting to compare the estimate in Eq. (12) with the theoretical result Kθ=2γow obtained at the beginning of this section. For both to be true, the interfacial hydrocarbon/water interaction, represented as γow, would have to be larger for DMPC than for DOPC because Sx-ray is larger for DMPC. Nevertheless, it is encouraging that estimates from three theoretical models, from simulations, and from experimental order parameters give similar values of Kθ.</p><p>There is also a general issue that arises because there are two fundamentally different tilt moduli that are not necessarily equal. One is the theoretical modulus Kθ in Eq. (5) that is input into a statistical mechanical model. This is the modulus that is evaluated by fitting simulation data, using either the spectra such as those in Eqs. (6) and (7) or fitting the angular distribution to a potential of mean force (Khelashvili and Harries, 2013). The other modulus is the thermodynamic modulus, let us call it Kθ+, that emerges from a statistical mechanical calculation of the model free energy. (Generally in statistical mechanics an input model parameter is not necessarily equal to a similarly named output thermodynamic quantity; the B modulus that appears in Eq. (9) is another example in membrane mechanics (Petrache et al., 1998).) For the thermodynamic modulus (May et al., 2004) calculated Kθ+ = Kθ + 3kT/AC and both terms have similar magnitudes so their Kθ+ is about a factor of two greater than Kθ. In supplementary material we show an alternative model that gives Kθ+/Kθ ≈ 1.5. Because Eq. (12) uses experimental order parameters, it would seem that the ensuing tilt modulus should be compared to Kθ+. However, all the corrections and approximations mentioned above for Eq. (12) would only decrease those values further, so there is an unresolved discrepancy regarding values of the thermodynamic tilt modulus Kθ+.</p><p>Finally, x-ray data of the sort shown in Fig. 2 may provide another independent experimental measure of the same Kθ that appears in Eq. (9) if the analysis can be extended to include a tilt modulus. Alternatively, even if the low angle x-ray scattering data are insufficient to obtain Kθ independently, estimates of Kθ may help improve the fit to the low angle x-ray data and provide better values of KC.</p><!><p>Our group has found it very encouraging in the past that x-ray values of KC agree so well with the MM values of (Rawicz et al., 2000). However, one has to consider that both sets of values may be artifactually too small. The reason that the x-ray values may be too small is that they are applied to data that are determined at a smaller length scale; in simulations, addition of a tilt degree of freedom increases the obtained value of KC (May et al., 2007, Watson et al., 2012) and something similar might occur for the x-ray values. The reason that the MM values of (Rawicz et al., 2000) might be too small is that other studies have suggested that use of 200 mM sugar reduces KC by a factor of two or more. Our new data show that glucose does not reduce the x-ray value of KC. Our sucrose result also does not indicate a reduction, although more data should be obtained at higher concentration. Although the x-ray method does not yet include the tilt degree of freedom, nevertheless, it seems unlikely that a tilt dependent implementation would alter the observance or non-observance of a sugar effect.</p><p>It was recently suggested that the MM values of KC may have been too small because the analysis did not include a tilt degree of freedom (Nagle, 2013). Importantly, we have here shown that adding a tilt degree of freedom is very unlikely to change either the MM or the SA values of KC, so other causes must be hypothesized for the larger KC values often obtained by the SA method compared to those obtained by the MM method. Nevertheless, it has been becoming increasingly clear that a tilt modulus Kθ, while secondary to the bending modulus KC, is warranted in membrane mechanics, even though it is difficult to measure experimentally. This paper strengthens the case for including a tilt modulus by showing that several distinct approaches, theoretical models, simulations, and experimental order parameters obtain similar values of Kθ.</p><p>In conclusion, there are still unexplained differences in the experimental values of the bending modulus KC, but the number of possible reasons has been reduced.</p><!><p>This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.</p>