Two intertwined facets of adherent membranes: membrane roughness and correlations between ligand–receptors bonds

We study equilibrium fluctuations of adherent membranes by means of Langevin simulations in the case when the interaction of the membrane with the substrate is twofold: a non-specific homogeneous harmonic potential is placed at large distances, whereas discrete ligand–receptor interactions occur at short distances from the flat substrate. We analyze the correlations between neighboring ligand–receptor bonds in a regime of relatively strong membrane fluctuations. By comparison with the random distribution of bonds, we find that the correlations between the bonds are always positive, suggesting spontaneous formation of domains. The equilibrium roughness of the membrane is then determined by fluctuations in the number density of bonds within the domains. Furthermore, we show that the excess number of bonds arising due to correlations and the instantaneous roughness of the membrane both follow master curves that depend only on the instantaneous bond density and not on the intrinsic binding strength of the ligand–receptor pair. The master curves show identical trends, further corroborating the link between membrane roughness and bond correlations.


Introduction
The adhesion between fluid membranes and flat substrates mediated by specific ligand-receptor bonds has been extensively studied in recent decades (for a recent review see [1]). One driving force underlying such an extensive scientific activity is the relevance of this system to the adhesion of living cells [2]. A typical experimental realization of an adhesion assay involves a protein-carrying giant unilamellar vesicle that binds to a functionalized substrate [3]- [5]. The first step of adhesion is associated with the sedimentation of the vesicle into a non-specific potential minimum at relatively large distances from the substrate [6,7]. The fluctuations of the membrane allow then for the formation of first ligand-receptor bonds [8] upon which tight contacts close to the substrate are established and specific adhesion takes place. Such behavior is well described as adhesion in a potential acting between the membrane and the substrate that contains two wells: a non-specific one further away from the substrate and a specific one in its vicinity [8,10,11].
If bare glass substrates were involved in the adhesion assay, the non-specific interactions driving the first step of adhesion would be dominated by the relatively strong van der Waals attraction of the glass. Similarly, substrates with strong Coulomb attraction could be envisaged. In these cases, the Helfrich steric repulsion imposed by fluctuations of the membrane would be overcome and the vesicle membrane would strongly adhere on the substrate. The non-specific minimum would be positioned very close to the substrate and its depth considerable. For the control of the generic non-specific substrate-vesicle interactions, a variety of experimental techniques have been developed [9], often involving soft polymers incorporated either into the vesicle membrane and/or deposited on the substrate. As a result, the non-specific minimum is pushed to far distances from the substrate and is characterized by a considerably smaller strength and curvature than the specific minimum associated with the formation of ligand-receptor bonds. Consequently, the two minima are well separated, which is a prerequisite for a successful specific adhesion assay as it allows for proper differentiation between non-specific and specific adhesion [1].
From a theoretical point of view, adhesion is often studied in the presence of membrane fluctuations but without the non-specific potential, which ignores subtle effects induced by the presence of a secondary minimum in the free energy. Even so, for the case of an effective membrane-substrate potential with a single minimum, valuable information on cooperativity A membrane made to adhere on a substrate by receptor-ligand bonds. The rest length of the tether l 0 representing a bond is marked. The local instantaneous height h(r ), which deviates from the minimum of the non-specific potential placed at h 0 , is shown. The potential of strength γ is depicted with the red curve.
of bonds [12,13], membrane roughness [14,15] and the development of domains has been obtained [16,17]. The analysis of a Hamiltonian containing a non-specific interaction potential [18]- [20] shows that a large section of the phase space is actually associated with a single well effective free energy. However, deeper understanding of the relevant mechanisms driving the adhesion in the regime relevant to experimental realizations of the system must involve a secondary minimum. Then a study where microscopic degrees of freedom are monitored directly, and their impact on the macroscopic behavior of the system analyzed in detail, should be performed. To address that issue, we concentrate on the part of the phase space where the free energy has features of a double-well potential. We study the equilibrium of the bonded-membrane system that is propagated in time by Langevin dynamics simulations in which the membrane can experience significant fluctuations. Our goal is to understand the interplay between specific and non-specific interactions and show that by means of finite elasticity, membrane roughness couples to correlations between bonds formed on the neighboring sites.

Model
The system (figure 1) consists of a fluctuating membrane that interacts with a substrate by a weak homogeneous harmonic potential whose strength and distance from the substrate are derived from the position and curvature of the non-specific minimum found in experiments [6,8]. Further, specific ligand-receptor interactions are imposed on a regular square lattice and are modeled as harmonic tethers. Ligands embedded in the fluctuating membrane are arranged on a regular lattice and move with the membrane only in the direction perpendicular to the surface. Receptors are immobilized at the same positions r i on the flat substrate. Each bond has a statistical probability to be open or closed, depending on the local temporal position of the membrane (see the appendix for details). The appropriate energy functional for a system including both, the membrane and the bonds, is The first term expresses the bending energy of a membrane with projected area A and bending rigidity κ. The membrane height profile h(r) is given in the Monge gauge as a function of its projected position r on the substrate plane. The non-specific membrane-substrate interaction, given by the second term, is modeled by a harmonic potential with strength γ placed at a distance h 0 from the wall (l 0 < h 0 ). The last term describes the total interaction between ligands and receptors positioned opposite to one another at r j , j indexing each of N t binding sites. Since the formation of a bond is taken to be a statistical event, only formed bonds with b j = 1 contribute to the free energy. Otherwise, if the bond is open, b j = 0. Since both the binding and the unbinding rates (k on and k off , respectively) are associated with distributions, the energy that each bond contributes to the free energy is not a unique value, but is a distribution too. Consequently, the binding affinity on the level of a bond [8] is defined as the average of its respective distribution and is intimately related to the detailed balance condition (equation (A.1)) that governs the formation of bonds. For a system described by equation (1), the free energy of a binding site f (h) can be calculated explicitly for the case of an infinitely stiff membrane (κ → ∞) from the partition function [18,20] where ρ ≡ N t /A is the density of binding sites, and β ≡ k B T , k B being the Boltzmann constant and T the temperature. As shown in figure 2, this free energy can have two minima. If an infinitely large and stiff membrane is placed in such a potential, it will reside in the minimum with the lower energy. As the relative strength of the minima changes in response to the variation of some parameter (e.g. intrinsic binding strength or the spring constant), the membrane jumps between two minima and exhibits a first-order transition. At finite temperatures and in a system of finite size, the stiff membrane as a whole spends a fraction of time in one minimum and the remainder of time in the other minimum. The occupancy of the two minima is given by the respective Boltzmann factors, the latter themselves depending on the membrane size.
For the simulations of a membrane that experiences both a finite bending stiffness and thermal fluctuations (βκ = 10), the parameters are chosen in such a way that the underlying free energy possesses the two minima (figure 2). Specifically, the system is parameterized by a square 64 × 64 lattice with a lattice constant a 10 nm. The system has N t = 64 binding sites arranged on a square 8 × 8 lattice such that the distance between neighboring binding sites is 8a. The time scale of the system is given by the viscosity of the water surrounding the membrane, which, at a temperature of T = 300 K, is η = 10 −9 Js cm −3 or 2.4 × 10 −7 k B T sa −3 in simulation units. The reaction rate of the unstressed tether is set to k 0 = 2000 s −1 and the discrete time step in the simulation is chosen as t = 10 −9 s, which is well under the smallest membrane time scale τ = 4η/(βκk 3 max ) 1.1 × 10 −9 s. This reaction rate is somewhat faster than the typical biological rates. However, the choice of k 0 has no effect on the equilibrium averages of the system. For dynamic fluctuations around the equilibrium, due to the clear separation of bond and membrane time scales, k 0 simply sets the time unit, a conclusion that has been tested by varying k 0 over several orders of magnitude (data not shown). The strength and for a 2 β K = 1.25 for the case of an infinitely stiff membrane. The nonspecific potential is placed at h 0 = 12a, while the rest length of the tethers is l 0 = 8a. Two minima are observed for ln (1.39) < β b < 5.4. The so-called non-specific minimum at h 12a is the signature of the bare non-specific membrane-substrate potential. A second minimum with a higher density of ligand-receptor bonds resides at lower heights h 8a and is, because it originates from the ligand-receptor interaction potential, named the specific minimum. At the critical e b 2.7, the two minima have the same free energy, and a first-order phase transition takes place for an infinitely stiff and infinitely large membrane.
γ of the bare non-specific potential is γ = 3.125 × 10 −3 k B T a −4 , which equals 1.25 × 10 9 Jm −4 , and the potential minimum is placed at h 0 = 12a. The rest length of the springs is l 0 = 8a. The typical intrinsic binding energy b of a ligand-receptor pair is a few k B T , while the spring stiffness is set to K = 1.25 k B T a −2 (equivalent to K 5 × 10 −5 Nm −1 ). For each set of parameters, the equilibrium dynamics simulations were performed following the protocol described in the appendix.
Fluctuations of the membrane and the membrane elasticity change the potentials shown in figure 2, and the membrane can reside in both minima simultaneously. The phase diagram for the simulated system can be constructed by evaluating the overall mean heighth or the mean bond number densityφ as a function of the intrinsic binding strength. Hereh andφ are defined ash with the sums accounting for the average over all instances of discredited time t i , and i 1 equals the total number of frames available to the analysis. The brackets indicate averaging over spacial coordinates of the membrane profile h i (r) or bond realizations b j,i at each bond site r j , the latter giving rise to the instantaneous number density of formed bonds φ i . A bar indicates a property that has been averaged over both the temporal and spatial coordinates. This notation for averaging will be used as a convention throughout this paper. As can be seen from the simulation data shown in figure 3, the transition occurs between β b = 3.45 and β b = 3.50. At binding strengths lower than β b = 3.45 (lower branch), configurations with no bonds occur in equilibrium (specific adhesion is unstable), whereas at binding strengths larger than β b = 3.5, a finite number of bonds are always present in the system, and the specific adhesion is stable (upper branch). This behavior is reminiscent of the first-order transition discussed above. From the simulated time evolution of the membrane profile it is possible to extract the membrane height distribution function (data not shown). An effective free energy can then be defined as the logarithm of that height distribution. As such, the effective potential has also two minima, although in a somewhat reduced range of intrinsic binding strengths β b . Interestingly, the critical binding strength at which the transition occurs is larger than both the binding strength (β b 3.1) at which the effective potential obtained from simulations has two degenerate minima and the binding strength (β b 2.7) at which the effective potential for a stiff infinite membrane has two degenerate minima (figure 2).
Typically, a soft fluctuating membrane patch has, on average, a smaller number of bonds than the stiff membrane segment of the same size, which can be calculated by the use of the partition function Z. The only exception is very low intrinsic binding strengths at which a stiff membrane makes no bonds. In this regime, a fluctuating membrane establishes some unstable bonds contributing to the non-zero average bond density, due to random excursions to the close proximity of the substrate. Furthermore, the observation is that even in such cases, specific 7 adhesion often appears not as individual bonds but as few bonds together. This indicates that in adhesion of fluctuating soft membranes, correlations between bonds may arise due to the deformability of the bilayer.

Correlations between bonds
Correlations between neighboring bonds emerge from membrane deformations taking place in the vicinity of a bond in a region whose size is typified by the lateral correlation length. If there is another binding site within this region, the probability for forming a bond at this second site increases with respect to the part of the membrane fluctuating in the non-specific minimum [8], simply because the ligands and receptors are brought closer together. Since the membrane above the neighboring site is on average, due to deformations, closer to the substrate, the formation of bonds is locally promoted, or said in another way, correlations between bonds are induced. There are no correlations between bonds in an infinitely soft or infinitely rigid membrane, for which the lateral correlation length is either zero or infinite because the local environment where the membrane is deformed does not exist. On the other hand, the fluctuations of the membrane effectively decrease the strength of the bare non-specific potential and increase the distance between the membrane that is not specifically bound and the substrate [21], which in turn affects the deformation of the membrane between the bonds. Correlations hence arise from the balance between the deformation and fluctuations of the membrane with finite bending elasticity on the one hand, and the non-specific potential that attracts the unbound membrane and induces deformations of the bonds on the other hand.
We focus on correlations between the bonds at nearest neighboring sites by calculating the relative mean coordination number. This is achieved by determining first the mean coordination numberc from calculating the number of bonds with 0-4 bonded neighbors in 10 3 equilibrium realizations of the system. The average number of bonded neighbors is then determined by calculating the mean over all bonds and frames. However, bonds placed on fully random positions at the density φ r will also have a finite coordination number c r with c r = 4φ r , for a large lattice ignoring finite size effects. For φ r fixed at the value corresponding to the equilibrium bond density in a given simulation (φ r =φ), positive correlations occur when the relative coordination number of a bond (c − c r )/c r is larger than zero. Here (c − c r ) is the mean excess coordination number, which is the average surplus of bonded neighbors relative to the random (uncorrelated) distribution.
Inspection of the left panel in figure 4 shows clearly that positive correlations typify the parameter space investigated herein. At low densities, in the unstable branch, the correlations between the bonds are strong. The equilibrium life time of a bond is extended in comparison to the life time of isolated bonds signifying coupling effects. At high densities of bonds, the correlations saturate to zero. Ultimately, when all sites are occupied by bonds (φ = 1, β b → ∞) the system does not differ from a random distribution, and all correlations are lost.
The probability for the appearance of a certain coordination number n for various parameter sets is shown for both the simulation data and random distributions (top right panels in figure 4). For the random distribution and a large lattice, the expectation value for the probability that a bond will have a coordination number n = 0, . . . , 4 is simply given by the binomial distribution In the simulations of soft fluctuating membranes, it is most likely that a bond will have one neighbor bound, at low binding affinities or bond densities (φ = 0.06) . However, the probability that the bond will have more than one neighbor is considerably larger than in the system with randomly distributed bonds, and the likelihood of appearance of a bond with four neighbors is more than 1000 times larger (bottom right panel in figure 4). As the binding affinity and the density of bonds increase, the most likely number of bound neighbors increases gradually. However, the difference in the distribution of neighbors between the simulated and the random system gradually decreases, consistent with the decrease of correlations between the bonds at larger binding affinities and mean bond densities.

Instantaneous membrane roughness and the coordination of bonds
For a deformable membrane adhering in a potential with two wells, the non-specific potential will attract the membrane between the bonds inducing deformations in the vicinity of the bonds. Furthermore, membrane fluctuations produce random forces acting on each bond. These forces are of larger amplitudes if the bonds are further apart [8], which affects the overall stability of the bonds. Obviously, if the bonds are confined into a domain, the membrane fluctuations are minimized within the domain, which may in turn stabilize the agglomerate. It is thus pertinent to study the roughness of the membrane in order to understand the correlations between the bonds. We first calculate the instantaneous membrane roughness as a function of the instantaneous density of bonds given by with the instantaneous roughness ξ ⊥ i calculated according to Likewise, in order to understand how the domain formation is related to the instantaneous density of bonds and the roughness of the membrane, the excess coordination number for a certain density is determined as The delta function in the above definitions takes care that the time average is made only over those i's (i.e. frames) for which at t i the number density of bonds is equal to the selected φ.
In figure 5, 4 the roughness ξ ⊥ (φ) and the excess coordination c(φ) − c r (φ) are shown for several intrinsic binding strengths b in the top and the bottom panel, respectively. The most striking feature in these graphs is that all data follow master curves, and the dependence on the intrinsic binding strength falls within the spread of data points belonging to a curve associated with a particular b . Furthermore, the master curves show very similar trends. Specifically, the roughness of the membrane increases until φ 0.2, upon which a decrease in roughness as a function of the bond density takes place. This clearly suggests that correlations between bonds and the instantaneous roughness are interrelated through the instantaneous number of bonds. Such an analogy of behavior in the excess coordination of bonds and the roughness ξ ⊥ shows that there is a certain connection between these two properties of the system. It is consistent with the observation that at low binding strengths the probability that a bond has a large number of neighbors is significantly bigger than in the random system, as shown in figure 4. This suggests that even at low affinity the bonds correlate to form unstable, short-lived domains. Furthermore, the deformations of the tethers are also insensitive to the intrinsic binding strength (data not shown). Consequently, if one fixes the density of bonds, the morphology of the bound membrane is dependent on b to a very little extent. The binding strength acts only to set the statistical probability for the appearance of certain bond densities.

Equilibrium membrane roughness and fluctuations in bond density
The equilibrium membrane roughness¯ ⊥ (φ), which should be distinguished from the instantaneous roughness ξ ⊥ (φ) discussed above, arises from the long-time behavior of the system, as a property averaged over both time and spatial coordinates of the membrane profile. To explore how this roughness¯ ⊥ depends on the system parameters, we choose an experimentally accessible definition Hereh is the overall mean height of the membrane (figure 3) and¯ ⊥ is evaluated directly from the simulations. In order to compare the equilibrium mean roughness with the random roughness r (φ) (figure 6), the latter first had to be evaluated. Therefore, for each bond density φ, 1000 membrane configurations with randomly distributed bonds have been generated. We calculate the mean roughness from this set of 1000 realizations according to equation (9), except that now index i enumerates different random configurations.
We investigate the ratio¯ ⊥ / r . Inspection of the top panel in figure 7 shows that this ratio is always larger than unity, signifying that, like the correlations, the equilibrium roughness of the membrane is always larger than in the system with randomly distributed bonds. From the point of view of the argument that the bonds have a tendency to form domains, such increased roughness is intuitively an unexpected result. Within domains, the membrane fluctuations (hence roughness) are supposedly suppressed, which in turn could lead to a decreased roughness with respect to the random distribution of bonds. The understanding of the obtained results, however, emerges from re-inspecting figure 3. Generally, for a fixed binding strength, the instantaneous bond density φ i adopts values in a certain range around the mean bond densityφ, whereφ( b ) was shown previously (left panel in figure 3). As the mean density of bonds changes (due to the variation of b ), the attained range of φ i , e.g. the intensity of deviations φ from the mean equilibrium density of bonds, defined as differs significantly (right panels in figure 3). However, every φ i contributes to¯ ⊥ through a membrane profile of a characteristic roughness ξ ⊥ (figure 5). Consequently, large deviations from the mean equilibrium density of bonds lead to large variations in the instantaneous roughness ξ ⊥ , which in turn affect the average equilibrium roughness¯ ⊥ . It is thus evident that the intensity of fluctuations in the equilibrium density has an impact on the overall equilibrium roughness¯ ⊥ . Indeed, a remarkable resemblance in behavior of¯ ⊥ / r and φ can be seen by comparison of the two panels in figure 7.
A more detailed analysis of ξ ⊥ in figure 5 confirms the above reasoning. At small mean bond density or binding strength b (the unstable branch in figure 3), the ξ ⊥ (φ) curves do not extend far from φ = 0, which allows only for small deviations φ. Nevertheless, the fluctuations are large relatively toφ ( φ/φ is considerable) andφ badly represents the instantaneous φ i (see lower panels on the right in figure 3). As b increases to intermediate values, φ i assumes almost all values between zero and unity. It is in this regime that¯ ⊥ as well as φ are maximal. For b at which the unstable branch approaches the transition point, the conformations with a large number of bonds somewhat occur more frequently, leading to a slight decrease in both φ and ⊥ / r . In the stable branch,φ well characterizes the system at any instance of time (see right upper panels in figure 3). The fluctuations in the density of bonds are less intense than in the unstable branch, as φ never reaches zero. As the bonds forming the adhesion domains are strengthened by the increasing binding strength of the ligand-receptor pair, φ decreases even more. The equilibrium mean roughness reflects this behavior very well. Ultimately, when the mean densitȳ φ reaches unity, the distribution of bonds does not differ from a random one and¯ ⊥ / r = 1. Similarly, when there is on average one or fewer bonds formed (φ → 0), a random configuration cannot be distinguished from any other configuration and¯ ⊥ / r reaches unity as well. In this regime, the deviations from the mean density of bonds reduce to zero.

Conclusions
The work presented herein shows that the role of fluctuations in membrane adhesion is twofold. On the one hand, the fluctuations decrease the overall adhesiveness of membranes and a larger effective binding affinity must be involved to obtain the same member of bonds as in non-fluctuating membranes. On the other hand, membrane roughness is intimately associated with positive correlations between the bonds simply because the latter drive the formation of domains. If the effective binding affinity (the mean contribution of a bond to the free energy) is such that the formation of a bond decreases the total energy, domains are stable. In this regime, the equilibrium roughness is small, and the correlations between the bonds, although positive, decay toward zero with increasing the mean density of bonds. However, the correlations between the bonds are very important in the regime of unstable specific adhesion where the number of formed bonds within the domain fluctuates significantly and the entire domain occasionally unbinds. In this case, equilibrium roughness of the membrane is significantly larger than the roughness of a membrane with no bonds or with randomly distributed bonds.
Since the instantaneous roughness reflects the instantaneous size of the domains, it should be, moreover, possible to use fluctuations as an imaging tool. This would be particularly useful in optical techniques where the lateral resolution of the order of 100 nm is still too large to observe individual bonds. However, local fluctuations of the membrane can be measured accurately. Corroborated by our predictions, membrane roughness in equilibrium could be, in principle, mapped onto the density of bonds within the domain.
It was recently shown experimentally that mobility of ligands and receptors furthermore affects the distribution of bonds within domains [27][28][29], where an increase of fluctuations have been observed in the process of the formation of the nucleation center [29]. Mobility certainly provides an additional entropic contribution to the free energy of the system. The current results are relevant to the situation in which one of the binding partners is immobilized (typically the receptors on the substrate), and the other binding partner is abundant in the opposing membrane (as in [4])-e.g. a high-density regime. Furthermore, the results discussed herein should be relevant for systems in which the diffusion time scale is much larger than the reaction time scale-a very-low-density regime that also has biological implications. However, the true impact of the diffusion remains to be elucidated quantitatively, in simulations and by other theoretical methods, particularly from the point of view of the relation between membrane roughness and correlations between bonds, revealed here as an essential feature of specifically adhering membranes.
14 The equation of motion for the membrane shape configuration associated with the energy given in equation (1) is given by [22] ∂h(k, t) ∂t = ξ(k) − (k) It was shown previously [23] that the Onsager coefficient depends on the mean distance between the membrane and the substrate. However, for the parameters used herein, this distance is sufficiently large and only the four smallest k-modes are affected appreciably, a result emerging from the comparison of simulations with the full treatment of hydrodynamic interactions with the current setup. For simplicity, only the hydrodynamic interactions between the solvent (of viscosity η) and the membrane are taken into account: (k) = 1/(4ηk), where k ≡ |k| and k > 0, whereas the interaction with the wall has been neglected. On the same level of approximation, for the k = 0 mode, which describes the mobility of the center of mass of the membrane, the Onsager coefficient is set to (k = 0) = 3 √ A/(8πη) [24]. In any case, averaged properties of the system do not depend on the treatment of hydrodynamic interactions.
The membrane shape is evolved by discretely integrating the Langevin equation (A.4) in time with a constant time step t on a square lattice [8,22]. Following each membrane step, the configuration of all receptor-ligand pairs is updated where an existing bond is broken with the probability k off t, while a bond is formed with the probability k on t. As can be seen from equations (A.2) and (A.3), the typical time scale for these events is set by the reaction rate k 0 . The initial configuration for each run is a membrane with no bonds placed at h = 12a and at least the first 10 ms were used for the equilibration and a production run was carried out for another 50 ms. The data are written out every 0.05 ms, which in total provides 1000 configurations for analysis.