**2. Theoretical prediction of the CV shape**

Theoretical and experimental studies on (artificial) phospholipid vesicles have shown that the shape of the bilayer membrane-enclosed compartments can be theoretically well described by the properties of the membrane [34]. By considering that membrane is composed of many constituents, methods of statistical physics were used for its description [35, 36]. The key feature in the expression of the energy of a single constituent is the mismatch of the local curvature and the constituent intrinsic curvature (the one fitting the shape of the membrane constituent) [36]. In order to compose the membrane, the constituent attains the local membrane shape that usually differs from its intrinsic shape. Moreover, if the constituent is not symmetric with respect to the axis perpendicular to the membrane surface, the principal axes of the membrane and the constituent are in general rotated by an in-plane angle, meaning that the constituents can attain different in-plane orientations in the membrane which correspond to different energies. The consequence of the mismatch in curvature and orientation is that certain energy is required to insert the constituent at the site, this energy being higher if the mismatch is greater. The thermal motion opposes the complete orientational ordering in the direction with the lowest energy but the constituents will spend on average more time in the orientation with lower energy. As constituents are more or less free to move laterally in the membrane, they can redistribute in a way to minimize the mismatch between the intrinsic and the actual shape. Summing up the contributions of all the constituents and considering the entropic effects due to lateral and orientational ordering yields the expression for the free energy of the membrane *F* in terms of the mean curvature of the membrane surface *H* and the curvature deviator *D*, both composed of the two principal curvatures *C*<sup>1</sup> and *C*<sup>2</sup> [34],

$$H = (\mathbf{C}\_1 + \mathbf{C}\_2)/2,\tag{1}$$

$$D = (\mathbf{C}\_1 - \mathbf{C}\_2) / \mathbf{2},\tag{2}$$

$$F = -kT \left[ \sum\_{i} m\_i \ln \left( q\_i^0 \, \text{2} \cosh \left( d\_{i, \text{eff}} \right) \right) \, \text{d}A + k\_{\text{B}} T \left[ \sum\_{i} m\_i \ln \left( m\_i / m \right) \text{d}A, \right. \tag{3} \right]$$

$$q\_i^0 = \exp\left(-\frac{\xi\_i (H - H\_{i,\text{m}})^2}{2k\_\text{B}T} - \frac{\left(\xi\_i + \xi\_i^\*\right)\left(D^2 + D\_{i,\text{m}}^2\right)}{4k\_\text{B}T}\right),\tag{4}$$

$$d\_{i, \rm eff} = \left(\sharp\_i + \sharp\_i^\*\right) D D\_{i, \rm m} / k\_{\rm B} T,\tag{5}$$

where ξ*<sup>i</sup>* and ξ <sup>∗</sup> *<sup>i</sup>* are constants, *k*<sup>B</sup> is the Boltzmann constant,*T* is temperature, *mi* is the local area number density of the *i*-th kind of constituents and *m* is the number density taking all constituents over all membrane area. Integration is performed over membrane surface *A.* The summation accounts for all types of constituents that are characterized by index *i.* The intrinsic mean curvature of the membrane surface *H*<sup>m</sup> and the intrinsic curvature deviator *D*<sup>m</sup> refer to respective principal curvatures intrinsic to the shape of the constituent [34].

Free energy given by Eq. (3) consists of two terms—the first one deriving from the single-constituent energy and entropy of orientational ordering, and the second one deriving from lateral distribution of constituents. Eqs. (3) and (4) can be rewritten in the form [34].

$$F = W\_\text{B} - \int 2mk\_\text{B}T \cosh\left(d\_{i,\text{eff}}\right) \,\text{d}A + k\_\text{B}T \left[\sum\_i m\_i \ln\left(m\_i/m\right) \text{d}A\right] \tag{6}$$

where *W*<sup>B</sup> has the form of the bending energy of a laterally isotropic membrane [34].

A basic physical principle that the system will attain the shape corresponding to minimal free energy [34]

$$F = \min,\tag{7}$$

at relevant constraints to the system is taken into account such as the requirement of fixed membrane area *A* and fixed enclosed volume *V* [37]. In the absence of net external forces acting on the membrane, a convenient geometrical parameter is the relative volume *v* which represents the volume to area ratio largely determined by the osmotic equilibrium [34]

$$v = \sqrt{V^2 / 36\pi A^3}.\tag{8}$$

Also, other constraints can be imposed upon the system, for example, constant average mean curvature [34]

$$ = \int H \, \mathrm{d}A // \int \mathrm{d}A. \tag{9}$$

To find the free energy minimum, the above formalism [Eqs. (1)–(5)] is used to state and solve the so-called variational problem in which the principal curvatures are expressed in terms of convenient coordinates. Dimensionless parameters are used for clarity: *cj* ¼ *CjR*, *j* = 1,2 are the principal curvatures normalized with respect to the radius of the sphere with the surface area *<sup>A</sup>*, *<sup>R</sup>* <sup>¼</sup> ffiffiffiffiffiffiffiffiffiffiffi *<sup>A</sup>=*4<sup>π</sup> <sup>p</sup> , *<sup>h</sup>* <sup>¼</sup> *HR,* and *<sup>d</sup>* <sup>¼</sup> *DR* are the normalized mean curvature and curvature deviator of the membrane, respectively, and *hi*,m ¼ *Hi*,m*R* and *di*,m ¼ *Di*,m*R* are the intrinsic mean curvature and the intrinsic curvature deviator of the *i*-th type of constituents, respectively. The area element is normalized with respect to the area *<sup>A</sup>*, d*<sup>a</sup>* <sup>¼</sup> <sup>d</sup>*A=*4π*R*<sup>2</sup> . The constants ξ*<sup>i</sup>* and ξ <sup>∗</sup> *<sup>i</sup>* are taken to be equal for simplicity. The free energy is normalized with respect to the free energy of the sphere composed of chosen constituents, 8π*m*ξ*i*.

Consistently related distributions of the constituents, their in-plane orientations (if relevant), and the membrane shape are determined simultaneously in a mathematical procedure. Several methods have been developed for this purpose—for example, ansatz [37], numeric solution of differential equation [38], surface evolver [39], or finite element method [40]. The simplest method is the ansatz approach in which the space of possible solutions is assumed within a family of shapes with adjustable parameters. Such a solution can be analytical throughout (depending on the sophistication of the ansatz) and therefore transparent. The advantage of transparent methods is that with some basic mathematical skills the procedures can be repeated and used accordingly. Differential equations expressing minimization of free energy are derived by applying the Euler–Lagrange method [36] in a convenient parametrization. Considering a multicomponent system, the constituents are free to move laterally over the membrane which may present singularities that have not yet been fully explored [36]. Such solutions require numerical procedures that are implemented in respective customized software which normally requires manipulation by a skilled researcher. However, the set of possible solutions is considerably larger than in the ansatz case, as for some classes of shapes the relevant ansatz does not exist. With more freedom in finding a solution, the achieved energy attains lower values. The advantage of the rigorous solution of the system of differential equations is that the lowest possible energy can in principle be found, however, the method is rather demanding for shapes that are not axisymmetric which limits the set of possible solutions.

The theory could be refined by considering particular experimental evidence. For example, it was found that cellular vesicles may go through a process of active and passive solutes' permeation which may cause cyclical expansion and contraction [41].
