INTRODUCTION

TOP ABSTRACT INTRODUCTION THEORY AND METHODS RESULTS AND DISCUSSION CONCLUSION REFERENCES

There is considerable evidence in recent literature for the active role of phospholipid molecules in lipid bilayer membranes (LBMs). This includes protein activity as modulated by the physical properties of lipid molecules (Brown, 1994; Burack et al., 1994), lipid-mediated protein-protein interactions due to hydrophobic matching (Mouritsen,1993; May and Ben-Shaul, 1999; Harroun et al., 1999), and lipid-mediated protein-protein interactions related to lipid-packing effects caused by hydrophobic interactions between proteins and lipids (Sintes and Baumbartner, 1997; Lague et al., 1998). These effects clearly show that the phospholipid bilayer does not act as a passive structureless non-polar solvent phase as suggested by the fluid mosaic model of Singer and Nicolson (1972). In a previous study, we examined nonspecific lipid-mediated protein-protein interactions arising from perturbations of the lipid structure by the proteins themselves (Lague et al., 1998). In this study the LBM was composed of dipalmitoyl phosphatidylcholine (DPPC) molecules. Here we extend the previous study by examining the change in such nonspecific lipid-mediated protein-protein interactions in pure lipid bilayers composed of different phospholipids.

The present study is motivated by the fact that a great variety of phospholipid molecules is found in biological membranes. Furthermore, several authors (Mouritsen et al., 1993a,b; Gil et al., 1998; Crane et al., 1999; Marsh, 1995) have shown that this leads to lateral membrane heterogeneity, which clearly has an impact on protein-lipid interactions. Lipid molecules with one or two unsaturated aliphatic chains are ubiquitous in biological membranes, and their physical properties are different from those of saturated lipids. Specific differences include area per lipid molecule, hydrophobic thickness (Nagle, 1998), and molecular order and dynamics (Mitchell and Litman, 1998) in LBM and lipid-protein interactions in biological membranes (Castuma et al., 1993). It is important to characterize the effect of these different physical properties on lipid-mediated protein-protein interactions.

In order to investigate the lipid-mediated protein-protein interactions, an approach based on statistical mechanical integral equations theories developed for the study of liquids (Hansen and McDonald, 1986; Chandler et al., 1986) was recently developped by the authors of the present work to describe the average structure of the hydrocarbon chains of LBMs (Lague et al., 1998, 2000). The approach has the advantage of combining aspects of both mean-field theory and results from fully detailed atomic simulations, much in the spirit of the Pratt-Chandler theory of the hydrophobic effect (Pratt and Chandler, 1977). The theory was derived in terms of a hypernetted chain (HNC) integral equation projected onto the two-dimensional Cartesian space of the lipid bilayer plane. The only input required for the application of this theory is the exact lateral density-density response function of the hydrocarbon core, computed from the configurations of a MD simulation of a specific single component lipid bilayer (without protein inclusions). The output of the theory is the perturbed density of the hydrocarbon chains around protein inclusions, and the lipid-mediated potential of mean force (PMF) between two protein inclusions. The calculation of the corresponding PMF using only MD simulations would be prohibitively expensive and is currently impossible. The present approach, based on a combination of MD and integral equations, provides a unique way to extend the results from MD to gain new information about lipid-mediated protein-protein interactions.

In an earlier paper, the approach was used to calculate both the protein-induced lateral perturbations of the DPPC bilayer structure and the related lipid-mediated protein-protein interactions (Lague et al., 2000). The proteins were modeled by hard repulsive cylinders for which three size were used: a small cylinder of 2.5 Å radius, corresponding approximately to an aliphatic chain, a medium cylinder of 5 Å radius, corresponding to a poly-alanine -helix, and a larger cylinder of 9 Å radius, representing a small protein such as the gramicidin channel.

The method and the related numerical calculations are applied to several different single component LBMs. The lipids chosen were dioleoyl phosphatidylcholine (DOPC) (18:1;18:1) with two aliphatic chains composed of 18 carbon atoms with a single unsaturated bond, palmitoyloleoyl phosphatidylcholine (POPC)(16:0;18:1), DPPC (16:0;16:0) and finally dimyristoyl phosphatidylcholine (DMPC)(14:0;14:0). The calculation proceeds as follows. First, the lateral density-density response function for DOPC, POPC, DPPC, and DMPC is calculated from MD trajectories. Second, the lateral density-density response function is then used as an input to calculate the lipid density around protein inclusions and the lipid-mediated potential of mean force (PMF) between two identical protein inclusions for both protein sizes and for each single component LBM. In the next section, we briefly present the formulation of the approach. The results for the various LBMs are described in the third section. The paper is concluded with a brief summary of the results and a discussion of future work.

	THEORY AND METHODS

TOP ABSTRACT INTRODUCTION THEORY AND METHODS RESULTS AND DISCUSSION CONCLUSION REFERENCES

2D-HNC integral equation

The theoretical method employed in this work was previously presented by Lague et al. (1998, 2000). Briefly, protein inclusions embedded in a uniform lipid bilayer in the liquid-crystalline phase are considered here. It is assumed that their dominant effect on the LBM can be described through the lateral perturbation of the average structure of the hydrocarbon chains. For the sake of simplicity, protein inclusions are modeled as hard uniform repulsive cylinders of radius  which interact only with the aliphatic chains, whereas the polar head groups are assumed not to be directly affected by the protein.

The theory corresponds to a HNC integral equation projected onto the 2D plane of the lipid bilayer (2D-HNC). The 2D-HNC equation can be written in terms of a pair of coupled integral equations:

(1)

and

(2)

where r  (x, y) is a 2D vector in the membrane plane, U(r) is the repulsive potential between the hard cylinders and the aliphatic chains, is the average hydrocarbon density in 2D of an unperturbed membrane, c(r) is the direct protein-lipid correlation function, h(r)(r)/ is the protein-lipid correlation function and _m(r) is a response function related to the density fluctuations of carbon pairs in the unperturbed LBM. Eq. 2 is the Ornstein-Zernicke (OZ) equation for an isolated impurity in an infinite bulk system (Hansen and McDonald, 1986).

The response function _m is related to the density-density fluctuations of lipid carbon pairs in the unperturbed LBM at equilibrium (Lague et al., 2000). The density-density fluctuations of the carbon pairs can also be expressed in terms of the radial intramolecular, S_m(r) and intermolecular, H_m(r), pair correlation functions of the unperturbed single component LBM (Lague et al., 1998, 2000).

The lipid-mediated potential of mean force (PMF) between two protein inclusions was computed using the OZ route (Lague et al., 2000):

(3)

where u(r) is the direct repulsive potential between protein inclusions and c(r) is the direct protein-lipid correlation function for a single protein inclusion. The PMF has the character of a Helmholtz free energy because the density response function _m was extracted from MD simulations at constant cross-sectional area (Armen et al., 1998; Feller et al., 1997; E. Dolan, R.M. Venable, R.W. Pastor, and B.R. Brooks, in preparation). As explained in Lague et al. (2000), other routes yield essentially the same PMF.

Computational details

In the present study, the pair correlation functions were calculated from configurations generated by MD simulations of detailed atomic models of single component LBMs. Previously published simulations of POPC (Armen et al., 1998) and DPPC (Feller et al., 1997) as well as unpublished simulations of DOPC (E. Dolan et al., in preparation) were used to compute pair correlation functions of the single component LBMs. In addition, data from an MD simulation of DMPC specially performed for this work was used (see below for details). These phospholipids were chosen in order to investigate the effects of chain length as well as unsaturation on lipid-mediated protein-protein interactions. The parameters used in the simulations for the different LBMs used in this study are presented in Table 1. All the heavy atoms from the aliphatic chains and the glycerol backbone were counted in the calculation of the pair correlation functions. The number of atoms per lipid species as well as the average carbon density per unit area are also given in Table 1. The polar head group was not included. The response function _m was calculated as described previously (Lague et al., 1998, 2000).

A 2D discrete grid of N = 1024 × 1024 nodes with a spacing d of 0.12 Å was used to solve numerically the 2D-HNC equations for a system with a single protein inclusion. As previously, an iterative scheme with simple mixing was used to solve the 2D-HNC integral equations self-consistently. Less than 200 iterations were necessary for convergence. The double convolution in Eq. 2 was calculated by using a 2D Fast Fourier transform. All parameters involved in the computation of the PMF using the OZ route were obtained from a single calculation with a system having a single protein inclusion, thus saving computational time.

Molecular dynamics simulations

The MD simulations for DMPC bilayers were performed as follows. The model bilayer consists of 112 DMPC molecules with 56 molecules in each leaflet of the bilayer plus 4590 water molecules for a total of 26,986 atoms (40.98 water molecules per lipid). This constitutes the central unit of a periodic system and the spatial dimension for this system is 60 × 60 × 73ų, giving a cross-sectional area of 64 Å² per DMPC molecule (Nagle, 1993). The membrane normal is oriented along the Z-axis and the center of the bilayer is placed at Z = 0. Periodic rectangular boundary conditions were applied in the XY directions to simulate an infinite planar layer and in the Z direction in order to represent a multilayer system. The MD trajectory was calculated in the microcanonical ensemble with constant energy and volume. The average temperature of the system was set to 330 K, above the gel-liquid crystal phase transition of DMPC (Gennis, 1989). The potential energy function used for the calculations was the all-hydrogen PARAM 22 force field (MacKerrel et al., 1998) of the biomolecular CHARMM program (Brooks et al., 1983), which includes phospholipids (Schlenkrich et al., 1996) and the TIP3P water potential (Jorgensen et al., 1983). The initial configuration was constructed from pre-equilibrated and pre-hydrated lipids as described previously (Woolf and Roux, 1996; Bernèche et al., 1998). The average density profile of the water molecules, the aliphatic chains, and the polar head groups were calculated in order to characterize the structure of the lipid bilayer. These density profiles are almost identical to previous results (Gambu and Roux, 1997; Damodaran and Merz, 1994), indicating that the present simulations are adequate.

	RESULTS AND DISCUSSION

TOP ABSTRACT INTRODUCTION THEORY AND METHODS RESULTS AND DISCUSSION CONCLUSION REFERENCES

Pair correlation functions for unperturbed LBMs

The first step in the analysis consists of calculating the lateral density-density response functions from the MD trajectories. The results are given in Fig. 1. The top panel of Fig. 1 shows that the carbon-carbon intramolecular correlation functions, S_m(r), have a similar shape for each of the single component LBMs. The correlation functions all exhibit a large peak at about 3 Å and then slowly decay over a distance of 10 to 15 Å. The main difference between the various LBMs occurs from r = 0 Å to r = 3 Å while the remainder of the curves are almost identical. As discussed previously (Lague et al., 2000), this short range contribution to the intramolecular correlation arises mainly from nearest neighbor carbons along the aliphatic chains (i.e., carbon i with carbon i  1 and i + 1), though there are also contributions from intermediate range (second neighbor) and long range correlations along the chains themselves. These include correlations between the final carbon atom of the aliphatic chains and the carbon atoms of the glycerol backbone. The large peak is thus an indication of the significant amount of short and long range order in the lipid chains perpendicular to the plane of the bilayer. The long range contribution to the intramolecular correlation functions shown in Fig. 1 extends to 15 Å and is due to carbon atoms located on different aliphatic chains of a single lipid molecule.

View larger version (20K): [in this window] [in a new window]   FIGURE 1   Carbon-carbon density-density response functions m(r) (bottom) as extracted from molecular dynamics simulations of different model LBMs and used as input to the integral equations. The intramolecular Sm(r) (top) and intermolecular, Hm(r) (middle) pair correlation functions are also shown. (in order to match the dimensions of Å2 of Sm(r), the function m/ and the function Hm(r)× are also shown).

The middle panel of Fig. 1 shows that the intermolecular correlation functions, H_m(r), which involve carbons from two different lipid molecules also have a similar shape for each single component LBM. The strong negative contribution at short distances is due to the lipid-lipid core repulsion. As stated in a previous paper (Lague et al., 2000), it is interesting to note that for r = 0 these intermolecular correlation functions have a value higher than the corresponding value for a normal liquid, i.e., a value of 0.44 to 0.53 as compared to the normal value of 1.0 for a bulk liquid (Hansen and McDonald, 1986). This mostly reflects the fact that the lipid in the two leaflets are weakly correlated and that the carbons of different molecules can overlap when their 3D configuations are projected onto the plane of the LBM, resulting roughly in a little less than half the average hydrocarbon density at r = 0. The first peak, around 5 Å, is in good agreement with the carbon-carbon intermolecular pair correlation function of butane (Tobias et al., 1997), at which a characteristic double peak appears in the region between 4 and 6 Å. The contributions to this peak include intramolecular correlations between terminal methyl groups for molecules in the trans conformation as well as intermolecular correlations. This result suggests a 5-Å distance between carbon atoms of two different aliphatic chains, corresponding to a radius of 2.5 Å for a single aliphatic chain.

These considerations show that the response function, _m(r), which is the sum of the intramolecular and intermolecular correlation functions, exhibits the following general features for the various LBMs; a large peak for distances up to 3 Å arising from the intramolecular correlations and a second peak at a distance of 4-5 Å due to the intermolecular correlations, followed by an oscillatory decay with small positive peaks appearing around 9 and 14 Å. The main differences between response functions for different LBMs occur before the first peak at 3 Å and are due to the intramolecular correlations. Only total correlation functions are used as input in the integral equations whereas intra- and intermolecular correlation functions are just used to analyze differences in the structures of the various LBMs.

Perturbed lipid density around protein inclusions

The structure of the correlation functions for the different single component LBMs shown in Fig. 1 suggests that the lateral response to perturbations is similar for each LBM. We examined this response by calculating the protein-lipid radial distribution function, g(r) = h(r) + 1, around protein inclusions using the 2D-HNC formalism of Eqs. 1 and 2 of the previous section. Three different radii of protein inclusion were considered: 2.5 Å, 5 Å, and 9 Å. g(r) is shown as a function of r for the different single component LBMs in Fig. 2. This figure shows that for all cases, the perturbation of the bilayer structure extends to 20 Å away from the edge of the cylinder with strong oscillations separated by ~5 Å. The complex structure reflects the oscillatory behavior of _m(r) shown in the bottom panel of Fig. 1. There are, however, significant differences depending on the size of the protein inclusion and the lipid composition. In the case of the small cylinder, the values of g(r) from the cylinder edge to about 2 Å are greater than the bulk unperturbed value of unity for every model LBM considered here except DMPC. In contrast g(r) for DMPC bilayers is less than unity in this region. For r >3 Å, g(r) oscillates around unity up to a distance of 25 to 30Å for every model LBM considered here. The region next to the small cylinder can therefore be regarded as an enriched layer with a lipid density higher than the uniform bulk value, except for DMPC bilayers where a depletion layer with respect to the average lipid density is predicted.

View larger version (24K): [in this window] [in a new window]   FIGURE 2   Radially averaged lipid density normalized by the respective mean density of the unperturbed lipid memnbrane g(r) around the cylinders as calculated from HNC integral equation (Eq. 6) for 2.5 Å (top), 5.0 Å (middle), and 9.0 Å (bottom) radius cylinders. The distance is calculated from the edge of the cylinders. The mean lipid density of the unperturbed membrane is represented by horizontal lines.

In the case of the medium cylinder, a depletion layer is predicted for DMPC and DPPC bilayers, while there is a small enriched region next to the cylinder for DOPC bilayers. For POPC bilayers, there is a depletion layer with a small peak at 1.3 Å, where g(r) is higher than unity. The enriched region of the DOPC bilayers extends from 0 to 1.7 Å and is followed by a depletion region, until 5 Å. This depletion region is followed by an enriched region, extending from 5 to 25 Å, where g(r) relaxes to unity. For the other LBMs, the depletion region next to the cylinder edge is followed by an enriched region, extending from 5 to 13 Å, depending on the type of LBM, to 22 to 30 Å, where g(r) relaxes to unity. Exactly the same qualitative trend is observed for the lipid density around the larger cylinder of 9 Å radius for every model LBM considered here. In particular, enriched regions are more enriched and depletion regions are more depleted for the case of the larger cylinder. It is important to emphasize the strikingly behavior of DOPC bilayers, where, in contrast to the other LBMs, the average lipid density next to the protein cylinder was always predicted to be higher than its bulk unperturbed value. Note that DOPC is the only lipid molecule with two unsaturated aliphatic chains considered in this work. In addition, the behavior of POPC bilayers in the presence of the medium or the large cylinder is intermediate to DPPC and DOPC bilayer behaviors, i.e., it has a depletion layer next to the cylinder for DPPC bilayers and an enriched region for the DOPC bilayers while it has both depletion and enriched regions next to these cylinders for POPC bilayers.

The depleted or enriched layers next to the protein edge modify the area available for lipid molecules in this region. In order to quantify this effect, the average area per lipid molecule within a layer of a given size next to the cylinder edge was calculated as follows. First, the number of carbon atoms was found by integration of the density of the carbon atoms over the layer surface using the graphs of Fig. 2. Next, the number of lipid molecules is estimated by dividing the number of carbon atoms in the layer by twice the number of carbon atoms per single lipid molecule as given in Table 1 in order to find the number of lipids in one leaflet of the LBM. This is an approximation since carbon atoms can belong to different lipid molecules. Finally, the surface area of the given layer around the protein inclusion is divided by the number of lipid molecules to give the area per lipid molecule. Results for different layer sizes and different protein inclusions are given in Table 2. This table shows that, in general the lipid molecules next to the protein edge have a greater area per lipid molecule than the lipid molecules in the unperturbed LBM. This observation is in agreement with the results of Husslein et al. (1998), who performed MD simulations of a hydrated diphytanol phosphatidylcholine lipid bilayer containing an -helical bundle of four trans-membrane domains of the influenza virus M2 protein. It was observed that the area per lipid molecule in the vicinity of the protein increased to 85 Å² per molecule, as compared to 74.6 Å² per molecule for an unperturbed LBM under the same conditions.

The presence of a depletion layer of lipid molecules around a protein inclusion, and the corresponding increase in cross-sectional area per lipid molecule, suggest that there is an effective long range repulsion between the lipids and the protein. According to this view, a lipid molecule must reduce its disorder significantly to be able to come in close contact with an embedded protein. An effective lipid-protein repulsion arises because it is entropically unfavorable. This view appears, however, to contradict the conclusion of Marsh and Horvath (1998) based on spin labels electron paramagnetic resonance measurements that the mean ordering of protein-associated lipid chains is very similar to that in the bulk liquid-crystalline regions of the LBM. On the other hand, the conjecture that the origin of the protein-lipid repulsion is configurational entropy is supported by MD simulations of the gramicidin channel in which the lipid chains were observed to adopt ordered configurations with a higher carbon-deuterium order parameter in the vicinity of the channel (Woolf and Roux, 1996; Chui et al., 1999). It would be of interest in the future to examine the dynamics of spin-labels with MD simulations to clarify the interpretation of the electron paramagnetic resonance data.

Lipid-mediated interaction between protein inclusions

It has been demonstrated that the lipid-mediated interactions between two protein inclusions arise directly from the perturbation of the average hydrocarbon density around the proteins (Lague et al., 2000). The results for the free energy which describe this interaction were calculated using Eq. 5 for the lipid-mediated PMF and are given in Fig. 3 as a function of the distance between the two cylinders.

View larger version (21K): [in this window] [in a new window]   FIGURE 3   Lipid-mediated interaction energy, in kBT units, between two hard repulsive cylinders of 2.5 Å radius (top), 5 Å radius (middle), and 9 Å radius (bottom) as a function of their displacement between cylinder edges.

In the case of two small cylinders, there is a free energy barrier at a distance of 15 Å between the cylinder edges for DPPC and POPC bilayers followed by an attractive free energy well for distances of <8 to 10Å. For DMPC bilayers, the free energy barrier is at 22 Å, and the attractive free energy well begins at 12 Å. The repulsive barrier is ~0.5 kBT and extends from about 9 Å to 20 Å for DPPC and POPC bilayers and from 12 Å to 27 Å for DMPC bilayers. At separations >20 Å (27 Å for DMPC), the lipid-mediated protein-protein interaction is negligible. Finally, at protein-protein contact, the magnitude of the lipid-mediated potential is in the range of 1.4 to 2.2 kBT. For the case of DOPC, the lipid-mediated interaction oscillates around zero, with an increasing amplitude when the separation distance is decreased, and there is a free energy well of 0.6 kBT at protein-protein contact.

Similar trends are observed in the case of the 5-Å radius cylinders embedded in the POPC and DPPC bilayers. For these bilayers, there is a free energy barrier at a distance of 10 Å followed by an attractive free energy well for distances <3-5 Å. The repulsive barrier is ~4 kBT and extends from 5 to 20 Å. Again, at separations greater than 20 Å, the protein-protein interaction is very small and oscillatory and the magnitude of the lipid-mediated potential is ~6 kBT at protein-protein contact, for POPC bilayers, 10 kBT for DPPC bilayers. Again as in the case of the small cylinders, the results for DMPC bilayers extend on a larger range than for DPPC and POPC bilayers, i.e., the free energy barrier is at 15 Å, and the attractive free energy well begins at 9 Å. The repulsive barrier is approximately 3 kBT and extends from 9 to 25 Å. Finally, at protein-protein contact, the magnitude of the lipid-mediated interaction is 7.5 kBT for DMPC bilayers. For DOPC bilayers, as for other LBMs, a free energy barrier is observed, and is followed by a free energy well. The magnitude of the free energy barrier is 6.5 kBT for a separation distance of 13 Å, and the free energy well at protein-protein contact is 5 kBT. The values of the lipid-mediated interaction at protein contact, which lie between 13 to 5 kBT, are in good agreement with, though somewhat more negative than, previous theoretical estimates in the literature (Marcelja, 1976; Schroder, 1977; Owicki et al., 1978; Owicki and Bloom, 1979; Pearson et al., 1984).

For the case of two large cylinders of 9 Å radius, there is a free energy barrier at a distance of 3 to 5 Å between the the cylinder edges for DOPC, POPC, and DPPC bilayers and, in contrast to the case of two smaller and the two medium cylinders, there is no attractive free energy well for these LBMs. For these cases, the repulsive barrier is ~29 kBT for DOPC bilayers, 20 kBT for POPC bilayers, and 10 kBT for DPPC bilayers, and extends from protein-protein contact to 20 Å. For DMPC bilayers, in contrast to other LBMs, there is a very small attractive free energy well for distances below 2 Å, with a lipid-mediated interaction of 6.0 kBT at contact. The free energy barrier between two large cylinders in DMPC bilayers extends from 25 Å to 2 Å between the the cylinder edges, with a maximum of 9 kBT at 10 Å.

In conclusion, the general response of POPC bilayers tends to be intermediate to DOPC bilayers and DPPC bilayers. In all the figures, the POPC curves are effectively always somewhere between DOPC and DPPC curves. For DMPC bilayers, the general behavior is the same as for DPPC bilayers except that the curves are smoother, i.e., the height of the curves is smaller and they extend to a greater distance than those of DPPC bilayers. Finally, the lipid-mediated interaction between two protein inclusions in DOPC bilayers is always less favorable than in DPPC bilayers. On a speculative note, it seems that an important consequence of the weaker association of two helices in DOPC relative to DPPC could result effectively in an increased flexibility of proteins formed as an assembly of -helices in such environment. This observation is consistent with the familiar idea that LBM composed of unsaturated lipids are more fluid than those composed of saturated lipids (Gennis, 1989).

	CONCLUSION

TOP ABSTRACT INTRODUCTION THEORY AND METHODS RESULTS AND DISCUSSION CONCLUSION REFERENCES

The dependence of lipid-mediated interactions between protein inclusions on protein size and membrane composition was investigated using a theory for examining the structure of the aliphatic chains around protein inclusions embedded in a lipid bilayer. This theory, which is based on the HNC integral equation theory for liquids, was recently developed (Lague et al., 1998, 2000) and uses the exact lateral density-density response function of the hydrocarbon core of LBMs, computed from configurations taken from the MD simulation of different single component lipid bilayers, as an input to the calculations. We first examined the density variation for the different model LBMs in the neighborhood of a single protein modeled as a hard cylinder for three different sizes: a small cylinder of 2.5 Å radius corresponding an aliphatic chain, a medium cylinder of 5 Å radius, corresponding an -helical poly-alanine protein, and a larger cylinder of 9 Å radius, representing a small protein such as the gramicidin channel (Woolf and Roux, 1996). The results showed that the average lipid order is perturbed over a distance of 20 to 30 Å from the edge of the protein. For small cylinders, the average density is higher than its bulk unperturbed value next to the cylinder for every model LBM considered here except DMPC where the average density is lower for this region. In addition, for every model LBM considered here, the average lipid density oscillates about its bulk unperturbed value up to a displacement of 25 to 30 Å. In the case of medium and large cylinders, the depletion layer is observed for every model LBM considered here except the DOPC, where again an enriched region is again present next to the cylinder edge. The enriched region of DOPC is followed by a depleted region, after which the density slowly relaxes to its unperturbed bulk value. Exactly the opposite trend was observed for other bilayers, where the depletion layer next to cylinder edge was followed by the enriched region and the density then relaxed slowly to the unperturbed bulk value at ~22-30 Å.

Next, the lipid-mediated interaction between two identical protein inclusions for the three inclusion sizes of inclusions was computed for every model LBM considered here. The lipid-mediated interaction was found to be repulsive at intermediate range and attractive at short range for two 2.5-Å cylinders in POPC, DPPC, and DMPC bilayers, while it oscillated between attractive and repulsive behavior in DOPC bilayers. For the case of two medium cylinders and for every model LBM considered here, the lipid-mediated interaction was repulsive at an intermediate range but attractive at short range. For the case of the larger cylinder, exactly the same trend was observed for DMPC bilayers, while the lipid-mediated interaction was predicted to be repulsive at all distance for the DOPC, POPC, and DPPC bilayers. The behavior of DMPC bilayers is qualitatively similar to that of DPPC bilayers. In contrast, POPC bilayers have an intermediate behavior between DOPC and DPPC bilayers.

A detailed comparison of the present results with experimental data is difficult. Although the highly specific interactions that drive the formation of helix dimers have been extensively studied (Lemmon et al., 1992; Lemmon and Elgelman, 1994), there is very little information about nonspecific helix-helix interactions in planar bilayers (Marsh and Horvath, 1998; Watts, 1998). One of the few direct measurements of the magnitude of the nonspecific interaction has been obtained recently by Matsuzaki and co-workers (Yano et al., 2001). They designed a 21-residue hydrophobic peptide composed of leucines and alanines and demonstrated that it forms a transmembrane helix. No specific sidechain interactions are expected to favor helix association for this peptide. Using a fluorescence resonance energy transfer technique, they showed that the free energy for helix-helix association was 4.8 kBT in POPC bilayers. This value compares very well with our result of 6.5 kBT for a cylinder of 5 Å (Fig. 3). In the future, it should be particularly interesting to compare theory and experiment for different lipid composition.

Extension of the approach will be aimed at investigating the influence of lipid composition and cholesterol on the stability of proteins formed as bundle of transmembrane -helices. In addition, several developments such as the inclusion of the coupling between lateral and transversal response of the membrane in a three-dimensional form of the the integral equation theory are currently in progress.