On The Critical Casimir Interaction Between Anisotropic Inclusions On A Membrane

TL;DR

This study investigates the critical Casimir interactions between anisotropic inclusions in membrane-like binary mixtures, revealing long-range, affinity-dependent forces that influence clustering and patterning of membrane proteins.

Contribution

The paper introduces a lattice model and thermodynamic integration scheme to analyze how anisotropic inclusions interact at criticality, highlighting the effects of domain composition and orientation.

Findings

Interactions are long-range and depend on inclusion affinity.

Opposing affinities lead to repulsive interactions, same affinities to attraction.

Multiple domains cause anisotropic interactions that influence clustering.

Abstract

Using a lattice model and a versatile thermodynamic integration scheme, we study the critical Casimir interactions between inclusions embedded in a two-dimensional critical binary mixtures. For single-domain inclusions we demonstrate that the interactions are very long range, and their magnitudes strongly depend on the affinity of the inclusions with the species in the binary mixtures, ranging from repulsive when two inclusions have opposing affinities to attractive when they have the same affinities. When one of the inclusions has no preference for either of the species, we find negligible critical Casimir interactions. For multiple-domain inclusions, mimicking the observations that membrane proteins often have several domains with varying affinities to the surrounding lipid species, the presence of domains with opposing affinities does not cancel the interactions altogether. Instead we can observe both attractive and repulsive interactions depending on their relative orientations. With increasing number of domains per inclusion, the range and magnitude of the effective interactions decrease in a similar fashion to those of electrostatic multipoles. Finally, clusters formed by multiple-domain inclusions can result in an effective affinity patterning due to the anisotropic character of the Casimir interactions between the building blocks.

Equations

$$H_{\alpha }=J\sum _{i=1}^N\sum _{⟨i⟩}|s_i-s_j|.$$

$$H_{\beta }=J\sum _{i=1}^N\sum _{⟨i⟩}|s_i-s_j|+h\sum _{i=1}^{n_{p1}}(s_i-s_{t1}{)}^2.$$

$$H_{\gamma }=J\sum _{i=1}^N\sum _{⟨i⟩}|s_i-s_j|+h\sum _{i=1}^{n_{p2}}(s_i-s_{t2}{)}^2.$$

$$H_{\delta }(d/r)=J\sum _{i=1}^N\sum _{⟨i⟩}|s_i-s_j|+h\sum _{i=1}^{n_{p1}}(s_i-s_{t1}{)}^2+h\sum _{i=1}^{n_{p2}}(s_i-s_{t2}{)}^2.$$

$$\Delta F_{fin/ini}\equiv F_{\mathrm{fin}}-F_{\mathrm{ini}}={\int }_0^1{⟨\frac{\partial H(\lambda )}{\partial \lambda }⟩}_{\lambda }d\lambda .$$

$$H(\lambda )=H_{\mathrm{ini}}+\lambda (H_{\mathrm{fin}}-H_{\mathrm{ini}}).$$

$$\Delta F_{fin/ini}={\int }_0^1{⟨H_{\mathrm{fin}}-H_{\mathrm{ini}}⟩}_{\lambda }d\lambda .$$

$$F_{CC}(d/r)=F_{\delta }(d/r)-(F_{\beta }+F_{\gamma })+F_{\alpha },$$

$$F_{CC}(d/r)=\Delta F_{\gamma /\delta }-\Delta F_{\alpha /\beta }={\int }_0^1{⟨H_{\delta }-H_{\gamma }⟩}_{\lambda }d\lambda -{\int }_0^1{⟨H_{\beta }-H_{\alpha }⟩}_{\lambda }d\lambda$$

Full text

I Introduction

Originally, the cell membrane was considered just as a physical barrier that kept cell components together. However, in the last decades, advances in experimental techniques, including atomic force and fluorescence microscopies, have enabled the probing of membranes’ inner structure and composition to such extents that it has been necessary to rethink our understanding of its physics, chemistry and its role in biology Shan and Wang (2015). Far from being a simple barrier, cell membranes are in fact very complex mixtures in which lipids and proteins meet, interact and self-assemble. A number of studies have shown that lipids are organized in the lateral dimension Kusumi and Sako (1996) and that most membrane proteins organize in clusters Uhles et al. (2003); Kai et al. (2006); Low et al. (2006); Sieber et al. (2006a). Further, there is now a growing consensus that such lateral organization and compartmentalisation play important roles in biological processes, such as in cell signalling and membrane trafficking Alonso and Millán (2001); Ikonen (2008); Martinez-Outschoorn et al. (2015); and it has also been suggested to be involved in a number of diseases from HIV Simons and Ehehalt (2002) to liver Holthuis and Menon (2014) and prion diseases like Alzheimer Critchley et al. (2004).

Biologically-motivated scenarios, such as the fence and pickets model Kusumi et al. (2005) or the lipid raft hypothesis Simons and Ikonen (1997), have been put forward to explain lateral membrane organizations. However, to date it remains unclear whether such scenarios are feasible from a chemical physics standpoint and, if possible, to what extent they are harnessed by biological processes in vivo. In this context, our aim here is to contribute towards understanding how inclusions (e.g. membrane proteins) may assemble into clusters within a cell membrane. Assuming the aggregation process occurs in thermodynamic equilibrium conditions, two possibilities come to mind. Firstly, it may occur by the free diffusion, collision and irreversible chemical bonding of the inclusions in a process reminiscent of a Diffusion Limited Aggregation Saxton (1993a, b); Tian et al. (2016). Secondly, the inclusions may instead self-organise into metastable clusters, with a given equilibrium size distribution, in order to minimise the overall free energy. In this work we will concentrate on the latter.

There are various candidate mechanisms for interactions between somewhat large protein-like inclusions surrounded by smaller lipid species. They range from direct, specific interactions between the inclusions, such as via van der Waals or electrostatic forces Honig et al. (1986); Lee (2003), to indirect interactions mediated by the lipid membranes. For the latter, one possible mechanism arises from minimising membrane deformations. Both theory and experiments have demonstrated that, by inducing curvature, proteins and colloidal particles can self-assemble into complex structures such as lines, rings and lattices Reynwar et al. (2007); Briggs et al. (2009); Saric and Cacciuto (2013); Vahid and Idema (2016). More recently, ideas for a new class of indirect mechanisms have emerged, namely Casimir-like forces Machta et al. (2012); Bimonte et al. (2013); Bitbol and Fournier (2011); Woods et al. (2016); Bitbol et al. (2010); Golestanian et al. (1996); Yolcu et al. (2011); Lin et al. (2011); Weikl (2001), which could arise both from shape and composition fluctuations of the membrane. Here we will focus on composition fluctuations. Initially envisioned for bulk critical mixtures in three dimensionsFisher and de Gennes (1987), Casimir-like forces are expected to arise between any inclusions embedded in a fluctuating fluid when separated by a distance shorter than the correlation length of the fluid; and for fluid mixtures near a critical point, this correlation length diverges as we approach the critical temperature. Interestingly, close to room temperature, lipid membranes with critical composition pass through a miscibility critical point whose nature closely follows that of the two-dimensional Ising universality class Honerkamp-Smith et al. (2008). This behaviour has been observed in giant plasma membrane vesicles isolated directly from living cells Veatch et al. (2008), as well as in synthetic lipid mixtures Connell et al. (2013). It thus becomes plausible that Critical Casimir (CC) interactions between membrane proteins embedded in critical lipid mixtures could pave ways for their aggregation.

Machta et al. Machta et al. (2012) recently demonstrated that two “like” inclusions immersed in a 2-dimensional critical binary lipid mixture have an attractive interaction, while “unlike” ones have a repulsive interaction. In their model each protein inclusion has a set uniform affinity with one of the lipid components. These findings are in line with those previously reported both theoretically and experimentally for CC interactions in 3 dimensions Krech (1997); Hertlein et al. (2008); Gambassi et al. (2009); Trondle et al. (2009); Nellen et al. (2009). Recent extensive simulations of colloids immersed in a critical solvent further suggest strong similarities between the 2 and 3-dimensional cases whereby the observed variations in the phase behaviour appear to be accountable to differences in critical exponent and dimension of the solvent-inclusion interface Edison et al. (2015); Tasios and Dijkstra (2017).

Besides uniform inclusions, inhomogeneities and anisotropies at the fluid-inclusion interface have been theoretically shown to affect substantially CC interactions in 3 dimensions Trondle et al. (2009); Parisen Toldin et al. (2013); Nowakowski and Napiórkowski (2014); Labbe-Laurent and Dietrich (2016) to the extent that patterning can now be tuned and exploited for targeted colloidal self-assembly Nguyen et al. (2016); Labbe-Laurent and Dietrich (2016). This observation is particularly interesting and relevant for CC forces between membrane proteins as proteins often have multiple domains with varying affinities with respect to the lipid speciesCohen et al. (1995); Pawson and Nash (2003); Mayer (2015). Thus, if CC interactions are to play a role in membrane protein aggregations, it is key to understand how their sign and magnitude can be tailored by heterogeneities in the protein domains and different affinities between the inclusions and the lipid species.

To address this issue, in this work we employ Monte Carlo (MC) simulations in a square lattice to simulate a binary mixture of lipids with two inclusions embedded in it. The free energy corresponding to the CC interactions between these inclusions is explicitly computed by using a thermodynamic integration scheme. The flexibility of our scheme allows us to (i) explore the effect of gradually modifying the inclusion boundary conditions, and (ii) study how the presence of different domains (anisotropy) can affect the interactions.

II Methods

The Critical Casimir (CC) interaction is characterised by the free energy $F_{CC}(d)$ of a critical fluid mixture within which lie two inclusions separated by a fixed distance $d$ Fisher and de Gennes (1987). There are two common approaches to compute this CC free energy. Firstly, the free energy can be defined as minus the work done by the ensemble average of the mechanical force exerted by the fluid on the inclusions when brought from infinity to a finite separation $d$ Frenkel and Smit (2002); Machta et al. (2012). The equilibrium fluctuations of this mechanical force also inform experiments relying on mechanical probes to evidence CC interactions; however, some care must be taken when interpreting these fluctuations, depending on how the inclusions interact with the fluid mixtures Bartolo et al. (2002); Bitbol and Fournier (2011).

Secondly, from an equilibrium statistical thermodynamics viewpoint, the statistics of configurations will be governed by a Boltzmann weight that depends on $F_{CC}(d)$. Thus, a fruitful strategy is to compute the probability weights for inclusions at varying distances, which are readily available from simulations, and invert them to recover the CC free energy as a function of distance Edison et al. (2015).

However, both the force and relative probability approaches usually need to be supplemented by a theoretical estimate at a reference distance because 2D CC interactions are very long range and it can be impractical to extend the calculations to get an effectively zero interaction at large distances Machta et al. (2012); Edison et al. (2015). In this work, we devise a non-distance-based thermodynamic integration scheme precisely to avoid the need for a reference point, which can be difficult to obtain for complex inclusion geometries. We will detail our approach in the following subsections.

II.1 The Model

Our system consists of a critical lipid binary mixture in which two proteins are embedded. We use a lattice gas model of this system with a square lattice of $N$=150x150 cells. Each cell is occupied either by a species of the binary mixture or part of a protein. This lattice size is chosen as it gives a good trade-off between reducing computational costs and finite size effects. By this we mean that up to $N=$500x500, larger lattice sizes yield results whose difference with our $N$=150x150 lattice size are not statistically significant.

We simulate the lipid mixture by making use of a two-dimensional Ising model, notorious for its critical behaviour, in which the spin variables can take values $s=\pm 1$ corresponding to lipid species A and B. The protein inclusions are modelled as two blocky patches occupying $n_{p}=12$ cells each, whose size is determined by their effective radius, $r$, as shown in Fig. 1. For the inclusions, the average value of their spin variable is set to a target value $s_{t}$ which can be different for the two proteins. To avoid a large standard deviation from this mean, spins belonging to a protein region are modelled with Potts-like spin variables $s={k/2}$, with $k$ an integer such that $k\in[-2,2]$.

To get the free energy of the system with two protein inclusions we make use of the 4 states represented in Fig. 1. The distance between proteins, $d$, is measured as the distance between their centers. These states are the following:

• •

State $\alpha$ comprising only the binary mixture at criticality (top left of Fig. 1) with hamiltonian:

[TABLE]

• •

State $\beta$ comprising the critical mixture and protein 1 (top right of Fig. 1) with hamiltonian:

[TABLE]

• •

State $\gamma$ comprising the critical mixture and protein 2 (bottom left of Fig. 1) with hamiltonian:

[TABLE]

• •

State $\delta$ comprising the two proteins 1 and 2 separated by a distance $d/r$ within the critical mixture (bottom right of Fig. 1) and with hamiltonian:

[TABLE]

The first term in Eqs. (1)–(4) is equivalent to a 2D Ising model, while the second and third terms are quadratic terms whose role is to impose the average spin values inside the regions where the proteins are. Here $J>0$ is the coupling parameter which characterises the energy cost for having two different species next to each other, $\sum_{\langle i\rangle}$ stands for a sum over the 4 closest neighbours of cell $i$, $h$ is a parameter setting the strength of the external potential imposing the average spin values inside the proteins, and $s_{t1}$ and $s_{t2}$ are the spin target values for proteins 1 and 2, respectively.

II.2 Thermodynamic Integration

Thermodynamic integration is a versatile method to compute free energy differences between two thermodynamic states. Let us consider an initial state with hamiltonian $H_{\rm ini}$, a final state with hamiltonian $H_{\rm fin}$, and a parameter-dependent hamiltonian $H(\lambda)$ such that $H(\lambda=0)=H_{\rm ini}$ and $H(\lambda=1)=H_{\rm fin}$. It can be shown that the free energy difference between the final and initial states can be computed by Frenkel and Smit (2002):

[TABLE]

The equilibrium averages $\left<\frac{\partial H(\lambda)}{\partial\lambda}\right>_{\lambda}$ can readily be obtained from standard MC methods. It is important to note that, as a thermodynamic quantity, the free energy difference is only a function of the final and initial states, but not the path we have taken between them. Thus, for simplicity, we have taken a linear interpolation with crossover Hamiltonian

[TABLE]

This then yields

[TABLE]

Specific to the problem at hand, here we have four thermodynamic states $\alpha$, $\beta$, $\gamma$ and $\delta$ (see Fig. 1), and we denote their free energies as $F_{\alpha}$, $F_{\beta}$, $F_{\gamma}$ and $F_{\delta}(d/r)$. We finally denote $F_{CC}(d/r)$ as the CC free energy and define it as the work done to bring the two proteins to within a distance $d/r$ from infinity. It follows from this definition that:

[TABLE]

where $(F_{\beta}+F_{\gamma})-F_{\alpha}$ represents the free energy of the system with the two proteins being infinitely far apart in the critical mixture. From Eq. (7), $F_{CC}(d/r)$ can then be expressed as a function of thermodynamic integrals only:

[TABLE]

which can be readily estimated from our MC simulations. We note that, for a given set of protein inclusions, the second thermodynamic integral only needs to be computed once, while the first integral has to be repeated for various values of $d/r$. The advantage of Eq. (9), which resemble to some extent the parameter variation method used in Ref. Parisen Toldin et al. (2013), is that it avoids the problem of having to supply a theoretical estimate of the reference free energy, which tends to be very much system-dependent and difficult to compute analytically.

II.3 Simulation Details

We perform standard Metropolis Monte Carlo simulations at the critical temperature of the Ising model with non-conserved order parameter. In our simulations, it corresponds to $J/k_{B}T=2.27$, in agreement with analytical predictions by Onsager and previous simulation results Onsager (1944); Mon and Binder (1992). For each system we employ three random different initial configurations, which are equilibrated for $5\times 10^{5}$ cycles, and sampled for $5\times 10^{6}$ cycles. Each cycle consists of $N$ random single spin flips, which ensures that all particles in the system can be selected in each cycle, and configurations are saved every 1000 cycles. In order to determine the error in our calculations we use block analysis Flyvbjerg and Petersen (1989) to determine the number of independent measurements of the free energy. The error is obtained from the standard deviation with a confidence interval of 95%.

When performing the thermodynamic integration, there are two technical aspects worth commenting. The first aspect is the choice for the parameter $h$ in Eqs. (1)–(4). On the one hand, its value must be large enough to effectively pin the spin values inside an inclusion such that the associated energy variation $\langle H_{\rm{fin}}-H_{\rm{ini}}\rangle_{\lambda}\rightarrow 0$ as $\lambda\rightarrow 1$. On the other hand, the accuracy of our estimate of the CC free energy, $F_{CC}$, is better when we have more points substantially contributing to the numerical integration in Eq. (7). As a consequence, $h$ must be in a window of values guaranteeing both fidelity to the model we want to implement and reliability of the integration method.

We further illustrate this issue in Fig. 2, where we show a typical result of our thermodynamic integration method for three values of $h=50k_{B}T$, $250k_{B}T$ and $500k_{B}T$. In this example, we compute the CC free energy difference between states $\gamma$ and $\delta$. For the largest value of $h$, the contribution to the integral in Eq. (7) primarily comes from a small window of $\lambda$ close to $\lambda\rightarrow 0$ and only involves a couple of integration points, leading to a poor estimate of the integral. From this perspective, smaller values of $h$ are better. However, for small values of $h$, we notice that the integrand in Eq. (7) does not converge to zero as $\lambda\rightarrow 1$ (see the inset of Fig. 2). Physically this means the spin value of protein inclusions is not correctly set to the target value $s_{t}$. In this work, we find that we obtain equivalent results for the Critical Casimir free energy if we use $h=100-250\,k_{B}T$. For the rest of this paper, we use $h=250k_{B}T$.

The second aspect concerns the difference in degrees of freedom when a given lattice point represents a lipid species or part of a protein inclusion. For instance, when computing $\Delta F_{\alpha/\beta}$, the protein sites in state $\beta$ have five possible spin values (Potts model with $s={k/2}$, $k\in[-2,2]$), whereas the equivalent sites in state $\alpha$ only have two possible spin values (Ising model with $s=\pm 1$), because in state $\alpha$ they represent lipid species. Since we are computing the CC free energy at varying distances between the protein inclusions, in principle we must account for corrections due to variations in degrees of freedom explicitly in the computations for $F_{CC}(d/r)$. To do this, we carry out two sets of thermodynamic integration calculations. In the first set, see sketch in Fig. 3(b), the initial state corresponds to the case where there is no interaction neither between the cells where the inclusions are located nor with their surrounding cells. The final state is where all cells are interacting and they all have two possible spin values (Ising model). We denote this free energy difference as $\Delta F_{\rm Ising}(d/r)$, where $d$ once again is the separation between the two inclusions and $r$ is the radius of the inclusions. For the second set of calculations, see sketch in Fig. 3(b), the initial state is the same as before. However, for the final state, the cells where the inclusions are located can now have five possible spin values (Potts model). The surrounding sites still have $s=\pm 1$. The free energy difference for this case is denoted as $\Delta F_{\rm Potts}(d/r)$. The difference $\Delta F_{\rm Ising}(d/r)-\Delta F_{\rm Potts}(d/r)$ is therefore the correction in the CC free energy due to changes in the spin degrees of freedom as a function of the normalised separation $d/r$, and the typical results are shown in Fig. 2(a) for a 150x150 square lattice. In practise, these corrections are small and, in fact, within the uncertainty of our typical thermodynamic integration results.

III Results and Discussion

We study the effective interaction between two protein inclusions embedded in a lipid membrane depending on the binding affinity that the proteins have for the lipid species A and B present in the membrane. For instance, a protein (or domain in a protein) labelled as -1 strongly prefers being surrounded by lipid species A, while a protein (or domain in a protein) labelled as 1 strongly prefers lipid species B. A value of 0 means that it has no preference for either of the species.

III.1 Single-domain proteins

We first illustrate the approach to the critical temperature $T_{c}$ by plotting two-dimensional maps of the ensemble average spin value denoted $\varphi$ for various temperatures in the neighbourhood of $T_{c}$. The presence of uniform inclusions with set spin value $s_{t}=-1$ introduces inhomogeneities into the spin system which persist over a distance of the order of the correlation length. Plotting the $\varphi-$maps enables us to have a direct visualisation of the correlation length. Indeed, the divergence of the correlation length at $T_{c}$ gives rise to long-range CC forces to exist.

As shown in Fig. 4, approaching $T_{c}$ from above, the distance over which the inclusions influence the rest of the spins in the system is initially rather short (panel (a)), but increases more and more (panels (b) and (c)) until the correlation length covers the whole box at $T_{c}$.

We next test our thermodynamic integration scheme by computing the CC free energy for proteins that consist of one single domain. We explore the variations of the protein-protein interactions as a function of their separation for different protein binding affinities to the lipid species. In our approach, we gradually change the binding affinity of one of the proteins to the lipid species, while keeping the other one constant. More specifically, we consider systems in which protein 1 has a spin target value of -1, while protein 2 can get spin target values of 1, 0.5, 0, -0.5 or -1 (see Fig. 5 for a description). The results for these systems are summarised in Fig. 5.

We start by looking at the limit in which both proteins have opposite binding affinities for the two lipid species, system $(-1,1)$. Here, protein 1 has strong preference for lipid A and protein 2 has strong preference for lipid B, which we call the ”unlike” limit. In this limit the effective interaction is clearly repulsive (orange diamonds in Fig. 5) and can be quite strong with a free energy barrier rising up to $3.5k_{B}T$. Then, as the binding affinity of protein 2 is gradually varied, and its preference becomes weaker for lipid species B and stronger for A, the strength of the repulsive interaction decreases (blue triangles in Fig. 5). At the point in which protein 2 has no preference for either of the lipid species, system $(-1,0)$, we find a crossover between the two regimes, and CC interactions become negligible with respect to $k_{B}T$ (green stars in Fig. 5). Finally, as protein 2 keeps on increasing its affinity for lipid A we get into the attractive regime (red squares in Fig. 5). In the ”like” limit in which both proteins have the same strong affinity to one of the lipid species, system $(-1,-1)$, we find that the magnitude of attraction is the strongest with a well depth of about $1k_{B}T$ (black circles in Fig. 5). However, this attraction is much weaker than the magnitude of the repulsion in the ”unlike” limit, by more than a 3 fold difference. That like proteins attract supports the idea that proteins of the same kind tend to coalesce in critical lipid mixtures in order to minimise the total free energy of the system.

We note that the free energy curves in Fig. 5 appear to have an offset with respect to the zero free energy baseline. However, this apparent offset is just a consequence of the very slow decay of these curves to zero as $d/r\rightarrow\infty$ and thus of the very long range nature of CC interactions in two dimensions. We further find that the free energy dependence with distance can be well fitted to a power law $F_{CC}(d/r)/k_{B}T=\zeta(d/r)^{-\nu}$, shown as plain lines in Fig. 5. The values for $\zeta$ and $\nu$ are tabulated in panel (c) in Fig. 5. Smaller values of $\nu$ signify longer-range of interactions.

Our results for the like (black circles) and unlike (orange diamonds) limits in Fig. 5 corroborate the findings of Machta et al. Machta et al. (2012) where we find excellent agreement without any fitting parameter. In their simulation work using Bennett algorithm Bennett (1976); Jarzynski (1997), they take advantage of asymptotic results from conformal field theory for the reference free energy at large distances. Our thermodynamic integration approach here circumvents the need for supplying the reference free energy at large distances, and thus it allows more complex scenarios to be computed readily, as we shall show in the next section for inclusions with multiple domains. Our results are also in good qualitative agreement with the experimental findings of Ref. Nellen et al. (2009), where the CC forces between colloids and a wall were tuned by varying the preferential absorption properties of the wall for the critical mixture’s components.

Focusing on the attractive regime which is the strongest for like inclusions, we next look at the magnitude of the CC interactions between like proteins, as we vary their binding affinity with the surrounding matrix from strong preference for lipid A (system $(-1,-1)$) to no preference for either lipid (system $(0,0)$). Note that we do not investigate explicitly the case where the two proteins prefer lipid B since it is equivalent to the system $(-1,-1)$. Our results are summarised in Fig. 5(b). It is found that the greater the affinity of the two proteins for one of the lipid species, the greater the magnitude of the resulting CC interactions. In fact, for system (0,0) where the inclusions have no preference for either of the lipid species, we cannot measure any net interaction within the limits of our computational accuracy. For system $(-0.5,-0.5)$ we clearly observe an attractive interaction (negative free energy). However, the CC force (the slope of the free energy with distance) is very small. It is akin to the case where the system sits in an effective negative square well potential.

III.2 Multiple-domain proteins

In the previous section we noted the absence of CC interaction at any distance when one or both of the proteins do not have any preference for either of the lipid species. In this section we investigate whether this still remains the case for identical proteins with multiple domains each of which being given a target affinity of either -1 or 1 — i.e. with strong binding affinity for either lipids A or B — but such that the overall net affinity of each protein is zero. To some extent, this approach with multiple-domain proteins represents a more realistic model of proteins which usually comprise several domains with different binding affinities to the surrounding lipids Cohen et al. (1995); Pawson and Nash (2003); Mayer (2015), thus making the present study more relevant to biological systems.

We start by studying interactions between two proteins with two different domains each. These proteins can have different relative orientations depending on which domains of the proteins are facing each other (see Fig. 7(a) for a description). Our aim is to study how orientation can affect these interactions and if there are preferential orientations that could lead to patterning in protein clusters. Given that these proteins do not have a net preference for any lipid type, one could naively expect no interaction between them as in the cases of the $(0,0)$ and the $(-1,0)$ pairs of single-domain proteins. However, our results as shown in Fig. 7(b) show again the presence of two different regimes and a crossover between them. We find that systems 1 and 2 appear in the repulsive regime while systems 4 and 5 fall in the attractive one. System 3, corresponding to a $90^{\rm o}$ relative orientation, is slightly attractive, though its magnitude is very small compared to the other systems.

A closer look at systems 1 and 2 shows that, in both cases, proteins are facing each other with domains that have opposite preferences for the lipid species, giving rise to predominantly ”unlike” interactions. Therefore, this situation is qualitatively analogous to the system $(-1,1)$ for single-domain proteins. On the other hand, systems 4 and 5 are confronting domains of the same kind, leading to predominantly ”like” interactions as in system $(-1,-1)$ for single-domain proteins. Lastly, in system 3, proteins are orientated in such a way that half of protein 1 has ”like” interactions with protein 2, while the other half has ”unlike” interactions. This case is a good example of the non-additive nature of CC interactions. If they were simply additive, taking advantage of the results in Fig. 5, we would have expected the net interaction to be repulsive. The magnitude of repulsive interactions between unlike uniform inclusions is stronger than attractions between like ones. However, as shown in Fig. 7, for system 3, we observe a very weak attractive interaction instead.

When compared with single-domain proteins, the range of the emerging CC interaction is clearly shorter in all cases, getting close to zero at distances of about ten times the radius of the protein; we interpret this as being a partial screening effect analogous to electrostatic screening in dielectric systems and owing to proteins comprising domains with opposing affinities. Panel (c) in Fig. 7 tabulates the power law fit of the CC free energy with distance, $F_{CC}(d/r)/k_{B}T=\zeta(d/r)^{-\nu}$. Here, the exponents $\nu$ are clearly larger for the values for single-domain inclusions in Fig. 5.

In order to further investigate the effect of the domains we next turn to proteins with four distinct domains. In this case only two possible orientations are of interest, see Fig. 8(a). We present our results in Fig. 8(b). Again, in spite of having no net preference for either of the species, proteins show either attractive or repulsive interactions. A similar analysis to that done with two-domain proteins show that system 6 confronts domains with opposite preferences for the lipid species, falling, therefore, in the repulsive regime. On the other hand system 7 has domains with the same preference facing each other, leading to attractive interactions. It is worth pointing out that, in these cases, both the range and the magnitude of the interactions are smaller than those observed for single- or two-domain inclusions. Fig. 8(c) tabulates the power law fit for the CC free energy with distance between two four-domain inclusions. Eventually in the limit of the domain size going to zero, we expect to recover the uniform case of $s_{t}=0$. Comparing the results in Figs. 5, 7 and 8, we also observe that the magnitude of the repulsive interactions falls faster than attraction interactions with increasing number of domains.

The observed behaviours for multiple-domain inclusions are in good qualitative agreement with the findings of Toldin et al. for CC forces between chemically striped surfaces immersed in critical liquid mixtures Parisen Toldin et al. (2013). In their study, they confront a surface with a strong preference for one component of the binary liquid mixtures with a striped surface in which the preference of the stripes for the liquid components is alternating. They observe a crossover between a regime with noticeable CC interaction and one with no noticeable CC interaction (similar to our $(0,0)$ and $(-1,0)$ cases for single-domain proteins) as the strip width is decreased from large values to zero.

Finally, we note that the range and strength of the CC interactions depend only on the lattice points which are at the boundary of the inclusions, and not the inner ones. This is as expected because we only have nearest neighbour interactions, and thus the inner sites do not interact directly with the surrounding cells representing the lipid mixtures.

IV Conclusions and Future Work

In this work, we have numerically computed the CC free energy between inclusions in two-dimensional critical binary mixtures. To this end, we have introduced a versatile thermodynamic integration scheme whose calculations for single-domain inclusions compare very favourably with the existing literature on the subjectMachta et al. (2012). This has enabled us to identify repulsive and attractive CC regimes respectively for unlike inclusions — with opposing affinities for the surrounding species comprising the binary mixture — or like inclusions — with the same affinity for the surrounding species. We noted that the CC interactions for single-domain inclusions are very long range and their magnitude is increasing with the strength of the affinity of the domain with either of the surrounding species.

A key advantage of our thermodynamic integration approach is that we do not need to supply a reference free energy at large distances, which can be notoriously difficult to compute analytically for complex scenarios. In turn this allows us to study inclusions comprising multiple domains such that, overall, they do not have any net preference for either of the surrounding species in the mixture.

The presence of domains with opposing affinities does not cancel the CC interactions altogether but instead leads to a decrease in the range and magnitude of the effective interactions in a manner reminiscent of electrostatic multipoles. For such multiple-domain inclusions, both repulsive and attractive CC interactions were observed depending on their relative orientations.

If biologically relevant lipid membranes are close to criticality as suggested by some studies Honerkamp-Smith et al. (2008); Veatch et al. (2008), then these findings may prove useful to better understand lateral spatial organisation within cellular membranes. As a matter of fact, attractive CC interactions between the single-domain inclusions of the same kind provide both sufficient magnitude and range for protein clustering to occur in equilibrium conditionsEdison et al. (2015); Nguyen et al. (2016). Utilising multiple-domain inclusions as more realistic models of membrane proteins lead us to the conclusion that the emerging CC interactions between them also suffice for them to aggregate, even where there is no overall net affinity to the surrounding species in the mixture.

There are several avenues for future work. Firstly, it would be interesting to study the clustering, and more generally the phase diagram, of the multiple-domain inclusions. For instance, in Fig. 9(a) we show the relative stability of the different orientations for two-domain proteins. Since some orientations are energetically favourable, we would expect proteins to rotate in order to minimise the overall free energy of the system. Our preliminary results suggest that, as these inclusions start to cluster, they will form an alternating patterning in affinities with respect to the surrounding species (Fig. 9(b)). In this context, it would be exciting to explore how CC interactions could be tailored to allow components in biological membranes achieve specificity as suggested in Ref. Sieber et al. (2006b). Such strategies could also be exploited for the self-assembly of anisotropic (e.g. Janus) colloidal particles Labbe-Laurent and Dietrich (2016); Noruzifar and Oettel (2009). Secondly, Casimir interaction due to shape fluctuations of the membrane has also been proposed as a non-specific mechanism for membrane proteins to interact Golestanian et al. (1996); Bitbol et al. (2010); Yolcu et al. (2011); Lin et al. (2011); Weikl (2001). Here it will be instructive to compare the strengths of the two possible Casimir interactions, and map regions in parameter space (temperature, composition and bending rigidity) where one mechanism dominates over the other, and where they interfere either constructively or destructively. For instance, Dean et al.Dean et al. (2015) suggests that taking into account membrane shape fluctuations can result in a shift in the critical temperature at which phase separation occurs for the lipid mixtures. Stress tensor-based methods following Ref. Bitbol et al. (2011) could help inform on these competing effects.

Conflicts of interest

There are no conflicts of interest to declare.

Acknowledments

We thank Pietro Cicuta, Mark Miller and Lorenzo di Michele for helpful discussions. We acknowledge funding from EPSRC, grant EP/J017566/1.