Article Information

• PDF (610 kb)
• Full Text with Large Figures
• Cited by in Scopus (4)

PubMed

• Articles by Qiang Zhu
• Articles by L. Amy Sung

Related Articles

• Spectrin Folding versus Unfolding Reactions and RBC Membrane...

  Spectrin Folding versus Unfolding Reactions and RBC Membrane Stiffness Biophysical Journal, Volume 94, Issue 7, 1 April 2008, Pages 2529-2545 Qiang Zhu and Robert J. Asaro Abstract Spectrin (Sp), a key component of the erythrocyte membrane, is routinely stretched to near its fully folded contour length during cell deformations. Such dynamic loading may induce domain unfolding as suggested by recent experiments. Herein we develop a model to describe the folding/unfolding of spectrin during equilibrium or nonequilibrium extensions. In both cases, our model indicates that there exists a critical extension beyond which unfolding occurs. We further deploy this model, together with a three-dimensional model of the junctional complex in the erythrocyte membrane, to explore the effect of Sp unfolding on the membrane's mechanical properties, and on the thermal fluctuation of membrane-attached beads. At large deformations our results show a distinctive strain-induced unstiffening behavior, manifested in the slow decrease of the shear modulus, and accompanied by an increase in bead fluctuation. Bead fluctuation is also found to be influenced by mode switching, a phenomenon predicted by our three-dimensional model. The amount of stiffness reduction, however, is modest compared with that reported in experiments. A possible explanation for the discrepancy is the occurrence of spectrin head-to-head disassociation which is also included within our modeling framework and used to analyze bead motion as observed via experiment. Abstract | Full Text | PDF (346 kb)

• Spectrin-Level Modeling of the Cytoskeleton and Optical Twee...

  Spectrin-Level Modeling of the Cytoskeleton and Optical Tweezers Stretching of the Erythrocyte Biophysical Journal, Volume 88, Issue 5, 1 May 2005, Pages 3707-3719 J. Li, M. Dao, C.T. Lim and S. Suresh Abstract We present a three-dimensional computational study of whole-cell equilibrium shape and deformation of human red blood cell (RBC) using spectrin-level energetics. Random network models consisting of degree-2, 3, …, 9 junction complexes and spectrin links are used to populate spherical and biconcave surfaces and intermediate shapes, and coarse-grained molecular dynamics simulations are then performed with spectrin connectivities fixed. A sphere is first filled with cytosol and gradually deflated while preserving its total surface area, until cytosol volume consistent with the real RBC is reached. The equilibrium shape is determined through energy minimization by assuming that the spectrin tetramer links satisfy the worm-like chain free-energy model. Subsequently, direct stretching by optical tweezers of the initial equilibrium shape is simulated to extract the variation of axial and transverse diameters with the stretch force. At persistence length =7.5nm for the spectrin tetramer molecule and corresponding in-plane shear modulus ≈8.3N/m, our models show reasonable agreement with recent experimental measurements on the large deformation of RBC with optical tweezers. We find that the choice of the reference state used for the in-plane elastic energy is critical for determining the equilibrium shape. If a position-independent material reference state such as a full sphere is used in defining the in-plane energy, then the bending modulus needs to be at least a decade larger than the widely accepted value of 2×10J to stabilize the biconcave shape against the cup shape. We demonstrate through detailed computations that this paradox can be avoided by invoking the physical hypothesis that the spectrin network undergoes constant remodeling to always relax the in-plane shear elastic energy to zero at any macroscopic shape, at some slow characteristic timescale. We have devised and implemented a liquefied network structure evolution algorithm that relaxes shear stress everywhere in the network and generates cytoskeleton structures that mimic experimental observations. Abstract | Full Text | PDF (907 kb)

• Ribosome Recycling, Diffusion, and mRNA Loop Formation in Tr...

  Ribosome Recycling, Diffusion, and mRNA Loop Formation in Translational Regulation Biophysical Journal, Volume 85, Issue 2, 1 August 2003, Pages 755-773 Tom Chou Abstract We explore and quantify the physical and biochemical mechanisms that may be relevant in the regulation of translation. After elongation and detachment from the 3′ termination site of mRNA, parts of the ribosome machinery can diffuse back to the initiation site, especially if it is held nearby, enhancing overall translation rates. The elongation steps of the mRNA-bound ribosomes are modeled using exact and asymptotic results of the totally asymmetric exclusion process. Since the ribosome injection rates of the totally asymmetric exclusion process depend on the local concentrations at the initiation site, a source of ribosomes emanating from the termination end can feed back to the initiation site, leading to a self-consistent set of equations for the steady-state ribosome throughput. Additional mRNA binding factors can also promote loop formation, or cyclization, bringing the initiation and termination sites into close proximity. The probability distribution of the distance between the initiation and termination sites is described using simple noninteracting polymer models. We find that the initiation, or initial ribosome adsorption binding required for maximal throughput, can vary dramatically depending on certain values of the bulk ribosome concentration and diffusion constant. If cooperative interactions among the loop-promoting proteins and the initiation/termination sites are considered, the throughput can be further regulated in a nonmonotonic manner. Experiments that can potentially test the hypothesized physical mechanisms are discussed. Abstract | Full Text | PDF (496 kb)

• …more

Copyright © 2007 The Biophysical Society. All rights reserved.
Biophysical Journal, Volume 93, Issue 2, 386-400, 15 July 2007

doi:10.1529/biophysj.106.094383

Biophysical Theory and Modeling

A Hybrid Model for Erythrocyte Membrane: A Single Unit of Protein Network Coupled with Lipid Bilayer

Qiang Zhu^*, Carlos Vera^†, Robert J. Asaro^*, Paul Sche^† and L. Amy Sung^†, ,  

^* Department of Structural Engineering, University of California, San Diego, La Jolla, California
^† Department of Bioengineering, University of California, San Diego, La Jolla, California

Address reprint requests to L. A. Sung, Tel.: 858-534-5250.

Although the mean size of erythrocytes (red blood cells, or RBC) in their natural state is ∼8μm, they routinely flow in deformed states through capillaries whose minimum lumen measures 5–7μm 1 or less and thereafter recover their biconcave shape in low stress states. Essential to this remarkable deformability is a composite membrane consisting of a lipid bilayer coupled to an actin/spectrin-based protein network, the membrane skeleton, which imparts flexibility and integrity to the overall structure. The interest in the mechanical properties of the erythrocyte membrane stems from multiple sources. For example, diseases or disorders such as hereditary spherocytosis, ovalocytosis, and elliptocytosis 2,3 have been shown to be associated with structural alterations in the erythrocyte membrane and the resulting changes in mechanical response 2,3,4,5. Another reason for the intense interest is that the erythrocyte membrane has been well characterized, so that it has served as a model system for developing a framework for analyzing the mechanical behavior of a broad range of biological membrane structures 5,6,7. In addition, useful properties of the composite erythrocyte membrane (strength, lightness, deployability, flexible connectivity, etc.) suggest potential synthetic applications in biomimetics.

Unlike the complicated three-dimensional cytoskeletons in other cells, a mature human erythrocyte has only a thin (yet three-dimensional) skeletal protein network beneath the lipid bilayer (Fig. 1). The network is composed of several major proteins: α- and β-spectrins, ankyrin, band 3, protein 4.1, protein 4.2, and actin, as well as some minor proteins such as myosin, tropomyosin, and tropomodulin (E-Tmod). Structurally, the network is organized into ∼33,000 repeating units. As revealed by transmission electron microscopy 8,9,10, each basic repeating unit is viewed as a small spoked but edge-free hexagon (Fig. 2) formed by a junctional complex (JC), with six long αβ spectrin dimers (Sp) radiating from a central short actin protofilament, and up to six suspension complexes (SC) 11. The radius of a spoked basic unit depends on the degree of extension of the six Sp. These repeating units connect with each other through head-to-head association of Sp from two neighboring units 12. It is important to note that not every Sp dimer is connected to a SC; however, most spectrin tetramers are associated with one (or two) SC. Furthermore, there are ∼3% pentagons and 8% heptagons in the network, in addition to 85% hexagons 10. The SC functions as the primary connection between the protein network and the lipid bilayer. A SC consists mainly of band 3, ankyrin, and protein 4.2 13. While band 3 is a major transmembrane protein (with multiple transmembrane domains) and an anion exchanger, protein 4.2 is myristylated (covalently associated with a completely saturated C14 fatty acid, and therefore tightly associated with the lipid bilayer) 14. The secondary site to link the protein network to the lipid bilayer is at JC, where protein 4.1 is associated with glycophorin C, another transmembrane protein with a single transmembrane domain 15,16.

Display large version of this figure Display high quality version of this figure

Figure 1

Schematic of the erythrocyte membrane showing the connectivity among the major proteins and with the lipid bilayer (modified from 51). SC, suspension complex. The inset illustrates two possible types of junctions for the Sp/actin attachment, i.e., point-attachment and wrap-around. Small closed circles indicate β-spectrin binding sites on G-actin.

Display large version of this figure Display high quality version of this figure

Figure 2

(a) Topology of a spoked hexagon (JC) in the erythrocyte membrane skeleton—a three-dimensional view showing a protofilament with six pairs of G actin associated with six Sps. Each Sp may connect to a SC in the lipid bilayer, forming a small hexagon without physical edges. The numbers in the SC indicate the order of G actin pairs, Sp, and SC. The angle between the actin protofilament and the lipid bilayer is not specified and molecules are not to scale. (b) A top view of a large hexagon, formed by six basic hexagonal units surrounding a center unit. The head-to-head interaction of Sp forms a spectrin tetramer between adjacent units. (Both a and b are from 11.)

Because the red cell membrane skeleton has a remarkably simple regular geometry compared to other cell types it has therefore motivated dozens of molecularly inspired modeling efforts, including our three-dimensional network model for a single unit 17. Other recent attempts to describe the mechanical response of the erythrocyte membrane include those of Boey et al. 18, Discher et al. 19, and Li et al. 20. These models are well suited for carrying out analysis of entire cell deformation but do not include hydrodynamic effects and are based on an idealized planar picture for the membrane. This makes it difficult to incorporate detailed molecular structure and structural changes that result from, for example, disease or mutation. The model by Hansen et al. 21 was concerned with describing the behavior of planar triangulated meshes of Sp segments modeled as a linear elastic spring. The lipid bilayer was not included, but their modeling did reveal some important effects of network connectivity and topology on the skeleton response and the predicted constitutive properties. The modeling of the lipid bilayer was achieved in Lin and Brown 22, using a Fourier space Brownian dynamics (FSBD) method introduced by Granek 23 to study the mobility of band 3 inside the lipid bilayer. This model did not consider the three-dimensionality and mobility of the protein skeleton and was not intended to probe membrane mechanical properties. It does, however, provide a viable framework for analyzing the dynamics of the lipid bilayer that we use herein.

A complete model for the membrane, i.e., one based on sufficient molecular detail and which would allow detailed correlations between the mechanical response and alterations in the skeletal network, the lipid bilayer, and/or their interconnectivity, would be of considerable value in understanding the dynamic properties of the erythrocyte membrane in normal and diseased states. For example, in disease states such as Southeast Asia ovalocytosis, mutations generate nearly rigid arrays of band 3 cross-linking Sp segments, changing the interconnectivity of the skeleton as well as the skeleton/bilayer connectivity, and leading to dysfunctional mechanical response of the cell 3. It is, accordingly, the intention of this work to develop a molecular-based model of the complete erythrocyte membrane (protein skeleton/lipid bilayer) that may be used to analyze the deformation behavior of this composite structure.

We propose a hybrid model wherein the three-dimensional dynamic response of a JC and its interaction with the lipid bilayer are accounted for. This model is based on a static model as reported in Vera et al. 17, which features a detailed three-dimensional description of a single JC unit as described by Sung and Vera 11. Numerical results based on this early model revealed that the protofilaments may be more tangent than normal to the lipid bilayer, in agreement with experimental observations made with the micropipette aspiration technique combined with fluorescent polarization microscopy 24,25. Our previous model, however, simplified the lipid bilayer by modeling it as a flat plate connected to the JC via a linear spring, and has not considered important effects such as hydrodynamic drag and thermal fluctuations. In addition, this model depicted the Sp-protofilament junction as a simple attachment point (ball-joint), even though the α- and β-spectrin may wrap around a protofilament 11. The proposed wrap-around junction is based on a three-dimensional model 11 and biochemical interactions between α- and β-spectrin 26,27,28,29. The wrap-around junction, adopted herein, features the extension of α- and/or β-spectrin around the protofilament, providing a mechanism to minimize the torque applied to the protofilament during deformation, affecting its roll and possibly yaw motions. It is also understood that the influence of lateral association on forced unfolding of antiparallel spectrin heterodimers plays an important role in the biomechanics of spectrin and the network 30, and this feature is to be incorporated into our models in the future. Currently, we consider the dynamic response of a fully coupled single unit of JC with the lipid bilayer. Careful attention is paid to describing the dynamic properties of the JC, the junction between Sp and the protofilament, and the interaction between the JC and the lipid bilayer. To achieve this, we model the protofilament as a rigid body and solve its motion by employing the Langevin equation, while the tethering force in Sp is obtained using the wormlike-chain (WLC) paradigm 31. A stochastic term representing thermal effects is incorporated. Domain unfolding in the Sp, as discussed by Rief et al. 32 and Lee and Discher 33 is not included in our present Sp model but can be incorporated as discussed below. The WLC model has been used in most other recent studies (e.g., 18,19,20) and thus is chosen here to provide linkage to these studies. The JC is connected to the lipid bilayer through the SCs. In the present study, because this model deals with one single unit of the network, the lateral mobility of these SCs is prohibited, while they are allowed to follow the heaving motion of the bilayer. Two types of Sp-protofilament junction, one corresponding to the point-attachment connection used in the previous study 17 and the other incorporating the torque-minimizing effect of the wrap-around connection are considered and compared with each other. Finally, we model the bilayer via the above mentioned FSBD method. This model is unique in its molecular detail, in its accounting for the vital three-dimensionality of the structure, in its inherently time-dependent response, and in its capacity of including fluid effects.

The rest of the article is organized as follows. In the next section we describe the physical model for the JC. The mathematical formulation and numerical algorithms are then detailed for the lipid bilayer and its connectivity to the skeleton. Numerical results, including both quasi-static response of the JC-bilayer structure and effects of thermal fluctuations, are provided. Preliminary calculations of important constitutive parameters such as the shear modulus are made and described. Finally, we summarize the primary findings. As noted, the present work is concerned with the development and implementation of a basic unit of JC/lipid bilayer model. Extensions of the work will include simulations of multiple JC units for the purpose of modeling the effects of topology and connectivity on mechanical response, including mutations of the skeleton/lipid bilayer structure on dynamic and quasi-static response.

Following our first report on the 3-D nanomechanics of an erythrocyte JC in deformation 17, we use Sung and Vera’s 3-D model for the JC 11. An actin protofilament contains 12–13 G-actin arranged in a helical fashion with a dihedral angle (angle between subsequent G actin) of 166.2° and with an axial distance of 5.5nm. The protofilament is reinforced by two dimeric tropomyosin, one in each of the two grooves, and is capped by an E-Tmod at one end (the pointed end). The positions of β-spectrin binding sites are based on the regularity of G-actin arranged in F-actin, and the predicted occupation of six out of twelve β-spectrin binding sites on the protofilament 11. As shown in Figure 2a, in a JC, the protofilament functions as a mechanical axis to anchor three (top, middle, and bottom) pairs of Sp. Each Sp pair (i.e., Sp₁/Sp₂, Sp₃/Sp₄, or Sp₅/Sp₆) is arranged in a back-to-back fashion and connects with the protofilament. The Sp pairs spiral down the protofilament in a right-handed manner from the pointed end with a dihedral angle of 55.4° between them. The detailed coordinates of the β-spectrin binding sites are listed in Table 1.

Our present model depicts a single JC consisting of an actin protofilament (modeled as a rigid circular cylinder with length b and radius a) and six Sp, coupled to the lipid bilayer. Since the length of the actin protofilament is ∼100 times smaller than the persistence length of the actin filament 34, it is reasonable to consider it as a rigid cylinder. The JC is connected to the lipid bilayer through SC. In this model, the actual volume and geometry of the SC are ignored, and thus the Sp is assumed to be connected to the bilayer at point junctions. Note that since we do not allow band 3 motion in the present simulations to move within the lipid bilayer, the Sp are pretensioned as described below.

For convenience, the problem is formulated using a Euler-Lagrangian dual-coordinate system (see Fig. 1), including a global space-fixed Euler coordinate system (X, Y, Z) and a local actin-fixed Lagrangian reference frame (x, y, z). The Euler coordinate system is defined so that the undisturbed lipid bilayer lies within the X-Y plane, and Z points toward the JC. The origin of the Lagrangian coordinate system (x, y, z) is at the center of the protofilament, x is tangential to the protofilament aligned toward the pointed end capped by E-Tmod, and both y and z are perpendicular to x. In this dual-coordinate system, an arbitrary vector g is represented either as (g_X, g_Y, g_Z) in the Euler system, or as (g_x, g_y, g_z) in the Lagrange system. The transformation between these two representations is

(1)

(2)

We note that the Euler parameters are not independent due to the relation . The rotation matrix, C, thus becomes 36

(3)

One of the structural characteristics determining the mechanical response of the JC is the exact manner in which a spectrin dimer is connected to the protofilament. In the previous numerical model of a JC 17, this important linkage was simplified as a point-attachment, i.e., an Sp was connected to the protofilament at a single point corresponding to the proposed β-spectrin binding site on the G actin 11 (see the inset in Fig. 1). According to this model, the forces applied through Sp will directly affect this binding site. The unfolding of domains would be the primary protective measure, relieving an overly large tensile force on the binding site. Corresponding to this model, the six Sps are pinned at the outer surface of the circular cylinder representing the protofilament.

Compared with the point-attachment junction, the proposed wrap-around junction 11 possesses many advantages in terms of structural stability. As displayed in the inset in Fig. 1, the α- and β-spectrins of an Sp are seen to split and wrap around the protofilament and reconnect at the back side. Thus, the tension inside Sp does not impose directly on the binding site; instead, it is distributed around the protofilament. The binding between Sp and protofilament is reversible in the event that dissociation occurs, since Sp is still wrapped around the protofilament so the configurational entropy increase of a released Sp is small. Furthermore, during deformation the α- and β-spectrin wrapping around the protofilament may stretch so that the roll moment applied to the protofilament is reduced, diminishing the coupling effect between the Sp and the protofilament in the roll mode. Although at this time it is difficult to estimate the exact magnitude of the reduction in coupling, it is reasonable to study an extreme case of complete decoupling in the roll mode. We mathematically simulated this condition by pinning the six Sps along the central axis of the protofilament so that the tension inside an Sp generates no roll moment on the protofilament. The exact locations of these pinning points, as well as the locations of the points where the Sps are attached to the lipid bilayer (corresponding to the locations of the SC), are listed in Table 1,Table 2 (see 11 for additional detail regarding the connection locations). In the following, we address the mathematical formulation and numerical issues in this coupled model.

In the following, we summarize the key steps involved in using the FSBD method 23 to model the dynamics of the lipid bilayer. Additional detail on the FSBD method is found in Lin and Brown 22. Considering a bilayer with uniform tension σ, the layer’s elevation, h(X, t) where X=(X, Y) is the horizontal position vector, is determined via the nonlocal Langevin equation. This takes the form

(4)

In our present analysis of a single JC unit, the interaction force (acting in the vertical Z direction) is composed of three parts, i.e.,

(5)

Let (i=1, 6) be the locations of the SCs (Table 1). Then, following Lin and Brown 22, we set

(6)

(7)

The repulsive force imposed by the Sp on the bilayer is evaluated in a similar manner. Consider the i^th Sp (i=1, …6); its projection, P_i, on the mean plane of the bilayer (Z=0), is expressed as a straight line, A_iX+B_iY+C_i=0, bounded by the projected points of the SC and the Sp-protofilament attachment point. From an arbitrary point on the bilayer, X, we draw a line that is perpendicular to P_i. The two lines intersect at X_Pi. The Z-coordinate of the Sp right above X_Pi is Z_Pi. The repulsive force is then assumed to be zero if X_Pi lies outside of the boundaries of P_i. Otherwise, this force is assumed to be a nonzero repulsive force that decays with the distance between X and X_Pi, i.e., and increases exponentially with −(Z_Pi-h(X, t)). Summing the contribution from all the Sp, we have

(8)

The random thermal fluctuation force, ζ(X, t), is modeled as a Gaussian variable with zero mean, i.e., 〈ζ(X, t)〉=0, and no temporal correlations, whose variance is

(9)

Among the parameters that quantify the interaction forces among the Sp, the bilayer, and the protofilament, some are chosen to be consistent with previous work 22, including c₁=7nm and c₂=0.2nm. Others (γ and η ) are chosen through our numerical tests. We find that when these parameters are not sufficiently large, their respective force expressions will not prevent the protofilament or the Sp from penetrating the bilayer; on the other hand, when their values are too large, extremely small time steps are required for numerical stability. Best results are obtained with γ=10²²Nm⁻⁴ and η=0.1Nm⁻².

We apply an efficient Fourier spectral algorithm to resolve the nonlocal Langevin equation 23. Assuming that the problem is doubly periodic (with period L) in both X and Y directions, we consider a piece of lipid bilayer whose projection on the X-Y plane is a square with length L in both dimensions. L is chosen to be much larger than the scale of the JC so that disturbances with wavelengths larger than L (representing global deformations) will not have a significant effect on our analysis, which is focused on local responses. Our numerical tests demonstrated that the responses of the JC were not sensitive to the value of L. For any doubly-periodic function f(X), the Fourier transform and reverse Fourier transform are, respectively,

(10)

(11)

(12)

When k=(0, 0), (− Nπ/L, 0), (0,−Nπ/L), or (− Nπ/l,−Nπ/L), the Brownian random term R(k, t) is a real Gaussian random variable with zero mean and variance 2L²k_BTΛ(k)Δt; for other k-values, R(k, t) is a complex function whose real and imaginary parts are both Gaussian random variables with zero mean and variance L²k_BTΛ(k)Δt. Once we obtain the updated elevation, the inverse transform is computed via a FFT so that the computational effort is only proportional to N² log N².

Owing to the small scale of the problem, the inertia of the protofilament is negligible. Its motion is thus determined by the balance among the hydrodynamic viscous drag exerted by the surrounding fluid, the Brownian force/moment, f^(ζ), M^(ζ), the force and moment exerted by Sp, f^(s), M^(s), and the repulsive force/moment due to its interaction with the lipid bilayer, f^(b), M^(b). Measured in the Lagrange coordinate system, let the translational velocity of the protofilament be (u_x, u_y, u_z) and its rotational velocity be (ω_x, ω_y, ω_z). Then according to the Langevin equation we have

(13)

The Reynolds number of the flow surrounding the protofilament is usually within the creeping/Stokes regime. Therefore the hydrodynamic effect on the protofilament is modeled by distributing a set of Stokeslets, as demonstrated by Lighthill 39. The normal and tangential drag coefficients per unit length are estimated as, respectively,

(14)

(15)

Following the fluctuation-dissipation theorem, the components of the thermal fluctuation force/moment , and are Gaussian random variables with zero mean. The variances for the forces are 2k_BTD_x/Δt, 2k_BTD_y/Δt, and 2k_BTD_z/Δt and for the moments we have 2k_BTD_Mx/Δt, 2k_BTD_My/Δt, and 2k_BTD_Mz/Δt, respectively. The force caused by interaction of the protofilament with the lipid bilayer, f^(b), is obtained by summing the contributions from the repulsive force acting on each element around its surface. Due to the principle of action and reaction, , in turn, is evaluated by integrating its corresponding part in the interaction force over the bilayer with a sign change. This yields

(16)

With translational and rotational velocities (u_x, u_y, u_z) and (ω_x, ω_y, ω_z) obtained through Eq. (13), the location of the center of the protofilament, at χ=(X₀, Y₀, Z₀) and the center’s Euler parameters are determined, via the procedure outlined in Tjavaras et al. 41, as

(17)

The strain-stress relation of the Sp is modeled using the WLC model. This model depicts a polymer such as a Sp as a continuous elastic cable possessing constant bending stiffness characterized by a persistence length, a measure of the average curvature along the contour 31. By considering thermal equilibrium, the WLC model provides a closed-form description of the stretch-tension relation of the polymer that is easy to implement. However, it ignores the transient behavior of the polymer and the randomness of the tension. To include the randomness, we use a modified form of the WLC and express the tension-stretching relation as

(18)

In their recent model studies, Discher et al. 19 and Li et al. 20 have modeled steric repulsion between the Sps of a JC by adding a term to the free energy of the form C/A, where A represents the area associated with the JC. C is a constant evaluated by setting the net pressure within the skeleton (assumed not attached to a bilayer) to zero once the JC attains a prescribed rest area; this also sets the length of the Sp at which steric repulsion, so modeled, balances the tension within the Sp. The steric contribution contributes to the area modulus but not to the area preserving shear modulus. Since in the present study the SCs are fixed to the bilayer, the bilayer supports pretension in the Sp segments and our skeleton does not collapse. Future modeling will include a detailed account of steric interaction such as we have presented here for steric interactions between the protofilament and the Sp with the bilayer itself. This will be necessary since we allow for SC mobility. For now we note that the effect of inter-Sp repulsion as modeled in the literature 18,19,20 is essentially on the computation of the area modulus, κ_s, of the skeleton itself (i.e., without the dominant effect of the bilayer on area stiffness)—and which can be analytically shown to amount to κ_s≈2μ (μ being the predicted skeleton shear modulus when computed at the rest area of the JC). As noted, the Sp WLC model used here does not include the extension softening that would come about due to either intramolecular transitions such as Sp domain unfolding, or interprotein events such as spectrin head-to-head disassociation. Models for such events are further considered in Conclusions and Discussion.

The elastic energy of the system contains three parts, the stretching energy stored in the bilayer owing to pretension, the elastic energy of the bilayer due to elevation by bending and via pretension (the combination of these is called Φ^(b)), and the strain energy in the Sp, called Φ^(s). In this study, the process that induces pretension is not considered—so that for the total elastic energy, we only consider Φ^(b) and Φ^(s). The elastic energy due to bending and pretension is

(19)

(20)

The total energy stored in the six Sp is called Φ^(s).

With the hybrid numerical algorithm described in the previous section, we conducted numerical simulations to illustrate the static and dynamic response of the JC coupled to the lipid bilayer. In our numerical simulations the size of the computational domain was chosen to be L=400nm. 64×64 Fourier modes were used to represent deformation of the bilayer. For numerical stability in the time integration, a small time step of 0.1ns was employed.

We first analyze the quasi-static response of the coupled JC-bilayer system without considering thermal fluctuations. Starting from an arbitrary initial state, the motion of this system was simulated until a steady state was reached. The steady values of the orientation angles of the protofilament, including the yaw angle and the pitch angle as defined in the inset in Fig. 4 and the deformation of the bilayer, as well as the strain energy Φ^(b) and Φ^(s), are summarized below. Moreover, to examine the dependence of the quasi-static response of the JC with respect to the lipid bilayer pretension (σ), we applied the point-attachment junction and simulated the yaw and pitch angles of the protofilament in the natural state as functions of σ. It was found that the pretension had no significant effect upon the yaw angle of the protofilament. What was found (data not shown) was that at values of σ >10⁻³ Nm the pitch angle decreased from 7° to ∼5°. This suggests a stiffening of the bilayer at sufficiently high levels of σ. Bilayer stiffening is manifested in a reduction in the ratio between energy associated with bending (the first term in Eq. (19)) and that associated with the pretension (the second term). At low levels of σ, the energy due to bending is much larger than that due to the pretension, yet when σ≳4×10⁻⁴Nm⁻¹ the two contributions are equal and reverse their order of magnitude at larger pretension. In the following examples we have chosen σ=10⁻⁴Nm⁻¹, where the effects of bending and pretension are comparable.

As shown in Figure 2a and Figure 3a, we assume that in its natural state the six SCs in a basic unit form a hexagon. A general in-plane stretch of the JC, corresponding to a topology-conserving redeployment of the SCs within the lipid bilayer, can be represented as a stretching along an arbitrary axis within the X-Y plane (denoted as axis 1) by a ratio λ₁, plus a stretching by a ratio λ₂ along an second axis which also lies within the X-Y plane and is normal to axis 1. The angle between axis 1 and the X axis, measured counterclockwise, is defined as the deformation angle τ. An equibiaxial deformation refers to an in-plane stretching in which λ₁=λ₂. As illustrated in Figure 3b, an equibiaxial deformation results in a similar hexagon enlarged by a factor of λ². Therefore, the equibiaxial deformation is τ-independent, and is related to an area change of the membrane skeleton. During the deformation, the lipids may flow into or out of the changing area spanned by the SCs, because the SCs may be dragged through the bilayer with little resistance, as noted above. Although the natural viscosity that this fluidic motion of lipid molecules contributes is not considered in this study, one can imagine that it will have significant effects on the dynamics of the coupled system. Fig. 4 shows the equilibrium orientation of the protofilament during a equibiaxial deformation from λ=1 to λ=3.5. The two Sp-protofilament junctions predict dramatically different behaviors of the protofilament. With the point-attachment junction, the yaw angle varies from 18° at the natural state to ∼30° at λ=3.5. The change in pitch angle is more profound, varying from 6° at λ=1 to 46° at λ=3.5, indicating a bow-up tendency with increasing dilation. With the wrap-around junction, on the other hand, the orientation of the protofilament is much less sensitive to the deformation. Within the whole range that we consider, the yaw angle remains at 30°, and the pitch angle is close to 0°, so that the axis of the protofilament is almost parallel to the bilayer.

Display large version of this figure Display high quality version of this figure

Figure 3

Schematic illustrating the natural state (a), a state corresponding to equibiaxial extension (b), and anisotropic deformations (c and d) as defined by the deformation angle τ. The pointed end of the protofilament is marked with an asterisk. Arrows point to the direction of deformation.

Display large version of this figure Display high quality version of this figure

Figure 4

Dependence of the attitude of the protofilament on the extension ratio (λ) in equibiaxial deformation. The inset defines the yaw and pitch angles. The pointed end of the protofilament is marked with an asterisk.

An anisotropic deformation is an in-plane stretching of the JC in which the area of the hexagon is preserved, i.e., λ₁λ₂=1. As shown in Figure 3cd, an anisotropic deformation depends on the deformation angle, τ, as well as the stretching ratio. In this case, as the area is conserved, we imagine that the lipid bilayer undergoes a viscous shear. This adds a viscous strain rate to the elastic shear deformability to be described. Figure 5ab, show the orientation of the protofilament at τ=0° and 135° and λ=1 to 3.5 obtained by using the point-attachment junction. The most noteworthy phenomenon is the bifurcation in the orientation of the protofilament. At both deformation angles, two equilibrium states are identified for each extension ratio λ provided that it is sufficiently large. For example, we consider the case with deformation angle τ=0°. The threshold λ for bifurcation is ∼1.05. For any λ larger than that value, two sets of yaw/pitch angles (each representing an equilibrium state) are shown. One of these states occurs at yaw angle ∼70°, the other at ∼−70°. The two states have similar pitch angles and are almost identical in terms of the total elastic energy, Φ^(s)+Φ^(b). This bifurcation can be explained by the close symmetry of the JC configuration with respect to the X axis, or an axis that is 60°, 90°, or 120°, from the X axis. The symmetry is not exact because the Sps have different lengths during deformation. With the point-attachment junction, the model also predicts dramatic variation of the yaw angle with respect to the extension ratio. For instance, at τ=0° near λ=1.05 the yaw angle quickly increases from 18° to 70°. This behavior, together with the bifurcation, indicates that the JC configuration described by using the point-attachment bilayer membrane is represented by its elastic moduli, whichjunction is quite unstable during anisotropic deformations. In contrast, using the wrap-around junction the model predicts a much more stable JC. As shown in Figure 5cd, in most cases the attitude varies slowly with the extension ratio. Among all the cases we studied, the only bifurcation occurs at τ=135°, and at a relatively high extension ratio λ=2. The pitch angle remains small.

Display large version of this figure Display high quality version of this figure

Figure 5

Dependence of pitch and yaw angles on stretch during anisotropic deformation with (a) τ=0° and using the point-attachment model, (b) τ=135° and using the point-attachment model, (c) τ=0° and using the wrap-around model, and (d) τ=135° and using the wrap-around model.

The predictions using the wrap-around junction (Figure 4 and Figure 5) are more consistent with experimental observations. For example, according to Picart et al. 24,25, during deformation the protofilaments remain tangent to the bilayer. The point-attachment junction, on the other hand, predicts large pitch angles in equibiaxial deformations (Fig. 4). Clearly, additional attention needs to be paid to the molecular structure of such junctions as the results show that the two, albeit extreme cases in terms of induced axial moment, impart quite different behavior to the JC.

To illustrate the effect of the lipid bilayer on the equilibrium state of the JC, simulations were also performed using a decoupled system without the bilayer. With the wrap-around junctions, the results are close to the predictions via the coupled JC-bilayer simulation mentioned above. In contrast, with the point-attachment junction, the results are different with and without the lipid bilayer. Without the lipid bilayer (and the repulsive force imposed) the attitude of the protofilament is not sensitive to deformations. In both equibiaxial and anisotropic deformations, yaw angles remain at ∼30° and pitch angles at ∼70°.

The disturbance on the bilayer exerted by the JC is highly localized and concentrated in the vicinities of the SCs and the contacting point between the protofilament and the lipid bilayer. The magnitude of this disturbance is significantly affected by the pitch angle of the protofilament. With the point-attachment junction, and, for example, for an equibiaxial extension of λ=3.5, the pitch angle of the protofilament is ∼46°, leading to a maximum elevation of 3nm (occurring near the SCs connected with Sp₁ and Sp₂). In contrast, the wrap-around junction leads to a pitch angle close to 0°, suggesting that the protofilament is almost parallel to the lipid bilayer. The corresponding maximum elevation is only 1nm (occurring right beneath the protofilament). In general, the pitch angles predicted by using the wrap-around junction are much smaller than those predicted by using the point-attachment junction (Fig. 4). These results indicate that by using the wrap-around junction, the disturbance by the protofilament on the bilayer is significantly reduced. This is clearly demonstrated in Fig. 6, which shows that during equibiaxial deformation, a JC with the wrap-around junction generates a much smaller elastic energy inside the bilayer than a JC with the point-attachment junction.

Display large version of this figure Display high quality version of this figure

Figure 6

The elastic energy in the lipid bilayer as a function of the extension ratio (λ) in equibiaxial deformation.

The elastic response of the hybrid protein skeleton/lipid bilayer membrane is represented by its elastic moduli, which are obtained herein via quasi-static deformations. For in-plane deformations, the relevant elastic moduli are the area and shear moduli. Since the lipid bilayer possesses very limited compressibility, it dominates the area modulus of the composite structure. The protein skeleton, on the other hand, contributes significantly to the shear modulus of the structure and thus is our focus here. Following Evans and Skalak 6, we calculate the shear modulus of the JC undergoing finite deformations using a method based on the strain energy, Φ^(s), stored in the Sp. In this approach, an arbitrary deformation is characterized by two independent deformation parameters α=λ₁λ₂−1 and β=()/(2λ₁λ₂)−1. Here α represents area change, and β a change of aspect ratio or eccentricity. The shear modulus, μ, is then given as

(21)

With the hybrid model, we find that the shear modulus is heavily influenced by the persistence length p of the Sp. Indeed, with the configuration of the JC prescribed and the extension ratio fixed, we are able to show that μ is roughly proportional to 1/p as would be expected by the 1/p scaling within the WLC model. Consequently, accurate determination of the persistence length p is critical to the evaluation of shear modulus via our model. Owing to the inconclusiveness concerning the exact value of p, in Fig. 7 we plot μ normalized by k_BT/pL_c. The point-attachment junction and the wrap-around junction provide similar predictions for μ, while the variation over different deformation angle τ is ∼20% with the maximum value occurring near τ=90°. We note that by using p=8nm, our model of the single JC unit predicts that at a small deformation, i.e., λ=1.25, μ=7×10⁻⁶N/m. This is close to predictions obtained in other reports 2,4,5,6,18,19,20,42,43. Fig. 7 indicates strain stiffening of the network consistent with other modeling results (e.g., 20), whereas evidence for increasing network compliance with extension has been obtained by measuring the fluctuations of fluorescent beads attached to the network during micropipette aspirations 33. As noted, this may be attributed to effects such as domain unfolding within the Sp or possibly interprotein disassociation.

Display large version of this figure Display high quality version of this figure

Figure 7

Shear modulus of a single JC unit coupled with the lipid bilayer normalized by k_BT/pL_c obtained with the point-attachment junction and wraparound junction.

We further studied the response of the hybrid JC/lipid bilayer subject to thermal fluctuations. As illustrated in Fig. 8, corresponding to the Brownian motion, a stochastic wavelike motion is seen to occur within the lipid bilayer. In addition, the JC itself is subject to a random vibration. The mean values and standard deviations of the attitude of the actin protofilament in an equibiaxial deformation obtained with both the point-attachment and the wrap-around junctions are extracted from the simulations. What is found is that the mean values of both the yaw and the pitch angles are close to the quasi-static predictions. With the point-attachment junction, the deviation of the yaw angle decreases from 64° at τ=1 to 44° at τ=3.5. Similarly, with the wrap-around junctions, this deviation decreases from 49° to 21°. Such oscillation reduction is due to an increased tethering force from the Sp. The deviation of the pitch angle, caused primarily by the thermal fluctuation of the bilayer, remains at ∼25° and is not sensitive to τ.

Display large version of this figure Display high quality version of this figure

Figure 8

Snapshot of an equibiaxially extended JC (λ=2.5) during thermal fluctuation obtained with the point-attachment junction. The protofilament is shown as a cylinder, the SCs are shown as spheres, and the elevation of lipid bilayer is shown by the contour (bright areas represent high elevation and dark areas are low elevation, with the contours in nanometers). The pointed end of the protofilament is marked with an asterisk.

Corresponding to the bifurcation discussed before, a mode-switching phenomenon is observed in anisotropic deformations. A mode-switching event is shown in Fig. 9, with τ=0° and λ=3.5, and using the point-attachment junction. In this simulation, the yaw angle jumps from ∼−70° (corresponding to state 2 in Figure 5a) to 70° (state 1) within a short time of ∼5μs. The transition process is often accompanied by an abrupt increase of pitch angle up to 80°, indicating a bow-up motion (see Fig. 9). This motion is expected to affect the structural stability of the network and increase the disturbance on the bilayer. The actual occurrence of mode switching, however, is expected to be affected by more subtle molecular structures of the membrane not yet included in the model, e.g., the secondary linkages between the protofilament and the lipid bilayer through protein 4.1 and glycophorin C.

Display large version of this figure Display high quality version of this figure

Figure 9

A mode-switching event during thermal oscillation of the JC in anisotropic deformation with λ=3.5 and deformation angle τ=0°. The result was obtained with the point-attachment junction. The pointed end of the protofilament is marked with an asterisk.

The mechanical properties of the membrane are important for erythrocyte deformation. The actin/spectrin-based network has been studied for more than three decades, and the wealth of information on the system continues to motivate modeling efforts. When modeled at the molecular level the investigation may shed new light into how elements in the network may respond to the deformation and what may be the physiological consequences. Among the key insights provided by this model is that how spectrin associates with actin may influence the dynamics of a protofilament in response to deformation.

Here we present the first computational model that incorporates the static and the dynamic responses of a single JC coupled to a lipid bilayer subjected to external loads. These studies are perhaps the first to examine the hydrodynamic aspect of a model with certain degrees of molecular detail. The coupling simulated in this report is through the primary linkage at the SCs. New features considered in this model include two types of junctions between the actin protofilament and the Sp (i.e., point-attachment and wrap-around). The conjectured molecular detail of the wrap-around model is shown in Fig. 10. In both equibiaxial (area-changing) and anisotropic (area-conserving) deformations, the structural responses with the wrap-around junction are much more stable than those with the point-attachment junction. In contrast to the point-attachment junction, in which the orientation of the protofilament is much affected by the amount of deformation and the presence of the lipid bilayer, the wrap-around junction diminishes the sensitivity to deformations and the JC-bilayer interactions. The wrap-around junction also leads to a reduced pitch angle for the protofilament, consistent with experiments 24,25, and the motion of the protofilament exerts much less disturbance on the lipid bilayer. In addition, the wrap-around junction suppresses the occurrence of bifurcations in anisotropic deformations, which is the existence of more than one equilibrium state at a certain extension ratio and deformation angle. The bifurcation originates from the symmetric distribution of the SCs (forming a regular hexagon), the three-dimensional geometry of the protofilament, the protofilament-bilayer interaction, and the manner of Sp-protofilament connection, which allow the existence of multiple minima in terms of the total strain energy stored in the Sp. In the dynamic response of the system with thermal fluctuations, a bifurcation may cause mode switching, in which the system switches from one equilibrium state to another within a short time interval. A mode-switching event is often accompanied by a bow-up/bow-down motion of the protofilament, i.e., a sudden increase/decrease of the pitch angle. This type of motion imposes a large transient load on the lipid bilayer, which may affect the structural stability of the membrane and its physiological performance. These simulations suggest the importance of three-dimensional molecular junctions and the lipid bilayer/protein network coupling on cell membrane mechanics.

Display large version of this figure Display high quality version of this figure

Figure 10

A proposed Sp wrap around the actin protofilament annotated with known spectrin domains. The tail end of Sp begins with the nonhomologous domains in the C-terminus of α-spectrin (C-α) and N-terminus of β-spectrin (N-β) 26,52. Nonhomologous domains are involved in the interchain binding at the tail end 29, but not in the dimer nucleation 27. The dimer nucleation sites spanning homologous domains of β1–4 and α18–21 are proximal 27. The actin-binding domains (ABD) are on β-spectrin, comprised of two Calponin Homology (CH) subdomains 53. There are two EF hand domains at the C-terminal α-spectrin 54. Some of Cys residues (*) at the tail end may be active, capable of forming disulfide bonds 27. Protein 4.1R stabilizes the binding between β-spectrin and F-actin, with a split β-spectrin binding site binding to CH1 and CH2, and an actin binding site positioned in between 28,45,55.

This hybrid model forms a framework upon which a sequence of more sophisticated models, which will include additional molecular detail as well as other factors such as the viscoelasticity of the lipid bilayer and that of the Sp (see e.g., 44), will be developed. The resulting multiple-scale model may provide a method of understanding the exact role each molecular component may play in terms of the overall mechanical performance of the membrane. A near-term extension of our model will include the secondary linkage between the protein network and the lipid bilayer at JC involving actin, protein 4.1, and the transmembrane glycophorin C 45,46. The mechanical function of these secondary linkages and their physiological importance to the membrane will be simulated with sufficient three-dimensional molecular detail to describe the interactions among them. Our preliminary simulations show that these secondary linkages may diminish the thermal fluctuation of the protofilament.

Future studies will involve simulations of multiple units of JCs (see Figure 2b) attached to a dynamic bilayer as described herein. The goal will be to model the effects of SC mobility, the connectivity between JC and SC, and molecular detail of the linkages between the protein network and the lipid bilayer on overall membrane mechanical response. Based on published values for the short range diffusion coefficient of band 3 22, we have estimated the mobility for the SC. Band 3 is mobile within the plane of the lipid bilayer and this mobility would give rise to a time-dependent viscoelastic response of the network. The short range diffusion coefficient, D, of band 3 in the lipid bilayer at 37°C is 0.53μm²s⁻¹22, which corresponds to a mobility (via Brownian dynamics) of D/k_BT=1.28×10⁸m(Ns)⁻¹. Typical levels of tension inside an Sp during network deformation, based on our previous model 17 and the simulations reported on here, is estimated in the range of 2–10 pN. Under this forcing, the time for the SC to drift a distance comparable to the size of a JC (∼70nm) is ∼55–270μs. It is thus expected that the SC will be seen to be dragged through the lipid bilayer in a timescale associated with cell deformations, e.g., ∼0.1s in micropipette aspirations. In the future analysis of ensembles of JCs we will allow SC mobility within the bilayer, and we note that this imparts an entirely natural viscosity to the protein network/lipid bilayer response as yet unaccounted for. On the other hand, during the deformation of the entire cell the lipid bilayer must follow the global deformation and thus undergo shear in parallel with the skeleton. Later in describing the computed shear modulus we will add a viscous term accounting for this viscous contribution to the shear resistance. Multiple JC units, connected through periodic boundary conditions, will allow for the free motion of SCs—i.e., via a force-mobility relation as prescribed with the mobility described above. In this way lipid molecules will be seen to flow in and out of the areas occupied by individual JCs. Particular studies will include topological changes that are characteristic of certain diseases, where mutations have been identified in membrane skeletal components, such as those in hereditary spherocytosis, ovalocytosis, or elliptocytosis. Simulations involving multiple JCs, as depicted in Figure 2b, will include descriptions of Sp disassociation from actin or from the head-to-head dimer connections, or Sp unfolding/folding, as influenced by Sp tension 32,47. The formalism for such descriptions can follow an analysis based on thermally activated kinetics whereby one expresses rate of disassociation, or transitions, in the form

(22)

Our hybrid model is capable of predicting the dynamics of individual components of the erythrocyte membrane, and has demonstrated good agreement with other models (e.g., 18,19,20) and experimental results (e.g., 24,25). It also could be applied in the interpretation of experiments involving localized behavior of the membrane, e.g., the thermal fluctuation of membrane attached beads or quantum-dots. By simulating the thermal oscillation of a bead attached to the protofilament using our model incorporating unfolding of Sp repeats in large stretching, our preliminary results suggest an increase in bead fluctuation accompanying large shear deformation, a trend observed in experiments 33. This increase in fluctuation may be attributed to both the unstiffening effect due to unfolding of repeats, the onset of mode switching, or other mechanisms. Further investigations to quantitatively correlate our numerical predictions with experimental measurements are underway and will be reported in a forthcoming report. There are recently developed technologies that make measuring the three-dimensional orientation of biomolecules possible. For example, single-molecule fluorescence polarization microscopy has been used to determine the three-dimensional orientation of a single F-actin in vitro 48; and single-pair fluorescence resonance energy transfer using a wide-field total-internal-reflection fluorescence microscopy has been used to image single molecules in living cells, providing a quantitative measuring of inter- and intramolecular distances and angles 49. These innovative techniques may allow the attitude changes of a protofilament to be recorded during membrane deformation, thus providing a way to test our model.