$$\require{cancel}$$

# mathematical continuum descriptions of membranes

Membranes play an important role in cellular functions. They act as barriers to the outside environment, but also must allow exchange of nutrients and waster. They must also be physically flexible to allow for cellular growth and movement. Membranes are composed of phospholipid bilayers. The phospholipids are amphipathic, containing a hydrophilic head and a hydrophobic tail. The membrane is composed of two adjacent layers of lipids, with the tails facing each other, and the hydrophilic heads facing the outside aqueous environment (Figure 1). These membranes have a high aspect ratio with a low thickness of about 10 nm, as the membrane is only two lipid molecules thick. However, the hydrophobic interior of the bilayer is a strong barrier against the transport of hydrophilic molecules across the membrane. While diffusion across the membrane is very slow for hydrophilic molecules, lateral diffusion can be quite quick for lipids and even transmembrane proteins. Due to the described properties of the cellular membranes, they are commonly modeled as a 2-dimensional planar, elastic sheet in a 3-dimensional space (Figure 2). Biological structures like vesicles and even entire cells can be described as a closed form made by these sheets.

Figure 1: Structure of Biological Membranes (Adapted from [1])

To understand the shape and dynamics of membranes, we will need to describe the geometry of the membrane surface itself. In addition, we will describe how the energy changes for different membrane configurations. In particular, we will look at how the energy changes when a membrane is stretched, bent, sheared, and changed in the thickness (Figure 2).

Figure 2: Biological membranes can be modeled as a two-dimensional sheet that stretches, bends, and changes in thickness. Adapted from [1].

Area Changes:

The area deformation can easily be parameterized by a single parameter $$\Delta a$$ or a function $$\Delta a(x, y)$$ if the area deformation is dependent on position. We change treat the lipids in the membrane as system of connected elastic springs, with increasing the area equivalent to stretching the springs (Figure 3).

As such, the energy of stretching a piece of membrane is:

$E_{a}=\frac{1}{2} \frac{K_{a} \Delta a^{2}}{a_{0}}$

where $$K_a$$ is the area modulus and $$a_0$$ is some reference area from which change in area is calculated (i.e. $$\Delta a=a-a_{0}$$). This equation is equivalent to a Hookean spring. Of course, when looking at a whole system, this term has to be integrated over the total area:

$E_{a}=\frac{1}{2} \frac{K_{a} \Delta a^{2}}{a_{0}}$

If the area change is constant over the surface of the membrane, then the integral reduces to Equation (1).

The area modulus has values of about 55-70 kBT or 230-290 mN/m.

Figure 3: The area change of a membrane can be modeled by treating the membrane as a network of springs. Adapted from [1].

Thickness Changes

Membrane compression can simply be described by a function $$w(x,y)$$ that outputs the thickness at point $$(x,y)$$. Again, we can model the thickness change as a system of springs (Figure 4). If the half-width of the membrane (essentially the lipid length in the membrane at equilibrium) is defined as $$w_0$$, and the thickness is $$2w_0$$ (Figure 4), then we get a similar equation for the energy penalty of changing the thickness:

$E_{\text {thickness}}=\frac{K_{t}}{2} \int\left(\frac{\left(w(x, y)-w_{0}\right)^{2}}{w_{0}}\right) d a$

Figure 4: The thickness change of a membrane can be modeled by treating the membrane as a network of springs. Adapted from [1].

Shear Changes:

Pure phospholipid membranes are essentially two-dimensional fluids. Due to the fluid nature of phospholipid bilayer, the membrane cannot support a shear strain. Therefore, a fluid membrane does not have a shear modulus. Therefore, static shear does not deform the membrane and does not cost any energy. However, many biological systems are not pure fluid membranes. For example, red blood cells has an underlying cytoskeletal network, which can support shear. Therefore, to model the equilibrium shapes and the dynamics of such cells, the shear must be accounted for. These details will not be described here, but an example of this is described below for the case of red blood cells.

An Interlude on Membrane Geometry

To understand how the energy costs of bending is determined, it is important to describe the geometry of the membrane.

If we take a reference two-dimensional plane, we can define a function $$h(x,y)$$ that characterizes the height of the membrane above the reference plane at point $$(x,y)$$. In most situations, this can accurately represent the shape of the membrane. This representation is known as the Monge parameterization (Figure 5). From the Monge parameterization we can now calculate the curvature of a membrane. This is important as the free energy cost of bending the membrane is dependent on the curvature.

Figure 5: Monge representation. Adapted from [1].

In a single dimension, if we have a smooth function $$y=f(x)$$, then at a point $$x$$ we can define the curvature $$\kappa (x)$$ to be the inverse of the radius of a circle that closely approximates the curve at point $$x$$ (Figure 6A):

$\kappa=\frac{1}{R}$

It might be expected as straight line as a curvature of 0. Indeed, since the circle that approximates the straight line is that of an infinite radius, the curvature is $$\kappa=\frac{1}{\infty}=0$$ (Figure 6).

Figure 5: A) The curvature is defined as the inverse of the radius of the circle that best approximates the curve at a given point. B) For a surface, the curvature is defined by two values, which are the curvatures of the curve when the surface is cut by orthogonal planes.

Figure 6: A line has an infinite radius of curvature and therefore a curvature of zero

To generalize this to the two-dimensional case, we could cut through the two-dimensional surface with a plane, obtaining a curve, from which we can calculate the curvature. However, depending on the orientation of the plane used to cut the surface, the value of the curvature. It can be shown that the curvatures in two directions describe a single point on a surface. Specifically, there is one particular choice of orthogonal planes that the curvatures take extreme values. These values $$c_1 , c_2$$ are known as the principal curvatures (Figure 5B). In fact, one can calculate the curvature in any other direction using the formula:

$c_{X}=c_{1} \cos ^{2} \alpha_{X}+c_{2} \sin ^{2} \alpha_{X}$

where $$c_{x}$$ is the curvature in the direction X, and $$\alpha_{X}$$ is the angle between the direction X and the direction belonging to curvature $$c_{1}$$. These statements are known as Euler’s theorem.

As an example of principal curvatures, the principal curvatures of a cylinder would be $$R$$ around the circumference of the cylinder, but $$0$$ along the length of the cylinder.

Other curvature quantities calculated from the principal curvatures are also commonly defined:

$\begin{array}{c}{K=c_{1}+c_{2}} \\ {K_{G}=c_{1} c_{2}} \\ {H=\frac{K}{2}=\frac{c_{1}+c_{2}}{2}}\end{array}$

where $$K$$ is the is the extrinsic curvature, $$G$$ is the Gaussian curvature, and $$H$$ is the mean curvature.

With the Monge parameterization, the curvatures can be expressed in terms of the height function. It can be shown that the extrinsic curvature is:

$K=\nabla \cdot\left(\frac{\nabla h(\boldsymbol{r})}{\sqrt{1+(\nabla h(\boldsymbol{r}))^{2}}}\right)$

For a small gradient, the extrinsic curvature can be approximated as $$K \approx \Delta h(\boldsymbol{r})$$.

Why is the Monge parameterization and expressing the curvature using the parameterization important? These are important as the bending energy is dependent on the         curvature, and to determine the equilibrium shapes of a membrane system, the height function that minimize the total free energy of the system will be determined to represent the final membrane shape.

This is only an introduction membrane/surface geometry, which itself is only a small part of the field of differential geometry. For more information on differential geometry in relation to membranes, the reader is referred to [2]. For a more general and advanced introduction to differential geometry, the reader is referred to [3].

Bending Changes – Helfrich theory:

The bending energy of a membrane is described by the well-known Helfrich theory. The Helfrich bending energy is defined in terms of the extrinsic and Gaussian curvature:

$E_{\text {bend}}=\int\left(\frac{1}{2} \kappa_{b}\left(c_{1}+c_{2}-c_{0}\right)^{2}+\overline{\kappa}_{b} c_{1} c_{2}\right) d A=\int\left(\frac{1}{2} \kappa_{b}\left(K-c_{0}\right)+\overline{\kappa}_{b} G\right) d A$

where $$\kappa_b$$ is the bending rigidity and $$\overline{\kappa}_{b}$$ is the Gaussian bending rigidity. This equation also contains the spontaneous curvature $$c_0$$. It is nonzero if there is some curvature in the membrane at an equilibrium state, which occurs is there is an asymmetric bilayer (different lipid composition for each monolayer). Typical values of $$\kappa_b$$ is 10-20 kBT.

It is important to note that that Gaussian bending term can be ignored in most cases due to the Gauss-Bonnet theorem. The Gauss-Bonnet theorem says that if the topology of the system remains unchanged, the integral over the Gaussian curvature will be a constant. The topology (genus) can be determined by the number of holes/handles the surface has. For example, a sphere has no holes, and therefore has a genus of 0, but a donut has one hole and a genus of 1. Therefore, unless the genus of the system changes, the Gaussian term vanishes. This is the case for many biological membranes. However, for complicated morphologies of the endoplasmic reticulum and mitochondria, Gaussian curvature could have a significant contribution as membrane holes, tubes, and tori are created and destroyed.

We can simplify the above equation for Monge parameterization, with no spontaneous curvature and with the Monge parameterization. It can be shown for the Monge parameterization, the area element is $$d A=d x d y \sqrt{1+(\nabla h)^{2}}$$. For the small gradient approximation, this is equivalent to $$d A=\left(1+\frac{1}{2}(\nabla h)^{2}\right) d x d y$$ (based on a simple Taylor expansion). We can also penalize area stretching with a term $$\frac{1}{2} \sum(\nabla h)^{2}$$ where $$\Sigma$$ is membrane tension. Therefore, a Helfrich energy equation for a Monge parameterization, under the small gradient approximation, and includes bending and tension is:

$E_{b e n d}=\frac{1}{2} \int d x d y\left[\kappa_{b}(\Delta h)^{2}+\Sigma(\nabla h)^{2}\right]$

Vesicles and Free Energy Cost

Vesicles are spherical lipid bilayer structures. Vesicles are used for many cellular processes for transporting large amounts of cargo within cells or even between cells. An example of the latter is the transport of neurotransmitters between neurons at the synaptic cleft.

We can calculate the free energy cost of the vesicles. The principal curvatures are the same: $$1/R$$ Therefore, the free energy cost is:

$E_{v e s i c l e}=\frac{\kappa_{b}}{2} \int\left(\frac{2}{R}\right)^{2} d A=\frac{\kappa_{b}}{2}\left(\frac{2}{R}\right)^{2} 4 \pi R^{2}=8 \pi \kappa_{b}$

Interestingly, this means there is a fixed energy cost to create a vesicle regardless of the size of the vesicle. If $$\kappa_{b} \approx 10 k_{B} T$$ then $$E_{\text {vesicle}} \approx 250 \mathrm{k}_{B} T$$, which means there must be an active mechanism for the formation of vesicles in cells. Indeed vesicle-coating proteins like clarithin decrease this energy barrier such that vesicle formation is favorable.

Applications:

Membrane tethers:

When a pulling force is applied to a single point on a vesicle, it will generate an extruded membrane tubule known as membrane tethers. Membrane tethers have many biological implications. For example neutrophils can use membrane tethers to attach to the endothelial wall when rolling. Membrane tethers also have an important role in the endomembrane system. In addition, the energetics and dynamics of membrane tethers is similar to that of other biological processes, such as endocytosis, and actin projections during cell motility.

Experiments use optical tweezers or molecular motors to apply a force at a single point on the membrane. We will use a simplified model (Figure 8) of the process where we have a sphere of radius $$R$$ attached to a cylinder of radius $$r$$ and length $$L$$. At the end of the cylinder is a hemispherical cap.

Figure 8: A simplified model of membrane tethers. Adapted from [1].

We can write the free energy of the vesicle+tether system by calculating the contributions of bending,area stretching, and compression. The bending energy of the system is:

$\begin{array}{c}{E_{\text {bending}}=8 \pi \kappa_{b}+\pi \kappa_{b} \frac{L}{r}+4 \pi \kappa_{b}} \\ {=12 \pi \kappa_{b}+\pi \kappa_{b} \frac{L}{r}}\end{array}$

The first term corresponds to the spherical vesicle (as demonstrated previously). The second term corresponds to the cylinder, which is calculated by squaring the sum of the principal curvatures $$\left(0+\frac{1}{r}\right)^{2}$$and multiplying by the area. The last term corresponds to the hemispherical cap.

Area stretching also contributes to the free energy of the system:

$E_{\text {area}}=\frac{K_{a}}{2} \frac{\left(a-a_{0}\right)^{2}}{a_{0}}$

where $$a_0$$ is the reference area and $$a=4 \pi R^{2}+2 \pi r L$$ which is the area of the system. Note that the area of the hemispherical cap $$4 \pi r^2$$ is ignored as it is assumed $$r \ll R, L$$.

Note that this system has no thickness change of the membrane and therefore does no contribute to the free energy of the system.

In the case of a closed bilayer like a vesicle, the pressure difference across the membrane $$\Delta p$$ also contributes to the energy of the system:

$E_{p V}=-\Delta p \int d V=-\Delta p\left(\frac{4}{3} \pi R^{3}+r^{2} \pi L\right)$

Again, the contribution of the hemispherical cap is ignored.

Finally, the force pulling the membrane must also contribute to the free energy of the system:

$E_{l o a d}=-f L$

where $$f$$ is the force applied.

Putting this all together, we get:

$E_{t o t}=12 \pi \kappa_{b}+\pi \kappa_{b} \frac{L}{r}+\frac{K_{a}}{2} \frac{\left(a-a_{0}\right)^{2}}{a_{0}}-\Delta p\left(\frac{4}{3} \pi R^{3}+r^{2} \pi L\right)-f L$

At equilibrium, the total energy is minimized. Therefore, the equilibrium shape can be determined by minimizing the free energy.

In this case, we can take the derivatives of the equations with respect to each of the variables to minimize the energy.

Note that

$\frac{\partial E_{\text {area}}}{\partial R}=K_{a} \frac{a-a_{0}}{a_{0}} \frac{\partial a}{\partial R}=\tau 8 \pi R$

where $$\tau$$ is the membrane tension which can be related to area stretching by the formula $$\tau = \frac{K_{a}\left(a-a_{0}\right)}{a_{0}}$$. Similarly, it can be shown that  $$\partial E_{\text {area}} / \partial r=\tau 2 \pi L$$ and $$\partial E_{\text {area}} / \partial L=\tau 2 \pi r$$.  Therefore:

$\begin{array}{l}{\frac{\partial E_{t o t}}{\partial R}=0 \rightarrow 8 \pi \tau R-4 \pi \Delta p R^{2}=0 \rightarrow \Delta p=\frac{2 \tau}{R}} \\ {\frac{\partial E_{t o t}}{\partial r}=0 \rightarrow-\pi \kappa_{b} \frac{L}{r^{2}}+2 \pi \tau L-2 \pi \Delta p r L=0} \\ {\frac{\partial E_{t o t}}{\partial L}=0 \rightarrow \pi \kappa_{b} \frac{1}{r}+2 \pi \tau r-\pi \Delta p r^{2}-f=0}\end{array}$

The first minimization leads to an equation known as the Laplace-Young relation, which is a well-known equation that relates the pressure difference $$\Delta p$$ across a fluid interface (here the membrane) and the surface tension $$\tau$$. A more general form takes into account both radii of curvature:

$\Delta p=\tau\left(\frac{1}{R_{1}}+\frac{1}{R_{2}}\right)$

We can plug in the Laplace-Young relation into the other minimized equations:

$\begin{array}{l}{-\pi \kappa_{b} \frac{L}{r^{2}}+2 \pi \tau L-\frac{4 \pi \tau r L}{R}=0} \\ {\pi \kappa_{b} \frac{1}{r}+2 \pi \tau r-f-\frac{2 \pi \tau r^{2}}{R}=0}\end{array}$

Since $$r \ll R$$ the last terms can be ignored:

$2 \pi \tau L=\pi \kappa_{b} \frac{L}{r^{2}}$

$r=\sqrt{\frac{\kappa_{b}}{2 \tau}}$

$\pi \kappa_{b} \frac{1}{r}+2 \pi \tau r-f=0 \rightarrow \pi \kappa_{b} \sqrt{\frac{2 \tau}{\kappa_{b}}}+2 \pi \tau \sqrt{\frac{\kappa_{b}}{2 \tau}}=f$

$f=2 \pi \sqrt{2 K_{b} \tau}$

We now get a relation between the radius of tether, bending energy, and membrane tension. In addition, we get a relation between the force applied, bending energy and tension.

General framework of determining equilibrium shapes (variational calculus):

The previous section looked at how relationships between meaningful physical quantities can be derived if the shape of the system is known. This was done by calculating the free energy cost for the membrane shape and minimizing it with respect to the variables of the system.  However, a different approach can also be taken. What if the physical quantities of the system are known, but the equilibrium shapes need to be determined? This can be done by minimizing an energy functional of the system, which is the basis of variational calculus.

We can define a functional as a function that assigns a scalar value for a given function. In this case, we can think of the free energy of the system $$E_{tot}$$ as a functional that is dependent on the height function $$h(x,y)$$. Variational calculus provides tools to find functions that minimize the functionals. This is exactly what we need to determine the equilibrium shapes! Variational calculus can be used to find the shape function $$h(x,y)$$ that minimizes the energy $$E_{tot}$$

In variational calculus, the minimization of the functional:

$J[y]=\int L\left(x, y, y^{\prime}\right) d x$

is solved by the solution of the differential equation:

$\frac{\partial L}{\partial f}-\frac{d}{d x} \frac{\partial L}{\partial f^{\prime}}=0$

where $$y=f(x)$$ minimizes the functional. This equation is known as the Euler-Lagrange equation. Indeed, the equation for $$E_{tot}$$ has this form.

As a side note, many equations in physics has this form and the Euler-Lagrange equation is used. In fact the dynamics of a physical system can be represented as a functional known as the Lagrangian, and solving the Euler-Lagrange equations for this Lagrangian lead to the equations of motion. This approach to classical mechanics is known as Lagrangian mechanics.

Going back to membranes, if we use Equation (9), we can get the following Euler-Lagrange equation:

$\kappa \Delta \Delta h-\Sigma \Delta h=0$

This partial differential equation can be solved for the height function that describes the shape of the membrane at equilibrium. This is called the shape equation.

Minimization of red blood cells:

While solving the shape equations are non-trivial and won’t be shown here, some results in the case of red blood cells will instead be shown.

It is important to note that there are often two additional energy terms added to the Helfrich bending energy to obtain the total free energy for red blood cells.

When an initially flat membrane is bent, the outer monolayer is stretched, and the inner monolayer is compressed and the number of lipids per unit area is higher for the inner surface and lower for the outer surface. The energy cost of this is added to the model as the following model:

$E_{A D E}=\frac{\kappa_{A D E}}{2} \frac{\pi}{A D^{2}}\left(\Delta A-\Delta A_{0}\right)^{2}$

where $$\kappa_{ADE}$$ is the associated bending energy, $$A$$ is the membrane area, and $$D$$ is the membrane thickness. This is known as the area-difference-elasticity (ADE) model.

$E_{e}=\frac{K_{A}}{2} \oint \alpha^{2} d A+\mu \oint \beta d A$

When the total free energy $$E_{b e n d}+E_{A D E}+E_{e}$$ is minimized, the equilibrium shapes of the red blood cells can be determined. The below figure shows a comparison of the theoretically determined shapes, and those that have been found experimentally.

Figure 9: A comparison of experimentally observed red blood cell shapes (left side) and theoretically determined equilibrium shapes (right side). Adapted from [6]

References:

1. Phillips, R.; Kondev, J.; Theriot, J.; Garcia, H. G. Physical Biology of the Cell2013, Garland Science, Taylor & Francis Group, LLC: Abingdon, UK.
2. Boal, D. Mechanics of the Cell. Cambridge, UK: Cambridge University Press. 2002
3. Deserno, Markus. Fluid Lipid Membranes – A Primer. Max-Planck-Institut fu ̈r Polymerforschung. Mainz, Germany. 2006
4. Deserno, Markus. Notes on Differential Geometry. Max-Planck-Institut fu ̈r Polymerforschung. Mainz, Germany. 2004
5. Li, Xuejin, et al. Continuum- and Particle-Based Modeling of Shapes and Dynamics of Red Blood Cells in Health and Disease. Soft Matter, vol. 9, no. 1, 2013, pp. 28–37., doi:10.1039/c2sm26891d.

6. W., G. Lim H., et al. “Stomatocyte-Discocyte-Echinocyte Sequence of the Human Red Blood Cell: Evidence for the Bilayer- Couple Hypothesis from Membrane Mechanics.” Proceedings of the National Academy of Sciences, vol. 99, no. 26, 2002, pp. 16766–16769., doi:10.1073/pnas.202617299.