Active Forces of Myosin Motors May Control Endovesiculation of Red Blood Cells

By using Monte Carlo (MC) simulations, we have shown that the active forces generated by (NMIIA) motor domains bound to F-actin may partially control the endovesiculation of the red blood cell (RBC) membrane. The myosin generated active forces favor pancake-like (torocyte) RBC endovesicles with a large flat central membrane region and a bulby periphery. We suggest that the myosin generated active forces acting on the RBC membrane in the direction perpendicular to the membrane surface towards the interior of the RBC may influence also other RBC shape transformations and the stability of different types of RBC shapes and should be therefore considered in the future theoretical studies of the RBC vesiculation and shape transformations.


Introduction
A biological membrane forms a physical boundary between the inner volume of a biological cell and the external medium, as well as, within the cell, between the lumens of intracellular organelles and cytosol. The main building block of the biological membrane lipid bilayer consists of two monolayers of phospholipid molecules. 1 In the lipid bilayer, other constituents like membrane proteins are intercalated. 2,3 The transmembrane proteins are the pining points for membrane skeleton or/and cytoskeleton. 4 Furthermore, some of the inclusions in the lipid bilayer, like proteins or nanoparticles, can impose local curvature of the membrane, [5][6][7][8][9] forcing the membrane to locally and/or globally adapt its shape. [10][11][12][13][14][15][16][17] The phase separation of membrane inclusions (nanodomains) is an important mechanism that may induce the local changes of membrane curvature and is therefore the driving force for transformations of the cell shape. [17][18][19] The spontaneous phase separation of membrane proteins, driven by the forces of actin polymerization and cell-substrate adhesion, was in a quantitative manner predicted for the first time by Veksler and Gov. 20 Besides the membrane inclusions' driven membrane shape changes, other mechanisms also determine the membrane cell shape. Among them, constant energy consuming forces are acting in the living cell. 17 In experiments, the cells or lipid bilayer vesicles (as model systems) may also be deformed by extracellular forces due to pressure differences, 21 stretching 22 and fluid flow. [23][24][25] An external force to the membrane surface can also be experimentally generated, for example by a cantilever of the atomic force microscope. 26 The cell membrane resistance to the shape changes depends on its mechanical properties 27 mainly governed by the composition of the membrane. The dynamic response of the cell membrane shape to the force 17 is important for different cell functions that range from adhesion and migration to division, differentiation, and cell death. 28 Therefore, the membrane mechanical properties are of great scientific interest. Consequently, new theoretical approaches for deeper understanding of the mechanisms that rule the cell's functions are constantly being developed.
In red blood cells (RBC), besides the lipid bilayer also the membrane skeleton (spectrin-F-actin network attached to the inner lipid layer) 29 plays an important role in cell shape determination and cell transformations. [30][31][32] It has been shown that equilibrium RBC shapes that correspond to the minimal membrane elastic energy, should take into account the local and non-local bending energy 27,33,34 and the elastic energy of the membrane skeleton. [30][31][32]35 The shear deformation of the membrane skeleton plays an essential role in the stability of echinocytic RBC shapes. 30,31 Additionally, cytoskeleton induced protrusions in lipid membranes has been theoretically explained by the interplay between the elasticity of the membrane lipid bilayer and the membrane skeleton -firstly for axisymmetric shapes 36,37 and later also for non-axisymmetric shapes. 38 Non-local bending energy 27,33,39,40 depends on the change in the relaxed areas of the two lipid layers of the membrane bilayer. 41 Decreasing the difference between the relaxed areas of the outer layer and the inner lipid layer, favors the concave local membrane shape, i.e. the inward bending of the membrane. 32,35,41 On the other hand, increasing the difference between the relaxed areas of the outer and inner membrane layers promotes the formation of convex local shape, i.e. the outward bending. 31,35,42 It was shown that exogenously added amphiphiles predominantly bound in the outer lipid layer induce the transformation of the discoid RBC into the spiculated echinocytic RBC, while amphiphiles predominantly bound in the inner lipid layer induce the transformation into invaginated stomatocytic shapes. 43 When RBCs are incubated with high sublytic concentrations of amphiphiles, the microexovesiculation starts. 42,44,45 The amphiphile induced RBC microexovesicles are highly depleted in the membrane skeleton 37 suggesting that a local disruption of the interactions between the membrane skeleton and the membrane bilayer occurred prior to micro-or nano-exovesiculation. 45 In RBCs, most types of amphiphiles induce spherical microexovesicles that are formed from sphere-like membrane buds and are free of the membrane skeleton. On the other hand, there are also anisotropic amphiphiles which induce the growing of tubular membrane buds and the release of stable tubular microexovesicles. 44,45 The tubular budding and vesiculation of the RBC membrane can be theoretically explained by deviatoric membrane properties due to the in-plane orientational ordering of the anisotropic membrane inclusions induced by intercalated anisotropic amphiphiles. 44,45 The deviatoric properties of the membrane 12,15,42,[45][46][47][48][49][50][51][52] may explain the experimentally ob-served tubular membrane protrusion without the application of the local force. 16,42,45,51,53,54 For isotropic membrane, the application of the local force is necessary to theoretically predict the tubular membrane protrusions. 16,54,55 It was reported 44,56 that some amphiphilic molecules can induce large membrane invaginations (i.e. stomatocyte RBC shape) as shown in Figure 1C and endovesicles in RBCs. By means of transmission electron microscopy (TEM) and confocal laser scanning microscopy using fluorescent markers, it was also observed that many stomatocytogenic amphiphiles (for example chlorpromazine hydrochloride) can induce in RBCs with large concentrations of amphiphiles small spherical endovesicles ( Figure 1D). 44,56,57 On the other hand, the stomatocytogenic detergent polyethyleneglycol dodecylether (C 12 E 8 ) induces large endovesicles with a unique ring-like toroidal shape joined with a central flat membrane segment, i.e. flat membrane structures with a bulby periphery called also torocyte endovesicles (Figure 1 A,B and E) 57 with the shape very similar to the shape of Golgi bodies. 58 It was proposed that the observed RBC torocyte endovesicles were formed in a process in which an initially stomatocytic RBC invagination loses volume while maintaining a large surface area. 57,59 The results of theoretical modeling indicated that torocyte RBC endovesicles can be mechanically stabilized by non-homogeneous lateral distribution of laterally mobile anisotropic membrane inclusions, like for example by anisotropic detergent-membrane component complexes. 58,59 It was further shown in 2019 that the mechanical stability of torocyte (pancake-like) closed membrane vesicle shapes can be explained also by the coupling of the curved isotropic membrane inclusions and active (cytoskeletal) forces. 55 Until recently, 60 it was believed that the active forces are absent in the mechanisms of the determination of the RBC shape and vesiculation, as discussed above. It has been shown in 2018 that nonmuscle myosin IIA (NMIIA) motor domains may generate tension in spectrin-F-actin also in a 2-dimensional RBC membrane skeleton and in this way partially control the RBC shape. 60 The role of NMIIA contractility in generating tension in the RBC network and partially controlling the RBC shape was thus in the past completely neglected, until recently. 60 The length of NMIIA filament is around 200 nm. 60 The influence of NMIIA bipolar filaments, associated with a 2-dimensional RBC membrane skeleton can be thus described by the local force acting on the RBC membrane in the direction to the interior of the RBC. In the present paper, we shall study the interplay between the laterally mobile membrane inclusions and the inward-oriented local forces of the NMIIA bipolar filaments in the formation of torocyte endovesicles (flat membrane structures with a bulby periphery). The present study was motivated by the results presented in references 44, 56, 57, 59, 60, but it is not entirely/directly connected to the experimental and theoretical results presented in these articles. of simplicity, we shall consider the closed lipid bilayer membrane with one type of inclusions only. In our model, the single inclusion can induce local membrane bending due to its negative intrinsic curvature and also because of inward-oriented active forces. This means that we have joined the effect of the intrinsic curvature of the membrane inclusions (nanodomains) and the effect of the local active forces of NMIIA bipolar filament domains in a single type of membrane inclusions (nanodomains).

Monte Carlo Simulations
Monte Carlo (MC) triangulated mesh was used to numerically model and investigate the vesicle shape and the lateral distribution of membrane inclusions by means of computer simulation. 55 Phospholipid bilayer membranes can be treated due to their small thickness in the first approximation as a two dimensional surface, allowing the continuum approach in the theoretical description of membrane surfaces. 61 In the model we discretize the membrane into patches consisting of many molecules ( Figure  2). A single patch is represented by a vertex in a triangulated surface model. The main model parameter that defines mechanical bilayer properties is bending stiffness.
The vesicle is represented by a set of N vertices that are linked by bonds of flexible length d to form a closed, randomly triangulated, self-avoiding network. 62,63 The lengths of the tethers can vary between a minimal (d min ) and a maximal (d max ) value. The self-avoidance of the network can be implemented by ensuring that no vertex can penetrate through the triangular network. The maximal Figure 1: TEM and SEM micrographs showing stomatocytic human RBC with plate-like torocytic invaginations/endovesicles with a bulby periphery (panels A and B) formed in a process in which initially stomatocytic RBC invaginations (panel C) loose volume while maintaining a large surface area. RBCs were incubated by amphiphilic molecules octaethyleneglycol dodecylether (C 12 E 8 ). Panels D and E show confocal laser scanning microscopy of RBCs incubated with amphiphilic molecules chlorpromazine inducing small spherical endovesicles (panel D) and with C 12 E 8 inducing toroidal endovesicles (panel E) as presented also in panels A and B. Adapted from. 57 possible random displacement of the vertex in a single step (s) should be small enough so that the fourth vertex cannot move through the plane of the other three to the minimal allowed distance, d min , from the three vertices. In our scenario, we use s=0.15 d min and d max =1.7 d min . For details about the expressions to calculate self-avoidance constraint d max , see. 64 The initial state of triangulated surface is a pentagonal bipyramid with all the edges divided into equilateral bonds so that the network is composed of 3(N-2) bonds forming 2(N-2) triangles. N c randomly selected vertices are given non-zero isotropic intrinsic curvature of c 0 , thus they become the model of the membrane inclusion. The rest of the vertices have zero intrinsic curvature. Positive curvature means that the membrane will locally bulge towards the exterior, negative curvature will force the membrane to bulge towards the interior compartment of the vesicle.
The system is developed into the thermal equilibrium state. The evolution of the system is measured in Monte Carlo sweeps (mcs). One mcs consists of individual attempts to displace each of the N vertices by a random increment in the sphere with radius s -the action we will refer to as a vertex move. Membrane fluidity is maintained by flipping bonds within the triangulated network. In each mcs, the vertex move attempts are followed by a 3N attempts to flip a randomly chosen bond. A single bond flip involves the four vertices of two neighboring triangles. The tether between the two vertices is cut and re-established between the other two, previously unconnected vertices (a detailed description was published elsewhere). 65 Each individual Monte Carlo step (vertex move or bond flip) is accepted with probability according to Metropolis Hastings algorithm, based on free energy change due to the Monte Carlo step. The energy W in simulation consists of three parts: (1) where W b is the bending energy of the membrane, W d is the energy of the direct interaction between vertices with intrinsic curvature and W F is the energy due to myosin forces acting on the membrane.
For the bending energy W b of the membrane, we use the standard Helfrich expression for a tensionless membrane with a term that represents intrinsic curvature. 66 The membrane keeps fixed topology, thus the contribution of the Gaussian curvature to the change of bending energy is cancelled out. (2) where κ is the bending stiffness of the membrane, c 1 , c 2 and c 0 are the two principal curvatures and the intrinsic curvature of the vesicle membrane at the point under consideration. Note that only points where inclusions are located have non-zero intrinsic curvature c 0 , whereas all other points on the membrane have intrinsic curvature c 0 set to zero. The integration is performed over the membrane area A.
For modelling attraction force between the vertices with intrinsic curvature energy term: 55 (3) where w is a direct interaction constant. The energy is summed over all inclusion pairs with their in-plane distance r ij , r 0 is the range of direct interaction and H is a Heaviside step function with (r 0 -r ij ) being function's argument.
The energy contribution of the local protrusive forces due to the myosin motor domains/inclusions: 55 (4) where F is the magnitude of the force, the sum runs over all proteins, is the outwards facing normal to the membrane at the location of the vertex with the inclusion i and x i is the position vector of the vertex with the inclusion i. Figure 3 shows the MC simulations of RBC shape transformations induced by laterally mobile membrane inclusions (nanodomains) with negative intrinsic curvature. It can be seen in Figure 3 that the non-homogeneous lateral redistribution of the mobile inclusions with negative intrinsic shape causes the RBC membrane to locally change the curvature resulting in the global transformation of the RBC shape. The MC predicted RBC shapes presented in Figure 3 depend on the inclusion concentration, the inclusion intrinsic curvature, the strength of the direct attractive interaction between inclusions and on the active forces exerted by the inclusions. , while in the lower panels C and D the myosin generated actin active forces F = -1 kT/l min are taken into account. The negative sign of F denotes that the myosin generated forces points into the vesicle interior (perpendicular to the local membrane surface). The triangulated membrane surface is drawn semitransparent to uncover its interior shape. Red vertices on the mesh represents the locations of the membrane inclusions. In the panels C and D black and white outlines are added to better visualize the shape of endovesicles and membrane necks. The values of other model parameters are: bending rigidity κ = 25 kT and direct interaction parameter w = 1.25 kT. Figures 3A and B show that the lateral accumulation of membrane inclusions (nanodomains) can induce the formation of long undulated thin inward membrane protrusions (buds). Long undulated membrane protrusions may be further transformed into small independent spherical endovesicles, as observed in Figure 1D, due to the frustrations in the orientational ordering of membrane components in the highly curved membrane necks. 67 The same mechanism can also be responsible for the detachment of the complete inward membrane protrusion from the parent membrane 67 and the consequent formation of the endovesicles as shown in experimental Figure 1.

Results and Discussion
In calculations presented in Figures 3C and D, it is taken into account that the membrane inclusions (nanodomains) exert also active forces in the direction perpendicular to the membrane surface towards the interior of the RBC. As already discussed above, the active forces in the RBC membrane are generated by myosin (NMIIA) motor domains (inclusions) bound to F-actin of the RBC membrane skeleton. 60 It can be seen in panels of Figure 3 that at a smaller concentration of membrane inclusions exerting force on the membrane, the MC predicted RBC shape has one large invagination (Fig. 3C) as can be observed in some experiments (see refs. 44,56 and the references therein). Large invaginations can be separated from the parent cell due to the frustrations in the orientational ordering of membrane components in the highly curved membrane necks connecting the invagination and the parent cell, 67 resulting in the formation of a large endovesicle.
Furthermore, Figure 3 shows that a larger concentration of membrane inclusions exerting force on the membrane, the MC predicted RBC shape has one small and two large pancake-like torocyte membrane invaginations (Figure 3D) as can also be observed in the experiments ( Figure  1, panels A, B and E). Again, the necks connecting the torocyte structures as well as the neck connecting the complete invagination to the parents cell are supposed to be ruptured due to the frustrations in the orientational ordering of membrane components in the highly curved membrane necks. 67 Note that in the torocytic membrane invaginations ( Figure 3D), the myosin motor domains generated active forces that acted in the outward direction with respect to the torocyte invagination, i.e. in the direction towards the inner RBC solution. It is also important to point out that the myosin motor domains/inclusions are mainly accumulated at the bulby rim of the torocyte membrane invaginations ( Figure 3D).
The stability of torocyte RBC endovesicles can be theoretically explained also by anisotropic membrane inclusions which exhibit orientational ordering in a highly curved bulby periphery. 57,59 In this work we have shown that the stability of torocyte endovesicles may be additionally favored also by active forces on the RBC membrane, generated by myosin motor nanodomains (i.e. membrane inclusions/nanodomains composed of myosin-actin-spectrin-lipids complex, see also references 7, 19, 68, 69).
Note that in the present work the total number of myosin motor nanodomains (inclusions) in Figures 3C and D is larger than the actual number of myosin motor nanodomains found in the RBC membrane. 60 Therefore, the proposed mechanisms of myosin inclusions/nanodomains driven formation and stabilization of torocyte invaginations and torocyte endovesicles in RBCs can be considered as an additional and complementary mechanism to the non-homogeneous lateral distribution and the orientational ordering of anisotropic membrane inclusions/ nanodomains in the RBCs membrane. 57,59 In accordance with experimental observations (Figure 1) we have predicted in this work the invaginated stomatocyte RBC shapes having different shapes of invaginations, like torocytic, spherical, undulated necklace-like, etc. This is an extension of the previously theoretically predicted shape classes of the invaginated stomatocytic shapes which were mostly limited to the simple (axisymmetric) stomatocytic shape with only one large invagination ( Figure 3C) (see for example ref. 32), experimentally observed also in a giant unilamellar lipid vesicle. 70

Conclusions
Numerical computer modelling of the cell membrane shapes and shape transformations is a widely used method to investigate physical models of various phenomena found experimentally in cellular systems. In the present work, we used MC simulations to theoretically study the influence of active forces on the red blood cell (RBC) shape and vesiculation. We have shown that the active forces, generated by myosin motor domains, may partially control endovesiculation of the RBC membrane and also the RBC shape changes in general. Among others the myosin domains generated forces on the RBC membrane favor experimentally observed pancake-like (torocyte) RBC endovesicles with a large flat central membrane and a bulby periphery which were also experimentally observed. 57 These theoretically predicted shapes ( Figure 3D) are very similar to the shapes of Golgi bodies, so the theoretical study presented in this work may be relevant also for better understanding of the physical mechanisms determining the shape of Golgi body. 58 In conclusion, by using MC simulations, it is shown in this work that the recently discovered myosin motor domains generated active forces on the RBC membrane 60 may partially control the endovesiculation of the red blood cells and the RBC shape changes in general. We conclude that the myosin generated active forces acting on the RBC membrane should be therefore necessarily considered in the future relevant theoretical studies of the RBC vesiculation and shape transformations. Povzetek Z uporabo Monte Carlo (MC) simulacij stomatocitnih oblik eritrocita smo pokazali, da lahko aktivne sile miozinskih domen, pripetih na F-aktinske molekule na notranji strani membrane eritrocita, delno kontrolirajo endovezikulacijo membrane eritrocita. Miozinsko-aktinsko generirane aktivne sile, ki delujejo na membrano eritrocita, tako med drugim vplivajo tudi na nastanek ploščatih torocitnih endoveziklov eritrocitov. Na osnovi dobljenih rezultatov MC simulacij vezikulazije eritrocitov zaključujemo, da lahko miozinsko-aktinsko generirane aktivne sile, ki delujejo v smeri pravokotno na površino membrane eritrocita proti notranjosti eritrocita, pomembno vplivajo tudi na ostale transformacije oblik eritrocita ter na stabilnost različnih oblik eritrocita in bi jih bilo zato potrebno upoštevati pri bodočih teoretičnih študijah transformacije oblik in vezikulacije eritrocitov.