Membrane remodeling due to mixture of multiple curvature proteins: Insights from mesoscopic modeling

We present an extension of the Monte Carlo based mesoscopic membrane model, where the membrane is represented as a dynamically triangulated surface and the proteins are modeled as anisotropic inclusions formulated as in-plane nematic field variables adhering to the deformable elastic sheet. The local orientation of the nematic field lies in the local tangent plane of the membrane and is free to rotate in this plane. Protein-membrane interactions are modeled as anisotropic spontaneous curvatures of the membrane and protein-protein interactions are modeled by the splay and bend terms of Frank’s free energy for nematic liquid crystals. In the extended model, we have augmented the Hamiltonian to study membrane deformation due to a mixture of multiple types of curvature generating proteins. This feature opens the door for understanding how multiple kinds of curvature-generating proteins may be working in a coordinated manner to induce desired membrane morphologies. For example, among other things, we study membrane deformations and tubulation due to a mixture of positive and negative curvature proteins as mimics of various proteins from BAR domain family working together for curvature formation and stabilization. We also study the effect of membrane anisotropy, which manifests as membrane localization and differential binding affinity of a given curvature protein, leading to insights into the tightly regulated cargo sorting and transport processes. Our simulation results show different morphologies of deformed vesicles that depend on the curvatures and densities of the participating proteins as well as on the protein-protein and membrane-proteins interactions.


INTRODUCTION
Bio-molecular structures dynamically assemble, perform highly orchestrated biochemical functions, and disassemble seamlessly before reassembling to continue the cycle. Since the molecular structure and the biochemical pathways in biological processes are optimized over evolutionary timescales to perform elaborate functions, the "biological circuitry" is extremely sophisticated. Despite being under strong thermal noise and embedded in a crowded molecular milieu, the molecular machines and motors encrypt and process complex set of information in a precise manner. This information is embedded at multiple length and time scales that range from structures at atomic resolution to formation of molecular complexes and leading to subcellular and cellular processes such as ion transport across channels and pumps, fusion and fission in endocytic recycling, presynaptic neurotransmission, morphogenesis and metastasis to name a few [1][2][3].
A paradigmatic example of this "biological circuitry" is the Clathrin-mediated endocytosis (CME) in vesicular trafficking. After five decades since its discovery in 1964 and countless seminal subsequent discoveries, it is now well established that more than 50 proteins make up this molecular machinery [4][5][6]. Proteins, either individually or in complex are "called upon" from the cytosolic reservoir in a highly coordinated manner at different stages that include (i) initiation/identification of the endocytic site, (ii) recruitment of the cargo to the site, (iii) invagination/remodelling of the membrane, (iv) leakage-free scission (v) disassembly/uncoating of proteins and (vi) vesicle release for trafficking. These functions are carried out by disparate proteins that have been studied in great detail. However, we still lack a good understanding of how all these different components work together in a highly coordinated manner to drive vesicle formation. There are several examples like CME in the cell, which take place due to a large number of proteins and other biomolecules working as a 'biological circuitry'.
BAR domain class of proteins are the quintessential curvature inducing scaffold proteins and are ubiquitously present in the endocytic recycling pathways and play a crucial role in vesicular transport processes in the cell [5,[7][8][9][10][11][12]. BAR family of proteins have a wide range of curvatures with proteins such as F-Bar and N-BAR exhibiting positive curvatures [5,5,13] and proteins such as I-BAR exhibiting negative curvature [14]. When these proteins bind to the membrane surface, they impose their local curvature on the membrane and also work in concert with each others at longer length scales to induce non-local curvatures on the membrane [15][16][17][18][19]. Though the sub-cellular processes of vesicular transports in mediated by multiple curvature proteins interacting in concert with each other, this aspect has not been studied very much in reconstitution experiments [20][21][22][23][24]. Also, a handful of theoretical and computational mesoscopic modeling work have provided important insights into the synergy between multiple scaffold proteins during membrane remodeling processes [25][26][27][28] but the studies are quite few in numbers and a more extensive exploration in this area is needed.
In this work, we develop a new framework to study multiple curvature proteins at mesoscopic scales. We use this framework to address the evolution of membrane morphologies as a function of mixture of differently curved proteins with different compositions and membrane binding strength as well as with a variety of proteinportein interactions behaviour. In our mesoscopic membrane mechanics based approach, we model the membrane as a triangulated surface and proteins as nematics. We use an augmented form of Helfrich type free energy function and perform Monte-Carlo simulations to model this process [29,30]. Membrane deformation with a single fixed scaffold protein is extensively studied by this model [29][30][31][32]. Here we have evolved the underlying Hamiltonian so that the method can be used to study membrane remodeling due to mixture of proteins with different curvatures working together.

MODEL
We model proteins as nematics and deformable fluid membrane as a triangulated surface. In this model [29,30], nematic lies on the local tangent plane of the vertex of the triangle. These nematics (proteins) are free to rotate in the local tangent plane and are denoted byn( r). Here, the interactions between protein and membrane are modelled as anisotropic spontaneous curvatures of the membrane in the vicinity of the protein filament. Protein-protein interactions are modelled by the splay and bend terms of the Frank's free energy for nematic liquid crystals. The total energy is the sum of membrane energy, energy due to interaction between protein and membrane and the protein-protein interaction energy. The energy of the bare membrane, which is the first of the total energy is written as: where E 1 is the Canham-Helfrich elastic energy for membranes [33,34], κ is bending rigidity and H is membrane mean curvature with H = (c 1 + c 2 )/2. c 1 and c 2 are the local principal curvatures on the membrane surface along the orthogonal principal directionst 1 and t 2 .
The coupling energy for the proteins and membrane interactions is written as E 2 below.
Here, N is the number of different curvature values of proteins. κ and κ ⊥ are the induced membrane bending rigidities and c and c ⊥ are the induced intrinsic curvatures alongn and t.n is the orientation of the protein andt is its perpendicular direction in the local tangent plane. H and H ⊥ are the membrane curvatures in the direction ofn andt, where H = c 1 cos 2 φ+c 2 sin 2 φ and H ⊥ = c 1 sin 2 φ+c 2 cos 2 φ. Here φ is the angle between the direction of nematic orientationn and the principal direction t 1 .
The protein-protein interactions is formulated as below: Here, E 3 is the energy for the nematic-nematic interactions. It is taken from Frank's free energy for nematic liquid crystals. Here∇ is the covariant derivative on the curved surface. K 1 and K 3 are the splay and bending elastic constants for the in plane nematic interactions. Here, a discrete form of this energy is used [29] that makes the implementation amenable for the Monte-Carlo simulations and is expressed as: The above formulation is inspired from the Lebwohl-Lasher model [35]. Here, LL is strength of the nematic interaction with a constant approximation (K 1 = K 3 ). The sum i>j is over all the nearest neighbour (i, j) vertices on the triangulated grid, promoting alignment among the neighbouring orientation vectors. Total energy is the sum of these three energy.

PARAMETERS
In this work, we use multiple curvature proteins with different densities and interactions. Below, we elucidate the input parameters used in the simulation for multiple curvature proteins expressed as matrices. If total number of different curvature proteins is N , then the nematic density is given as: where ρ 1 , ρ 2 , ......ρ N are the densities of different curvature proteins(nematics) and n v is the total number of the vertices on the triangulated surface. Here we have one nematic on every vertex. So the total number of nematics cannot exceed the total number of vertices. Different proteins can have different curvatures, induced membrane bending rigidities and interactions.
LL is the interaction between same curvature of nematics (same interaction) and i,j LL is the interaction between different curvature of nematics (cross interaction). We have considered both cases i,i LL > i,j LL and i,i LL < i,j LL and show that how they affect the membrane deformation.

RESULT
In this work, our main focus is on the membrane remodeling due to mixture of different curvature proteins or due to a single protein with variable curvature. Curvature of the protein can vary due to its linker regions. For this we have used nematic based model and developed it for the multiple curvature proteins. Fig.1(a-d) show simulation results for the membrane remodeling due to two, three, four and five different curvature proteins. Here the total protein coverage is 100%, which is equally divided into different curvature proteins. The bending rigidity of the membrane κ = 20K b T , induced bending rigidity κ = 5K b T and nematic-nematic interaction LL = 2K b T . Initially the proteins are randomly distributed on the spherical vesi-  Here the lower panel of Fig.2 is reproduced with different value of membrane protein interaction κ = 10 K b T . Other parameters are same as Fig.2. cle. Fig.1a shows the converged solution of the Monte-Carlo simulation for the two different proteins with curvatures ±0.5. Similarly, Fig.1b-d show the membrane remodeling results for three, four and five different proteins with curvature values (0, ±0.5), (±0.5, ±0.8) and (0, ±0.5, ±0.8) respectively. Fig.2 shows the vesicle deformation due to two different proteins with opposite curvatures ±0.5. The total protein coverage is 100% but densities of both proteins are changing regularly. On the X axis, density of positively curved protein is decreasing 90% to 10% while the density of negatively curved proteins in increasing 10% to 90%. Here the different shapes are obtained due to different value of nematic-nematic interaction and different proteins density. As nematicnematic interaction in increasing, proteins those have higher density are coming together. Fig.3 shows the vesicle deformation due to two different proteins with curvature ±0.5. Here the induced membrane bending rigidity in higher then Fig.2. These results show if the proteins bind to membrane more strongly, more deformation occurs on the vesicle. Here the induced membrane bending rigidity κ is 10 K b T and the other parameters are same as Fig.2.
In Endoplasmic reticulum, tubes are stabilized by reticulon and reep proteins while the junction between them is stabilized by protein lunapark. Fig.4a shows the membrane tubulation with two different proteins with curvature 0.5 and−0.3. Positive curvature proteins form tube while the negative curvature proteins form the junction between the tubes. Here the protein density for positively and negatively curved proteins are 90% and 10% respectively. Other parameters are κ = 20, κ = 15and LL = 2 in K b T unit. Fig.4b shows the membrane tubulation for different curvature values of two different proteins. Here the triple junction is forming like Endoplasmic reticulum. Positively curved proteins are forming the tubes while the negatively curved proteins are stabilizing at the junction. Fig.5 shows the simulation results for two different curvature of proteins with different densities. Here the interaction between the different curvature proteins is also included and denoted by Cross LL . Protein coverage in upper and lower panels are 40% and 80% respectively. Here the curvature of the two different proteins are ±0.5. Proteins with positive curvature are shown by blue color while the negatively curved proteins are shown by red color. Fig.5a shows the simulation result with Cross LL = 0. Although there is no interaction between positively and negatively curved proteins still they have arranged close to each other because positively curved proteins generate a negatively curved region on the bare membrane at the boundary of there arrangement and negatively curved proteins prefer that region due to there own curvature. Proteins have arranged parallely due to interactions between them. In Fig.5b   Here the total protein coverage is 40% which is equally divided for two different proteins. Fig.6c,d shows the similar results like 6a,b with 80% protein coverage.
Some proteins have variable curvature due to its intrinsically disordered linker regions. In these proteins, curvature can vary around a fixed mean value. Fig.7 shows the vesicle deformation due to proteins with three different curvatures. Here protein curvatures are 0.5, 0.75 and 1.0 with densities 30%, 40% and 30% respectively. In Fig.7 different shapes are obtained due to different values of protein membrane interaction and protein protein interaction. As expected, we observe that the deformation is increasing with the membrane protein interaction. Fig.8 shows the vesicle deformation results due to three different curvature proteins where x-axis represents the densities of the different curvatures of proteins while y-axis represents proteinprotein interaction. Different shapes are obtained due to different coverage ratio of curvature proteins and protein-protein interactions. In these results, the curvature of the proteins are 0.3, 0.4 and 0.5 and the membrane-protein interaction κ is 10 K b T . Membrane tubes are quite common in cell functions and are formed by the curved proteins and bundling tendency of the proteins. Tube formation due to bundling tendency of protein is briefly discussed in ref [30]. In this work, we show membrane tubulation due to combination of two different curvature proteins. Fig.9a shows the membrane tubulation due to two curved proteins with curvature c = 1.0 and -0.05. We show the positively curved proteins by Blue color those form the tube while the negatively curved proteins form valley and we show that by Red color. Proteins that don't have curvature but can form bundles like microtubule are also capable of creating membrane tubes. Fig.9b shows the membrane tubulation due to bundling of the proteins. Here one protein is has zero curvature (Red) c ⊥ = 1.0, c = 0 and another one have a positive curvature (Blue) c = −0.05. Arrangement of the protein on the tube is different in Fig.9a and b. In Fig.9a, protein arrange themselves in the transverse direction to the length of the tube while in Fig.5a, they are along the length of the tube. This happens because in first case, tubulation occurs due to curved protein while n second case, bundling tendency of the proteins is the driver of tube formation.
Experimentally, it is observed that two different proteins can interact with each other during membrane remodeling processes [20][21][22][23][24]. In this section, we show some simulation results in which membrane is deforming with two different proteins having opposite curvatures. Also, two different proteins with same curvature but different membrane protein interactions, also deform the membrane. Fig10 shows the membrane remodeling with two different proteins shown by blue(A) and red(B). We fixe the membrane-protein interaction of protein A, which is 10 K B T and membrane protein interaction of protein B is varying on the Y Axis. The density of both proteins are changing on the X Axis but the total density is fixed and equal to 50%. Here, the protein-protein interaction between same proteins is 2 K B T and in different proteins is 1.5 K B T . We find that tubulation is possible is certain cases where the density of tightly membrane bound proteins is 10% and loosely bound proteins is 40%. The membrane protein interaction of both proteins are 10 K B T and 20 K B T , respectively. are same as Fig10. The protein density is varying on the X axis and the membrane pro-tein interaction is varying on the Y axis. It is clear from this figure, Around 70% density of protein is required for the tubulation but Fig10 shows that when two different proteins interact with the membrane, less density is required for tubulation (around 50%). After comparison of both results, we find that less protein density is required for membrane remodeling when two different proteins interact the membrane with different membrane protein interactions. Tubulation is also possible with two different proteins when they are not interacting with each other as shown in Fig12. Here the total protein density is same as Fig10 which is equal to 50% but the ratio of both protein's density is different. The membrane protein interaction for tightly and loosely bound proteins are 40 K B T and 20 K B T and the other parameters are same as Fig10.

CONCLUSION
Membrane deformation in cell occurs due to multiple proteins working in synergy with each other. Some processes like clathrin mediated endocytosis [4] and formation of three way junctions in the endoplasmic reticulum [36,37] are examples of membrane deformation and stabilisation of induced curvature due to multiple proteins. In this work, we have shown the deformation of membrane due to multiple proteins. These proteins can have different curvatures and interactions. We have included different values of the curvature and interaction and we also include the interaction between different curved proteins. From our Monte Carlo simulation results. we show that tubes are forming by positively curved proteins and valley (junction) are formed and stabilized by negatively curved proteins. We also show the effect of interaction between different curved proteins, where we report two particular cases: (i) Cross LL < LL and (ii) Cross LL > LL . We show that the deformation occurs when Cross LL < LL between the two oppositely curved proteins. If the cross interaction is high ( Cross LL > LL ), then no deformation occurs because in this case two oppositely curved proteins arrange at random position parallel manner. In this case, proteins are unable to deform the membrane. Further, we show that the membrane tubulation by the mixture of different proteins having a variety of curvatures leading to formation of different morphologies.

AUTHOR CONTRIBUTIONS
Anand Srivastava conceived the idea and formulated the Hamiltonian with help of Gaurav Kumar. Gaurav Kumar wrote the software and carried out the simulations and all the analyses. Gaurav Kumar and Anand Srivastava wrote the paper.