PHYSICAL REVIEW E 94, 032406 (2016)

Budding transition of asymmetric two-component lipid domains
Jean Wolff*
Institut Charles Sadron, UPR22-CNRS 23, rue du Loess BP 84047, 67034 Strasbourg Cedex, France

Shigeyuki Komura†
Department of Chemistry, Graduate School of Science and Engineering, Tokyo Metropolitan University, Tokyo 192-0397, JapanDavid Andelman‡
Raymond and Beverly Sackler School of Physics and Astronomy, Tel Aviv University, Ramat Aviv, Tel Aviv 69978, Israel
and CAS Key Laboratory of Theoretical Physics, Institute of Theoretical Physics, Chinese Academy of Sciences, Beijing 100190, China(Received 17 July 2016; published 14 September 2016)
We propose a model that accounts for the budding transition of asymmetric two-component lipid domains,
where the two monolayers (leaflets) have different average compositions controlled by independent chemical
potentials. Assuming a coupling between the local curvature and local lipid composition in each of the leaflets, we
discuss the morphology and thermodynamic behavior of asymmetric lipid domains. The membrane free-energy
contains three contributions: the bending energy, the line tension, and a Landau free-energy for a lateral phaseseparation. Within a mean-field treatment, we obtain various phase diagrams containing fully budded, dimpled,and flat states as a function of the two leaflet compositions. The global phase behavior is analyzed, and dependingon system parameters, the phase diagrams include one-phase, two-phase, and three-phase regions. In particular,
we predict various phase coexistence regions between different morphologies of domains, which may be observedin multicomponent membranes or vesicles.
DOI: 10.1103/PhysRevE.94.032406
I. INTRODUCTION

The cytoplasmic membrane separates the living cell fromits extracellular surroundings, while other intracellular membranes compartmentalize cellular organelles. Biomembranes
are constructed from two monolayers (leaflets) in a backto-back arrangement and are in general asymmetric in theirlipids composition [1,2]. For example, in human red blood
cells, the inner cytoplasmic leaflet is composed mostly ofphosphatidylethanolamine (PE) and phosphatidylserine (PS),while the outer cytoplasmic leaflet is composed of phosphatidylcholine (PC), sphingomyelin (SM), and a varietyof glycolipids [3,4]. The asymmetric nature of the cell
membrane plays a key role in a variety of cellular processes
such as endocytosis [5], vesicle budding, and trafficking [6].
Furthermore, in living cells, the composition asymmetry isan active and energy-consuming process. It is maintained byseveral membrane proteins such as flippase and floppase thatallow lipids to exchange between the two leaflets with the aidof adenosine triphosphate (ATP) [7].
In artificial multicomponent lipid bilayers, the two membrane leaflets can undergo a lateral phase separation. Several
authors made the connection between such a phase separationin artificial membranes and existence of small dynamic
domains (“rafts”) in biological membranes [8]. It should be
noted, however, the size of rafts in biological membranes are
expected to be in the range of 10–100 nm [9,10]. Raft are
believed to be enriched mixtures of cholesterol and SM ina liquid-ordered phase (Lo ), embedded in a background of a*jean.wolff@ics-cnrs.unistra.fr
komura@tmu.ac.jp‡
andelman@post.tau.ac.il†

2470-0045/2016/94(3)/032406(9)liquid-disordered phase (Ld ). Despite the lack of an ultimate
proof for the existence of rafts, they have been advocated inrelation with their potential influence on biological cellular
processes. It has been suggested that rafts act as organizing
centers for the assembly of signaling molecules, influencing
membrane fluidity, and regulating receptor trafficking [9,10].
A theoretical model for domain-induced budding of planar
membranes was proposed by Lipowsky [11] some years
ago, and later was extended for closed vesicles [12,13]. In
the model, the competition between the membrane bending
energy and domain line tension leads to a budding transition.
More recently, a model describing domain-induced budding inbilayers composed of a binary mixture of lipids was proposedby us [14]. In particular, we have shown that dimpled domains
are formed and remain stable due to an asymmetry between thetwo compositions of the corresponding domain leaflets, giventhat the line tension along the domain rim is not too large. The
calculations in Ref. [14] were done for a specific case where
the relative concentration between the two lipids on the bilayer
domain stays constant, while the lipids are allowed to freely
exchange between the two leaflets.
In the past decade, however, techniques such as LangmuirBlodgett or Langmuir-Schäefer have enabled a control over theasymmetric composition of artificial membranes [15,16]. For
example, unsupported bilayers via the Montal-Müller method
have been used to form asymmetric bilayers in which thecomposition of each leaflet was independently controlled [17].
A key ingredient in understanding those experiments is thefact that the flip-flop process that exchanges lipids acrossthe leaflets is slower than experimental times. Given theseexperimental findings, it is worthwhile to consider bilayers
where each of the leaflet domain composition (rather than theoverall bilayer composition) can be controlled in a separate
and independent way.

032406-1©2016 American Physical Society

JEAN WOLFF, SHIGEYUKI KOMURA, AND DAVID ANDELMAN

FIG. 1. Schematic vertical cut through a membrane consistingof several domains. The domains are embedded in a lipid matrix(yellow) and have a line tension acting along the domain rim. Three
domain types can be seen and are further explained in the caption ofFig. 2: flat (F), budded (B), and dimpled (D). Each domain is formedby two lipids (red and blue) that partition themselves differently on
the two domain leaflets.

In this paper, we generalize our previous budding
model [14] and extend it to asymmetric two-component lipid
domains. Each domain leaflet has a conserved lipid composition that is independent from that of the other leaflet. We consider the possibility of domains curved in the third dimension,
which can produce buds as shown in Fig. 1. In our model,
the composition-dependent spontaneous curvature leads tocoupling between curvature and lipid composition in each of
the domain leaflets [18–20]. Such a coupling leads to richphase behavior including various phase coexistence regions.
For the sake of clarity, we do not take into account any directinteraction between two domains that face each other, althoughvarious possibilities have been previously proposed [21].
In addition, we would like to mention that in a separate setof studies for completely planar membranes, asymmetricbilayers composed of two modulated monolayers (leaflets)
were considered [22,23]. In these works, the static anddynamic properties of concentration fluctuations have beenpresented together with the related microphase separation.
Our model may have several experimental implications.
We predict that bilayer domains can exist in three states having
different equilibrium shapes: fully budded, dimpled, and flatstates. Their relative stability depends on controlled systemparameters: temperature, degree of compositional asymmetry
between the two leaflets, and domain size. We find that the
dimpled state is the most stable one in some of the parameter
range and this is in accord with recent experiments [24].
Based on the calculated phase diagrams, we anticipate that
membranes should exhibit in some parameter range two-phaseand three-phase coexistence between different domain states.
The outline of this paper is as follows. Section II generalizes
our previous budding model [14], and the free-energy describing asymmetric two-component lipid domains is discussed. In
Sec. III, we explain the conditions for various phase equilibria
and how to calculate the phase diagrams. We then proceedby presenting the phase diagrams in Sec. IV, and discuss theresulting global phase behavior. Finally, a more qualitative
discussion is provided in Sec. V.
II. MODEL

We consider a membrane consisting of two monolayers (theterms “monolayer” and “leaflet” will be used interchangeably

PHYSICAL REVIEW E 94, 032406 (2016)

in this paper), each composed of an A/B mixture of lipids,
which partition themselves asymmetrically between the twoleaflets. We assume that the membrane can undergo a lateralphase separation creating domains of different lipid composition. The domains are taken to fully span the two monolayers,
but the leaflet compositions in these domains can be different.
Hence, the formed lipid domains are, in general, asymmetric.
We further assume that the two leaflet compositions areconserved and can be taken as independent from each other.
In the following, we discuss the thermodynamic behavior ofa system consisting of a large number of such asymmetricdomains as shown in Fig. 1. For this purpose, we explain below
the different terms that contribute to the domain free-energy.
In Fig. 2, we show a vertical cut through three possible
domain states: flat (F), dimpled (D), and fully-budded (B). The
flat circular domain (F) in (a) has an area, S = π L2 , whichis assumed to remain constant during the budding process.
For simplicity, we consider in (b) dimpled buds (D) whose
shape is a spherical cap of radius 1/C, and in (c) the extreme
case of a completely detached spherical bud (B). The totalbending energy of the budded domain is given by the curvature
contributions from its two monolayers [25,26],
Ebend = 2π L2 κ (C − C0 )2 + (C + C0 )2 ,(1)
where κ is the bending rigidity modulus (assumed to be aconstant that is independent of lipid composition) and C0 themonolayer spontaneous curvature. As shown in Fig. 2, the twomonolayer deform with curvatures +C and −C, respectively.

FIG. 2. A schematic vertical cut through three possible domain
shapes: (a) a flat bilayer domain composed of a mixture of A andB lipids (red and blue, respectively), and embedded in an otherwise
flat membrane (yellow). The circular flat domain (F) has a radius
L and area S = π L2 . (b) A partial bud (dimpled domain) curved in
the out-of-plane direction. The bud (D) of the same area S forms aspherical cap of radius 1/C, where C is the curvature. (c) A fully
budded domain (B) has a spherical shape of total area S and justtouches the flat membrane. The line tension γ acts along the line
boundary between the domain (red/blue lipids) and the flat membrane
matrix (yellow lipids).

032406-2

BUDDING TRANSITION OF ASYMMETRIC TWO- . . .

PHYSICAL REVIEW E 94, 032406 (2016)

The second energy contribution is the domain edge energythat is proportional to the perimeter length and to the line
tension, γ [11]:
Eedge = 2π Lγ 1 − (LC/2)2 .
(2)Note that in the extreme case, when the domain buds into acomplete spherical shape (B) as in Fig. 2(c), C = ±2/L and
Eedge = 0. In the above, the strong variation in composition
between the domain and its surrounding matrix is effectively
taken into account through the line tension, γ , which is treatedas an external control parameter. This situation can be justified
for a strong segregation that results in a sharp boundary
between the domain and its flat matrix surroundings.
As each domain is composed of an A/B mixture, we defineφA (φB ) as the area fraction (assumed to be equal to the molar
fraction) of the A lipid (B lipid) in the upper leaflet domain, andsimilarly, ψA (ψB ) for the lower leaflet domain. We assume
that each monolayer is incompressible so that φA + φB = 1
and ψA + ψB = 1. Hence, the two relevant order parameters
are the relative composition in the upper leaflet,
φ = φA − φB ,(3)dimensionless form, ε = Etot /2π κ, is expressed asε(φ,ψ,c) = (c + βφ)2 + (c − βψ)2 +
+

 

t 2
1 U
1 4
2
4
(φ
(φ
+
ψ
)
+
+
ψ
)
.
ξ 2 2κ2
42(c + βφ) + 2(c − βψ) −

c

= 0.
4ξ 1 − c2 /4

As in any A/B mixture, the possibility of a phase separation
can be described by a phenomenological Landau expansion of
the free-energy in powers of φ and ψ. This expansion is done
separately for each monolayer, and the total contribution to thefree-energy is the sum over the two monolayers,


U t
1
(5)
Ephase = π L2 2 (φ 2 + ψ 2 ) + (φ 4 + ψ 4 ) , 2
4where  ≡ κ/γ is the invagination length, U is a parameterthat sets the energy scale, and t ∼ (T − Tc )/Tc is the reduced
temperature (Tc being the critical temperature). In Eq. (5)above, we multiply by the domain area, π L2 , to obtain thedomain free-energy.
Hereafter, we will use several dimensionless variables: arescaled curvature c ≡ LC, rescaled spontaneous curvaturec0 ≡ LC0 , and rescaled invagination length ξ ≡ /L. The
coupling between curvature and composition is taken into
account by assuming a linear dependence of the spontaneous
curvature c0 on the relative composition in each of theleaflets [18–20],
c0 (φ) = c̄0 − βφ,(6)c0 (ψ) = c̄0 − βψ,(7)where c̄0 is the monolayer spontaneous curvature for thesymmetric 1:1 composition, φ = ψ = 0, and β is a coupling
parameter that has the same value for the two monolayers.
Since c̄0 is a constant, it merely shifts the origin of the chemical
potential, and will be dropped out without loss of generality.
The total free-energy per domain is then given by the sumof Eqs. (1), (2), and (5): Etot = Ebend + Eedge + Ephase , and its(9)Then, by substituting the curvature c = c(φ,ψ) into Eq. (8),results in a partially minimized free-energy, ε∗ε∗ (φ,ψ) = ε(φ,ψ,c(φ,ψ)),

(4)

(8)

We note that Eq. (8) depends on three dimensionless parameters: β, ξ , and U/(2κ), while the thermodynamic variables
are the reduced temperature t and the three order parameters:
φ, ψ, and c. In the calculations presented hereafter, we set
U/(2κ) = 1 and vary the values of β and ξ .
Within mean-field theory, the equilibrium states and phasetransitions are determined by minimizing ε with respect toφ, ψ, and c, under the condition that φ and ψ are conservedorder parameters while c is not. From the minimization of εwith respect to c, we obtain the condition

and in the lower leaflet,ψ = ψA − ψB .

1
1 − c2 /4
ξ(10)

as a function of φ and ψ.
Typical experimental values of domain sizes are in the rangeof L  50–500 nm [9], the bending rigidity κ  10−19 J [26],and the line tension γ  0.2–6.2 × 10−12 J/m [27,28]. Hence,the scaled invagination length is estimated to be in the range,ξ  0.01–10. These values will be used in the next section
where we calculate numerically the phase diagrams.
III. PHASE EQUILIBRIA CONDITIONSIn order to obtain various phase coexistence regions, ε∗in Eq. (10) should be further minimized with respect to theconserved order parameters, φ and ψ. Hence, we consider thefollowing thermodynamical potential:g(φ,ψ) = ε∗ (φ,ψ) − μφ φ − μψ ψ,(11)where μφ and μψ are the chemical potentials coupled with the
A/B relative compositions in the upper and lower domains,respectively. They act as Lagrange multipliers that take into
account the conserved φ and ψ compositions. In general, thesetwo chemical potentials have different values, μφ = μψ . The
special case of μφ = μψ , for which only the total relative
composition, φ + ψ, is conserved was investigated in ourprevious work [14], while here we deal with a general situation
where each of the compositions, φ and ψ, are conserved
independently.
The thermodynamic equilibrium between the two coexisting phases denoted as “1” and “2” and characterized by (φ1 ,ψ1 )and (φ2 ,ψ2 ), satisfies the conditions [29]032406-3∂φ g(φ,ψ) 1 = ∂φ g(φ,ψ) 2 = 0,∂ψ g(φ,ψ) 1 = ∂ψ g(φ,ψ) 2 = 0,g(φ1 ,ψ1 ) = g(φ2 ,ψ2 ).

(12)

JEAN WOLFF, SHIGEYUKI KOMURA, AND DAVID ANDELMAN

FIG. 3. (a) Plot of the partially minimized free-energy ε∗ (φ,ψ)
as a function of φ and ψ for ξ = 0.25, β = 1, and t = −0.2 (below
Tc ). (b) A cut through the free-energy landscape ε∗ in the particular
direction, φ + ψ = 0, plotted as a function of= (φ − ψ)/2.
The two P points are cusps at which the slope of ε∗ changes
discontinuously, although the Hessian H remains positive (seeAppendix). The two Q points correspond to the location at which
H vanishes.

Similarly, for a three-phase coexistence between phases “1,”“2,” and “3,” the following set of conditions should be satisfied:∂φ g(φ,ψ) 1 = ∂φ g(φ,ψ) 2 = ∂φ g(φ,ψ) 3 = 0,∂ψ g(φ,ψ) 1 = ∂ψ g(φ,ψ) 2 = ∂ψ g(φ,ψ) 3 = 0, (13)g(φ1 ,ψ1 ) = g(φ2 ,ψ2 ) = g(φ3 ,ψ3 ).
In Fig. 3(a), we show an example of the partially minimized
free-energy, ε∗ (φ,ψ), as a function of φ and ψ at a fixed
temperature t = −0.2 (below Tc ), and for given values of ξand β. In order to have a better view of the free-energy surface,
we show in Fig. 3(b) a cross-section cut of the free-energy
surface in the direction of= (φ − ψ)/2, while keeping
φ + ψ = 0. Here, we see two singular cusps (points P) where
the determinant of the Hessian matrix, H (see Appendix)
does not vanish. At these cusps, H changes discontinuously
although it remains positive. At points Q, on the other hand,the Hessian H vanishes. More details on the Hessian matrix
and determinant, and their relation to the phase stability and

PHYSICAL REVIEW E 94, 032406 (2016)

spinodal lines, are presented in the Appendix. By calculation
the Hessian H , it is possible to derive the spinodal lines andcritical points. Note that in most cases, the critical points andspinodal will not be shown on the phase diagrams, because
they are preempted by the first-order phase transition lines andcoexistence regions.
The phase diagrams are obtained by further minimizing ε∗ ,with respect to the two independent variables, φ and ψ. Convex
regions of the free-energy correspond to single thermodynamical phases. Two-phase coexistence regions correspond tononconvex regions, where we can construct a common tangent
plane. The plane touches the free-energy surface at two points
that determine the two phases in coexistence. A more special
three-phase coexistence region corresponds to a plane thattouches the free-energy surface at three points. More details
on the numerical procedure of finding the phase diagrams are
given below.
The numerical computation of the phase diagram is performed using a public-domain software called “Qhull” [30].
The Qhull software generates initially a fine grid of triangulation, which approximates the free-energy surface, ε∗ . The
three-phase coexistence regions correspond to facets with allsides being much larger than the initial discretization. The
two-phase coexistence regions are associated with elongatedtriangles having one short side that is much smaller thanthe other two longer sides that approximate the tie-lines.
Finally, small triangular facets of the free-energy surface areassociated with stable one-phase regions. The projection of thetriangulated free-energy surface onto the composition plane
provides a systematic approximation for the phase diagram.
The Qhull results are then used as an initial condition in
calculating more precisely the equilibrium phase diagrams,
including the various phase coexistence, Eqs. (11)–(13).

IV. RESULTS
A. Phase diagrams

Four representative types of phase diagrams are shown in
Figs. 4–7. In the first two figures, we show the calculated phase
diagram for β = 1 (Fig. 4) and β = 3 (Fig. 5), while the otherparameters are fixed to ξ = 0.25 and t = −0.45. The latter
two phase diagrams are for t = −0.2 (Fig. 6) and t = −0.02(Fig. 7), while keeping β = 1 and ξ = 0.25. For presentation
purposes, in (b) of each figure, we have superimposed thecolor plot for the curvature c on the phase diagram presentedin (a). As denoted above, the various stable phases are thefully budded phase (B) with curvature |c| = 2, the dimpled or
partial budded phase (D) with curvature 0 < |c| < 2, and theflat membrane (F) with zero curvature c = 0. The subscripts± denote whether the bud is curved positively or negatively
with respect to the positive normal direction of the planar
membrane.
Figure 4 includes five homogeneous phases: B± , D± , and F.
The B+ and B− occupy the regions around the φ = −ψ = 1
and φ = −ψ = −1 corners, respectively. On the other hand,
the D± occupy the two remaining corners: φ = ψ = ±1. The
difference between D+ and D− is only associated with theircurvature (c ≶ 0), which changes continuously. This is further
clarified in part (b) of the figure. For all t values, the flat

032406-4

BUDDING TRANSITION OF ASYMMETRIC TWO- . . .

PHYSICAL REVIEW E 94, 032406 (2016)

FIG. 4. (a) Phase diagram in the (φ,ψ) plane for β = 1, ξ = 0.25,
and t = −0.45. The corners of the phase diagram indicates the four
one-phases: D± and B± , while the flat F phase (not drawn) strictly lies
only on the diagonal φ = ψ. The black lines represent tie-lines in the
two-phase regions, and the two triangles are the three-phase regions(see text for more details). (b) The phase diagram is plotted as in (a)but with a superimposed colored plot for the curvature, c. As one
crosses the major diagonal, φ = ψ, there is a smooth change fromD+ (c > 0) through the flat F (c = 0) to the D− (c < 0). Furthermore,
the curvature also changes smoothly inside the D± phases, but thegradient in orange (D− ) and light blue (D+ ) colors is not shown forclarity. However, the curvature has a jump between B− (c = −2, red)and D− (−2 < c < 0, orange) regions, as well as between B+ (c = 2,
dark blue) and D+ (0 < c < 2, light blue) ones.

F phase (not seen on the figure) strictly exists only on theφ = ψ diagonal, but is important for the multiphasecoexistence regions (as discussed below). As one approaches
this line from above or below, the D± phases change oversmoothly into the F phase (with c = 0) on the diagonal line.
Moreover, five two-phase coexistence regions are shown in
Fig. 4 together with their calculated tie lines: two B+ /D+ and
two B− /D− along the boundaries of the phase diagram, andone B+ /B− along the major diagonal, φ = −ψ. In addition,

FIG. 5. (a) Phase diagram in the (φ,ψ) plane as described in thecaption of Fig. 4 but with β = 3, ξ = 0.25, and t = −0.45. (b) The
phase diagram is plotted as in (a) but with a superimposed colored
plot for the curvature, c.

two three-phase coexistence regions: B+ /B− /F can be seen.
These are the two triangular regions lying above and below theφ = −ψ diagonal. Note that the F corner of the three-phaseregion lies close to the D± corners, but since it lies on theφ = ψ diagonal, it is identified as the F phase with c = 0.
The three-phase coexistence region between F and B± phasesmeans that each point inside the triangular region is composedof three relative area fractions of the three coexisting phases:
the flat (F) and budded (D± ) phases. The nearly horizontal(D− /B− ) or vertical (D+ /B+ ) tie lines on the boundaries ofphase diagram indicate that the two φ and ψ monolayers arealmost decoupled, because either φ or ψ do not vary along thetie line. On the other hand, tie lines that lie along the major
diagonal, φ = −ψ, indicate a strong coupling between the φand ψ monolayers in the B+ /B− coexistence region.
In Fig. 5 with β = 3 and the same temperature as in Fig. 4,
we see that the central binary coexistence region (B+ /B− )becomes much larger. On the other hand, the four two-phasecoexistence regions of Fig. 4 have shrunken because the extent

032406-5

JEAN WOLFF, SHIGEYUKI KOMURA, AND DAVID ANDELMANPHYSICAL REVIEW E 94, 032406 (2016)

FIG. 6. (a) Phase diagram in the (φ,ψ) plane as described in thecaption of Fig. 4 but with t = −0.2, ξ = 0.25, and β = 1. (b) The
phase diagram is plotted as in (a) but with a superimposed colored
plot for the curvature, c.

FIG. 7. (a) Phase diagram in the (φ,ψ) plane as described in thecaption of Fig. 4 but with t = −0.02, ξ = 0.25, and β = 1. The
red circles correspond to the critical points. (b) The phase diagramis plotted as in (a) but with a superimposed colored plot for thecurvature, c.

of the dimpled phase becomes smaller for larger coupling
parameter, β. When β gets large values in Eqs. (6) and (7),the composition of each monolayer induces higher curvaturethat promotes budding. In this case, the central two-phaseregion, B+ /B− , occupies a large fraction of the phase diagram,and most of its tie-lines are parallel to the diagonal φ = −ψ,suggesting a strong coupling between the two monolayers.
On the other hand, the two three-phase coexistence regions,
B+ /B− /F, become smaller in Fig. 5. Furthermore, we note thatin Fig. 5(b), the regions of the dimpled phases, D± , represented
by the light blue and orange regions are narrower as comparedwith those of Fig. 4(b).
It is of interest to explain how one of the monolayers
induces a phase transition in the second monolayer. As an
illustrative example, let us consider a point in Figs. 4 and 5with average leaflet composition, (φ,ψ) = (0, −0.9). At thesecompositions, the bilayer separates in Fig. 4 into a D−phase with compositions (−0.75, −0.9) and a B− phase with
(0.75, −0.9). In this weak coupling case, β = 1, the phaseseparation in one monolayer does not induce any instability
leading toward phase separation in the second monolayer,
because the tie-line is nearly parallel to the horizontal φ axis.
On the other hand, in Fig. 5 with a larger value of β = 3,
the bilayer separates into a B− phase with (0.1, −0.9) anda B+ phase with (−0.9,0.1). As the tie-line in this case lies
along the major diagonal, φ = −ψ, the two monolayers are
influencing each other. Such a situation results from a strong
composition-curvature coupling when the parameter β is large
enough. A more general dependence of the phase diagram on
the parameter β will be further discussed in the next subsection.
Figure 6 is plotted for a higher temperature t = −0.2 thanin Figs. 4 and 5, while we fix β = 1 as in Fig. 4. The dimpled
region expands both toward the corners and the middle of thephase diagram. Moreover, there are two new one-phase regionsof the dimpled phase (D± ) resulting in four additional twophase coexistence regions: two D+ /D+ and two D− /D− . The
system exhibits a first-order phase transition in composition,
while the transition is second-order in curvature.

032406-6

BUDDING TRANSITION OF ASYMMETRIC TWO- . . .

PHYSICAL REVIEW E 94, 032406 (2016)

the (1/ξ,β) plane for t = −0.25 in (a) and t = −0.02 in (b).
We show how the four different types of phase diagrams of
Figs. 4, 5, 6, 7, labeled as I, II, III, and IV, respectively, evolve
as function of ξ and β.
Figure 8(a) summarizes the result for t = −0.25 and isvalid even for lower t values. For β > 2 (strong coupling), thediagram is that of type II in which the dimpled phase (D) region
shrinks while the budded phase (B) expands. For β < 2 (weakcoupling), on the other hand, the diagram is mostly of type I for
which the wide three-phase coexistence occurs. For β < 1 and1/ξ < 4, type III diagram is found, and contains coexistence
regions of the dimpled phases (D+ /D+ and D− /D− ).
In Fig. 8(b), for t = −0.02, which is closer to the criticaltemperature, type III is replaced by type IV for β < 1 and1/ξ < 5, while type II also extends to smaller β when 1/ξ
is large enough. When the temperature is higher, the phasetransition in curvature from the flat (F) phase to the dimpled(D) one, becomes continuous [14].
V. DISCUSSION

FIG. 8. Schematic phase-behavior plot in the (1/ξ,β) plane fortemperatures (a) t = −0.25 and (b) t = −0.02. The diagram types, I,
II, III, and IV, correspond to representative examples as in Figs. 4, 5, 6,
and 7, respectively.

In Fig. 7, the temperature is increased to t = −0.02, while
β and ξ stay as in Fig. 4. The chosen t value is higher than in
the previous figures, and approaches the critical temperature,
tc = 0.04 [14]. The central region of the phase diagram isdominated by the dimpled phase (D± ), while the budded
regions (B± ) exist only close to the two corners, with twophase coexistence regions, B+ /D+ and B− /D− . As seen in(b), the regions of the dimpled phases (D± ) are determinedby β, and do not depend on the temperature t as long as ξ is
fixed. Furthermore, there are four critical points appearing in
the central region (marked by red circles), for which the system
exhibits a second-order phase transition both in curvature andcomposition.
B. Global phase behavior

Next, we investigate how the different phase-diagram types
appear and change as we adjust the system parameters in a
global way. In Fig. 8, we present the global phase behavior in

In this paper, we have discussed a phenomenologicalmodel that accounts for the budding transition of asymmetrictwo-component lipid domains, where the two domains havein general different average compositions, represented by φand ψ. Assuming a linear composition dependence of the
spontaneous curvature, we have taken into account a coupling
between the local curvature and local lipid composition ineach of the two leaflets. We then explored the morphologicalchanges between the flat and budded domains by using athermodynamic argument. Our free-energy model contains
three contributions: bending energy, accounting for domain
deformation in the normal direction; line tension along therim of the budded or flat domain; and a Landau free-energy
expansion, which accounts for a phase separation of thetwo-component lipid domains. We have assumed, in addition,that the domain area remains constant during the budding
process.
Our model predicts three different states for the domains:
fully budded (B), dimpled (D), and flat (F) states. In particular,
in some ranges of parameters, the D state is found to be the moststable one, as observed in the experiment [24]. Within meanfield theory, we have calculated various phase diagrams interms of two compositions for different temperatures t, domain
size 1/ξ , and coupling parameter β. The resulting phase behavior is very rich. The calculated phase diagrams include various
one-phase, two-phase (e.g., B+ /D+ , B− /D− , or B+ /B− ), andthree-phase coexistence regions (e.g., B+ /B− /F), dependingon the curvature-composition parameter β as well as thetemperature. Finally, four different types of phase diagram
morphologies are found, and we have analyzed the global
phase behavior in terms of the coupling parameter β, and thedomain size ξ . The model analysis suggests that the asymmetry
in the lipid composition between the two leaflets can lead tocomplex morphological and thermodynamic behavior of lipid
domains.
The most important mechanism that leads to the interleaflet
correlated phase separation, such as the two-phase coexistence
seen in Fig. 5 between two budded phase (B+ /B− ), is
the composition-dependent spontaneous curvature, introduced

032406-7

JEAN WOLFF, SHIGEYUKI KOMURA, AND DAVID ANDELMAN

in Eqs. (6) and (7). With this mechanism, the coupling
between domain curvature and its composition is controlled
by the parameter β. Although we did not include any directinteraction between the domains occurring on the two leaflets,
the fact that the two domains are in full registry to eachother and have opposite curvatures, ±C, results in a strong
correlation between the opposing domains. We consider thatsuch an effect generally exists in biomembranes, even in theabsence of specific proteins that maintain the compositional
asymmetry between the two leaflets.
When the coupling parameter β is made larger at a fixed
temperature (compare Figs. 4 and 5), the two-phase coexistence region, B+ /B− , dramatically expands due to the coupling
between the two leaflets. Such a strongly correlated phaseseparation can be in general observed for lower temperatures.
We also note that the stable budded (B) phase occupies the two
asymmetric corners of the phase diagram, i.e., φ = −ψ = ±1.
Moreover, the region of the dimpled (D) phase increases as
the temperature is raised toward the critical temperature frombelow (compare Figs. 4, 6, and 7). Especially, in Fig. 7, thedimpled one-phase region dominates the central region of thephase diagram and four different critical points are expected
to appear. In the intermediate temperature, as in Fig. 6, thephase diagram contains several types of two-phase coexistence
regions.
The present work for asymmetric domains is a generalization of our previous model [14], where only the averagecomposition of the domains in the two leaflets is controlled bya single chemical potential. Here, we have considered a moregeneral situation, where each of the two domain compositions
is controlled in a separate and independent way. This is doneby introducing two independent chemical potentials coupled
to the two domain compositions, φ and ψ, as in Eq. (11).
Hence, each domain has a conserved lipid composition that is
independent from the other one. The resulting phase diagrams
show a much richer phase behavior, while the previous
results [14] can be recovered by considering the special caseof φ + ψ = const.
Our assumption that the two opposing domains in different
leaflets are correlated to each other is in accord with several experimental observations. For studies of planar but asymmetriccomposition in the two leaflets, Collins et al. [17] reported thatin some cases, one leaflet can induce a phase separation in theother leaflet, depending on local lipid composition. Whereas
in other cases, the two leaflets do not interact. A similar experimental phenomenon was observed for membranes composedof lipids extracted from biological cells [15]. Such situations
for asymmetric membranes can be partially described by thedifferent phase behaviors of asymmetric membrane dependingon the coupling parameter β as shown in Fig. 8. Moreover,both in our model and in experiments [17], domain-induced
processes take place in lipid membrane without any proteins.
These results may suggest a cellular mechanism for regulating
protein function by modulating the local lipid composition orinterleaflet interactions.
Although we have mainly discussed the domain-induced
budding, our model can be applied to describe the formationof vesicles in mixed amphiphilic systems [31]. It was observedin experiment that mixtures of anionic and cationic surfactantsin solution form disk-like bilayers in some range of the relative

PHYSICAL REVIEW E 94, 032406 (2016)

amphiphilic composition. As these disk-shaped bilayers growin size, they transform into spherical caps and eventually
become spherically closed vesicles. In such cases, the spontaneous curvature of bilayer membranes may be induced due tothe compositional asymmetry between the two monolayers.
Finally, our model suggests that the asymmetry in thelipid composition between the two leaflets leads to a complex
behavior of lipid domains even in the absence of any specificenzymes or proteins, which can induce additional coupling
between the two leaflets [32]. The importance of such apure physical mechanism can be verified in experiments onasymmetric model membranes involving only a lipid mixture.
They also can be of relevance to signal transduction [33],membrane fusion [34], or penetration of viruses into cells [35].
We hope that additional experiments will address these issues
in the future.
ACKNOWLEDGMENTS

J.W. acknowledges support from the Service deCoopération Scientifique et Universitaire de l’Ambassade deFrance en Israël, the French ORT association, and ORT schoolof Strasbourg. S.K. acknowledges support from the Grant-inAid for Scientific Research on Innovative Areas “Fluctuationand Structure” (Grant No. 25103010) from the Ministry of
Education, Culture, Sports, Science, and Technology of Japan,
and the Grant-in-Aid for Scientific Research (C) (Grant No.
15K05250) from the Japan Society for the Promotion of
Science (JSPS). D.A. acknowledges the hospitality of the
KITPC and ITP, Beijing, China, and partial support fromthe Israel Science Foundation (ISF) under Grant No. 438/12,
the United States–Israel Binational Science Foundation (BSF)
under Grant No. 2012/060, the ISF-NSFC joint research
program under Grant No. 885/15, and a CAS President’s
International Fellowship Initiative (PIFI, China).
APPENDIX: THE HESSIAN AND STABILITY ANALYSISIn order to discuss the stability of the free-energy, we
consider the 2 × 2 Hessian matrix of ε∗ (φ,ψ) given by∂φφ ε∗
H=∂ψφ ε∗

∂φψ ε∗,∂ψψ ε∗(A1)

where ∂ij ε∗ are the second-order partial derivatives of ε∗ .
We recall that by minimizing ε(φ,ψ,c) with respect to c, we
obtained ε∗ as in Eq. (10). This leads to ∂c ε = 0, ∂φ ε∗ = ∂φ ε,and ∂φ c = −∂cφ ε/∂cc ε, as well as similar expressions for theψ derivatives. Using these relations, one can show that thecomponents of H are given by

032406-8∂φφ ε∗ = ∂φφ ε −(∂cφ ε)2,∂cc ε∂ψψ ε∗ = ∂ψψ ε −(∂cψ ε)2,∂cc ε∂φψ ε∗ = ∂ψφ ε∗ = ∂φψ ε −(∂cφ ε)(∂cψ ε)
.
∂cc ε

(A2)

BUDDING TRANSITION OF ASYMMETRIC TWO- . . .

PHYSICAL REVIEW E 94, 032406 (2016)

At the spinodal condition, the Hessian defined as the determinant of the matrix, H = det H, vanishes. This condition canbe written as


1
4β 2
2β 2 + 2 (t + 3φ 2 ) −
ξ
∂cc ε


16β 4
1
4β 2
2
2
−
= 0, (A3)× 2β + 2 (t + 3ψ ) −ξ∂cc ε(∂cc ε)2

The critical point can be derived by considering another2 × 2 matrix
∂φφ ε∗ ∂φψ ε∗

,
(A5)
H =
∂φ H
∂ψ H
and its determinant H  = det H . One can explicitly show that∂φ H = (∂φφφ ε∗ )(∂ψψ ε∗ ) − (∂φψ ε∗ )(∂φφψ ε∗ )+ (∂φφ ε∗ )(∂φψψ ε∗ ) − (∂φφψ ε∗ )(∂φψ ε∗ ),∂ψ H = (∂φφψ ε∗ )(∂ψψ ε∗ ) − (∂φψ ε∗ )(∂φψψ ε∗ )

where

(A6)

+ (∂φφ ε∗ )(∂ψψψ ε∗ ) − (∂φψψ ε∗ )(∂φψ ε∗ ).
∂cc ε = 4 −1.
4ξ (1 − c2 /4)3/2(A4)Then, the conditions for the critical point are given by [29]H = 0 and

[1] P. F. Devaux and R. Morris, Traffic 5, 241 (2004).
[2] G. van Meer, F. R. Voelker, and G. W. Feigenson, Nat. Rev. Mol.
Cell. 9, 112 (2008).
[3] J. A. F. Op den Kamp, Annu. Rev. Biochem. 48, 47 (1979).
[4] P. F. Devaux, Biochemistry 30, 1163 (1991).
[5] T. Pomorsky, R. Lombardi, H. Riezman, P. F. Devaux, G. vanMeer, and J. C. Holthuis, Mol. Biol. Cell. 14, 1240 (2003).
[6] D. L. Daleke, J. Biol. Chem. 282, 821 (2007).
[7] D. L. Daleke, J. Lipid Res. 44, 233 (2003).
[8] See e.g., S. Komura and D. Andelman, Adv. Colloid Interface
Sci. 208, 34 (2014), and references thererin.
[9] K. Simons and E. Ikonen, Nature 387, 569 (1997).
[10] S. Munro, Cell 115, 377 (2003).
[11] R. Lipowsky, J. Phys. II (France) 2, 1825 (1992).
[12] F. Jülicher and R. Lipowsky, Phys. Rev. Lett. 70, 2964
(1993).
[13] F. Jülicher and R. Lipowsky, Phys. Rev. E 53, 2670 (1996).
[14] J. Wolff, S. Komura, and D. Andelman, Phys. Rev. E 91, 012708(2015).
[15] V. Kiessling, J. M. Crane, and L. K. Tamm, Biophys. J. 91, 3313
(2006).
[16] H.-T. Cheng, Megha, and E. London, J. Biol. Chem. 284, 6079
(2009).
[17] M. D. Collins and S. L. Keller, Proc. Natl. Acad. Sci. U.S.A.
105, 124 (2008).
[18] S. Leibler and D. Andelman, J. Phys. (France) 48, 2013 (1987).
[19] S. A. Safran, P. Pincus, and D. Andelman, Science 248, 354(1990).

H  = 0.

(A7)[20] F. C. MacKintosh and S. A. Safran, Phys. Rev. E 47, 1180
(1993).
[21] S. May, Soft Matter 5, 3148 (2009).
[22] Y. Hirose, S. Komura, and D. Andelman, ChemPhysChem 10,2839 (2009).
[23] Y. Hirose, S. Komura, and D. Andelman, Phys. Rev. E 86,
021916 (2012).
[24] J. E. Rim, T. S. Ursell, R. Phillips, and W. S. Klug, Phys. Rev.
Lett. 106, 057801 (2011).
[25] W. Helfrich, Z. Naturforsch. 28c, 693 (1973).
[26] S. A. Safran, Statistical Thermodynamics of Surfaces, Interfaces, and Membranes (Addision Wesley, Reading, 1994).
[27] T. Baumgart, S. T. Hess, and W. W. Webb, Nature 425, 821(2003).
[28] A. Tian, C. Johnson, W. Wang, and T. Baumgart, Phys. Rev.
Lett. 98, 208102 (2007).
[29] J. W. Gibbs, The Scientific Papers of J. Willard Gibbs, Vol. 1
Thermodynamics (Longmans, Green, and Co., London, 1993).
[30] C. Barber, D. Dobkin, and H. Huhdanpaa, ACM Trans. Math.
Soft. 22, 469 (1996).
[31] A. Shioi and T. A. Hatton, Langmuir 18, 7341 (2002).
[32] G. Gri, B. Molon, S. Manes, T. Pozzan, and A. Viola, Immunol.
Lett. 94, 247 (2004).
[33] K. Simons and D. Toomre, Nat. Rev. Mol. Cell. Biol. 1, 31(2000).
[34] N. Puri and P. A. Roche, Traffic 7, 1482 (2006).
[35] N. Chazal and D. Gerlier, Microbiol. Mol. Biol. Rev. 67, 226(2003).

032406-9