Role of Transmembrane Proteins for Phase Separation and Domain Registration in Asymmetric Lipid Bilayers

It is well known that the formation and spatial correlation of lipid domains in the two apposed leaflets of a bilayer are influenced by weak lipid–lipid interactions across the bilayer’s midplane. Transmembrane proteins span through both leaflets and thus offer an alternative domain coupling mechanism. Using a mean-field approximation of a simple bilayer-type lattice model, with two two-dimensional lattices stacked one on top of the other, we explore the role of this “structural” inter-leaflet coupling for the ability of a lipid membrane to phase separate and form spatially correlated domains. We present calculated phase diagrams for various effective lipid–lipid and lipid–protein interaction strengths in membranes that contain a binary lipid mixture in each leaflet plus a small amount of added transmembrane proteins. The influence of the transmembrane nature of the proteins is assessed by a comparison with “peripheral” proteins, which result from the separation of one single integral protein into two independent units that are no longer structurally connected across the bilayer. We demonstrate that the ability of membrane-spanning proteins to facilitate domain formation requires sufficiently strong lipid–protein interactions. Weak lipid–protein interactions generally tend to inhibit phase separation in a similar manner for transmembrane as for peripheral proteins.


Introduction
Lipids in membranes tend to mix nonideally [1]. Many lipid mixtures are known for their ability to phase separate or form domains [2]. Of special interest is domain formation in biomembranes because of its putative functional role associated with the membrane raft hypothesis [3][4][5]. The plasma membrane of mammalian cells is asymmetric and multicomponent, but its lipid composition has often been described-in a first-order approximation-as consisting of phosphatidylcholine (PC) and sphingomyelin (SM) in the outer leaflet, phosphatidylserine (PS) and phosphatidylethanolamine (PE) in the inner leaflet, and cholesterol as being able to populate both leaflets [6]. It is well known from experiments in model membranes that the lipids in the outer leaflet appear to represent a mixture of saturated and unsaturated lipids with cholesterol that forms liquid-ordered (lo) domains [7,8]. No such tendency is observed for the lipids in the inner leaflet [9]. However, there are some hints-concluded mostly from computer simulations [10,11]-that suggest lo domains could also exist in the inner leaflet of the plasma membrane and that they are spatially registered with those in the outer leaflet [6].
The raft hypothesis remains controversial [12,13], but it has sparked a large number of experimental [14,15], computational [16], and theoretical [17] studies about domain formation in model membranes, with an increasing focus on inter-leaflet domain coupling in asymmetric bilayers [7,[18][19][20]. Despite being in their fluid state, sufficiently large domains located in the apposed leaflets of a lipid bilayer tend to register due to a domain mismatch energy on the order of 0.1 − 0.2 k B T/nm 2 [21] (k B T is the thermal energy unit: Boltzmann constant times absolute temperature). There is experimental evidence that the mismatch energy is large enough to not only register preexisting domains, but to even induce domains in one leaflet by an existing domain in the apposed leaflet [22]. The origin of the mismatch energy has been suggested to be mostly entropic [23,24], stemming from a more efficient dynamic penetration of the bilayer's midplane by the lipid tails in the registered as compared to the unregistered domain arrangement. Sufficiently small domains may antiregister to minimize the line tension by hydrophobic domain matching [25,26]. Recent theoretical modeling on the mean-field level of a lattice gas has addressed the calculations of phase diagrams in asymmetric membranes [25,[27][28][29][30][31][32].
Here, the domain mismatch energy penalty drives domain registration, but domain formation itself is driven foremost by interactions of the lipids in the same leaflet. This can lead to a rich phase behavior depending on the lipid-lipid interaction strength within each leaflet and the strength of the inter-leaflet domain coupling.
As mentioned, lipid domains can be coupled across the membrane "thermodynamically" through a domain mismatch energy. Here, registered domains are energetically (but not structurally) connected across the bilayer. There is another possibility that has been suggested [33][34][35] but not further pursued: transmembrane proteins or peptides, or membrane-spanning lipids (such as bolalipids [36]), provide a "structural" domain coupling mechanism that may act in conjunction with the above-mentioned thermodynamic mechanism of energy penalties for mismatching domains. Obviously, membrane-spanning proteins are able to physically connect the domains they are associated with across the membrane, irrespective of the inter-leaflet domain interaction energy. In addition, one single transmembrane protein has a lower in-plane translational entropy in a membrane as compared to two equivalent "peripheral" proteins that result from the separation of the transmembrane protein into two independent units. The lower entropy too is expected to favor domain registration. On the other hand, transmembrane proteins of different hydrophobic lengths invoke hydrophobic mismatch penalties in membranes [37] that will affect their ability to induce phase separation. This was observed, for example, by Ackerman and Feigenson [38] in a coarse-grained molecular dynamic simulation of a four-component lipid membrane in the presence of additional transmembrane WALP peptides of varying lengths. Independently of their length, however, all WALP peptides were observed to increase domain alignment. The structural coupling mechanism is not confined to transmembrane proteins; it also applies to bolalipids [39] and even to lipids with long tails such as monosialotetrahexosylganglioside (GM1) [40] and other long saturated acyl chains [11] that interact with the lipids in the apposed leaflet.
The objective of the present work is to propose and analyze a minimal model for phase separation in a mixed lipid bilayer that is subject to the two distinct inter-leaflet coupling mechanisms: a thermodynamic one due to the presence of a compositional mismatch between the two leaflets and a structural one due to the presence of transmembrane proteins. The term "transmembrane protein" stands as a representative for any type of membrane-spanning molecule that is able to interact with the lipids in both leaflets, including integral proteins, transmembrane peptides, bolalipids, and even long-chain lipids. We propose a bilayer-type lattice model, with two two-dimensional lattices stacked one on top of the other. Lipids of two different types (referred to as A and B) occupy one lattice site each, whereas transmembrane proteins consist of two lattice sites that span the bilayer. Hence, each leaflet contains a ternary mixture consisting of two different lipid types and the protein. We introduce all relevant lipid-lipid and lipid-protein interactions and analyze the model on the mean-field level by calculating spinodal surfaces, critical points, tri-critical points, as well as coexistence regions and tie lines in some cases. We demonstrate the ability of transmembrane proteins to facilitate phase transition and to register domains across the bilayer. Our work represents a first attempt to approach an understanding of the three-dimensional phase diagram of a mixed protein-containing bilayer-the richness of the features in the phase diagram justifies the simplicity of our lattice approach, including the neglect of effects due to hydrophobic mismatch, membrane bending, and multi-body interactions.

Theory
We consider two-dimensional lattice models for the external ("ext") and internal ("int") leaflets of a lipid bilayer that contains a fixed number of transmembrane proteins. The two lattices have the same coordination number z (for example, z = 4 for a cubic and z = 6 for a hexagonal lattice) and reside on top of each other so that each lattice site on the external lattice contacts exactly one lattice site on the internal lattice. Each lattice has a total of M lattice sites; the external one hosts P transmembrane proteins, A ext lipids of type A, and B ext = M − P − A ext lipids of type B. Similarly, the internal lattice hosts P transmembrane proteins, A int lipids of type A, and B int = M − P − A int lipids of type B. Because transmembrane proteins span the entire bilayer, the protein positions in each leaflet are exactly the same for each microstate. The illustration of one specific microstate in Figure 1 shows the correlations of protein numbers and positions across the two lattices.

Free Energy of a Lipid Membrane That Contains Transmembrane Proteins
We consider a mean-field Helmholtz free energy F = U − TS of our lattice model at fixed temperature T. Its internal energy U = U ext + U int + U coupl reflects nearest-neighbor interactions within each lattice (U ext and U int ) and a coupling term across the two lattices (U coupl ). The entropy S = k B ln Ω (with k B denoting Boltzmann's constant) accounts for the number of available states, On the level of the random mixing approximation [41,42], the in-plane nearest-neighbor interaction energies U ext in the external layer and U int in the internal layer can be expressed as where are effective lipid A-lipid B, lipid A-protein, and lipid B-protein interaction strengths in the external and internal leaflets. These reflect the actual interaction strengths between lipid A-lipid A (ω ext AA and ω int AA ), lipid A-lipid B (ω ext AB and ω int AB ), lipid B-lipid B (ω ext BB and ω int BB ), lipid A-protein (ω ext AP and ω int AP ), lipid B-protein (ω ext BP and ω int BP ), and protein-protein (ω ext PP and ω int PP ), all expressed in units of The physical situation we are interested in is the presence of lipid A-lipid B interactions in each leaflet and preferential interactions of the proteins with only one lipid type. To this end, we may simply assume ω ext BP = ω ext AA = ω ext BB = ω ext PP = 0 for the external leaflet and ω int BP = ω int AA = ω int BB = ω int PP = 0 for the internal leaflet. This would leave us with χ ext L = zω ext AB /2, χ ext P = zω ext AP /2, χ int L = zω int AB /2, and χ int P = zω int AP /2, whereasχ ext P andχ int P both vanish. More generally, we assume everywhere in the present workχ ext P =χ int P = 0, whereas χ ext L , χ ext P , χ int L , χ int P may all be non-vanishing. If the interactions are symmetric across the bilayer, we are left with only two interaction strengths that we refer to as We also consider cases of asymmetric interactions where χ ext L = χ int L or χ ext P = χ int P . The case χ ext L = χ int L accounts for the different propensities of the lipids in the two leaflets of a plasma membrane to undergo phase separation [9].
Symmetry demands the lowest order term of the inter-leaflet coupling energy across the membrane to be quadratic; on the basis of the random mixing approximation, we obtain U coupl /k B T = Λ(M − P)[(A ext − A int )/M] 2 , where the coupling constant reflects the inter-leaflet interaction of lipid A with lipid B (ω AB ), of lipid A with lipid A (ω AA ), and of lipid B with lipid B (ω BB ). It is convenient to define the three mole fractions Using Stirling's approximation ln(x!) ≈ x ln x − x in the expressions for S as well as the definitions in Equation (5), we find for the total Helmholtz free energy f = F/(Mk B T), in units of k B T and per lattice site, where we again note that we assumeχ ext P =χ int P = 0. Without that assumption, the free energy in Equation (6) would contain the additional contributionχ ext . A specific goal of the present work is to quantify the phase behavior that follows from the function f = f (φ, ψ, α) specified in Equation (6), subject to fixing the interaction strengths χ ext L , χ int L , χ ext P , χ int P , Λ in the thermodynamic limits of fixed temperature T and an infinitely large membrane size, M → ∞. The compositional variables φ, ψ, and α can vary independently within the ranges 0 ≤ α ≤ 1, 0 ≤ φ ≤ 1 − α, and 0 ≤ ψ ≤ 1 − α. However, the interesting protein mole fraction-on which we focus in the present work-is that of small α. We generally consider 0 ≤ α < 0.05. This seems reasonable because only the most protein-rich biological membranes such as the purple membrane of Halobacterium halobium have protein-to-lipid molar ratios of less than 1:50 [43]. Larger mole fractions may occur locally [44], but their consideration would not only add another layer of complexity to the present work; they would also raise concerns about the structural stability of membranes that the present simple lattice model is not designed to address. The mole fractions φ, ψ, and α constitute three independent degrees of freedom. This renders the phase diagram three-dimensional, with a maximum of four phases that can coexist. We note two symmetries. The first one, f (φ, ψ, 0) = f (1 − φ, 1 − ψ, 0), applies to a lipid bilayer that does not contain proteins. The second one, f (φ, ψ, α) = f (ψ, φ, α), is valid for a membrane with symmetric interactions, χ ext L = χ int L and χ ext P = χ int P .

Spinodals and Critical Points
To investigate the phase behavior, we first consider the spinodal surface, which can be calculated from the vanishing of the determinant, detA = 0, of the stability matrix Points {φ, ψ, α} inside the spinodal surface, where A is negative definite, are locally unstable. Tie lines with end points {φ 1 , ψ 1 , α 1 } and {φ 2 , ψ 2 , α 2 } are determined by the familiar common tangent plane construction [42] where we have introduced the abbreviations The existence of three distinct points that satisfy the common tangent plane construction defines three-phase coexistence. We also note that the limit of tie lines with vanishingly small separation between the two coexisting compositions {φ 1 , ψ 1 , α 1 } and {φ 2 , ψ 2 , α 2 } defines a critical point. Critical points are located on the spinodal surface (as determined by Equation (7)); in addition, they fulfill the operator equation where the B φ , B ψ , B α are the cofactors of A along one arbitrarily chosen row or column. For example, taking the middle row implies We are not aware of previous approaches to express the critical point condition in the operator form of Equation (9). In Appendix A, we briefly discuss the derivation of Equation (9) and state equivalent criteria that appear elsewhere in the literature [45][46][47].

Example for Calculation of Spinodals and Critical Points
In Figure 2, we show two examples for the calculation of spinodal surfaces according to detA = 0 (with A specified in Equation (7)) and its critical points (Equation (9), with the cofactors specified in Equation (10)). To this end, we choose f (φ, ψ, α) according to Equation (6) with fixed Λ = 0.05. Note that the magnitude of Λ can be obtained by dividing the domain mismatch energy by the cross-sectional area per lipid (typically 0.7nm 2 ). However, the domain mismatch energy is not well known. In the Introduction, we referred to the range 0.1 − 0.2 k B T/nm 2 estimated by Risselada et al [21] through Molecular Dynamics simulations, but Putzel et al. [29] argued the domain mismatch energy could be an order of magnitude lower. Clearly, our value Λ = 0.05 should be regarded as a rough estimate of a quantity that remains poorly understood.
The left diagram of Figure 2 refers to the absence of proteins, α = 0. It displays four spinodal lines (curves in blue color labeled "a"-"d") and the solution of Equation (9) (the curve in red color), all calculated at fixed χ ext L = 2.1 and the four different choices χ int L = 2.1 (spinodal labeled "a"), χ int L = 2.045 ("b"), χ int L = 2.0272 ("c"), and χ int L = 1.95 ("d"). Note that there is only a single red curve, independent of χ int L because the third derivatives present in Equation (9) remove the quadratic dependencies of the lipid-lipid interaction strengths. Intersections of the blue and red curves mark the critical point locations. The spinodal marked "a" exhibits two critical points, "b" four critical points, "c" two tri-critical points, and "d" zero critical points.  Figure 2 will serve as useful reference in the Results section below.

Numerical Calculation of Coexisting Phases
In order to compute phase diagrams, we need to determine the location of coexisting phases. To this end, it is convenient to minimize the composite thermodynamic free energy of a potentially phase-separated membrane. In Equation (11), we allow for up to four coexisting phases, Hence, the minimization of f th = f th (φ, ψ, α) according to Equation (11) with respect to 12 independent variables and fixed interaction parameters χ ext L , χ int L , χ ext P , χ int P , Λ completely specifies the phase behavior and thus can be used to compute all coexisting phases that correspond to a given point {φ, ψ, α} in the phase diagram. We point out that the accurate calculation of complete three-dimensional phase diagrams as function of all parameters and their meaningful visualization is a formidable task beyond the scope of the present work. Instead, we focus on a few examples that illustrate the role of transmembrane proteins for domain registration across the bilayer.

Free Energy of a Lipid Membrane That Contains Peripheral Proteins
The expected mechanism of how transmembrane proteins couple domains across a lipid bilayer is a structural one, based on the ability of the proteins to protrude into (and interact with) both the external and internal membrane leaflets. To assess how effective this mechanism is, we compare transmembrane proteins with "peripheral" proteins, where we obtain 2P peripheral proteins by cutting each of the P transmembrane proteins in the middle. The two peripheral proteins that result from one single transmembrane protein are able to independently relocate in their host leaflet; see the illustration in Figure 3.
Two modifications of the free energy are associated with transitioning from transmembrane to peripheral proteins. The first is an increase in the number of available states due to the presence of twice as many proteins which will give rise to the additional free energy contribution α ln α + (1 − α) ln(1 − α) in the free energy per lattice site. The second is a modification of the inter-leaflet coupling term due to the presence of lipid-protein interactions across the bilayer. To deduce the latter, we first introduce the presence of inter-leaflet interaction strengths between lipid A and protein (ω AP ), between lipid B and protein (ω BP ), and between protein and protein (ω PP ). These interaction strengths are present in addition to the inter-leaflet interaction strengths of lipid A with lipid A (ω AA ), of lipid B with lipid B (ω BB ), and of lipid A with lipid B (ω AB ) as introduced in Equation (4). Recalling that we operate on the level of the random mixing approximation, we recognize that the total inter-leaflet interactions are given by U coupl /(Mk B T) = Λ(φ − ψ) 2 , where the coupling constant Λ, defined in Equation (4), is independent ofω AP ,ω BP , andω PP . The fact that the inter-leaflet interactions of the peripheral proteins are irrelevant is an immediate consequence of employing the random mixing approximation and the conservation of the number of proteins in the external and internal leaflets. As a consequence, the coupling term in Equation (6) for transmembrane proteins, Λ(1 − α)(φ − ψ) 2 , must be replaced by a coupling term Λ(φ − ψ) 2 for peripheral proteins. Hence, we can write for the mean-field free energy per lattice site in the presence of peripheral proteins where f (φ, ψ, α) has been inserted from Equation (6). In summary, separating all transmembrane proteins in the membrane into pairs of peripheral proteins that reside in opposite leaflets corresponds to adding the mixing entropy term α ln α + (1 − α) ln(1 − α) and an inter-leaflet-interaction contribution Λα(φ − ψ) 2 . We expect the latter to be small because we assume α < 1 throughout this work.

A Simple Example for Phase Stability in the Presence of Transmembrane versus Peripheral Proteins
In order to most clearly illustrate the different influence of transmembrane versus peripheral proteins on the phase behavior, we consider the specific situation χ int L = χ ext L = χ int P = χ ext P = χ, where the proteins interact the same with lipids of type A as lipids of type B interact with lipids of type A. The motivation behind this choice is to eliminate the difference between lipid B and peripheral protein, rendering each leaflet effectively a binary system. We obtain an especially simple free energy expression by demanding identical compositions of the two leaflets, φ = ψ. This removes the inter-leaflet coupling term and thus reduces the free energy difference between a membrane with transmembrane proteins and peripheral proteins to the mixing term α ln α + (1 − α) ln(1 − α).
More specifically, for the peripheral protein, this leads to the free energyf (φ, φ, α) , which amounts to two identical contributions from the two symmetric leaflets. According to Equations (7) and (9), the critical point is located at Phase separation can only be observed for χ > 2, independent of α. Indeed, α is merely a dummy variable because there is no difference anymore between lipid B and peripheral protein: they only differ by their name.
We contrast this to the presence of transmembrane proteins, for which the free energy amounts Here, the critical point is implying that phase separation can already be observed for χ < 2 if proteins are present. Let us focus on the case φ = 1/2 because that always passes exactly through the critical point as α is varied. We can choose α from its minimal value α = 0 to its maximal value α = 1 − φ = 1/2. In the former case, the membrane consists of 50% lipid A and 50% lipid B, which produces a critical nonideality parameter χ = 2. In the latter case, 50% lipid A and 50% transmembrane protein are present, implying a critical nonideality parameter χ = 1. Because every transmembrane protein interacts with two leaflets, the lipid-protein interaction strength is effectively doubled. This implies a reduction of the critical nonideality parameter from χ = 2 to χ = 1. In the intermediate case, for 0 < α < 1/2, the critical point is predicted by Equation (16) to decrease linearly with increasing mole fraction of transmembrane proteins. The presence of lipid-protein interactions is crucial for the ability of the transmembrane proteins to facilitate phase separation. If we repeat the same calculation as in the preceding paragraph, yet with χ int L = χ ext L = χ L and χ int P = χ ext P = 0, we obtain the same critical point irrespective of the proteins being transmembrane or peripheral. Here, the proteins do not exhibit any interactions with the lipid; they merely dilute the two-component lipid mixture, and they do so in the same way for transmembrane and peripheral proteins. This elevates the critical value for χ and thus suppresses phase separation.

Phase Behavior in the Absence of Proteins
We start our analysis by recalling the previously analyzed case α = 0, where no proteins are present in the membrane [24,27]. To this end, Figure 4 displays two phase diagrams, both calculated for α = 0 and fixed Λ = 0.05. They show spinodal lines in blue color with the location of critical points (if present), marked as blue bullets. They also show tie lines as straight solid lines in black, with the two coexisting phases indicated by black bullets at the two ends of each tie line. Regions enclosed by three connected tie lines (present in the left diagram) exhibit three-phase coexistence. The left diagram was calculated for χ ext L = χ int L = 2.1 and the right diagram for χ ext L = 2.1 and χ int L = 1.95. We first discuss the left diagram, which can be viewed as a specific example for a class of systems with χ ext L = χ int L = χ L . The existence of three-phase coexistence requires a sufficiently small coupling constant Λ < Λ with [27] For Λ > Λ and Λ < Λ , the maximum number of coexisting phases is two and three, respectively. The former may be referred to as the strong-coupling regime. At Λ = Λ , the phase diagram contains two tri-critical points φ = φ and ψ = ψ = 1 − φ with Second, for φ = ψ, the membrane is symmetric and the inter-leaflet coupling term vanishes. The phase behavior is then unaffected by Λ and is determined solely by the non-ideality parameter χ L . That is, the spinodal points for φ = ψ are defined by the quadratic equation 2φ(1 − φ) = 1/χ L . Third, the presence of the two three-phase regions for a lipid layer with intermediate asymmetry (where |φ − ψ| is neither very small nor large) reflects the competition between inter-leaflet coupling and the tendencies of each leaflet to phase separate. Indeed, two of the three coexisting phases exhibit the same compositional difference, whereas the remaining third phase (which has the largest compositional difference) serves as host for the "non-matching" lipids. We finally note that the critical points (marked by the blue bullets) are inside the three-phase coexistence regions, which render them irrelevant for the thermodynamically observed phase behavior.
Next, we discuss the right diagram of Figure 4. The parameters χ ext L = 2.1 and χ int L = 1.95 imply that the external but not the internal leaflet tends to phase separate on its own. Hence, phase separation is completely suppressed at small and large φ when the outer leaflet resides outside its binodal region. We still observe two-phase coexistence regions, but no three-phase coexistence. Even a different choice of Λ > 0 will not give rise to three-phase coexistence. Instead, when Λ grows (starting from Λ = 0.05 at fixed χ ext L = 2.1 and χ int L = 1.95), the spinodal detaches from the borders of the phase diagram (that happens at Λ = 0.1, with two critical points appearing at {φ, ψ} = {0.5, 0} and {φ, ψ} = {0.5, 1}) and then forms an increasingly more circular shape. In the limit Λ → ∞, the spinodal is a circle of radius 0.079 centered at φ = ψ = 0.5, with the two critical points {0.556, 0.444} and {0.444, 0.556} attached. The tilt of the tie lines in the right diagram of Figure 4 is a consequence of the inter-leaflet coupling. The observed positive slope of the tie lines emerges from Λ > 0 and the ensuing tendency to minimize the local compositional difference across the bilayer. Any non-vanishing tilt implies distinct compositions of the two coexisting phases in both leaflets. Hence, phase separation in the external leaflet induces phase separation in the internal leaflet. The maximal degree of this "enslaved" phase separation in the inner leaflet is adopted for the tie line that passes the point φ = ψ = 1/2.
How the two phase diagrams in Figure 4 transform into each other upon decreasing χ int L from 2.1 to 1.95 is revealed by the spinodals shown in the left diagram of Figure 2: a tri-critical point exists at χ int L = 2.027, which separates the presence of three-phase regions (for χ int L > 2.027) from its absence (for χ int L < 2.027). It is thus interesting to note that for, say, χ ext L > 2.1, χ int L > 2.02, and Λ = 0.05, there exists no three-phase coexistence despite the presence of a non-vanishing inter-leaflet coupling parameter and despite the tendency of both leaflets to phase separate.

Phase Behavior in the Presence of Proteins
Here, we investigate the influence of transmembrane proteins on the phase behavior and compare it with peripheral proteins. For transmembrane proteins, we use f (φ, ψ, α) according to Equation (6) and for peripheral proteinsf (φ, ψ, α) according to Equation (14).
Consider first the symmetric case with χ ext L = χ int L = χ L = 2.1, χ ext P = χ int P = χ P , and α = 0.04. We initially discuss the left two diagrams of Figure 5. They only differ in the protein type: transmembrane for the top left diagram and peripheral for the bottom left diagram. Each diagram shows eight spinodals (displayed in blue color) with the locations of the associated critical points (blue bullets), calculated for χ P = 0 (the innermost spinodal) and changing in increments of 0.3 until reaching χ P = 2.1 (the outermost spinodal). The light-blue bullets mark additional critical points without the corresponding spinodal lines being displayed. The innermost spinodal in the upper left diagram, calculated for χ P = 0 and α = 0.04, has already been displayed in the right diagram of Figure 2. We recall that two sets of three critical points, residing in close vicinity to each other, are located on that innermost spinodal. We have calculated a number of tie lines for the innermost spinodal and added them to the phase diagram (black lines): clearly, the phase diagram exhibits two three-phase regions, but their small size prevents them from being visible, given our choice of tie line positions. Note that, for χ P = 0, all coexisting phases have the same protein mole fraction, α = 0.04, thus preserving the two-dimensional nature of the phase diagram.
Starting from the innermost spinodal in the upper left diagram of Figure 5 (calculated for χ P = 0 and α = 0.04), upon slightly increasing χ P , the three critical points merge into a single one. Indeed, the second innermost spinodal has only two single critical points left, and so does the third one (which is calculated for χ P = 0.6). Immediately after that, for χ P slightly larger than 0.6, two tri-critical points appear and then give rise to a total of four additional critical points. To visualize this, we have added multiple light-blue bullets that mark critical point locations between the two spinodals for χ P = 0.6 and χ P = 0.9. Hence, the next spinodal, calculated for χ P = 0.9, contains two sets of three critical points. Two of the six critical points have moved outside the phase diagram boundaries for χ P = 1.2, and four of the six critical points have moved outside the phase diagram boundaries for all subsequent spinodals (the three outermost spinodals, calculated for χ P = 1.5, χ P = 1.8, and χ P = 2.1).  The lower left diagram of Figure 5 exhibits a similar scenario as the upper left diagram with two differences. First, the innermost spinodal (which has χ P = 0) is very similar in size and shape. All differences result entirely from the different coupling terms (that is, (1 − α)Λ(φ − ψ) 2 for transmembrane proteins versus Λ(φ − ψ) 2 for peripheral proteins) but not from the different mixing entropies of the proteins. However, the innermost spinodal in the lower left diagram exhibits only two (instead of six) critical points. That is, the phase diagram for χ P = 0 and α = 0.04 exhibits three phase regions for transmembrane proteins, but not for peripheral proteins. As in the upper diagram, we have added a number of tie lines to the innermost spinodal; here, no three-phase region is present in the phase diagram. Second, the spinodals for the peripheral proteins enclose smaller regions than the corresponding spinodals for transmembrane proteins. The latter is a major finding of the present work. It quantifies the ability of the transmembrane proteins-even when present at small mole fractions-to induce membrane phase separation, whereas, for equivalent peripheral proteins, the membrane remains uniform. In addition, transmembrane proteins tend to induce or widen three-phase coexistence regions. For example, all tie lines displayed in the left two diagrams of Figure 5 encode for two-phase coexistence where the two phases have the same degree of inter-leaflet mismatch (|φ − ψ| is the same in both phases). A third phase, one with larger inter-leaflet mismatch, does not form because of the prohibitively large inter-leaflet domain-coupling energy. However, transmembrane proteins are able to couple mismatching domains structurally, thus counteracting the unfavorable inter-leaflet domain-coupling mechanism.
The outermost spinodals in the upper and lower left diagrams of Figure 5 refer to χ P = 2.1, where the interactions of the transmembrane proteins with the A-lipids are the same as the interactions of the B-lipids with the A-lipids (that is, χ L = χ P = 2.1). Peripheral proteins then behave the same as B-lipids, reducing each leaflet effectively to a binary mixture of A-lipids with the union of indistinguishable B-lipids and proteins. Hence, all differences in the phase diagram arise from the presence (for transmembrane proteins) or absence (for peripheral proteins) of the inter-leaflet connectiveness between the proteins segments, without being further affected by the different lipid-lipid and lipid-protein interactions within each leaflet. We have discussed this for the special case of a symmetric membrane (φ = ψ) in Section 2.6. The outermost spinodals in the upper and lower left diagrams of Figure 5 clearly demonstrate the widening of the region where the membrane is unstable due to the inter-leaflet connectiveness of the transmembrane proteins.
We have selected one particular spinodal from the upper and lower left diagrams of Figure 5-the fourth one counted from outside, corresponding to χ P = 1.2-and calculated a number of tie lines; the results are shown in the two right diagrams of Figure 5. Note that the phase diagrams are three-dimensional, with tie lines and three-phase regions extending out of the plane of the displayed {φ, ψ, α = 0.04} section of the phase diagram. Hence, unlike in Figure 4, each unstable point {φ, ψ} has its own individual tie line (or three-phase region). This makes a meaningful visualization of the phase diagram more challenging. On the two right diagrams of Figure 5, we have simply selected a few points (marked by triangles) located on the spinodal (the solid blue lines with the critical points marked as blue bullets) and calculated the corresponding phase behavior. In all cases, we found two-phase coexistence as characterized by tie lines. The α-values of the end points of the tie lines are color-coded, with green hues for α < 0.04 and red hues for α > 0.04; see the legend on the bottom right diagram of Figure 5. From the effective lipid A-protein interactions terms in the free energy, χ ext P αφ + χ int P αψ, it follows that proteins preferentially reside in a phase rich in lipid B (that is, small φ and small ψ), given that χ ext P > 0 and χ int P > 0. This has two consequences in the phase diagrams. First, coexisting phases with smaller φ and ψ values tend to have larger protein content (i.e., red triangles on the left bottom end and green triangles on the right top end of each tie line). In addition, second, phase separation tends to become more pronounced in regions of larger φ and ψ. Regarding the latter, compare, for example, the cases φ = 0 and φ = 1 − α = 0.96 for transmembrane proteins (the upper right diagram in Figure 5). Clearly, there is no phase separation for φ = 0, but there is phase separation for φ = 0.96. Comparing the free energy in Equation (6) shows the presence of the additional term χ int P α(1 − α) for φ = 1. This term is thus responsible for the additional destabilization of the membrane. The physical origin of the difference is that, for φ = 0, the upper leaflet contains only lipid type B, whereas, for φ = 1, the upper leaflet contains only lipid type A. There is no preferential interaction of the proteins with lipid B, but a tendency for segregation when the protein resides in a matrix of lipid A. Hence, we observe a stronger tendency for phase separation when φ = 1 as compared to φ = 0.
Similar considerations are also valid for membranes with asymmetric interactions, χ ext L = χ int L or χ ext P = χ int P . We exemplify this in Figure 6 for a membrane with χ ext L = χ ext P = 2.1 and χ int L = χ int P = 1.95. As in the two right diagrams of Figure 5, we show a {φ, ψ, α = 0.04} section of the phase diagram with the spinodal line and a selected number of tie lines displayed. The left and right diagrams in Figure 6 correspond to transmembrane and peripheral proteins, respectively. In fact, the only difference when comparing the upper right diagram of Figure 5 with the left diagram of Figure 6 (for transmembrane proteins) and the lower right diagram of Figure 5 with the right diagram of Figure 6 (for peripheral proteins) is the parameter change from χ int L = χ int P = 2.1 to χ int L = χ int P = 1.95. Hence, because the latter is below the critical point, we no longer observe phase transitions in the small and large φ-regions of the phase diagram. By comparing the phase diagram in the right diagram of Figure 4 with those in Figure 6, we directly observe the influence of adding transmembrane or peripheral proteins of mole fraction α = 0.04 and interaction strengths χ ext P = 2.1 and χ int P = 1.95. We also note that the spinodal lines in Figure 6 reproduce the small and large ψ-regions of those in Figure 5 for χ P = 2.1 (the outermost spinodals in the two left diagrams of Figure 5). As a result of decreasing χ int P from 2.1 to 1.95, no critical points are present anymore and phase separation always produces exactly two coexisting phases. Most importantly and similar to our discussion of Figure 5, transmembrane proteins produce substantially larger unstable regions in the phase diagram as compared to peripheral proteins. We emphasize that the phase diagrams shown in Figures 5 and 6 only contain limited information-spinodals and tie lines at a few selected positions of a {φ, ψ, α = 0.04}-section. Complete phase diagrams contain information about all coexisting phases at every point {φ, ψ, α} in the three-dimensional phase diagram. Visualizing full phase diagrams in a meaningful way and analyzing them comprehensively for the set, χ int L , χ ext L , χ int P , χ ext P , Λ, of interaction parameters is a work in progress.

Conclusions
Our work is a first step into the direction of analyzing the interplay between "thermodynamic" and "structural" coupling of domains across a lipid membrane. "Thermodynamic" coupling results from the inter-leaflet interactions of the lipids, whereas the coupling becomes "structural" for transmembrane proteins (or other membrane-spanning components such as bolalipids) that stretch across the entire lipid bilayer. We have considered a mean-field, lattice-based model with nearest neighbor pair interactions. This type of model is highly simplistic and in many ways oversimplifies a protein-containing lipid membrane. However, it captures the principal aspect of introducing membrane-spanning molecules (which we refer to as transmembrane proteins) into a lipid bilayer and its influence on the membrane phase behavior. We find that this influence depends on the interaction strength of the transmembrane proteins with the lipids in the two leaflets. Weakly interacting proteins suppress phase separation by merely diluting the lipids. If the lipid-protein interaction strength resembles that among the lipids (which leads to domain formation in the first place), then transmembrane proteins are indeed able to couple domains and enhance or even induce their formation. Our present work may trigger a number of extensions. Perhaps one of the most interesting is related to the ability of immobilized pinning sites to limit the growth of lipid domains [48], which is one among a number of mechanisms [49,50] that have been suggested to rationalize the nanoscopic size of membrane rafts. Cytoskeletal coupling of the inner leaflet of a plasma membrane creates pinning sites, but it is the outer leaflet that has a high propensity to phase separate.
We emphasize the simplicity of our model. For example, all energy penalties due to hydrophobic mismatch (either between the lipids among each other or between the lipids and proteins) are lumped into effective interactions' parameters. Hence, our model is not capable of making predictions of how domain coupling depends on lipid chain length or on the hydrophobic protein thickness. Moreover, we represent lipids by the sites of two coupled two-dimensional lattices, thus ignoring any molecular specificities such as head group size, degree of lipid chain unsaturation, etc. The special role of cholesterol does not enter explicitly into our model, neither do membrane bending [51], line tension effects of domains [52,53], protein-mediated coupling of domains to its surroundings [54][55][56], or any multi-body interactions. Finally, we treat our lattice model on the mean-field level, which ignores any type of correlations between interacting membrane constituents [57]. However, even with all these approximations, our simple model leads to a surprisingly rich phase behavior. . We also acknowledge financial support through the Phospholipid Research Center, Heidelberg, Germany.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Appendix A. Thermodynamic Criterion for Critical Point
Equation (9) (together with Equation (10)) defines the location of critical points on the spinodal surface. It can be derived starting from the tangent plane construction in Equation (8) by expressing the two coexisting compositions as φ 1 = φ + φ, ψ 1 = ψ + ψ, α 1 = α + α, φ 2 = φ − φ, ψ 2 = ψ − ψ, α 2 = α − α. Upon carrying out a series expansion up to first order in x = { φ, ψ, α}, we can express the three equal-slope equations in the first line of Equation (8) as Ax = 0, thus reproducing the equation detA = 0 that defines the spinodal. With that, the first non-vanishing term in the expansion with respect to φ, ψ, α of the second line of Equation (8) is of third order. The two eigenvectors of x allow us to express two of its components by the remaining one; for example, φ = φ( ψ) and α = α( ψ). Inserting these into the third-order expansion of the second line of Equation (8) yields where our choice of the independent variable in expressing the two eigenvectors determines which row (or, equivalently, column) the cofactors B φ , B ψ , B α refer to. Equation (A1) is equivalent to Equation (9). We have not observed Equation (9) to be stated previously in the literature. However, multiple equivalent ways to express Equation (9) are well known; see the discussions by Akasaka [46] and by Bell and Jäger [47]. Among those is a set of equations that goes back to Heidemann and Khalil [45] and results from a third-order stability analysis of a Taylor series expansion of the Helmholtz free energy about the equilibrium state. Another method involves the vanishing of the determinant detB = 0 where the matrix B is produced from matrix A (see Equation (7)) by replacing any one of its rows by the row vector {∂(detA)/∂φ, ∂(detA)/∂ψ, ∂(detA)/∂α} [58] .