Curvature correction to the mobility of fluid membrane inclusions

Using rigorous low-Reynolds-number hydrodynamic theory on curved surfaces, we provide, via a Stokeslet-type approach, a general and concise expression for the leading-order curvature correction to the canonical, planar, Saffman-Delbrück value of the diffusion constant for a small inclusion embedded in an arbitrarily (albeit weakly) curved fluid membrane. In order to demonstrate the efficacy and utility of this general result, we apply our theory to the specific case of calculating the diffusion coefficient of a locally curvature inducing membrane inclusion. By including both the effects of inclusion and membrane elasticity, as well as their respective thermal shape fluctuations, excellent agreement is found with recently published experimental data on the surface tension dependent mobility of membrane bound inclusions.


Introduction
The lateral diffusion of membrane embedded proteins is vital to many physiological processes required in order to maintain life, such as regulating ion transport, maintaining cellular adhesion, and signal transduction [1]. The diffusion constant for planar fluid membranes was calculated using hydrodynamic theory some time ago by Saffman and Delbrück [2], and has been widely used and applied since, however, many physical membranes typically possess some degree of curvature. In this work we use covariant, classical, low-Reynolds number, hydrodynamic theory in order to rigorously elucidate the intrinsic curvature dependence of the diffusion coefficient of an embedded inclusion, in a generally (albeit weakly) curved fluid membrane, or surface. This is carried out using a Stokeslet-type approach, which has been successfully applied in previous work to calculate the diffusion constant for a small particle embedded in non-planar membranes [3,4].
Much modelling work has been carried out recently on trying to understand the mobility of membrane deforming inclusions, using various theoretical approaches [5][6][7][8][9][10]. While one would hope that there exists a unique physical mechanism and explanation for the experimentally observed results of e.g. [5], there exist numerous important and significant issues with theories such as that presented in [5], for example. Firstly, the effect of curvature on membrane hydrodynamics is not taken into account in [5], unlike our approach which explicitly and directly uses the governing low-Reynolds-number hydrodynamic equations. Secondly, there still exists a need for a classical hydrodynamic explanation, as provided by our work, since a e-mail: D.R.Daniels@swansea.ac.uk the Saffman-Delbrück result must be inserted by hand in the theory of [5]. Thirdly, in order for the theory of [5] to fit the experimental data, other additional dissipative mechanisms must be invoked, mainly due to membrane shear. Fourthly, the fit to the data in [5] requires the inclusion size to be much larger than is typically measured experimentally, such that it is postulated in [5] that the inclusion drags along with it a large patch of membrane lipids. Finally, moreover, the theory of [5] requires a membrane/inclusion coupling coefficient in order to fit the experimental data.
Indeed, the theory presented here is much closer in approach to the recent work of [10], which also used classical hydrodynamical theory to calculate membrane-bound inclusion mobility, though not via a Stokeselet-type approach as used in this work. Moreover, the work of [10] discusses a single, particular, geometry only, namely that of [5], unlike the general result presented here. Furthermore, the final result contained in [10] does not obviously display an explicit dependence of inclusion mobility on the intrinsically meaningful, covariant, quantity available for surfaces, namely the Gaussian curvature. A direct, term by term, comparison of the rather lengthy, numeric, expression for inclusion mobility contained in [10], with our more general result contained in eq. (8) of below, is therefore rather non-trivial. Additionally, in order to reproduce the experimental data, a variation in the mean contact angle of the membrane inclusion with membrane tension was found to be necessary, which required the introduction of an additional fit parameter, in the form of a proposed torsional stiffness constant in [10]. In this work we alternatively model inclusion elasticity via a harmonic potential for the membrane inclusion contact angle, with the stiffness of the harmonic potential playing the role of the torsional stiffness in [10]. Finally and importantly, given the relatively concise form obtained in this work for the diffusion constant, we are additionally able to investigate the role of both inclusion and membrane shape fluctuations on the mobility of membrane bound inclusions. Interestingly, and unlike the work presented in [10], by including such fluctuations we are able to account for variations in inclusion diffusivity with membrane tension even in the absence of a mean, or average, contact angle between the membrane and inclusion.
Our motivation in this work is to provide a simple, analytical, and general expression for inclusion mobility on arbitrarily (albeit weakly) curved membranes, including fluctuations. This result is valid for any small deformation membrane geometry, is not tied to a specific geometry, and hence does not require entire hydrodynamic calculations to be redone, on a case by case basis, thus avoiding substantial computational effort in future work. This result is given by eq. (8) of below, and constitutes the main result of this paper. Having obtained such a general expression for inclusion mobility, we naturally wish to demonstrate its efficacy and usefulness by applying it to the recently studied experimental data of [5] on inclusion mobility.

Low-Reynolds hydrodynamics in curved membranes
The low-Reynolds-number hydrodynamics for incompressible fluid flow, with velocity u α , in a curved membrane with viscosity η, is described by the following equations [4,10]: In eq. (1), Δ = D β D β , and D α is the covariant derivative given by D β u α = ∂ β u α + Γ α βγ u γ , where Γ α βγ = 1 2 g αδ (∂ γ g δβ +∂ β g δγ −∂ δ g βγ ) in terms of the metric g αβ [11]. Additionally, K is the local Gaussian curvature, p is the pressure, and we have also introduced a Stokeslet term [3,4,12,13], F α = 1 η f α δ 2 (x)/ √ g, located at the origin, which we require in order to enforce the no-slip condition at, and hence calculate the hydrodynamic drag on, a small, testparticle embedded in the membrane. In this work, we assume that the membrane fluid viscosity, μ m , is much greater than the viscosity of the embedding, or bulk, fluid, μ f , (typically μ m /μ f ∼ 10 2 [5]). One cannot however completely ignore the role of the bulk fluid, compared to that of the membrane, when calculating the viscous drag on a membrane inclusion, due to the socalled Stokes paradox. This leads to the well-known result that low-Reynolds-number hydrodynamics in planar, twodimensional, surfaces typically produce fluid flows that diverge logarithmically at large distances [10,14]. By solving the full boundary-value problem (incorporating a no-slip condition with the embedding fluid) it was found in [10,14] that the main influence of the viscous drag associated with the bulk fluid is to effectively introduce a long-distance cut-off to the two-dimensional surface hydrodynamics (for a nice discussion of this point see also [15]). This large distance cut-off is given by the Saffman-Delbrück length, l SD = η/μ f , where η = μ m h 0 , and h 0 is the membrane thickness. Moreover, following [10], we can also assume that the bulk fluid flow remains uninfluenced by any deformation of the membrane, provided that the characteristic length scale associated with such membrane deformations, r 0 = κ/σ, is much less than the Saffman-Delbrück length, l SD (where σ and κ are the membrane tension and rigidity respectively). Using the following values from [5,10] of η = 6 × 10 −10 kg s −1 , μ f = 10 −3 kg m −1 s −1 , and κ = 20 k B T , it can be shown [10] that this condition is indeed satisfied for membrane tensions σ 2 × 10 −7 N m −1 , which is always typically the case physiologically. Moreover, using the above values for η, μ, and κ, and assuming a typical membrane tension of σ = 10 −5 N m −1 , we obtain r 0 = 90 nm and l SD = 600 nm.
Due to the tensorial nature of the covariant derivatives involved, and the concomitant proliferation of indices, the solution of eq. (1) represents a formidable challenge, even in the perturbative regime. Hence we proceed as follows [4,10]. We can straightforwardly satisfy the incompressibility condition present in eq. (1) by introducing a stream function ψ, such that: u α = ε αβ D β ψ, where we use the epsilon notation throughout such that ε αβ = αβ / √ g and ε αβ = αβ √ g, and αβ is the ordinary two-dimensional, flat space, antisymmetric tensor 1 . Using this expression for u α , we can eliminate the pressure term from eq. (1) by taking the two-dimensional "curl" (see footnote 1 ) of eq. (1), in order to arrive at the governing equation for ψ A similar equation for ψ was found in [4,10], albeit without the Stokeslet term described by F α . Using the following fluid Green functions, defined via: Hence the fluid velocity, u α = ε αβ D β ψ, becomes simply 1 We use the following properties of ab throughout: ab Since we wish to work perturbatively to first order in the Gaussian curvature K(x), we can replace the full stream function ψ(x) appearing in eq. (4), with the zerothorder (in K(x)) term from eq. (3), to give consistently As discussed above, membrane fluid flows obtained using eq. (5), will typically diverge at large distances [10,14,15]. Moreover, following on from the discussion above (and as outlined in [10,14,15]), such large distance divergences can hence be remedied in a consistent manner (taking into account the viscous drag associated with the surrounding fluid [10,14,15]), via the use of an effective longdistance cut-off, given by the Saffman-Delbrück length: l SD = μ m h 0 /μ f . Furthermore, we reassuringly find below that the use of such an approximate treatment for the viscous resistance of the embedding fluid, will allow us to both re-derive the classical Saffman-Delbrück result for the diffusion constant on planar membranes [2], as well as enable us to accurately capture the experimental data on the mobility of curved membrane bound inclusions, as found for example in [5].
In the work of [5][6][7], the surrounding fluid was alternatively included via the use of an Oseen-type approximation, which gave rise to additional dissipative terms entirely due to viscous losses in the surrounding solvent. Our work, therefore, which does not use this type of approach, naturally does not contain such additional terms. Moreover, it was shown in [5] that such additional dissipative terms (as calculated in [5][6][7] and due entirely to viscous losses in the surrounding solvent) were unable to account for the experimental results of [5], and still required extra dissipation due to membrane shear, for example. Our work, similar to [10], provides a classical hydrodynamic explanation for the additional drag, and in doing so we find no additional dissipative mechanisms are required in order to fit the experimental data of [5].

General expression for the diffusion constant D in curved membranes
In order to calculate the diffusion constant, we need to find u a (0). Integrating eq. (5) by parts, inserting the expression for the Stokeslet force F α (x) given above, and evaluating at the location of the membrane embedded, point-like, test-particle (x = 0), we get Due to the symmetry of the problem, the integral over all x required in eq. (6) must give a contribution proportional to the symmetric tensor δ σ β , and we can also conveniently use the result: . By additionally using the properties of the antisymmetric tensor ε αβ (see footnote 1 ), we arrive at our final expression for u a (0) Using the relationship f a = ξu a (0) [12,13], that pertains between the applied force and the fluid velocity at the particle's location (where ξ is the coefficient of friction), we can straightforwardly read off the value of the diffusion coefficient D = k B T /ξ from eq. (7) as being The canonical Saffman-Delbrück result: for a planar membrane [2], with a longdistance cut-off given by the Saffman-Delbrück length l SD = η/μ f , and a short-distance cut-off a given by the size of the inclusion (with γ being Euler's constant), can be obtained from the first term in eq. (8), via careful consideration of the Green function G 0 (x) (as outlined below). Our expression, eq. (8), for D = D 0 + δD additionally contains within it the first-order correction (δD) due to membrane curvature. This general expression for the diffusion constant, as given in eq. (8), with its explicit dependence on the intrinsic curvature K(x), represents the main result of this work. We re-emphasise that this general result is valid for any arbitrary (weakly curved) surface geometry.

Membrane and inclusion elasticity including fluctuations
As a physically important and specific application of our general theory outlined above (and in particular eq. (8)), we now derive the mobility of an embedded, membrane deforming, inclusion. We use the following total Hamiltonian H in order to describe both the membrane and inclusion elasticity The membrane part [16,17] is written in the Monge representation, where the membrane shape is given perturbatively via a height function h(x), which contains surface tension (σ), and rigidity (κ) controlled terms. We have also included in eq. (9) a term governed by λ which constrains the shape of the membrane bound inclusion to be h (a) = ∂h(r) ∂r | r=a to be equal to β. The final term in eq. (9) describes the elastic response of the inclusion via a simple harmonic potential for β (with strength k p ) around some preferred shape β 0 . A similar, though not identical approach to modelling inclusion elasticity, and associated shape-change, was used in [10], via the introduction of a torsional stiffness force. However, fluctuations in the inclusion shape governed by β were ignored in [10], while they are included in the work presented here. Following [18], we split the membrane contribution into a ground-state part (h 0 ) and a fluctuating part (δh). Minimising H in eq. (9), subject to appropriate boundary conditions, it can be straightforwardly found that the dominant membrane ground-state contribution is given in terms of modified Bessel functions by: h 0 = −βr 0 K 0 (r/r 0 )/K 1 (a/r 0 ) [5,10], where r 0 = κ/σ. Note that in this approach [18], fluctuations of the inclusion (governed by β) are included in the equilibrium part of the membrane Hamiltonian (via h 0 ). Membrane fluctuations, given by δh, can thus be seen to be independent of the inclusion shape β in this approach [18], such that they must now satisfy the boundary condition: δh (a) = 0. Our total Hamiltonian H thus becomes where k m is given by: k m = 2πκ a r0 K 0 (a/r 0 )/K 1 (a/r 0 ), which can be seen to depend on the membrane tension σ via r 0 . Due to the underlying rotational symmetry assumed in this work, and by inspection of eq. (8), we can simplify matters considerably by observing that all the quantities we require in order to calculate the diffusion constant D must depend on the radial distance r only. Hence, in the Monge gauge, the Gaussian curvature K(r) is given perturbatively by: K(r) = 1 2r ∂ ∂r h 2 (r), and we can write the membrane fluctuation modes as: Performing the required thermal averages, using eq. (10) we find where K 1 (a/r 0 )I 1 (r/r 0 ) − I 1 (a/r 0 )K 1 (r/r 0 ) K 1 (r/r 0 )/K 1 (a/r 0 ).
The first term in eq. (11) includes both the elasticity and fluctuations of the membrane bound inclusion shape β. The second term in eq. (11) describes membrane shape fluctuations around the ground state solution.
In this work, we neglect in our membrane Hamiltonian a possible term given by 1 2 κ G √ gd 2 xK due to Gaussian curvature, which has been considered in related work on inclusion bound membrane elasticity [19][20][21]. However, it can be shown via a straightforward calculation that such a term can easily be included in our work, should we wish to do so, via a renormalisation of the parameter k m defined above as follows: k m → k m − πκ G , where κ G is the elastic modulus for Gaussian curvature.

Fluid Green functions G 0 and G 1
We now proceed to calculate perturbatively the fluid Green functions on a curved surface, as required in eq. (8) for the diffusion constant D, in terms of the membrane height function h(r). With a metric component: g rr = 1 + h 2 and determinant such that √ g = r(1 + h 2 ) 1/2 , the radial Green function G 0 (r) satisfies (for r > 0) Evaluating eq. (13) perturbatively in membrane curvature, and matching the solution thus obtained at large distance to the analogous Green function obtained by solving the full boundary-value problem (incorporating a noslip condition with the embedding fluid), as carried out similarly in [10,14], we obtain that where for convenience we have defined the large-distance length-scale l SD as l SD = l SD exp(−γ). Using this approximate treatment we guarantee, by design, that the leading-order contribution to the diffusion constant (D 0 ) ultimately agrees precisely with that calculated in [14], which rigorously includes the effect of the surrounding bulk fluid. Additionally, as outlined above, we can safely assume that the bulk fluid flow remains uninfluenced by membrane deformations, for all physically reasonable membrane tensions [10]. Similarly, from above, the Green function G 1 (r) satisfies to lowest order (as required consistently by eq. (8)) which, incorporating the appropriate boundary condition at large distances, such that: G 1 (r) → 0 as r → l SD , has the solution

First-order correction δD to the inclusion diffusion constant
Using the Green functions G 0 (r) and G 1 (r) in our general result given by eq. (8) of above, and integrating by parts, we can write the first-order correction to the inclusion diffusion constant, δD, due to membrane curvature, in the Monge gauge as with h 2 (r) given by eq. (11) and eq. (12). Via the variation of r 0 = κ/σ with σ, we can see that the magnitude of the first-order correction to the diffusion constant δD, as given by eq. (17), necessarily depends on the value of the membrane tension. It is interesting therefore to investigate the limits of eq. (17) at both low and high membrane tensions σ. Assuming l SD /a 1, in the low tension limit, σ → 0, we find that δD → − kB T 8πη (β 2 0 + 1 kp − 1 4πκ ) log 2 (l SD /a), while for the large tension limit, σ → ∞, we additionally find from eq. (17) that δD → − kB T 8πη 1 2πa √ κσ log 2 (l SD /a), to leading order.

Comparison with recent experimental data
Single-particle tracking experiments [5] on reconstituted membranes have found that membrane curvature, surface tension, and inclusion shape can have a significant effect on the mobility of membrane embedded proteins. As a specific application, and to probe the generic membrane tension dependence of eq. (17), we can compare our general theory with the recently obtained experimental data of [5] on membrane inclusion mobility. In [5] it was found that the diffusion constant for a voltage-gated potassium channel (KvAP) protein was significantly increased, as the surface tension was increased, whereas the mobility of a water channel aquaporin 0 (AQP0) protein was, relatively, fairly insensitive to membrane tension. This is thought to occur due to the fact that KvAP locally bends the membrane considerably, forming an effectively conical inclusion, whereas AQP0 is thought to negligibly deform the bilayer locally [5,22]. Shown in fig. 1 is the experimental data from [5] on the variation of the diffusion constant D with membrane tension for AQP0 (triangular data points) and KvAP (circular data points). Also shown in fig. 1 are the theoretical fits (solid black lines) using eq. (17) for δD, where we have also included the contribution from D 0 = kB T 4πη log(l SD /a), such that D = D 0 +δD. In this work we fix the following parameters, used for both AQP0 and KvAP, to those from [5,10]: η = 6 × 10 −10 kg s −1 , μ f = 10 −3 kg m −1 s −1 , κ = 20 k B T , a = 5 nm, and k B T = 4 × 10 −21 J. However, in order to  Fig. 1. Plot of the diffusion constant D versus the logarithm of the membrane tension log 10 σ. Experimental data from [5] is shown for AQP0 (triangular data points) and KvAP (circular data points). Also shown are the theoretical best fits (solid black lines) using eq. (17), and the model parameters quoted in the text. For KvAP we require β 0 = 0.5 rad and kp = 50kBT , while for AQP0 we find β 0 = 0 rad and kp = 32kBT .
best fit the mobility data for KvAP, it was found necessary to let β 0 = 0.5 rad and k p = 50k B T , whereas AQP0 required the best fit parameter of β 0 = 0 rad and k p = 32k B T . Reassuringly, the best fit model parameters used in this work compare favourably with those similarly reported elsewhere [5,10].

Discussion
This theoretical work addresses the general problem of calculating via classical hydrodynamic theory, the leading-order curvature correction to the canonical, planar, Saffman-Delbrück value of the diffusion constant for a small inclusion embedded in an arbitrarily (albeit weakly) curved membrane. This result was then applied to the specific example of a locally membrane deforming protein membrane inclusion, as contained in the recent experimental case of [5]. The application of our main result, given by eq. (8), to other possible membrane geometries of interest (such as the catenoid, for example) will be left to future work. We also look forward to further experimental investigations of our general result, including possible applications to "man-made", patterned (or "bumpy"), twodimensional thin fluid films, in addition to more familiar bio-membranes. As in [10], we found that the rigidity of the protein is an essential ingredient in order to fit the experimental data of [5]. However, the theory outlined in this work was also able to include the role of thermal fluctuations on inclusion mobility. Such fluctuations have an effect on the membrane as well as the inclusion shape. Indeed, unlike [5,10], we find above that by taking into account thermal fluctuations we are also able to accurately capture the membrane tension dependent mobility of an AQP0 inclusion even when its average contact angle ( β ∼ β 0 ) vanishes (β 0 ∼ 0) [22]. Furthermore, using our theory we can now be seen to be able to directly probe and investigate inclusion elasticity (governed by the elastic constant k p ) via measurements of the inclusion's mobility.
We can also compare our result for the diffusion coefficient given by eq. (8) to previous results obtained in non-planar membranes for inclusions which do not deform the membrane locally [4,3,21]. In the case of a straight, membrane tube, for example, the Gaussian curvature is exactly zero, and so our expression for the diffusion constant agrees with that given in [3], if we choose the longdistance cut-off in this case to be given by the radius of the membrane cylinder. Indeed, for any developable surface (K = 0) the intrinsic curvature vanishes identically, and so the mobility reverts to the usual form of the Saffman-Delbrück result. For spherical membranes, where the Gaussian curvature is constant, the integral required in eq. (8) gives a negligible contribution, which compares favourably with expressions found in [4,3], as long as we choose the long-distance cut-off in this case to be given by the radius of the membrane sphere.
Finally, the general result of the work presented here, as given by eq. (8) of above, provides a relatively simple expression for inclusion mobility on arbitrarily (albeit weakly) curved membranes, avoiding the use of intensive numerical computation, and can therefore be straightforwardly used for the practical analysis of future experimental data on membrane hydrodynamics and inclusion mobility. By utilising this result, therefore, substantial computational effort can be avoided in future when considering inclusion mobility on other, practically important, membrane geometries. The theoretical work presented here is also likely to be highly relevant to experiments pertaining to many biological processes involving the mobility of membrane embedded inclusions and their dynamic spatial positioning, such as membrane receptor clustering and ligand association [1], for example.
Useful discussions with M.S. Turner and R.G. Morris are gratefully acknowledged. We would also like to thank the anonymous referees for their helpful comments.

Open Access
This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.