M Misra

Basile Audoly

I G Kevrekidis

S Y Shvartsman

Shape Transformations of Epithelial Shells

Keywords: vertex models, epithelial morphogenesis, continuation, bistability Running Title: Shape transformations of epithelial shells

HAL is

Introduction

Regulated deformations of epithelial sheets play important roles during tissue morphogenesis and can be driven by a variety of mechanisms, the simplest of which rely on spatial patterns of apical cell contractility [START_REF] Lecuit | Cell surface mechanics and the control of cell shape, tissue patterns and morphogenesis[END_REF][START_REF] Sawyer | Apical constriction: a cell shape change that can drive morphogenesis[END_REF][START_REF] Zartman | Unit operations of tissue development: epithelial folding[END_REF]. These prepatterns are commonly established by upstream signaling processes, which upregulate the activity of the actomyosin networks on the apical surfaces or apical edges in a subset of cells within the epithelium [START_REF] Martin | Pulsed contractions of an actin--myosin network drive apical constriction[END_REF][START_REF] Wrighton | Forces of Change: Mechanics Underlying Formation of Functional 3D Organ Buds[END_REF]. Intracellular nonuniformities in contractility trigger shape changes of individual cells, leading to three-dimensional (3D) deformations, such as localized tissue invaginations. This sequence of processes, from patterning of apical contractility, to spatially restricted cell shape changes, to 3D tissue deformations has been documented in a wide range of experimental systems [START_REF] Sherrard | Sequential activation of apical and basolateral contractility drives ascidian endoderm invagination[END_REF][START_REF] Tamulonis | A cell--based model of Nematostella vectensis gastrulation including bottle cell formation, invagination and zippering[END_REF][START_REF] Osterfield | Three--dimensional epithelial morphogenesis in the developing Drosophila egg[END_REF][START_REF] Du | Computational analysis of three--dimensional epithelial morphogenesis using vertex models[END_REF][START_REF] Polyakov | Passive mechanical forces control cell--shape change during Drosophila ventral furrow formation[END_REF].

A canonical example is provided by the early stages of mesoderm invagination in Drosophila, where apical constriction of cells on the ventral side of the embryo transforms a convex epithelial shell into a more complex shape with an omega-like crosssection [START_REF] Leptin | Gastrulation movements: the logic and the nuts and bolts[END_REF]. In this system the apical surfaces of epithelial cells are facing towards the rigid membrane surrounding the embryo. In another well-studied experimental model, the developing Drosophila egg chamber, the epithelium has opposite polarity, with the apical surfaces oriented towards the oocyte, which is enclosed by the epithelial sheet [START_REF] Osterfield | Three--dimensional epithelial morphogenesis in the developing Drosophila egg[END_REF][START_REF] Berg | The Drosophila shell game: patterning genes and morphological change[END_REF]. In this case, a two-dimensional (2D) patch of the follicle cells evaginates, bending in the direction of the membrane surrounding the egg chamber. This deformation is thought to be guided, at least in part, by the embedded contour of apically constricting epithelial cells.

The general principles governing 3D deformations induced by patterns of apical contractility remain poorly understood. It is in general unclear what determines the direction in which the patterned sheet is going to bend and whether the resulting deformation will be discontinuous or smooth. Several mathematical and computational models have been proposed to explore these questions, but none of them are sufficiently versatile for analyzing the interplay of multiple physical and geometrical factors [START_REF] Sherrard | Sequential activation of apical and basolateral contractility drives ascidian endoderm invagination[END_REF][START_REF] Tamulonis | A cell--based model of Nematostella vectensis gastrulation including bottle cell formation, invagination and zippering[END_REF][START_REF] Polyakov | Passive mechanical forces control cell--shape change during Drosophila ventral furrow formation[END_REF][START_REF] Odell | The mechanical basis of morphogenesis. I. Epithelial folding and invagination[END_REF][START_REF] Davidson | How Do Sea--Urchins Invaginate --Using Biomechanics to Distinguish between Mechanisms of Primary Invagination[END_REF][START_REF] Conte | A 3D finite element model of ventral furrow invagination in the Drosophila melanogaster embryo[END_REF][START_REF] Pouille | Hydrodynamic simulation of multicellular embryo invagination[END_REF][START_REF] Jones | Modelling apical constriction in epithelia using elastic shell theory[END_REF][START_REF] Hocevar Brezavscek | A model of epithelial invagination driven by collective mechanics of identical cells[END_REF]. For instance, multiple models have been used to describe the invaginations in the early embryos [START_REF] Sherrard | Sequential activation of apical and basolateral contractility drives ascidian endoderm invagination[END_REF][START_REF] Tamulonis | A cell--based model of Nematostella vectensis gastrulation including bottle cell formation, invagination and zippering[END_REF][START_REF] Polyakov | Passive mechanical forces control cell--shape change during Drosophila ventral furrow formation[END_REF][START_REF] Hocevar Brezavscek | A model of epithelial invagination driven by collective mechanics of identical cells[END_REF]. However, most of these models consider only a 2D crosssection of the epithelial shell, which limits the class of shape transformations that can be analyzed.

Here we provide what appears to be a minimal, yet physically realistic computational framework that enables systematic exploration of epithelial deformations induced by spatial patterns of active and passive cell properties. Our approach is based on the recently proposed energy formulation in which every cell is modeled as a prism with apical, basal, and lateral surfaces that can have different properties [START_REF] Hannezo | Theory of epithelial sheet morphology in three dimensions[END_REF]. We use this model to construct epithelial shells enclosing a fluid-filled volume and surrounded by a hard membrane, mimicking a scenario encountered in developing tissues and tissue organoids. We then use numerical continuation algorithms to explore how these model epithelial shells deform in response to the prepatterns of apical contractility. We find that, for most cases, the deformations are smooth, i.e. they happen without bifurcations, when viewed as a function of the amplitude that characterizes the spatial pattern of apical contractility. We also demonstrate that a simpler 2D model, in which epithelial shells are constructed from polygonal cells, can describe the effects predicted by the 3D model and can be used to explore a wide range of morphogenetic processes.

Materials and Methods

3D vertex model

Vertex models constructed from 3D cells have been proposed to study epithelial morphogenesis in multiple experimental systems [START_REF] Hannezo | Theory of epithelial sheet morphology in three dimensions[END_REF][START_REF] Honda | A three--dimensional vertex dynamics cell model of space--filling polyhedra simulating cell behavior in a cell aggregate[END_REF][START_REF] Okuda | Reversible network reconnection model for simulating large deformation in dynamic tissue morphogenesis[END_REF][START_REF] Okuda | Vertex dynamics simulations of viscosity--dependent deformation during tissue morphogenesis[END_REF][START_REF] Bielmeier | Interface Contractility between Differently Fated Cells Drives Cell Elimination and Cyst Formation[END_REF]. Based on these previous studies, we used the following energy functional to model an epithelial monolayer constructed from cells with distinct apical and basal surfaces (Fig. 1(A,B)):

𝐸 !! = 𝜎 𝑙 ! + ! 𝛼 𝑆 ! + 𝛾 𝑆 ! + 𝐵 𝑉 ! -𝑉 ! ! ! ! ! ! . (1)
In this expression, the first term sums over all apical edges 𝑒, the second term sums over lateral surfaces 𝑙, the third term sums over basal surfaces 𝑏, and the last term sums over cells 𝑐 . The first term corresponds to the line tension along the apical edges of neighboring cells; 𝑙 ! is the edge length and 𝜎 is the apical line tension coefficient. The second and third terms correspond to the contributions from lateral and basal tensions, which are proportional to the lateral and basal surface areas (𝑆 ! , 𝑆 ! ) with coefficients 𝛼 and 𝛾, respectively. These coefficients model the resistance to deformation due to the cytoskeletal meshwork that underlies the surface. The last term penalizes the deviation of the cell volume, 𝑉 ! , from its target value, 𝑉 ! ! , with compression modulus 𝐵 to approximate cytoplasmic incompressibility. For more details see Section S1 in the Supporting Material.

In our 3D model, we define prepatterned apical edges to represent a contractile ring 𝑅 embedded in the apical surface (shown in red in Fig. 2(A1), 3(A1)) or an apically constricting patch of cells 𝑃 (shown in red in Fig. 2(B1), 3(B1)). The first pattern is motivated by the embedded 'compressive cable' observed in follicle cells during Drosophila appendage formation and early sea urchin embryo, which displays localization of the motor protein myosin [START_REF] Osterfield | Three--dimensional epithelial morphogenesis in the developing Drosophila egg[END_REF][START_REF] Nakajima | The initial phase of gastrulation in sea urchins is accompanied by the formation of bottle cells[END_REF]. The second pattern is motivated by biological settings in which a group of cells are influenced by uniform actomyosin constriction across the apical surface [START_REF] Sawyer | Apical constriction: a cell shape change that can drive morphogenesis[END_REF][START_REF] Mason | Apical domain polarization localizes actin--myosin activity to drive ratchet--like apical constriction[END_REF][START_REF] Martin | Apical constriction: themes and variations on a cellular mechanism driving morphogenesis[END_REF]. Apical edges belonging to the contour or the patch are assigned an apical line tension (𝜎 + 𝛤), which is larger than the apical line tension 𝜎 of the other non-patterned edges. This is implemented through an additional energy term for the prepatterned edges

𝐸 ! = 𝛤 𝑙 ! ! ∈ ! !" ! , (2)
where the index 𝑖 runs over all edges belonging to the apical ring 𝑅 or apical patch 𝑃.

2D vertex model

Most of the previously published 2D vertex models were used to describe the dynamics of epithelia constrained to two spatial dimensions [START_REF] Nagai | A dynamic cell model for the formation of epithelial tissues[END_REF][START_REF] Farhadifar | The influence of cell mechanics, cell--cell interactions, and proliferation on epithelial packing[END_REF][START_REF] Fletcher | Vertex models of epithelial morphogenesis[END_REF]. More recently, we have used these models to explore 3D deformations induced by prepatterning of cell properties [START_REF] Osterfield | Three--dimensional epithelial morphogenesis in the developing Drosophila egg[END_REF][START_REF] Du | Computational analysis of three--dimensional epithelial morphogenesis using vertex models[END_REF][START_REF] Murisic | From Discrete to Continuum Models of Three--Dimensional Deformations in Epithelial Sheets[END_REF]. We extended these models so that they can account for the natural curvature of the sheet, and, importantly, for the intrinsic apicobasal polarity of epithelial cells. The energy in this model is defined in the following way:

𝐸 !! = 𝜇 𝐴 ! -𝐴 ! ! ! + 𝜎 𝑙 ! + 𝛽 1 -𝑵 !! ! ! • 𝑵 !! ! ! + 𝑘 ! ! 𝑵 !! ! ! - ! ! ! ! 𝑵 !! ! ! • 𝒖 ! ! . (3)
The first term sums over all cells 𝑐, the second term sums over all edges 𝑗, and the last term sums over all interior edges 𝑗 ! , i.e. edges shared by two adjacent cells 𝑠1(𝑗 ! ) and 𝑠2(𝑗 ! ). The first term corresponds to sheet elasticity, where 𝐴 ! is the area of cell 𝑐, 𝐴 ! ! is its target area, and 𝜇 is the stretching modulus. The second term captures intercellular interactions in the form of tensile forces along cell-cell edges, where 𝑙 ! is the length of the edge 𝑗 and 𝜎 is the line tension coefficient. The last term represents the bending energy term; where 𝛽 is the bending elasticity coefficient and 𝑘 ! ! represents the local natural curvature of the shell (𝑘 ! ! is roughly the rest value of the dihedral angle between the faces adjacent to edge j' and it scales like the natural curvature of the midsurface times the typical cell size; see Supporting Material Section S2.1 for details). Here 𝑵 ! denotes the 'outward' unit normal to a cell 𝑠 (i.e. the normal vector of each cell points towards the 'outside' of the closed shell). 𝒖 ! ! is the unit vector joining the cell centers (Fig. 1 (C,D)) ,

𝒖 ! ! = 𝑪 !!(! ! ) !𝑪 !!(! ! ) |𝑪 !! ! ! !𝑪 !! ! ! | . ( 4
)
where 𝑪 !! ! ! and 𝑪 !! ! ! are cell centers, defined as the mean position of the nodes belonging to the cells 𝑠1 and 𝑠2, respectively (Fig. 1(C,D)). The bending energy term defines an "at rest" value of the angle between the normals of adjacent cells, which depends on 𝑘 ! ! (and is proportional to 𝑘 ! ! for small 𝑘 ! ! ) and hence implements the notion of natural curvature. For more details see Section S2 in the Supporting Material.

In our 2D model, we idealize the epithelial sheet as the midsurface of a 3D epithelial monolayer (Fig. 1(C)). We proceed to show how we mimic the effect of the prepatterns in our 2D model that are present on the apical surface of an epithelium (Fig. 1(D)). From now on, we will refer to this representation of prepatterns as being offset from the midsurface. Such prepatterns include a contractile ring 𝑅 (Fig. 5(A1)) and a constricting patch 𝑃 (Fig. 5(B1): a contractile ring 𝑅 is defined by a closed path joining the centers of a ring of adjacent cells (Fig. 5(A1)), while a constricting patch P represents a subset of the mesh and is implemented by adding a line tension along all segments joining the adjacent cell centers that are contained within this subset (Fig. 5(B1)). Segments belonging to either type of prepattern are assigned a line tension 𝛤 and an offset parameter h, and an additional energy term is considered:

𝐸 ! = 𝛤 (|𝑪 !! ! ! ∈! !" ! -𝑪 !! ! | -ℎ 𝑵 !! ! -𝑵 !! ! • 𝒖 ! ) . ( 5
)
where the index 𝑖 runs over all edges belonging to the apical ring 𝑅 or apical patch 𝑃, and h is the half the thickness of the real tissue.

The key feature of equation ( 5) is that it represents the line tension of segments that are offset by a distance h from the midsurface (see Supporting Material Section S2 for details), but does so by using solely the degrees of freedom of the 2D model, namely the positions of the nodes lying onto the midsurface. The line tension parameter 𝛤 acts in two ways: it brings the cell centers closer together and, when h ≠ 0, it bends the surface so that the end points of the normal vectors come close together (h > 0) or further apart (h < 0). Note that the line tension prepatterning parameter 𝛤 in the model with 3D cells has qualitatively similar effects.

In the special case 𝑘 ! ! = 0 and ℎ = 0, this model reduces to the 2D vertex model which has been used in the literature to understand epithelial deformations for flat epithelial sheets in different contexts [START_REF] Osterfield | Three--dimensional epithelial morphogenesis in the developing Drosophila egg[END_REF][START_REF] Du | Computational analysis of three--dimensional epithelial morphogenesis using vertex models[END_REF][START_REF] Farhadifar | The influence of cell mechanics, cell--cell interactions, and proliferation on epithelial packing[END_REF][START_REF] Murisic | From Discrete to Continuum Models of Three--Dimensional Deformations in Epithelial Sheets[END_REF]. Hereafter, we will refer to that particular case as the naturally planar 2D vertex model (since 𝑘 ! ! = 0, the bending energy term corresponds to a naturally flat sheet). Note that the new terms which we introduce in the present work, namely the one proportional to 𝑘 ! ! in Eq. ( 3) and the one proportional to ℎ in Eq. ( 5), both depend on the nodal positions through the same expression

𝑵 !! ! -𝑵 !! ! • 𝒖 ! .
As a result, both terms have similar numerical implementation.

Modeling the effects of constraints

To study the effects of epithelial shell curvature, a fluid-filled inner cavity and an outer stiff membrane during epithelial bending, we consider an initial homogeneous spherical configuration and include additional terms in the energy functions:

𝐸 = 𝐵 ! 𝑉 ! -𝑉 ! ! ! + 𝜖 ! ! ! !! ! ! + 𝐸 !!/!! ! + 𝐸 ! . (6)
Here, the first term penalizes deviations from the initial volume of the inner cavity. 𝑉 ! is the volume enclosed by the closed surface, 𝑉 ! ! is the initial volume, and 𝐵 ! is a compression modulus. The second term corresponds to the outer membrane stiffness [START_REF] Pouille | Hydrodynamic simulation of multicellular embryo invagination[END_REF][START_REF] Hocevar Brezavscek | A model of epithelial invagination driven by collective mechanics of identical cells[END_REF] that runs over all vertices on the outer surface of the shell and restricts the radial motion of vertices within a sphere of radius 𝑅 ! . This sphere is concentric with the homogeneous spherical configuration from which the system is initialized. The center of the initial spherical configuration acts as reference center to calculate the radial distance of the vertices. The term 𝑅 ! denotes the radial distance of the vertex 𝑘, and hence 𝑅 ! -𝑅 ! represents the membrane thickness near that vertex. 𝜖 is a membrane stiffness parameter, and 𝑛 is the exponent of the repulsive potential term that models the effects of the outer stiff membrane.

Numerical methods

To find an equilibrium shape for a typical configuration in the 3D vertex model, we solved the system of nonlinear algebraic equations that correspond to vanishing of the gradient of the energy with respect to node positions by Newton-Raphson iteration. The number of equations was halved for the corresponding configurations in the 2D vertex model. Initial guesses were obtained by direct forward Euler integration following overdamped gradient dynamics. In some cases, rotational and translational symmetries had to be factored out through appropriate pinning conditions, as described in Supporting Material Section S3. Pseudoarclength continuation [START_REF] Keller | Lectures on numerical methods in bifurcation problems[END_REF][START_REF] Nayfeh | Applied nonlinear dynamics: analytical, computational and experimental methods[END_REF] was used to follow solution branches in parameter space and to go around turning points. The eigenvalues of the Jacobian upon convergence quantify the stability of the computed equilibrium shapes. Branches were terminated when equilibrium solutions featured vertices too close to each other (edge length < 0.01). The merging of vertices and the rearrangement of cell neighbors were not considered here. We used the Armadillo, a C++ linear algebra library to solve the linear systems of equations and calculate leading eigenvalues [START_REF] Sanderson | Armadillo: An open source C++ linear algebra library for fast prototyping and computationally intensive experiments[END_REF].

Homogeneous configurations and model parameters

In all cases considered below, we start with a system with spatially uniform cell properties at a mechanical equilibrium. This equilibrium configuration is then used as a starting point for numerical continuation that analyzes the effects of spatial patterns of cell properties, with the amplitude of the pattern chosen as continuation parameter. To construct the initial equilibrium states for our model epithelial shells, we analyzed an idealized system sheet with uniform cell properties. For the 2D vertex model, we considered a shell constructed from regular hexagons (Fig. 1(D)). For the 3D model, we considered a shell constructed from a 'lampshade' shape with hexagonal apical/basal surfaces (Fig. 1(A)). In both these cases, one can analyze the explicit expressions for the overall energy as a function of the geometric and model parameters. The geometric parameter for the 2D model is the edge length of the unit hexagonal cell. For the 3D model, the geometric parameters are the edge length of the hexagonal apical and basal surfaces and the height of the unit cell. The minima of such functions provide a relationship between equilibrium cell shape and tissue properties.

Once the parameters were found for an idealized sheet at an equilibrium, we used DistMesh to realize an epithelial shell with a finite number of cells tiling a closed surface [START_REF] Persson | A simple mesh generator in MATLAB[END_REF]. This meshing package builds triangular tessellations on a sphere for the given edge length. As an input geometry for the 2D vertex model, we used these tessellations to construct polygonal cells by considering centers of triangular mesh elements as vertices of epithelial cells. To construct the 3D vertex model input geometry, we first constructed the inner surface of the epithelial shell using the 2D model input geometry approach and then extended the vertices radially to a given height to form 3D cells. This produces a globally curved epithelial shell with finite number of 3D 'lampshade'-shaped cells, i.e. cells with asymmetric apical and basal surface. Further details are provided in the Supporting Material Section S6.1.

Results

3D vertex model

3D vertex models have previously been used to study dynamic tissue morphogenesis for cell aggregates and spatially uniform cell sheets [START_REF] Okuda | Reversible network reconnection model for simulating large deformation in dynamic tissue morphogenesis[END_REF][START_REF] Okuda | Vertex dynamics simulations of viscosity--dependent deformation during tissue morphogenesis[END_REF][START_REF] Nagai | A dynamic cell model for the formation of epithelial tissues[END_REF]. Recently Hannezo et al. extended a similar framework to model different shapes of epithelial cells [START_REF] Hannezo | Theory of epithelial sheet morphology in three dimensions[END_REF]. Here, we use a simplified 3D vertex model (Eq. 1), along with the additional term for prepatterns (Eq. 2), to numerically find equilibrium shapes in an epithelial monolayer induced by the in-plane prepatterns of apical contractility. We also analyze the effects of a fluid-filled cavity and an outer membrane (Eq. 6) on out-of-plane deformations induced by such prepatterns. All of these factors are present in different developmental contexts including Drosophila appendage formation [START_REF] Berg | The Drosophila shell game: patterning genes and morphological change[END_REF][START_REF] Osterfield | Diversity of epithelial morphogenesis during eggshell formation in drosophilids[END_REF] and the blastula-to-gastrula transition in different species [START_REF] Sherrard | Sequential activation of apical and basolateral contractility drives ascidian endoderm invagination[END_REF][START_REF] Tamulonis | A cell--based model of Nematostella vectensis gastrulation including bottle cell formation, invagination and zippering[END_REF][START_REF] Polyakov | Passive mechanical forces control cell--shape change during Drosophila ventral furrow formation[END_REF][START_REF] Hocevar Brezavscek | A model of epithelial invagination driven by collective mechanics of identical cells[END_REF][START_REF] Solnica--Krezel | Gastrulation: making and shaping germ layers[END_REF][START_REF] Rauzi | Physical models of mesoderm invagination in Drosophila embryo[END_REF] To understand the effect of model parameters on cell morphology and to provide good initial guesses for the numerical calculation of equilibrium shapes, we calculated equilibrium shapes for flat, as well as curved homogeneous epithelial monolayers. The analysis was performed for ideal cases where apical and basal surfaces were approximated as hexagons. For the 3D vertex model with planar cells, increasing the lateral tension coefficient (𝛼) makes the equilibrium cell shape more columnar, whereas increasing the basal tension coefficient (𝛾) makes the cell more squamous (Fig S4 (C,D)). We discuss the effects of basal and lateral tension coefficients in greater detail in Supplementary Information Section S6.1. We use this information to construct the initial homogeneous configuration (𝛤 = 0) with given model parameters (see Materials and Methods) and then apply different prepatterns (𝛤 > 0) and solve for the equilibrium shapes.

We started with a closed shell monolayer of 3D cells, which had polarity similar to the cells in the follicular epithelium in the developing Drosophila egg, where the apical side forms the inner surface of the epithelium (apical-in case, Fig. S1(D),( 8)). Consistent with recent observations [START_REF] Gelbart | Volume conservation principle involved in cell lengthening and nucleus movement during tissue morphogenesis[END_REF], we kept the value of the cell-volume compression modulus B large, to make sure that the deviation from the initial cytoplasmic volume of each cell was small. We then introduced spatial patterns of apical contractility. We considered two prepatterns of apical constriction in different settings: (I) an embedded contractile ring (Fig. 2(A1)); and (II) a patch of apically constricting cells (Fig. 2(B1)). As mentioned before, the first prepattern is motivated by the study of Drosophila appendage formation [START_REF] Osterfield | Three--dimensional epithelial morphogenesis in the developing Drosophila egg[END_REF]. The second prepattern is found in diverse developmental settings in which a group of cells undergo apical constriction driving localized tissue bending [START_REF] Sawyer | Apical constriction: a cell shape change that can drive morphogenesis[END_REF][START_REF] Mason | Apical domain polarization localizes actin--myosin activity to drive ratchet--like apical constriction[END_REF][START_REF] Martin | Apical constriction: themes and variations on a cellular mechanism driving morphogenesis[END_REF].

We found one continuous branch of stable equilibria in which the epithelium bent in the basal direction (evagination) for both apical contractility prepatterns (Fig. 2(A4,B4)). Due to cytoplasmic incompressibility, a patch of cells undergoing apical constriction leads to the expansion of the basal surface and thus causes the sheet to bend in the apicalto-basal direction. Inclusion of an outer membrane term kept the configuration more spherical (Fig S8 ), and increasing the inner-fluid compression modulus (𝐵 ! ) decreased deflection without affecting the qualitative nature of the equilibrium shapes (Fig S9 ).

We then considered a scenario in which the apical side forms the outer surface of the shell (apical-out case, Fig. S1(C)) by reversing the polarity and keeping the energy formulation the same. A similar scenario is observed during the blastula-to-gastrula transition in many species, such as sea urchin [START_REF] Nakajima | The initial phase of gastrulation in sea urchins is accompanied by the formation of bottle cells[END_REF] and Drosophila [START_REF] Rauzi | Physical models of mesoderm invagination in Drosophila embryo[END_REF]. We introduced similar prepatterns of apical contractility as before (Fig. 3(A1,B1)) and found one continuous branch of stable steady states for the case of a patch of apically constricting cells in which the epithelial shell formed an invagination (Fig. 3(B4)). However, in the case of prepatterning with a contractile ring, we observed bistability in the shapes of shells (Fig 3 (A4)). This shape bistability was not present in configurations with smaller shell curvature (Fig. S12) and with a one-cell wide ring (Fig. S11). We also considered a prepattern of line tension in the form of a two-cell wide ring of apically constricting cells, which was inspired by the pattern of actomyosin localization during primary invagination in the sea urchin embryo [START_REF] Nakajima | The initial phase of gastrulation in sea urchins is accompanied by the formation of bottle cells[END_REF][START_REF] Kimberly | Bottle cells are required for the initiation of primary invagination in the sea urchin embryo[END_REF]. Here, as well, we found one continuous branch of stable invaginated shapes (Fig 3 (C4)).

To summarize, in response to most of the analyzed patterns of apical contractility, the 3D model predicts a smooth transition to the evaginated state when the apical surface forms the inner surface and to the invaginated states when the apical surface forms the outer surface of the shell.

Naturally planar 2D vertex model

Naturally planar 2D vertex models provide a computationally simpler framework to model epithelial deformations, while still incorporating the essential physical features of the apical surface of the epithelium that drive such transformations. In this section we determine what aspects of the 3D model can be captured with a naturally planar 2D vertex model. We use the model presented by Murisic et al. (30) as our starting point. This model can be viewed as a special case of the 2D model introduced here (ℎ = 0; 𝑘 ! ! = 0). Similar to our analysis of the model with 3D cells, we examined the case of a closed shell with different prepatterns of apical cell contractility. In contrast to the 3D model, the naturally planar model showed no out-of-plane deformations for a constricting patch prepattern. However, for a contractile ring, we found two disconnected branches of steady state solutions (Fig. 4). Equilibrium shapes with cells deflecting outward (evaginated state) formed a continuous branch of steady states from an initially homogeneous configuration, whereas solutions with cells deflecting inward (invaginated state) formed a disconnected branch (Fig. 4(E)). This steady state diagram can be traced back to the pitchfork bifurcation observed by Murisic et al. [START_REF] Murisic | From Discrete to Continuum Models of Three--Dimensional Deformations in Epithelial Sheets[END_REF] in the presence of contour forces for flat configurations. However, the pitchfork bifurcation here is imperfect, as the sheet curvature breaks the symmetry present in the flat configuration case (Fig S5 ). We observed that the turning point of the branch of invaginated states shifts to the right as the curvature of the initial homogeneous configuration is increased. This is consistent with the past analysis of buckling in Elastica (a thin strip of elastic material) with imperfections, such as preferred curvature or misaligned axial load [START_REF] Tauchert | Large deformation and postbuckling behavior of an initially deformed rod[END_REF][START_REF] Bažant | Stability of structures: elastic, inelastic, fracture and damage theories[END_REF][START_REF] Xu | Stability Analysis of Nonlocal Elastic Columns with Initial Imperfection[END_REF].

The presence of an outer membrane term tends to keep the configuration more spherical (Fig. S6A,B). As expected, increasing the inner-fluid compression modulus (𝐵 ! ) and bending elasticity coefficient (𝛽) made it more difficult to bend the epithelium. It decreased the deflection of evaginated shapes for a given value of 𝛤 and shifted the limit point to the right (Fig S6D ,E). We observed similar imperfect pitchfork bifurcations for different stretching moduli 𝜇.

In sum, we observed a typical imperfect pitchfork bifurcation for the closed shell configuration. The addition of an inner fluid-filled cavity and outer membrane term does not affect the qualitative response of the naturally planar 2D vertex model, influencing only the domain of shape bistability. The results of the naturally planar 2D vertex model suggest that the initially curved and homogeneous epithelium can smoothly bend, due to contour forces, only in the outward direction, which is similar to the response of the apical-in case of the 3D model to a contractile ring. Thus, this model does not capture all of the effects predicted by the 3D model, particularly the deformations driven by an apically constricting patch.

2D vertex model

As noted above, the naturally planar 2D vertex model cannot capture different aspects of the 3D model, including apicobasal polarity, natural curvature of the sheet, and a smooth transition from homogeneous configuration to invaginated states in the presence of an apical constriction pattern. Here we show that the 2D model with natural curvature, as introduced in eqs. 3-5, is fully consistent with the 3D model responses to different prepatterns, while maintaining the simpler framework of the 2D model. The model includes additional terms making the sheet naturally curved and incorporates apicobasal polarity by offsetting the prepatterns away from the sheet's midsurface such that they lie on the apical side of the shell.

Similar to previous sections, we considered two prepatterns of apical contractility: a contractile ring (Fig. 5(A1)) and a constricting patch (Fig. 5(B1)), now offset from the 2D sheet. We started with a homogeneous closed shell configuration at equilibrium, for which we set the value of 𝑘 ! ! equal to the average angle between the normals of the adjacent cells (see Supporting Material Section S2.1 for details). We then introduced the prepatterns of constriction in two settings corresponding to apical-in (ℎ < 0) and apicalout (ℎ > 0) cases (note that we defined the cell normals in both cases such that they always point 'outwards' with respect to the shell). For the apical-in case, we found that each prepattern produces a continuous branch of equilibrium shapes as we increased 𝛤 in which the epithelium bent in the 'outward' direction (evaginated state) (Fig. 5(A3), 5(B3)). For the apical-out case, the constricting patch prepattern produces one continuous branch of equilibrium shapes with the epithelium bending in the 'inward' direction (invaginated state) (Fig. 6(B3)). However, we observed bistability in the shapes of shells in the case of the contractile ring prepattern (Fig. 6(A3)). Increasing the inner-fluid compression modulus (𝐵 ! ) decreased deflection without affecting the qualitative nature of steady state diagram (Fig S10). Thus, the results obtained with the 2D vertex model are qualitatively similar to results from the 3D vertex model, suggesting that the additional terms can effectively capture the features of the complete 3D description of the epithelial sheet and its responses to different patterns of apical cell contractility.

Discussion

We used vertex models to explore how the spatial patterns of cell contractility lead to 3D deformations of epithelial shells. Previous studies of such deformations used appropriately modified elastic shell theories to make connections between biophysical descriptions of epithelia and a large body of tools and ideas from continuum mechanics [START_REF] Conte | A 3D finite element model of ventral furrow invagination in the Drosophila melanogaster embryo[END_REF][START_REF] Jones | Modelling apical constriction in epithelia using elastic shell theory[END_REF][START_REF] Hohn | Dynamics of a Volvox embryo turning itself inside out[END_REF][START_REF] Krajnc | Theory of epithelial elasticity[END_REF]. Although these studies can mimic a number of canonical processes, including vertebrate neurulation, they are difficult to adapt to different biological systems, such as morphogenesis driven by cell rearrangements, directed migration, or cell proliferation [START_REF] Solnica--Krezel | Gastrulation: making and shaping germ layers[END_REF]. Furthermore, the spatial patterns of apical contractility that drive 3D epithelial deformations are commonly fine grained. For instance, the early stages of sea urchin gastrulation are triggered by a contractile ring that is only two cells wide [START_REF] Nakajima | The initial phase of gastrulation in sea urchins is accompanied by the formation of bottle cells[END_REF]. These observations provide a clear motivation for models that can resolve individual cells.

A number of cell-based models, similar to the vertex models used in our work, have been used to describe deformations in 2D cross-sections of epithelial tissues [START_REF] Sherrard | Sequential activation of apical and basolateral contractility drives ascidian endoderm invagination[END_REF][START_REF] Tamulonis | A cell--based model of Nematostella vectensis gastrulation including bottle cell formation, invagination and zippering[END_REF][START_REF] Fletcher | Vertex models of epithelial morphogenesis[END_REF][START_REF] Murisic | From Discrete to Continuum Models of Three--Dimensional Deformations in Epithelial Sheets[END_REF]. The most extensive work has been done to model the formation of the ventral furrow in the early Drosophila embryo [START_REF] Polyakov | Passive mechanical forces control cell--shape change during Drosophila ventral furrow formation[END_REF][START_REF] Rauzi | Physical models of mesoderm invagination in Drosophila embryo[END_REF]. We extended these models to three dimensions, making it possible to investigate deformations induced by a broader class of spatial patterns of apical contractility. Our results can be summarized as follows. Model epithelia constructed from 3D cells, with distinct apical and basal properties, can readily generate both invaginated and evaginated shapes. The direction of bending is dictated by the apicobasal orientation of the epithelium: cells with increased apical contractility bend out of the sheet in the apical-to-basal direction. In contrast, naturally planar 2D vertex models predict that two different equilibrium shapes (invaginated and evaginated) coexist for the same set of parameters. The invaginated state exists as an isolated branch, which cannot be reached by continuous changes in the parameters of the spatial prepatterns; this is in contrast to the predictions of the 3D model. It also lacks the ability to include the apicobasal polarity or the natural curvature of the epithelial monolayer.

To address these issues, we proposed a modified 2D vertex model that better accounts for the 3D nature of the epithelium. This was accomplished by incorporating additional terms in the energy function. These terms effectively offset contractile patterns from the sheet and include the natural curvature of the epithelium. The modified model can readily capture the full effects of the 3D model. Most importantly, it can generate both evaginated and invaginated shapes depending on the apicobasal polarity of the epithelial shell. This suggests that the modified 2D vertex model can be used to describe shape transformations in a broader class of systems, especially those in which 3D deformations are driven by asymmetric apical contractility patterns, while maintaining the computational simplicity of the 2D framework. In our current work, we are using these models to explore shape transformations that include cell rearrangements and directed cell migration. S3 in the Supporting Material. S4 in the Supporting Material.

Figure 1 :

 1 Figure 1: Schematic representation of our vertex models. (A-B) Schematic representation of the 3D model. (A) A model 3D cell with distinct apical, basal, and lateral surfaces. (B) Schematic highlighting different tension terms in the 3D vertex model. (C) Schematic explanation of the modified 2D vertex model. The modified 2D model 'idealizes' the epithelial sheet as the midsurface (shown in black) of an epithelial monolayer with finite thickness '2h' and assumes that a contractile segment (shown in red) joins the centers of the apical faces. (D) The figure shows two adjacent 2D cells S1 and S2 with shared edge j, cell centers (C S1 , C S2 ) and unit 'outward' normals (N S1 , N S2 ), pointing towards the outer shell. u j is the unit vector joining the cell centers and h is the offset parameter. See text for more details. The red line highlights a contractile segment that is offset by a distance h from the surface along the cell normals..

Figure 2 :

 2 Figure 2: 3D deformations induced by prepatterns of line tension in a model with 3D cells for the apical-in case; see text for details. (A1-B1) Initial configurations, showing the inner surface (i.e. apical side and omitting the outer basal surface) of the shell, with different heterogeneities: (A1) Contractile ring. (B1) Patch of apically constricting cells. (A2-B2) Representative evaginated states, showing the inner surface, where the enclosed patch bends outward (positive deflection, 𝛿 > 0, as defined in Fig. S3). (A3-B3) Crosssection representation of equilibrium shapes (A2-B2) respectively. (A4-B4) Steady state diagrams showing deflection 𝛿 with increasing parameter 𝛤. Solid line: stable steady states. Cross-sections of representative steady states (A1-B1, A2-B2) are shown as insets. Parameter values are listed in Table S2 in the Supporting Material.

Figure 3 :

 3 Figure 3: 3D deformations induced by prepatterns of line tension in a model with 3D cells for the apical-out case; see text for details. (A1-C1) Initial configurations, showing the outer surface (i.e. apical side) of the shell, with different heterogeneities: (A1) Contractile ring. (B1) Patch of apically constricting cells. (C1) A ring of apically constricting cells, two cells wide. (A2-C2) Representative invaginated states, showing the outer surface, where the enclosed patch bends inward (negative deflection, 𝛿 < 0 as defined in Fig. S3). (A3-C3) Cross-section representation of equilibrium shapes (A2-C2) respectively. (A4-C4) Steady state diagrams showing deflection 𝛿 with increasing parameter 𝛤 . Solid (dashed) line: stable (unstable) steady states. Cross-sections of representative steady states (A1-C1, A2-C2) are shown as insets. Parameter values are listed in Table S2 in the Supporting Material.

Figure 4 :

 4 Figure 4: 3D bending of closed epithelial shells in the naturally planar 2D vertex model due to the embedded contractile ring. (A) Initial equilibrium configuration with the contractile ring (shown in red). (B) Representative invaginated state, where the enclosed patch bends inward (negative deflection, 𝛿 < 0). (C) Representative evaginated state, where the enclosed patch bends outward (positive deflection, 𝛿 > 0). Dashed red line represents patterned edges that lie behind the evaginated patch. (D) Representative unstable equilibrium configuration. (E) Steady state diagram showing deflection 𝛿 with increasing parameter 𝛤. Solid (dashed) line: stable (unstable) steady states. Insets (I) and (II) highlights the breakup of the pitchfork bifurcation in the naturally planar 2D vertex model due to curvature (discussed more in Fig. S1). Cross-sections of representative steady states (A-D) are shown as insets. Parameter values are listed in TableS3in the Supporting Material.

Figure 5 :

 5 Figure 5: 3D deformations induced by prepatterns of line tension in a model with 2D cells for the apical-in case; see text for details. (A1-B1) Initial configurations with different heterogeneities: (A1) Contractile ring. (B1) Constricting patch. Both are implemented as subsets of the dual mesh, which is used to offset contractile patterns; see text for details. (A2-B2) Representative evaginated states. (A3-B3) Steady state diagrams showing deflection 𝛿 with increasing parameter 𝛤 corresponding to the prepatterning defined in (A1-B1) respectively. Solid line: stable steady states. Cross-sections of representative steady states (A1-B1, A2-B2) are shown as insets. Parameter values are listed in TableS4in the Supporting Material.

Figure 6 :

 6 Figure 6: 3D deformations induced by prepatterns of line tension in a model with 2D cells for the apical-out case; see text for details. (A1-B1) Initial configurations with different heterogeneities: (A1) Contractile ring. (B1) Constricting patch. (A2-B2) Representative invaginated states. (A3-B3) Steady state diagrams showing deflection 𝛿 with increasing parameter 𝛤 . Solid (dashed) line: stable (unstable) steady states. Schematic cross-sections of representative steady states (A1-B1, A2-B2) are shown as insets. Parameter values are listed in TableS4in the Supporting Material.

Acknowledgements: We thank Miriam Osterfield, Alexander Fletcher, Satoru Okuda, XinXin Du, Nebojsa Murisic, and Vincent Hakim for helpful discussions during the course of this work and comments on the manuscript. S.Y.S. acknowledges support from the 1R01GM107103 grant from NIGMS. I.G.K acknowledges support from the CS&E program of the NSF.


H Turlier

D A Fedosov

B Audoly

T Auth

N S Gov

C Sykes

J.-F Joanny

G Gompper

T Betz

Equilibrium physics breakdown reveals the active nature of red blood cell flickering

Red blood cells, or erythrocytes, are seen to flicker under optical microscopy, a phenomenon initially described as thermal fluctuations of the cell membrane. But recent studies have suggested the involvement of non-equilibrium processes, without definitively ruling out equilibrium interpretations. Using active and passive microrheology to directly compare the membrane response and fluctuations on single erythrocytes, we report here a violation of the fluctuation-dissipation relation, which is a direct demonstration of the non-equilibrium nature of flickering. With an analytical model of the composite erythrocyte membrane and realistic stochastic simulations, we show that several molecular mechanisms may explain the active fluctuations, and we predict their kinetics. We demonstrate that tangential metabolic activity in the network formed by spectrin, can generate curvature-mediated active membrane motions. We also show that other active membrane processes represented by direct normal force dipoles may explain the observed membrane activity. Our findings provide solid experimental and theoretical frameworks for future investigations of the origin and function of active motion in cells.

Main Text

One of the first observations of red blood cells (RBC) by optical microscopy during the nineteenth century was the 'vibrating' fluctuation of its plasma membrane [START_REF] Browicz | Further observation of motion phenomena on red blood cells in pathological states[END_REF] . Important biological functions of this flickering have been suggested, such as decreased interactions between cells 2 or facilitated diffusion of membrane proteins [START_REF] Lin | Dynamics of Pinned Membranes with Application to Protein Diffusion on the Surface of Red Blood Cells[END_REF] . Today, flickering analysis is envisaged as a promising tool for high-rate malaria diagnosis [START_REF] Cojoc | Toward fast malaria detection by secondary speckle sensing microscopy[END_REF] . However, although equilibrium fluctuation models [START_REF] Strey | Measurement of erythrocyte membrane elasticity by flicker eigenmode decomposition[END_REF][START_REF] Betz | ATP-dependent mechanics of red blood cells[END_REF][START_REF] Yoon | Flickering Analysis of Erythrocyte Mechanical Properties: Dependence on Oxygenation Level, Cell Shape, and Hydration Level[END_REF] have been used to infer mechanical properties of RBCs since the seminal work of Brochard & Lennon [START_REF] Brochard | Frequency spectrum of the flicker phenomenon in erythrocytes[END_REF] , the fundamental question whether membrane fluctuations are driven by an active process or simply by thermal agitation remains controversial [START_REF] Evans | Fluctuations of the red blood cell membrane: relation to mechanical properties and lack of ATP dependence[END_REF][START_REF] Yoon | Red blood cell dynamics: from spontaneous fluctuations to non-linear response[END_REF][START_REF] Boss | Spatially-Resolved Eigenmode Decomposition of Red Blood Cells Membrane Fluctuations Questions the Role of ATP in Flickering[END_REF][START_REF] Puckeridge | Membrane flickering of the human erythrocyte: physical and chemical effectors[END_REF] .

Possible active metabolic origins were conjectured more than 60 years ago [START_REF] Blowers | Flicker phenomenon in human erythrocytes[END_REF] , but were readily questioned by the observation that flickering arises in reconstituted membranes as well. By comparing membrane fluctuations directly between healthy RBCs and cells depleted of ATP (adenosine triphosphate), recent experiments find a decrease of fluctuations upon ATP-depletion, hence suggesting that an active process may contribute to membrane undulations [START_REF] Betz | ATP-dependent mechanics of red blood cells[END_REF][START_REF] Tuvia | Cell membrane fluctuations are regulated by medium macroviscosity: evidence for a metabolic driving force[END_REF][START_REF] Park | Metabolic remodeling of the human red blood cell membrane[END_REF][START_REF] Rodríguez-García | Direct Cytoskeleton Forces Cause Membrane Softening in Red Blood Cells[END_REF] . However, the passive mechanical properties of the RBC membrane themselves depend strongly on ATP, as illustrated by the stiffening of the RBC membrane on starvation [START_REF] Weed | Metabolic dependence of red cell deformability[END_REF][START_REF] Chabanel | Increased resistance to membrane deformation of shape-transformed human red blood cells[END_REF][START_REF] Discher | Molecular maps of red cell deformation: hidden elasticity and in situ connectivity[END_REF][START_REF] Betz | ATP-dependent mechanics of red blood cells[END_REF][START_REF] Park | Metabolic remodeling of the human red blood cell membrane[END_REF] , that correlates with a sudden transition from discocyte to spiculated echinocyte shapes [START_REF] Nakao | Adenosine Triphosphate and Maintenance of Shape of the Human Red Cells[END_REF][START_REF] Sheetz | On the mechanism of ATP-induced shape changes in human erythrocyte membranes. I. The role of the spectrin complex[END_REF] . Therefore a direct and compelling explanation for the decrease of fluctuations observed in ATP-depleted cells is the stiffening of the membrane [START_REF] Manno | Modulation of Erythrocyte Membrane Mechanical Function by -Spectrin Phosphorylation and Dephosphorylation[END_REF][START_REF] Manno | Modulation of Erythrocyte Membrane Mechanical Function by Protein 4.1 Phosphorylation[END_REF] rather than the loss of putative active (that is, non-equilibrium) fluctuations. To make a different assessment of the involvement of metabolic activity in flickering, an attractive approach was used by Tuvia et al [START_REF] Tuvia | Cell membrane fluctuations are regulated by medium macroviscosity: evidence for a metabolic driving force[END_REF] : at equilibrium, the mean-square fluctuations amplitude is a thermodynamic variable, which, by definition, cannot depend on the dynamics, and in particular on medium viscosity. Tuvia et al [START_REF] Tuvia | Cell membrane fluctuations are regulated by medium macroviscosity: evidence for a metabolic driving force[END_REF] reported a viscosity dependence of the mean-square fluctuation amplitude in a frequency range of 0.3-30Hz on changing the suspending medium and therefore drew the conclusion of a metabolic origin of flickering. However, recent analysis of fluctuations over a larger bandwidth (0.1-1000Hz) shows no detectable dependence of fluctuation amplitudes on medium viscosity [START_REF] Betz | ATP-dependent mechanics of red blood cells[END_REF] , as also confirmed here (Supplementary Figure S1 and SI Appendix C), hence making the conclusions of Tuvia et al equivocal. Recent publications [START_REF] Park | Metabolic remodeling of the human red blood cell membrane[END_REF][START_REF] Monzel | Measuring fast stochastic displacements of bio-membranes with dynamic optical displacement spectroscopy[END_REF] reported non-Gaussian features in the RBC membrane fluctuations distribution and associated them to an active mechanism. However, an equilibrium fluctuations distribution is also likely to exhibit non-Gaussian tails, as soon as anharmonicities are present in the membrane [START_REF] Van Kampen | Fluctuation Phenomena in Solids[END_REF][START_REF] Schmidt | Signature of a Nonharmonic Potential as Revealed from a Consistent Shape and Fluctuation Analysis of an Adherent Membrane[END_REF] . On the other hand, the fluctuations distribution of a nonequilibrium system can be perfectly Gaussian [START_REF] Van Kampen | Fluctuation Phenomena in Solids[END_REF][START_REF] Ben-Isaac | Effective Temperature of Red-Blood-Cell Membrane Fluctuations[END_REF] . As a consequence, the observation of Gaussian fluctuation distributions [START_REF] Betz | ATP-dependent mechanics of red blood cells[END_REF][START_REF] Evans | Fluctuations of the red blood cell membrane: relation to mechanical properties and lack of ATP dependence[END_REF][START_REF] Boss | Spatially-Resolved Eigenmode Decomposition of Red Blood Cells Membrane Fluctuations Questions the Role of ATP in Flickering[END_REF][START_REF] Puckeridge | Membrane flickering of the human erythrocyte: physical and chemical effectors[END_REF] does not necessarily imply that RBC flickering is purely thermal. ATP-dependent temporal correlations, consistent with an active component, were recently identified [START_REF] Rodríguez-García | Direct Cytoskeleton Forces Cause Membrane Softening in Red Blood Cells[END_REF] , however their interpretation is model-dependent and does not provide unequivocal evidence for active forces.

An unquestionable hallmark of equilibrium systems is alternatively found in the fluctuationdissipation theorem (FDT), which relates directly mechanical properties to fluctuations in the linear regime [START_REF] Callen | Irreversibility and Generalized Noise[END_REF] . By measuring the response function and spontaneous fluctuations for the same system and for the same experimental conditions, it is possible to assert directly the presence of " 2 non-equilibrium processes, as recently demonstrated in vitro for hair bundles [START_REF] Martin | Comparison of a hair bundle's spontaneous oscillations with its response to mechanical stimulation reveals the underlying active process[END_REF] , active gels [START_REF] Mizuno | Nonequilibrium Mechanics of Active Cytoskeletal Networks[END_REF] and suspended bone cells [START_REF] Mizuno | High-resolution probing of cellular force transmission[END_REF] . A series of techniques allow the measurement of cellular mechanics [START_REF] Hoffman | The consensus mechanics of cultured mammalian cells[END_REF] , and the RBC membrane response was measured previously using micropipette aspiration [START_REF] Discher | Simulations of the erythrocyte cytoskeleton at large deformation. II. Micropipette aspiration[END_REF] , electric fields [START_REF] Engelhardt | Viscoelastic properties of erythrocyte membranes in high-frequency electric fields[END_REF] , or optical tweezers [START_REF] Yoon | Red blood cell dynamics: from spontaneous fluctuations to non-linear response[END_REF][START_REF] Lenormand | Direct Measurement of the Area Expansion and Shear Moduli of the Human Red Blood Cell Membrane Skeleton[END_REF] , but never compared directly to the fluctuations in the same cell. Here, we use active and passive microrheology to measure directly single erythrocyte membrane response and fluctuations for the same cell and over timescales from milliseconds up to tens of seconds, a range which remains unattained so far [START_REF] Yoon | Red blood cell dynamics: from spontaneous fluctuations to non-linear response[END_REF] .

The experimental approach relies on a time-shared multiplexed optical tweezer, generating four independent traps that control the position of four streptavidin-functionalized polystyrene beads attached to a biotinylated RBC at opposite sites along its rim (figure 1, Supplementary figure S2). The beads allow us to apply well-defined forces on the cell membrane, but they can also track local membrane shape deformation with sub nm-accuracy and avoid problems related to the convolution between the cell edge and the detection laser when directly tracking the membrane [START_REF] Betz | ATP-dependent mechanics of red blood cells[END_REF][START_REF] Betz | Time resolved membrane fluctuation spectroscopy[END_REF] . Using this setup, we can detect both the free fluctuations and mechanical response of the cell membrane with sub-millisecond temporal resolution, thanks to an interferometric detection system (figure 1a). Three beads (figure 1a inset) serve as handles to hold the cell 20µm above the coverslip surface in the observation chamber. The fourth bead is used as a probe (figure 1a inset) to detect free fluctuations, or to measure the mechanical response of the membrane under external driving.

In a first experimental step, the free membrane fluctuations are determined by reducing the laser trap power on the probe bead from 5mW to 0.1mW. This value is insufficient for stable trapping but it allows laser-tracking of the bead position. We checked carefully that the presence of this weak trap does not restrict the fluctuations (Supplementary figure S3). Therefore, the bead motion reflects directly membrane fluctuations (figure 1b), which are characterized by their power-spectrum-density (PSD) expressed as function of the frequency f. Furthermore, the influence of the bead on the PSD remains furthermore lower than measurement error in the frequency range probed here when compared to the free membrane (Supplementary figure S4).

Immediately after PSD measurement, the mechanical response is determined for the same cell in a second experiment, where the three handle traps remain unchanged but the average laser power on the probe bead is increased to 5mW. Sinusoidal oscillations of the probe trap at a driving frequency f (figure 1c, blue) displace it over a 200nm amplitude (figure 1c, green). The complex response function " is calculated as the ratio between the Fourier transforms of the position and force at the driving frequency f. The real and imaginary parts of the response function, and , correspond respectively to the elastic and the dissipative components of the linear response of the RBC membrane to the applied force. To obtain the frequencydependent response function, we repeat this active rheology measurement for driving frequencies ranging from 0.2Hz to 1kHz. We check the linearity of the response by investigating the forcedeformation relation as well as the absence of higher harmonics in the bead displacement and response spectra for each driving force (Supplementary figure S5). For both active and passive

C( f ) χ f ( ) x f ( ) F f ( ) ʹ χ ʹʹ χ " 3
microrheology measurements, the handle beads move around their equilibrium position. Hence, the probe bead movement combines the mechanics of the cell and the handle traps. Fluctuations and mechanical response function are measured under the exact same experimental conditions, which enables us to directly test the validity of the FDT 28 , which reads (1) where kBT is the thermal energy, with kB being Boltzmann's constant.

We check for violation of the FDT in fresh RBCs, as well as in cells that were ATP-depleted by incubation in glucose-free medium (GFM) at 37°C for 16h (partial ATP-depletion) and >48h (full ATP-depletion). For fully ATP-depleted cells (nine experiments on seven cells) figure 2a displays the dissipative part of the response function and the predictions from the FDT based on the PSD measurements (equation ( 1)). The observed collapse of both data sets demonstrates the validity of the FDT for a purely passive RBC, as expected. In strong contrast, fresh RBCs (21 experiments on 9 cells) show a clear violation of the FDT in the low frequency regime f<10Hz (figure 2b). This result is a direct and conclusive demonstration that an active mechanical process contributes to the flickering of the RBC membrane at timescales above 100ms. Conversely, the collapse of free fluctuations and response at high frequencies (f>10Hz) indicates the validity of the FDT for rapid fluctuations. Therefore, thermal agitation is the dominant source of fluctuations at short timescales (<100ms). We conclude that fluctuation measurements can be reliably used to infer mechanical properties of the RBC membrane only in this high frequency limit. The validity of previous studies, which used equilibrium theories to extract mechanical information on the RBC membrane, has therefore to be critically checked. The comparison of the dissipative response for the three conditions (zero, partial and full ATP depletion) in figure 2c confirms the strong metabolic dependence of passive membrane properties noted previously [START_REF] Betz | ATP-dependent mechanics of red blood cells[END_REF][START_REF] Puckeridge | Membrane flickering of the human erythrocyte: physical and chemical effectors[END_REF][START_REF] Chabanel | Increased resistance to membrane deformation of shape-transformed human red blood cells[END_REF][START_REF] Manno | Modulation of Erythrocyte Membrane Mechanical Function by -Spectrin Phosphorylation and Dephosphorylation[END_REF][START_REF] Manno | Modulation of Erythrocyte Membrane Mechanical Function by Protein 4.1 Phosphorylation[END_REF] .

To measure directly the metabolic energy input which drives active fluctuations, we use the concept of effective energy [START_REF] Betz | ATP-dependent mechanics of red blood cells[END_REF] , defined as the thermal energy at a frequency f that would be necessary to recover the FDT in the non-equilibrium system . The metabolic energy input driving active flickering is visualized directly above in figure 2d.

Interestingly, ATP-depletion of RBCs is not fully complete after 16h of incubation in glucosefree medium [START_REF] Weed | Metabolic dependence of red cell deformability[END_REF][START_REF] Chu | Identification of cytoskeletal elements enclosing the ATP pools that fuel human red blood cell membrane cation pumps[END_REF] . The effective energy coincides with the thermal energy at high frequency for all conditions, but it increases up to 23kBT at low frequencies for fresh RBCs (ref [START_REF] Ben-Isaac | Effective Temperature of Red-Blood-Cell Membrane Fluctuations[END_REF] ). Integrating this deviation over the full spectrum gives an estimate of the power Pf ~0.025 fW used for driving the membrane fluctuations. This is only a very small fraction of the total metabolic energy turnover in the RBC (ref [START_REF] Bandmann | Clinical Physiology: Microcalorimetric Measurements of Heat Production in Whole Blood and Blood Cells of Normal Persons[END_REF] ) of PRBC~10 fW.

Diverse metabolic processes have been identified in the RBC membrane. We use here analytical and simulation models to explore their relevance in explaining the observed active fluctuations.

C f ( ) = 2k B T ʹʹ χ f ( ) 2π f ( ) E eff f
( )

E eff f ( ) = C f ( ) π f ʹʹ χ f ( ) 1 k B T " 4
The membrane consists of a network of spectrin filaments anchored to a lipid bilayer, as illustrated in figure 3a. Biochemists identified several phosphorylation sites directly in the network: at the transmembrane protein complexes 4.1R (refs [START_REF] Manno | Modulation of Erythrocyte Membrane Mechanical Function by Protein 4.1 Phosphorylation[END_REF][START_REF] Subrahmanyam | Phosphorylation of protein 4.1 on tyrosine-418 modulates its function in vitro[END_REF][START_REF] Gauthier | Phosphorylation-Dependent Perturbations of the 4.1R-Associated Multiprotein Complex of the Erythrocyte Membrane[END_REF] ) and ankyrin [START_REF] Lu | Phosphorylation of ankyrin decreases its affinity for spectrin tetramer[END_REF] , which anchor the network in the bilayer, but also at the spectrin chains [START_REF] Manno | Modulation of Erythrocyte Membrane Mechanical Function by -Spectrin Phosphorylation and Dephosphorylation[END_REF] . These phosphorylations have been associated with a decreased mechanical strength of the membrane [START_REF] Manno | Modulation of erythrocyte membrane mechanical function by beta-spectrin phosphorylation and dephosphorylation[END_REF][START_REF] Manno | Modulation of Erythrocyte Membrane Mechanical Function by Protein 4.1 Phosphorylation[END_REF][START_REF] Gov | Red Blood Cell Membrane Fluctuations and Shape Controlled by ATP-Induced Cytoskeletal Defects[END_REF][START_REF] Betz | ATP-dependent mechanics of red blood cells[END_REF][START_REF] Park | Metabolic remodeling of the human red blood cell membrane[END_REF] , but the underlying molecular mechanism remains unclear. Instead of considering one specific process, we argue here that any relevant spectrin phosphorylation should lead to a local decrease of the network shear modulus, as this single parameter characterizes its mechanical properties. Our mean-field approach encompasses therefore the filament detachment hypothesis [START_REF] Gov | Red Blood Cell Membrane Fluctuations and Shape Controlled by ATP-Induced Cytoskeletal Defects[END_REF][START_REF] Gov | Active elastic network: Cytoskeleton of the red blood cell[END_REF] , but also other potential spectrin metabolic processes. More precisely, we assume a simple two-state dynamics, where phosphorylation sites can switch between active and inactive states with transition rates and respectively, defining hence a mean activity . The network shear modulus fluctuates accordingly around a mean value , that decreases with phosphorylation. Modeling the RBC membrane mechanics remains a challenging task, as multiple length and timescales are involved. Most previous composite models are hence numerical [START_REF] Boal | Dual network model for red blood cell membranes[END_REF][START_REF] Peng | Multiscale modelling of erythrocytes in Stokes flow[END_REF][START_REF] Li | Two-Component Coarse-Grained Molecular-Dynamics Model for the Human Erythrocyte Membrane[END_REF][START_REF] Peng | Lipid bilayer and cytoskeletal interactions in a red blood cell[END_REF] or limited to flat geometries [START_REF] Gov | Red Blood Cell Membrane Fluctuations and Shape Controlled by ATP-Induced Cytoskeletal Defects[END_REF][START_REF] Gov | Active elastic network: Cytoskeleton of the red blood cell[END_REF][START_REF] Fournier | Fluctuation spectrum of fluid membranes coupled to an elastic meshwork: jump of the effective surface tension at the mesh size[END_REF][START_REF] Dubus | A Gaussian model for the membrane of red blood cells with cytoskeletal defects[END_REF][START_REF] Zhang | Cytoskeleton mediated effective elastic properties of model red blood cell membranes[END_REF][START_REF] Rodríguez-García | Direct Cytoskeleton Forces Cause Membrane Softening in Red Blood Cells[END_REF] . We derive here an active composite membrane model in spherical geometry, while considering explicitly the finite elasticity of the spectrin network and its sliding relative to the bilayer.

The spectrin cytoskeleton is idealized as a perfect triangular network of springs. In line with previous studies [START_REF] Discher | Simulations of the erythrocyte cytoskeleton at large deformation. II. Micropipette aspiration[END_REF][START_REF] Gov | Red Blood Cell Membrane Fluctuations and Shape Controlled by ATP-Induced Cytoskeletal Defects[END_REF][START_REF] Paulose | Fluctuating shells under pressure[END_REF] , we assume it to be prestressed (see figure 3a). The filaments are homogeneously prestretched by a factor " , which preserves the sixfold symmetry of the network [START_REF] Kantor | Entropic elasticity of tethered solids[END_REF] . We can hence homogenize the discrete filament energy into an effectively isotropic continuous elastic membrane that we describe in the regime of finite elasticity (see SI p13-16). Its energy, expressed as a function of the incremental deformations from the prestressed state, is Hookean and includes a prestress term. Integrated over the membrane surface , it reads

, where is the constant energy density due to prestress and is the non-linear incremental membrane strain tensor. The prestressed network mechanics is characterized by an effective spectrin tension " and an incremental shear modulus , which are both functions of the prestretch ratio , the phosphorylation activity and the stress-free shear modulus . For vanishing prestress ( ) the spectrin tension vanishes and (Supplementary figure S6). Extending a previous approach [START_REF] Rochal | Viscoelastic dynamics of spherical composite vesicles[END_REF] , we derive the elastic forces of the spectrin network and expand them into vector spherical harmonics, which allows

k a k i n a = k a k a + k i µ = µ 0 1-n a ( ) ξ > 1 A E = dA ε i + Se αα + M 1 2 e αβ ( ) 2 + e αβ 2 ⎛ ⎝ ⎜ ⎞ ⎠ ⎟ ⎡ ⎣ ⎢ ⎤ ⎦ ⎥ ∫ ε i e αβ S M ξ n a µ 0 ξ = 1 M = µ " 5
one to identify the three elementary modes of deformations: bending, stretching and shear. As usual in thin shell theory, bending and stretching deformations are linearly coupled for curved membranes, as illustrated in figure 3c for the harmonic l=3; by contrast, shear is decoupled here and is therefore ignored. The inextensible lipid bilayer is classically characterized by a bending modulus " and a surface tension " , and its elastic forces can be derived from the Canham-Helfrich energy [START_REF] Canham | The minimum energy of bending as a possible explanation of the biconcave shape of the human red blood cell[END_REF][START_REF] Milner | Dynamical fluctuations of droplet microemulsions and vesicles[END_REF] , where C measure the mean local curvature.

The dissipative forces on the membrane have several sources. The inner and outer fluids have different viscosities and exert viscous stresses [START_REF] Rochal | Viscoelastic dynamics of spherical composite vesicles[END_REF] , while the relative movement of the spectrin with respect to the lipid bilayer leads to a tangential friction force. Surprisingly this friction force is found to largely overcome viscous damping for large wavelengths; it therefore controls the tangential relaxation dynamics (Supplementary figure S6). A Langevin membrane dynamics is obtained by combining elastic and dissipative forces with thermal noise. This leads, for a given spherical harmonic (l,m), to an effective dynamics coupling bending " and stretching " deformations (2) where is the relaxation frequency matrix for the membrane deformation .

The thermal noise naturally contributes to both bending and stretching fluctuations, while metabolic shear modulus fluctuations generates, in an isolated network, direct active forces both in the tangential and normal directions. When the network is coupled to the bilayer, however, we predict the local bilayer incompressibility to systematically cancel the direct active normal force and we end up with an additional active source of noise for the stretching modes only, . This is therefore in contrast with previous phenomenological active models based on the existence of a direct active normal force [START_REF] Rodríguez-García | Direct Cytoskeleton Forces Cause Membrane Softening in Red Blood Cells[END_REF][START_REF] Gov | Red Blood Cell Membrane Fluctuations and Shape Controlled by ATP-Induced Cytoskeletal Defects[END_REF] . Tangential metabolic noise nevertheless perturbs normal membrane motions via the curvaturemediated coupling between bending and stretching modes (non-diagonal components of the matrix ). This effect is dampened by the friction between the bilayer and the network, and depends therefore critically on the tangentially fluid character of the bilayer. The membrane shape fluctuation spectrum for each spherical harmonic (l, m) can be calculated in Fourier space as a sum of the dissipative part of the response function and this active contribution

κ σ 0 Η = κ 2 C 2 +σ 0 ⎡ ⎣ ⎤ ⎦ ∫ dA u lm (r ) u lm (1)
dU lm t ( ) dt = -W lm ⋅ U lm t ( ) + Z lm TH t ( ) + Z lm A t ( ) W lm U lm = u lm (r ) , u lm (1)
( )

Z lm TH t
( )

µ lm t ( ) Z lm A t ( ) = 0, Ω lm A µ lm t ( ) ( ) W lm ʹʹ χ lm f ( ) " 6 (3)
where captures the complex mode-and frequency-dependent propagation of tangential active noise into membrane shape fluctuations. The deviation from FDT is characterized by an active timescale , which controls the frequency, below which a departure from equilibrium is observed for a specific mode (l,m) : for active timescales much lower than the typical membrane relaxation frequency " , the membrane does not have time to move and the metabolic activity is ineffective, whereas for larger active timescales the FDT violation increases and finally saturates at a value set by the maximum energy input in the system (see Supplementary figure 6). We calculate the fluctuation spectrum as a sum over the relevant modes and predict a violation of FDT for non-zero metabolic activity and for a prestressed spectrin network , in agreement with previous hypotheses [START_REF] Gov | Red Blood Cell Membrane Fluctuations and Shape Controlled by ATP-Induced Cytoskeletal Defects[END_REF][START_REF] Gov | Active elastic network: Cytoskeleton of the red blood cell[END_REF] . In contrast, we find here that the violation vanishes for flat membranes ( ). The membrane curvature is therefore a key feature for propagating tangential spectrin activity into normal shape motions, and we predict an increase of active fluctuations with membrane curvature, which may explain the higher fluctuations recently measured along the rim of the RBC, compared to its central region [START_REF] Park | Metabolic remodeling of the human red blood cell membrane[END_REF][START_REF] Monzel | Measuring fast stochastic displacements of bio-membranes with dynamic optical displacement spectroscopy[END_REF] . To fit the experimental data for rheological tests, we use typical mechanical parameters for the RBC, that have been either experimentally measured or estimated in simulations, and four fitting parameters (see Supplementary Table ): the active rates and , leading to a mean fraction of active sites , the network prestretch , controlling both passive and active spectrin mechanics, and the lipid bilayer excess area which determines the bilayer tension [START_REF] Seifert | The concept of effective tension for fluctuating vesicles[END_REF][START_REF] Loubet | Effective tension and fluctuations in active membranes[END_REF] (see SI p29-31). Our model shows excellent agreement with the passive response of the membrane and can reproduce the active fluctuation spectrum (figure 3d). A broader distribution of active timescales may refine this agreement and could arise from a Bell-like dependence of the transition rates on membrane extension [START_REF] Bihr | Nucleation of Ligand-Receptor Domains in Membrane Adhesion[END_REF] . Interestingly, we find a slightly negative lipid bilayer tension , which reveals a tangential compression of the bilayer, in agreement with a previous independent estimation [START_REF] Fournier | Fluctuation spectrum of fluid membranes coupled to an elastic meshwork: jump of the effective surface tension at the mesh size[END_REF] . Writing the force balance for the harmonic , which reduces to the Laplace's law " , we find that the elastic contribution arising from compression of the bilayer

C lm f ( ) = 2k B T 2π f ʹʹ χ lm f ( ) + 2 n a 1-n a ( ) τ a 1+ 2π f τ a ( ) N lm f ( ) 2 N lm f ( ) τ a = k a + k i ( ) -1 W lm C f ( ) (l ≥ 2, -l ≤ m≤ l ) n a 1-n a ( ) > 0 ξ > 1⇔ S> 0 R→ ∞ k a ≈ 1.12 s -1 k i ≈ 5.6 s -1 n a ≈ 17% ξ ≈ 1.04 Δ ≈ 1.9% σ 0 ≈ -6.510 -7 J.m -2 l = 0 ΔP = 2 S+σ 0 ( ) R " 7
is large enough to counterbalance intrinsically the spectrin tension S>0 in the membrane, even in the absence of pressure jump (figure 3b). Our analytical model is therefore self-consistent and shows that the composite nature of the RBC membrane allows it to maintain an intrinsic prestress, a necessary feature for the emergence of spectrin-based active fluctuations.

To investigate the role of other potential metabolic processes in flickering, such as ion-pumps [START_REF] Prost | Shape fluctuations of active membranes[END_REF] or phospholipid translocase [START_REF] Graham | Flippases and vesicle-mediated protein transport[END_REF] activities, we employ a cell-scale in-silico model, which allows us to closely mimic the experimental situation, including the biconcave RBC geometry, the membrane curvature and stretching elasticity, the attached beads, and the optical tweezers. In our simulations, the hydrodynamics of the internal and external fluids and thermal fluctuations are taken into account by a particle-based mesoscale simulation technique, the dissipative particle dynamics method [START_REF] Hoogerbrugge | Simulating Microscopic Hydrodynamic Phenomena with Dissipative Particle Dynamics[END_REF][START_REF] Español | Statistical Mechanics of Dissipative Particle Dynamics[END_REF] . Our RBC model (figure 4a) is based on a two-dimensional triangulation of the composite membrane [START_REF] Fedosov | A Multiscale Red Blood Cell Model with Accurate Mechanics, Rheology, and Dynamics[END_REF] , which corresponds to a network of viscoelastic springs with the inclusion of membrane bending energy and constraints for the conservation of membrane area and volume (see Supplementary Note, Section 1). The spring elements mimic the stretching elasticity of the spectrin network, and simultaneously capture viscous dissipation within the membrane. The bending energy describes bending resistance of the lipid bilayer, while the area and volume conservation constraints mimic the incompressibility of the lipid bilayer and cytosol, respectively. Modeling the RBC membrane at the spectrin level is still very challenging for realistic physical conditions (i.e., appropriate ambient temperature, membrane elasticity, and long times). To obtain an accurate averaging of fluctuations at the time scales of several minutes, which is much longer than a typical RBC relaxation time, we need to employ a coarse-grained model of the membrane, restricting its triangulation to about a thousand vertices. This coarsegrained description for a RBC has proved to reproduce quantitatively and accurately the viscoelasticity of the RBC membrane as well as deformation of RBCs due external forces (laser tweezers) and under flow [START_REF] Fedosov | A Multiscale Red Blood Cell Model with Accurate Mechanics, Rheology, and Dynamics[END_REF][START_REF] Fedosov | Multiscale modeling of blood flow: from single cells to blood rheology[END_REF] .

To simulate active processes, we furthermore add random active forces acting normally on membrane vertices (figure 4c). For ion pumps or other membrane protein activity the force exerted by an active center has to be counteracted by an exact opposite force in the surrounding fluid. The force density field of active centers have therefore a zero monopole moment, and the first contribution expected is a force dipole [START_REF] Manneville | Active membrane fluctuations studied by micropipet aspiration[END_REF] .

As in experiments, the three handle beads in simulations are held by a harmonic potential mimicking optical tweezers, while the probe bead is moved sinusoidally to determine the mechanical response function (figure 4a-b). The free fluctuations of the probe bead are measured in a separate simulation. The simulations reproduce accurately the experimental response function (figure 4d), and lead to the estimation of mechanical properties of RBCs (Supplementary Table ). The found mechanical properties are fully consistent with previous experimental estimations. Simulations reproduce the experimental violation of the FDT very precisely with the force amplitude for active dipoles and with the transition rates and between active and inactive states (figure 4c). The order of

F 0 ≈ 650 fN k a ≈ 0.42 s -1 k i ≈ 1.3 s -1 " 8
magnitude of the dipole force amplitude is consistent with simple estimations proposed for ionpumps [START_REF] Manneville | Active membrane fluctuations studied by micropipet aspiration[END_REF] . The timescales are larger than those predicted analytically with the spectrin activity hypothesis, a difference, which may be explained by the large bilayer friction, dampening the propagation of tangential spectrin activity into normal membrane shape fluctuations. When active processes are turned off, the simulated membrane fluctuations correspond to the FDT predictions (figure 4c, dashed lines) as expected. The resulting stochastic process for the membrane activity leads to an average fraction of active sites (i.e. those with the force applied at any moment of time) on a RBC of about 12%. Thus, our RBC simulations show that simple uncorrelated active force dipoles can reproduce experimental results with excellent accuracy. The agreement between experiments and simulation result, using force dipoles, indicates that the spectrin network, as suggested by our and other theoretical models, may not be the only source of activity in the RBC membrane. In particular, the activity of ion pumps in the membrane might offer a compelling explanation for the active flickering we measure here directly [START_REF] Prost | Shape fluctuations of active membranes[END_REF][START_REF] Manneville | Active membrane fluctuations studied by micropipet aspiration[END_REF] . Simulations also allow efficient testing of hypothetical active processes, such as local force monopoles, that may conserve momentum only globally. Such active force centers fit the data equally well with the same active rates as force dipoles but with ten fold lower force amplitudes , in agreement with previous estimations [START_REF] Gov | Active elastic network: Cytoskeleton of the red blood cell[END_REF] (Supplementary figure S7). Therefore, our simulations show that several, non-excluding, underlying mechanisms of uncorrelated active forces may explain membrane activity and suggest that further experimental characterizations are required to clearly identify the molecular mechanisms involved in flickering. They provide clear kinetic and amplitude predictions to guide future experiments. ). When activity is turned off (passive), the FDT is recovered as in fully ATP-depleted RBCs (Figure 2a). Shaded areas denote s.d.

F 0 ≈ 50 fN " 9

C( f )

E eff f ( ) = C f ( ) π f ʹʹ χ f ( ) " 11 b d c

Methods:

Experimental preparation and effects on tension: Biotinylated red blood cells are prepared as described previously [START_REF] Pierrat | Enforced Detachment of Red Blood Cells Adhering to Surfaces: Statics and Dynamics[END_REF] and mixed with streptavidin coated beads (3.28 µm diameter). 4 beads are attached to the RBC in an opposing geometry (Figure 1) by trapping and manipulating the beads and the RBC using the multiplexed optical tweezers. The attachment area of the bead on the membrane was visually estimated to be about 1µm in diameter (Supplementary Methods, Supplementary Figure S2), which corresponds to about 2% of the RBC membrane. The resulting increase of tension can be estimated to be only about 0.2%, and will hence not have a significant impact on the membrane fluctuations.

Bead and membrane motions: In the experiments we use the probe bead to determine the free fluctuations of the membrane, but also to apply a given force to the cell. The handle beads are exclusively used to keep the system in a stable 3D situation 20µm above the substrate. To test to which extend the probe bead fluctuations actually represent the membrane motion we also prepared sample cells with only three handle beads and measured the membrane fluctuations on the position where the fourth bead usually is attached using the methods reported previously [START_REF] Betz | ATP-dependent mechanics of red blood cells[END_REF] . Supplementary Figure S4 shows an overlay of the free probe-bead fluctuations and the free membrane. Within the error the two remain typically equal, and only at the low frequencies there is a systematic, but not significant difference between the probe-bead and the membrane. Interestingly, the membrane fluctuations in this 3D situation are a factor 5 higher than the fluctuations measured on the sedimented RBC (Supplementary Figure S4). Hence, when modeling the RBC membrane fluctuations, surface effects such as gravity mediated interactions as well as friction need to be taken into account for a correct analysis. To this end, our results can explain the previously found overestimation of RBC viscosities in the fluctuation analysis of sedimented cells [START_REF] Betz | ATP-dependent mechanics of red blood cells[END_REF] .

Effect of the trapping laser on measurements: To verify that the weak laser does not suppress any fluctuation, we compare the PSD of a free bead to the PSD of a bead attached to the RBC in the described experimental situation, using 0.1 mW in both cases. The resulting plot (Supplementary Figure S3) shows that the motion of the bead is indeed limited by the cell, as the free bead fluctuations are larger over the full frequency range.

Furthermore, we checked that the response function measured by the active microrheology does not simply reflect the trap stiffness of the handle beads. This was verified by comparing the PSD of the handle beads with the PSD derived from the response function measured (Supplementary Figure S8). The result shows that the response function of the cell is different to the response function imposed by the handle traps over the relevant frequency regime of 0.1 -100 Hz. Hence, the response function does indeed reflect the RBC mechanics.

Possible heating effects of the laser on the cell cannot be excluded. We can see hemolysis of cells if trapped with 100mW for more than 30 minutes. However, with the applied laser intensities we do not observe any effect during the experiments. As both, the active and passive measurements " 14 are done in the same situation (handle beads trapped with 5mW), the typical heating should be similar in both situations and not lead to a systematic error of the measurements.

Finally, we can use the diffusive motion of the handle beads in the traps during the passive experiments to ensure that the determined PSD of the probe bead is not simply reflecting the fluctuations of the whole cell. Supplementary Figure S9 (a) shows the PSD of the probe bead and the handle beads in the passive experiments (0.1 mW on the probe bead), while panel (b) gives the time-evolution of these motions. The fluctuations of the free probe bead are significantly larger which confirms that we can use it to infer information about the cell membrane fluctuations.

QPD and force calibration: For the correct measurement of both, the bead motion and the forces applied to the bead, the calibration of the QPD as well as the force calibration are key elements. These calibrations are done in situ on the beads attached to the cells. QPD calibration is obtained by a rapid scanning protocol that ensures no relevant bead motion during the scan (Supplementary Figure S10) and by an in situ force calibration (Supplementary Methods, Figure S9). As presented in Supplementary Figure S11, the presence of the cell does influence the position calibration by about 10%. This also leads to a change in trap stiffness in the situation where the cell is attached to the beads (Supplementary Methods for details). Both effects are accounted for in the in situ calibration. It should be noted that the force calibration is possible even in the absence of a priori knowledge over the friction in the cell-bead system. The in situ calibration strategy becomes also important when the cells change their shape, for example under ATP depletion, where the used procedure ensures a correct interpretation of the measurements.

As the QPD calibration depends on a precise interpretation of the QPD signal upon rapid laser motion, artifacts from the repositioning of the laser are suppressed by exclusively taking the data acquired during the time the laser rests at a given position (Supplementary Figure S12).

More details about the scanning routines, force calibration and the experimental procedure can be found in the Supplementary Methods.

" 15

FiguresFigure 1 :Figure 2 :

 12 Figures

Figure 3 :Figure 4 :

 34 Figure 3: Analytic composite model of the red blood cell membrane. (a) Schematic description of the membrane geometry and configurations: the membrane is composed of the lipid bilayer and a regular triangular network of spring-like spectrin filaments. It fluctuates around an intermediate state where the spectrin network is prestressed and the bilayer compressed on average. (b) A large enough excess area of lipid bilayer creates a negative bilayer tension σ0, as a result of its tangential compression. This bilayer elastic contribution can counterbalance the effective spectrin tension S, which results from filament prestretch. (c) The membrane deformation is a combination of bending and stretching modes, illustrated here for the spherical harmonic l=3. (d) The analytic model (plain lines) can reproduce the experimental fluctuations and response data (crosses). Shaded areas denote s.d.

  fluct. (active) Sim rheo. (active) Sim fluct. (passive) Sim rheo. (passive)

Acknowledgments

We thank V.L. Lorman and J. Prost for careful reading of the analytical model. T.B. was founded by ANR-11-JSV5-0002 and by the Deutsche Forschungsgemeinschaft (DFG), Cells-in-Motion Cluster of Excellence (EXC 1003 -CiM), University of Münster, Germany. H.T. acknowledges support from the Young Researcher Prize of the Bettencourt-Schueller Foundation. D.A.F. acknowledges funding by the Alexander von Humboldt Foundation. D.A.F., T.A., and G.G. gratefully acknowledge a CPU time grant by the Jülich Supercomputing Center. N.S.G. acknowledges the Institut Curie's Mayent-Rothschild Visiting Professor fund and Labex CelTisPhyBio for their support during the stay at the Institut Curie. N.S.G. is the incumbent of the Lee and William Abramowitz Professorial Chair of Biophysics, and thanks ISF Grant 580/12 for support.

Additional information

The authors declare no competing financial interests. Supplementary Information accompanies this paper on www.nature.com/naturephysics. Reprints and permissions information is available online at www.nature.com/reprints. Correspondence and requests for materials should be addressed to T.B.