A minimal selforganisation model of the Golgi apparatus
Abstract
The design principles dictating the spatiotemporal organisation of eukaryotic cells, and in particular the mechanisms controlling the selforganisation and dynamics of membranebound organelles such as the Golgi apparatus, remain elusive. Although this organelle was discovered 120 years ago, such basic questions as whether vesicular transport through the Golgi occurs in an anterograde (from entry to exit) or retrograde fashion are still strongly debated. Here, we address these issues by studying a quantitative model of organelle dynamics that includes: denovo compartment generation, intercompartment vesicular exchange, and biochemical conversion of membrane components. We show that anterograde or retrograde vesicular transports are asymptotic behaviors of a much richer dynamical system. Indeed, the structure and composition of cellular compartments and the directionality of vesicular exchange are intimately linked. They are emergent properties that can be tuned by varying the relative rates of vesicle budding, fusion and biochemical conversion.
Introduction
The Golgi apparatus is an intracellular organelle at the crossroad of the secretory, lysosomal and endocytic pathways (Heald and CohenFix, 2014). One of its most documented functions is the sorting and processing of many proteins synthesized by eukaryotic cells (LippincottSchwartz et al., 2000). Proteins translated in the endoplasmic reticulum (ER) are addressed to the Golgi, where they undergo posttranslational maturation and sorting, before being exported to their final destination. The Golgi itself is composed of distinct subcompartments, called cisternae, of different biochemical identities (Stanley, 2011). From the ER, cargoproteins successively reside in cisternae of the cis, medial and transidentities, after which they exit the Golgi via the transGolgi network (TGN). In some organisms such as the yeast Saccharomyces cerevisiae cisternae are dispersed throughout the cell (Suda and Nakano, 2012) and each cisterna undergoes maturation from a cis to a transidentity (Losev et al., 2006; MatsuuraTokita et al., 2006). In most other eukaryotes, cisternae are stacked in a polarized fashion, with cargoproteins entering via the cisface and exiting via the transface (Boncompain and Perez, 2013). This polarity is robustly conserved over time, despite cisternae constantly exchanging vesicles with each other, the ER and the TGN (LippincottSchwartz et al., 2000). Pharmacological treatments that alter the structure and stacking of Golgi compartments in mammalian cells affect the processing of certain proteins, and in particular glycosylation (Hu et al., 2005; Shen et al., 2007; Xiang et al., 2013; Bekier et al., 2017).
The highly dynamical nature and compact structure of a stacked Golgi makes it difficult to determine how cargoproteins are transported in an anterograde fashion from the cis to transside of the stack and how this transport is coupled to processing. Several models have been proposed to explain the transport dynamics of Golgi cargo (See sketch in Figure 1A). They mostly belong to two classes: one involving stable compartments exchanging components through fusionbased mechanisms and one involving transient compartments undergoing enblock maturation and in which fusion mechanisms are dispensable for cargo transport. Historically (Malhotra and Mayor, 2006), the 'Vesicular transport’ model postulated that cisternae have fixed identities and cargo progresses from one cisterna to the next via anterograde vesicular transport, while the 'Cisternal Maturation’ model postulated that cisternae themselves undergo maturation from the cis to the transidentity and physically move through the stack, while Golgi resident proteins remain in the stack by moving toward the cisface by retrograde vesicular transport. Several extensions of these two models have been proposed, including cisternal maturation with tubular connections (Trucco et al., 2004; Rizzo et al., 2013), rapid partitioning within Golgi cisternae between processing and exporting subcompartments (Patterson et al., 2008), or the cisternal progenitor model (Pfeffer, 2010), in which stable cisternae exchange material by the fusion and fission of subcompartments defined by their composition of Rab GTPases, which evolve over time through exchange with the cytosol. The strengths and weaknesses of these different models are nicely reviewed in Glick and Luini, 2011. Variations of these models, such as the ‘rim progression’ model also exist (Lavieu et al., 2013). It is noteworthy that these models do not provide a quantitative analysis of the generation and maintenance of Golgi compartments, nor do they attempt to relate the Golgi structure (number and size of compartments) and transport dynamics. Therefore, much remains unknown regarding the mechanisms that dictate the organisation and dynamics of the Golgi.
Although there is a large diversity in Golgi structures and dynamics among different species, the physiological function of the organelle as a sorting and processing hub is common to all species, suggesting that important mechanisms controlling Golgi dynamics are conserved. Works in evolutionary biology and biophysics have attempted to describe these mechanisms (Klute et al., 2011; Sens and Rao, 2013). Different classes of mathematical models have been proposed, from models of vesicle budding and fusion based on rate equations (Binder et al., 2009; Ispolatov and Müsch, 2013; Mukherji and O'Shea, 2014)  some specifically including space (Sachdeva et al., 2016), to continuous, discrete and stochastic models of protein sorting in the Golgi (Gong et al., 2010; Vagne and Sens, 2018b) and computer simulations based on membrane mechanics (Tachikawa and Mochizuki, 2017). To account for the ability of the Golgi to reassemble after mitosis (Wei and Seemann, 2017), many of these studies have sought to describe the Golgi as a selforganized system (Binder et al., 2009; Kühnle et al., 2010). However, the manner in which the kinetics of the Golgi two main functions (namely the transport and biochemical conversion of its components) interplay with one another to yield a robust steadystate has so far received much less attention (Sachdeva et al., 2016; Vagne and Sens, 2018b).
Most existing models discuss transport through a preexisting compartmentalised Golgi structure. The main goal of our model is to explain how this compartmentalised structure and the flux within it can spontaneously emerge through selforganisation. We propose here a model in which both Golgi selforganisation and vesicular transport are solely directed by (local) molecular interactions resulting in compositiondependent budding and fusion. In particular, we neglect the spatial structure of the Golgi, and assume that any two components of the system can interact with one another (see the Discussion section for the relevance of this approximation). Thus, the directionality of vesicular transport is an emergent property of the selforganisation process that coevolves with the size and number of Golgi subcompartments, rather than being fixed by arbitrary rules. This approach brings new conceptual insights that highlight the link between the largescale steadystate organisation of the Golgi, the directionality of vesicular transport, and the kinetics of individual processes at the molecular scale. We show that a wide variety of organelle phenotypes can be observed upon varying the rates of the local processes. In particular, our model uncovers a correlation between the size and composition (purity) of Golgi subcompartments. Besides, it emphasizes that compositiondependent exchange between compartments yields complex vesicular flux patterns that are neither purely anterograde (as in the ‘vesicular transport’ model) nor purely retrograde (as in the ‘cisternal maturation’ model).
Model
The model, shown in Figure 1 and fully described in the Appendix 1 Detailed model is a coarsegrained representation of the Golgi Apparatus. The parameters used in the model are summarised in Appendix 1—table 1, and the quantities used to charaterise a compartment and the steadystate of the system are summarised in Appendix 1—table 2 and Appendix 1—table 3, respectively. The system is discretised at the scale of small vesicles, which define the unit size of the system. Vesicles can fuse together or with existing compartments to form larger compartments and can bud from compartments. Compartments consist of a number of membrane patches (of size unity) with distinct identities. The biochemical identity of a membrane patch is defined in a very broad sense by its composition in lipids and proteins that influence its interaction with other patches, such as Rab GTPases (Grosshans et al., 2006) or fusion proteins such as SNAREs (Cai et al., 2007). Consistently with the three main biochemical identities reported for Golgi membranes (Day et al., 2013), vesicles and membrane patches are of either cis, medial or transidentity, and can undergo irreversible biochemical conversion from the cis to the transidentity. The system is fed by a constant influx $J$ of cis vesicles coming from the ER, and vesicles and compartments can fuse homotypically with the ER (itself of cis identity) and the TGN (of trans identity).
Compartments are defined by their size $n$ (a number of patches) and their composition, the fractions ${\varphi}_{i}$ (with $i=cis,medial,trans$) of patches of the three different identities composing it. These quantities are dynamically controlled by three, composition dependent, microscopic mechanisms – budding, fusion and biochemical conversion – as described below.
Fusion: Homotypic fusion – the higher probability of fusion between compartments of similar identities – is a feature of cellular organelles in general (Antonin et al., 2000) and the Golgi apparatus in particular (Pfeffer, 2010; Bhave et al., 2014). This is controlled here by a fusion rate between compartments proportional to the probability that they present the same identity at the contact site. The total fusion rate for two compartments $(a)$ and $(b)$ with composition ${\varphi}_{i}^{(a)}$ and ${\varphi}_{i}^{(b)}$ for each identity $i$ (with $i$ equals cis, medial and trans) is then:
(1) ${K}_{\mathrm{f}}\times \sum _{i}{\varphi}_{i}^{(a)}{\varphi}_{i}^{(b)}$where ${K}_{\mathrm{f}}$ is the fusion rate between two identical compartments. Fusion with the boundaries follow the same rules, with the ER containing a fraction ${\alpha}_{\mathrm{ER}}$ of cis patches and the TGN a fraction ${\alpha}_{\mathrm{TGN}}$ of trans patches. Numerical results described in the main paper are for ${\alpha}_{\mathrm{ER}}={\alpha}_{\mathrm{TGN}}=1$, and the role of the composition of the boundaries is discussed in the Appendix 7.
Budding: Budding is the emission of a vesicle from a larger compartment. Compositiondependent vesicular budding is an important sorting mechanism in cellular traffic (Bonifacino and Glick, 2004). Owing to the high specificity of the budding machinery, each vesicle is assumed to be composed of a single membrane identity. The budding flux for any components $i$ (equals cis, medial or trans) present in a compartment is assumed not to depend on the number ${n}_{i}=n{\varphi}_{i}$ of patches of identity $i$, but on the total size of the compartment:
(2) ${J}_{\mathrm{b},i}={K}_{\mathrm{b}}\times n.$This budding kinetics corresponds to the saturated regime of a MichaelisMenten reaction kinetics, and is appropriate to the case where budding actors (e.g. coat proteins) bind nonspecifically to the membrane of the compartment and find their budding partners (a patch of a particular identity) by diffusion (Vagne and Sens, 2018b). An alternative budding scheme where the budding flux is linear with the number of patches of a given identity (${J}_{\mathrm{b},i}={K}_{\mathrm{b}}{n}_{i}$) is discussed in Appendix 8.
Biochemical conversion: Each membrane patch is assumed to undergo stochastic biochemical conversion from a cis to a medial to a transidentity. This local mechanism of identity conversion is motivated by the Rab cascade, during which the membrane identity evolves over time through molecular exchange with the cytosol, as was observed in endosomes (Rink et al., 2005) and in the Golgi of Yeast (Losev et al., 2006; MatsuuraTokita et al., 2006). The Rab conversion mechanism has not been demonstrated in mammalian Golgi, but is thought to be universally important for the specificity of membrane trafficking (Grosshans et al., 2006), and many of the proteins involved in the Yeast Rab cascade are conserved in mammalian cells (Klute et al., 2011). Rab conversion is at the basis of the cisternal progenitor model (Pfeffer, 2010), and it is thus of fundamental importance to understand its potential implication on Golgi structure and dynamics with the use of a quantitative model. We stress that our model of biochemical conversion is not limited to Rab GTPases, but can also correspond to the modification of lipid components by enzymes in the Golgi (McDermott and Mousley, 2016). This is a local mechanism of identity conversion distinct from the maturation of entire compartments, which is also affected by the dynamics of budding and fusion. To limit the number of parameters, the biochemical conversion rate ${K}_{\mathrm{m}}$ is the same for the two reactions (cis→medial and medial→trans) and is independent of the concentration. Introducing cooperativity in the conversion process is expected to increase the robustness of the selforganisation process (Vagne and Sens, 2018a).
In its simplest form, our model contains only four parameters: the rates of injection, fusion, budding, and biochemical conversion ($J,{K}_{\mathrm{f}},{K}_{\mathrm{b}},{K}_{\mathrm{m}}$). By normalizing time with the fusion rate, we are left with three parameters: $j=J/{K}_{\mathrm{f}}$, ${k}_{\mathrm{b}}={K}_{\mathrm{b}}/{K}_{\mathrm{f}}$ and ${k}_{\mathrm{m}}={K}_{\mathrm{m}}/{K}_{\mathrm{f}}$. The dynamics of the system is entirely governed by these stochastic transition rates, and can be simulated exactly using a Gillespie algorithm (Gillespie, 1977), described in Appendix 1.
Codes used for the simulation can be found in Source code 1.
Results
Meanfield description of the system
The complexity of the system prohibits rigorous analytical calculation. Nevertheless, analytical results can be obtained for a number of interesting quantities in certain asymptotic regimes where simplifying assumptions can be made. We present below some of these derivations.
Denovo formation and steadystate composition in the wellsorted limit
The composition of the system at steadystate is difficult to compute due to the fact that the exit of components from the system, through fusion with the boundaries, depends on the composition of the exiting compartments, which cannot be derived analytically in the general case. This calculation becomes straightforward if the system is well sorted and all compartments are pure. As explained below, this corresponds to the limit of high budding rate: ${k}_{\mathrm{b}}\gg {k}_{\mathrm{m}}$. In the wellsorted limit, only cis components exit through the cis boundary and only trans components exit through the trans boundary, and meanfield equations can be derived for the total amount ${N}_{\text{\mathit{c}\mathit{i}\mathit{s}}}$, ${N}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}}$ and ${N}_{\text{\mathit{t}\mathit{r}\mathit{a}\mathit{n}\mathit{s}}}$ of each species:
With $j$ the injection rate of new cisvesicles in the system and ${k}_{\mathrm{m}}$ the biochemical conversion rates, both normalized by the fusion rate, and ${\alpha}_{\text{ER}}$ (${\alpha}_{\text{TGN}}$) is the fraction of cis (trans) species promoting homotypic fusions with Golgi components in the ER (TGN). At steadystate, the total amounts of cis, medial and transspecies are fixed by the balance between influx, exchange, biochemical conversion and exit:
To estimate the typical size of compartments, we assume that for each species, a single large compartment of size $n$ interacts through budding and fusion with a number ${n}_{\mathrm{v}}$ of vesicles of the same identity (so that $N=n+{n}_{\mathrm{v}}$). The compartments size for each species then satisfies:
To limit the number of parameters, we fix the average system’s size $N={N}_{\text{\mathit{c}\mathit{i}\mathit{s}}}+{N}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}}+{N}_{\text{\mathit{t}\mathit{r}\mathit{a}\mathit{n}\mathit{s}}}$ to a set value $N=300$ in the main text, which is suitable for Golgi ministacks whose total area is of the order of 1 μm^{2} (the area of a mammalian Golgi ribbon is much larger) (Yelinek et al., 2009), corresponding to about a few hundreds of vesicles of diameter ~ 1050 nm. This is obtained by adjusting the influx to the variation of the other parameters according to:
We show in Appendix 1—figure 1 that the overall structure of the phase diagrams for the compartments size and composition is robust upon variation of the average system’s size.
The denovo formation of the system can be approached with Equation 3. Starting with an empty system with an influx $J\phantom{\rule{veryverythickmathspace}{0ex}}(={K}_{f}j)$ of cis vesicles, the system grows stochastically to reach the steadystate size after a time of order $1/{K}_{\mathrm{f}}$. The evolution of the size of the system with time obtained from numerical simulations, shown in Appendix 4—figure 1 for different sets of parameters, agrees with Equation 3. At steadystate, the influx of vesicles is balanced by an outflux of material. The exit kinetics can be quantified by looking at the fate of a pulse of cargo released from the ER at a given time. Appendix 4—figure 2 shows that the exit kinetics is exponential, in agreement with results of FRAP (Patterson et al., 2008) or pulsechase (Boncompain et al., 2012) experiments.
Although Equation 6 is only valid in the wellsorted limit, it permits a satisfactory control of the average system size over the entire range of parameters, with only a 10% variation for small budding rates (Appendix 5, Appendix 5—figure 1). At steadystate, the system size exhibits stochastic fluctuations around the mean value that depend on the parameters (Appendix 5, Appendix 5—figure 1). The large fluctuation (about 30%) for low budding rate ${k}_{\mathrm{b}}$ stem from the fact that compartments are large and transient, as explained below.
Compartment size and composition in the maturationdominated regime
For low values of the budding rate ${k}_{\mathrm{b}}$, compartments do not shed vesicles and evolve in time independently of one another in a way dictated by a balance between fusionmediated growth and biochemical maturation. Their size and composition can be approximatively calculated by assuming that compartments grow in size by fusion from a constant pool of vesicles containing the same number ${n}_{\mathrm{v}}$ of vesicle for all three identities. Calling ${n}_{\text{\mathit{c}\mathit{i}\mathit{s}}}(t)$, ${n}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}}(t)$ and ${n}_{\text{\mathit{t}\mathit{r}\mathit{a}\mathit{n}\mathit{s}}}(t)$ the amount of components of the three identities in a given compartment of total size $n={n}_{\text{\mathit{c}\mathit{i}\mathit{s}}}+{n}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}}+{n}_{\text{\mathit{t}\mathit{r}\mathit{a}\mathit{n}\mathit{s}}}$, and neglecting vesicle budding, the meanfield equations satisfied by these quantities are:
Starting with a vesicle of cis identity for $t=0:\text{}{n}_{\mathit{\text{cis}}}(0)=1$, ${n}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}}(0)={n}_{\text{\mathit{t}\mathit{r}\mathit{a}\mathit{n}\mathit{s}}}(0)=0$, the compartment size evolves linearly with time: $n(t)=1+{n}_{\mathrm{v}}t$, and the composition of each species satisfies:
The fraction of each species in the compartment is thus independent on ${n}_{\mathrm{v}}$, and reads:
These results show that in this regime, a given compartment smoothly evolves from a mostly cis to a mostly trans identity over a typical time $1/{K}_{\mathrm{m}}$. The maximum amount of medial identity is obtain for ${k}_{\mathrm{m}}t=1$ and is ${\varphi}_{\mathrm{medial}}=1/e\simeq 0.37$. Therefore, the mechanism does not lead to pure medial compartments. These analytical results are confirmed by the full numerical solution of the system in the low budding rate regime shown in Appendix 5 and discussed in detail below.
Steadystate organisation
The steadystate organisation and dynamics of the system is described in terms of the average size and purity of compartments (introduced here and detailed in the Appendix 2). The stationary size distribution of compartments is a decreasing function of the compartment size, which typically shows a power law decay at small size, with an exponential cutoff at large size (Figure 2A). This is expected for a nonequilibrium steadystate controlled by scission and aggregation (Turner et al., 2005) and means that there exist many small compartments and very few large ones. The typical compartment’s size is defined as the ratio of the second over the first moments of the size distribution. This corresponds to a size close to half the size of the exponential cutoff (see Appendix 2) beyond which it is unlikely to find a compartment (Vagne et al., 2015).
The typical compartment size varies with parameters as shown in Figure 2B. Increasing the ratio of budding to fusion rate ${k}_{\mathrm{b}}$ decreases the compartment size (Turner et al., 2005). The size depends on the biochemical conversion rate ${k}_{\mathrm{m}}$ in a nonmonotonic manner: for a given budding rate, the size is smallest for intermediate values of the biochemical conversion rate (${k}_{\mathrm{m}}\simeq 1$). The conversion rate controls the composition of the compartment and hence their fusion probability. This dependency is well explained by the simple analytical model discussed above (Equations 3,4) that approximates the system by three pure compartments interacting with a pool of vesicles (solid lines in Figure 2B).
The purity of a given compartment is defined as the (normalized) Euclidean distance of the composition of the compartment (the fractions ${\varphi}_{i}$ in the different iidentities, $i=$ cis, medial and trans) from a perfectly mixed composition. $P=0,\text{}1/2,\text{}1$ correspond respectively to compartments that contain the same amount of the three identities, the same amount of two identities, and a single identity (Appendix 2). The purity of the system is the average purity of each compartment weighted by its size, ignoring vesicles. The dependency of the average purity with the parameters is shown in Figure 2C. As for the size of the compartments, it is nonmonotonic in the biochemical conversion rate, and regions of lowest purity are found for intermediate biochemical conversion rates ${k}_{\mathrm{m}}\simeq 1$, when all three identities are in equal amount. The purity increases with the budding rate ${k}_{\mathrm{b}}$ in a sigmoidal fashion (see also Figure 3A). This can be understood by analyzing the processes involved in mixing and sorting of different identities. In a first approximation, mixing occurs by biochemical conversion and sorting occurs by budding of contaminating species, suggesting a transitions from low to high purity when the budding rate reaches the biochemical conversion rate (${k}_{\mathrm{b}}\simeq {k}_{\mathrm{m}}$). A ‘purity transition’ is indeed observed in an intermediate range of biochemical conversion rates ($0.1\phantom{\rule{thinmathspace}{0ex}}<\phantom{\rule{thinmathspace}{0ex}}{k}_{\mathrm{m}}\text{}\phantom{\rule{thinmathspace}{0ex}}10$, see Figure 2C). It is most pronounced for ${k}_{\mathrm{m}}\simeq 1$. Beyond this range, the system is dominated by one identity and compartments are always pure. The purity variation can be qualitatively reproduced by an analytical model (Appendix 3) that accounts for the competition between biochemical conversion and budding, but also includes fusion between compartments and with the exits (solid lines in Figure 2C).
In summary, three types of organisation can be found, mostly controlled by the ratio of budding over fusion rate ${k}_{\mathrm{b}}$: a mixed regime at low ${k}_{\mathrm{b}}$, where compartments are large and contain a mixture of identities, a vesicular regime at high ${k}_{\mathrm{b}}$, where compartments are made of only a few vesicles and are very pure, and a sorted regime for intermediate values of ${k}_{\mathrm{b}}$, where compartments are both rather large, and rather pure. Snapshots of these three steadystates are shown in Figure 2D. The relationship between average size and average purity of compartments is shown in Figure 2E. The existence of a welldefined intermediate regime with large sorted compartments is promoted by our assumption of ‘saturated budding’ where the budding flux of any identity only depends on the compartment size (Equation 2). If the budding flux is linear with the number of patches of a given identity, the purity transition occurs for larger values of the budding rate ${k}_{\mathrm{b}}$, where compartments are smaller (see Appendix 8).
Vesicular transport
Vesicular exchange between compartments is quantified by following the dynamics over time of passive cargo injected from the ER. Each cargo in a compartment of size $n$ has a probability $1/n$ to join a vesicle budding from this compartment. When this vesicle fuses with another compartment or a boundary, the composition difference $\mathrm{\Delta}{\varphi}_{i}$ between the receiving and emitting compartment (a number between −1 and 1) is recorded for the three iidentities ($i=$cis, medial, trans). Averaged over all budding/fusion events, this defines the enrichment vector $\overrightarrow{E}$, whose three components ${E}_{i}$ are normalized for readability so that ${\sum}_{i}{E}_{i}=1$ (see Appendix 6—figure 1 for nonnormalized enrichment). Vesicular flux is predominantly anterograde if ${E}_{\mathit{\text{cis}}}\phantom{\rule{thinmathspace}{0ex}}<\phantom{\rule{thinmathspace}{0ex}}0$ and ${E}_{\mathit{\text{trans}}}\phantom{\rule{thinmathspace}{0ex}}>\phantom{\rule{thinmathspace}{0ex}}0$, while it is predominantly retrograde if ${E}_{\mathit{\text{cis}}}\phantom{\rule{thinmathspace}{0ex}}>\phantom{\rule{thinmathspace}{0ex}}0$ and ${E}_{\mathit{\text{trans}}}\phantom{\rule{thinmathspace}{0ex}}<\phantom{\rule{thinmathspace}{0ex}}0$.
We focus on systems containing comparable amounts of each species at steadystate (${k}_{\mathrm{m}}=1$  see Equation 6). The enrichment in cis, medial and transidentities are shown in Figure 3B as a function of the budding rate ${k}_{\mathrm{b}}$. The most striking result is the high correlation between purity of the compartment and the directionality of vesicular transport. For low values of the budding rate (${k}_{\mathrm{b}}\ll 1$), compartments are well mixed (purity $\lesssim 1/2$) and the vesicular flux is retrograde, with a gain in cisidentities and a loss in transidentities. For high values of the budding rate (${k}_{\mathrm{b}}\gg 1$), compartments are pure (purity $\simeq 1$) and the vesicular flux is anterograde, with a gain in transidentities and a loss in cisidentities. The archetypal behaviors most often discussed in the literature are thus, within the limits of our model, asymptotic regimes for extreme values of the ratio of budding to fusion rates. Remarkably, the cross over between these two asymptotic behaviors is rather broad (${k}_{\mathrm{b}}\sim 0.11$) and displays a more complex vesicular transport dynamics, mostly oriented toward medial compartments.
In our model organelle, vesicle budding and fusion are biased solely by local compositions, which is subjected to irreversible biochemical conversions. The interplay between these microscopic processes gives rise to an irreversible flux of vesicles across the system. To better understand the directionality of vesicular exchanges, the net vesicular flux leaving compartments with a given composition is represented as a vector field in the compartments composition space in Figure 4. This is shown as a triangle plots in which each point corresponds to a fraction of cis, medial and trans components in a compartment, and the three vertices of the triangle correspond to pure compartments (see Appendix 2). The variation of the vesicular flux with the budding rate ${k}_{\mathrm{b}}$ (Figure 4) is consistent with that of the composition enrichment in the different identities (Figure 3B). As ${k}_{\mathrm{b}}$ increases, the flux evolves from being mostly retrograde at low ${k}_{\mathrm{b}}$ (from transrich toward less mature compartments) to being anterograde at high ${k}_{\mathrm{b}}$ (from cis to medial, and medial to transcompartments). In between, the vesicular flux is centripetal toward mixed compartments, leading to a net enrichment in medialidentity. The relationship between structure and transport is explained below, and a dissection of this fluxes with respect to the composition of the vesicles is shown in Appendix 6. The net flux leaving a particular region of the triangular phase space is proportional to the total system’s size (number of compartments times their size) with this composition, while the flux arriving at a particular region depends on the number of compartments with that composition. Both quantities are shown on Figure 4.
At low budding rate (${k}_{\mathrm{b}}\ll 1$  top row in Figure 4), the system is in the mixed regime, where compartments are large (hundreds of vesicles in size) but mixed (purity $\lesssim 1/2$). The system is dominated by compartment maturation, and the compartment distribution is spread around a path going from purely cis to purely transcompartments, consistent with a dynamics where each compartment maturates independently from the other and given by Equation 9.  see top left panel of Figure 4. The majority of the transport vesicles are emitted by transrich compartments, which concentrate most of the components of the system, leading to a retrograde vesicular flux. Note that although vesicular flux is clearly retrograde in the maturationdominated regime, as seen from the increase in cis identity and the decrease in trans identity between emitting and receiving compartment (Figure 3B), there is no net flux of vesicles from medialrich to cis compartments in this limit (Figure 4). This differs from the classical ‘cisternal progression and maturation’ mechanism sketched in Figure 1A, and is due to the fact that the compartments that contain the highest fraction of medial identity in this regime are well mixed (${\varphi}_{cis}\simeq {\varphi}_{medial}\simeq {\varphi}_{trans}\simeq 1/3$). The vesicular flux leaving such compartments is equally split into cis, medial and trans vesicles, which fuse homotypically with compartments enriched in their identity, yielding a vanishing net flux averaged over the three identities. There is nevertheless a strong retrograde flux of cis vesicles toward cis compartments, as shown in Appendix 6—figure 1. In the classical ‘cisternal maturation’ mechanism, retrograde transport is needed to recycle cisGolgi components from more mature compartments. Within the present model, this can be achieved by targeting such components to $cis$ vesicles budding from more mature compartments.
At high budding rate (${k}_{\mathrm{b}}\gg 1$  bottom row in Figure 4), the system is in the vesicular regime, where compartments are fairly pure, but rather small, with a size equivalent to a few vesicles. Compartments are exclusively distributed along the cis/medial and medial/transaxes of composition, and accumulate at the vertices of the triangle. Membrane patches that just underwent conversion bud quickly from the donor compartment and fuse with pure compartments with the appropriate medial or transidentity, explaining the clear anterograde vesicular transport. This feature is further reinforced under very high budding rate, where the vesicular flux is dominated by budding vesicles undergoing biochemical conversion before fusion, which prohibits their back fusion with the donor compartment.
In the intermediate regime (middle row in Figure 4), the distribution of mixed compartments is more homogeneous along the line describing the maturation of individual compartments  both in terms of size and number of compartments. Cisrich and transrich compartments emit comparable amount of transport vesicles, while large medialrich compartments are absent. This leads to a centripetal vesicular flux with an enrichment in medial identity (see Figure 3B). The lack of large medialrich compartment is due to the fact that they can be contaminated by either cis or transspecies and fuse with cisrich or transrich compartments. If intercompartment fusion is prohibited, the enrichment in medialidentity in the intermediate regime is abolished, with a direct transition from anterograde to retrograde vesicular transport upon increasing ${k}_{\mathrm{b}}$ (see Appendix 8 and the Discussion section).
A graphical sketch consistent with our numerical observations for the typical size and composition of different compartments and for the dominant vesicular fluxes between them in the different regimes is shown in Figure 3C.
Discussion
The genesis and maintenance of complex membranebound organelles such as the Golgi apparatus rely on selforganisation principles that are yet to be fully understood. Here we present a minimal model in which steadystate Golgilike structures spontaneously emerge from the interplay between three basic mechanisms: the biochemical conversion of membrane components, the compositiondependent vesicular budding, and intercompartment fusion. As our model does not include space, we do not refer to the spatial structure of a stacked Golgi, but rather to the selforganisation of Golgi components into compartments of distinct identities undergoing vesicular exchange. We show that directional vesicular exchanges between compartments of different identities spontaneously emerge from local biochemical interactions, despite the absence of spatial information. This is in line with recent observations that the Golgi functionality is preserved under major perturbation of its spatial structure resulting from landlocking Golgi cisternae to mitochondria (Dunlop et al., 2017). Indeed, cargo transport still proceed, and is merely slowed down by a factor of two (from 20 min to 40 min) in landlocked Golgi despite the absence of spatial proximity between cisternae. This suggests that spatial information is less crucial to the Golgi function than biochemical information. Our main result is that the steadystate organisation of the organelle, in terms of size and composition of its compartments, is intimately linked to the directionality of vesicular exchange between compartments through the mechanism of homotypic fusion. Both aspects are controlled by a kinetic competition between vesicular exchange and biochemical conversion.
The results presented in Figure 2 can be used to explain a number of experimental observations. The predicted role of the budding and fusion rates is in agreement with the phenotype observed upon deletion of Arf1 (a protein involved in vesicle budding), which decreases the number of compartments and increases their size, particularly that of transcompartments which seem to aggregate into one major cisterna (Bhave et al., 2014). A comparable phenotype has been observed upon mutation of NSF (a fusion protein), which produces extremely large, but transient, transcompartments (Tanabashi et al., 2018), thereby increasing the stochasticity of the system. This is consistent with an increase of the fusion rate according to our model (see Appendix 5 for a quantification of the fluctuations in our model). Modifying the expression of VPS74 (yeast homolog of GOLPH3) leads to an altered Golgi organisation, comparable to the ARF1 deletion phenotype (Iyer et al., 2018). Both $\mathrm{\Delta}$arf1 and $\mathrm{\Delta}$vps74 present an enlargement of the Golgi cisternae and a disruption of molecular gradients in the system, which we interpret, in the limits of our system, as a lower purity due to slower sorting kinetics following an impairment of the budding dynamics. Importantly, the correlation between the size and purity of the compartments predicted by our model is in agreement with the observation that decreasing the budding rate by altering the activity of COPI (a budding protein) leads to larger and less sorted compartments (Papanikou et al., 2015). In yeast again, overexpression of Ypt1 (a Rab protein) increases the transition rate from early to transitional Golgi, and increases the colocalization of early and late Golgi markers (Kim et al., 2016). We interpret this as a decrease of the purity of Golgi cisternae upon increasing the biochemical conversion rate by unbalancing the ratio of budding to conversion rates. This suggests that the wildtype Yeast Golgi is closed to the line ${k}_{\mathrm{m}}={k}_{\mathrm{b}}$ shown in Figure 2.
One crucial model prediction is that the size and purity of Golgi compartments should be affected in a correlated fashion by physiological perturbations. This could be tested experimentally by further exploring the phase diagrams of Figure 2B–C by simultaneously varying the ratios of biochemical conversion over fusion rate ${k}_{\mathrm{m}}$ and of budding over fusion rate ${k}_{\mathrm{b}}$. Compartment purity can be altered without changing their size by acting on ${k}_{\mathrm{m}}$, while it should be correlated with a change of size by acting on ${k}_{\mathrm{b}}$. We thus predict that the purity decrease concomitant with the size increase observed in Papanikou et al., 2015 upon impairing COPI activity (decreasing ${k}_{\mathrm{b}}$) could be reversed without a change of cisterna size upon decreasing the biochemical conversion rate ${k}_{\mathrm{m}}$, for example by impairing Ypt1 activity (as done in Kim et al., 2016). Along the same line, we predict that the decrease of purity observed by Kim et al., 2016 upon Ypt1 overexpression (which increases the early to transitional conversion rate) should not be associated to change of cisternae size (at steadystate) if altering Ypt1 does not modify the budding rate ${k}_{\mathrm{b}}$. We can predict further that the loss of purity phenotype could be reversed upon increasing the budding rate by overexpressing COPI, but that this would be accompanied by an decrease of cisternae size. This experimental proposal is represented graphically on the size and purity phase diagrams in Appendix 9—figure 1, with somewhat arbitrary arrows, to give a feel for the way our model may be used to design experimental strategies.
The directionality of vesicular transport is intimately linked to the steadystate organisation of the organelle; vesicular transport is anterograde when compartments are pure, and it is retrograde when they are mixed. This can be intuitively understood with the notion of ‘contaminating species’. If a compartment enriched in a particular molecular identity emits a vesicle of that identity, the vesicle will likely fuse back with the emitting compartment by homotypic fusion, yielding no vesicular exchange. One the other hand, an emitted vesicle that contains a minority (contaminating) identity will homotypically fuse with another compartment. If budding is faster than biochemical conversion, the contaminating species is the one that just underwent conversion, and its budding and fusion with more mature compartment leads to an anterograde transport. In such systems, compartments are rather small and pure. If biochemical conversion is faster than budding, the contaminating species is the less mature one, leading to a retrograde transport. Such systems are rather mixed with large compartments. Thus, the directionality of vesicular exchange is an emergent property intimately linked to the purity of the compartments. The medialrich compartments are special in that regard: for intermediate values of the purity (${k}_{\mathrm{b}}\sim {k}_{\mathrm{m}}$), they may be contaminated both by yet to be matured cisspecies and already matured transspecies and can fuse both with cisrich and transrich compartments. If intercompartment fusion is allowed, large medialrich compartments are relatively scarce, and emit few vesicles. On the other hand, medialvesicles emitted by cis/medial and medial/transcompartments may fuse together to form (small) medialrich compartments. For intermediate values of the budding rates, this vesicular flux dominates vesicular exchange and leads to a centripetal vesicular flux towards medialcompartments (details in Appendix 6). If intercompartment fusion is prohibited, which for instance corresponds to a situation where Golgi cisternae are immobilized (Dunlop et al., 2017), medialrich compartments are present at steadystate, the centripetal vesicular flux is less intense, and the purity transition is accompanied with a direct transition from retrograde to anterograde vesicular flux (Appendix 8). Note that vesicular transport of cargo through the Golgi is not bound to follow the net vesicular flux between compartments discussed above. Indeed, the net flux is the average of the flux of vesicles of the three identities. A given cargo follows this flux  on average  if it does not interact preferentially with membrane of particular identity. On the other hand, a cargo that is preferentially packaged into vesicles of trans identity (for example) will be transported toward transrich compartment even if the net vesicular flux is mostly retrograde.
Regardless of the details of the model, we find that systems showing a well defined polarity with well sorted cisternae exhibit anterograde vesicular fluxes, whereas systems with mixed compartments exhibit retrograde fluxes. The former dynamics is expected when biochemical conversion is the slowest kinetic process and compartments are longlived, while the latter is expected when vesicular exchange is slow and the system is composed of transient compartments undergoing individual maturation. One can relate this prediction to the difference in organisation and dynamics between the Golgi of S. cerevisiae and the more organized Golgi of higher organisms such as vertebrates. Maturation of Golgi cisternae has been directly observed in S. cerevisiae (Losev et al., 2006), with colocalization of different identity markers within single cisternae (low purity) during maturation, whereas vesicular transport phenotypes have been indirectly observed (Dunlop et al., 2017) or inferred through modeling (Dmitrieff et al., 2013) in mammalian cells. Consistent with our predictions, Golgi dynamics is one to two orders of magnitude faster in S. cerevisiae Golgi cisternae (as measured by the typical maturation rate of cisternae $\sim 1/\mathrm{min}.$ [Losev et al., 2006; MatsuuraTokita et al., 2006]) than in mammalian cells (as measured by the typical Golgi exit rate of cargo $\sim 1/20\mathrm{min}.$ to $\sim 1/40\mathrm{min}.$ Bonfanti et al., 1998; Hirschberg et al., 1998; Patterson et al., 2008; Boncompain et al., 2012]). At this point, one may speculate a link between structure and function through kinetics. A well sorted and polarized Golgi is presumably required to accurately process complex cargo. The glycosylation of secreted cargo is key to the interaction of a vertebrate cell with the immune system of the organism it belongs to Ryan and Cobb, 2012, and glycans appear to be more diverse in higher eukaryotes  which also possess a highly organized Golgi  than in unicellular eukaryotes like yeast (Wang et al., 2017). The Golgi organisation in S. cerevisiae could thus be the result of an adaptation that has favored fast transport over robustness of processing, leading to a less organized Golgi characterized by cisternal maturation and retrograde vesicular transport. Remarkably, the slowing down of Golgi transport when S. cerevisiae is starved in a glucosefree environment coincides with the Golgi becoming more polarized (Levi et al., 2010), strengthening the proposal derived from our model of a strong connection between transport kinetics and steadystate Golgi structure.
The present model is designed to account for crucial dynamical processes at play in the Golgi with a limited number of parameters. It is sufficient to yield robust selforganised structures which can reproduce a number of experimental observations, but cannot be expected to account for the full richness of Golgi dynamics, especially under severe perturbation. The model includes a single kind of exchange mechanism based on vesicles. Exchange could also proceed via the pinching and fusion of larger cisterna fragments (Pfeffer, 2010; Dmitrieff et al., 2013), or via tubular connection between compartments (Trucco et al., 2004), although the relevance of the latter to Golgi transport is still debated (Glick and Luini, 2011). Exchange mechanisms that do not establish stable connections able to relax concentration gradient between compartments by diffusion could be included as extension of our vesicle exchange mechanism (at the expense of introducing additional parameters to characterise their composition and sizedependent nucleation, scission and fusion rates). Stable Golgi tubulation without tubule detachment is observed under brefeldin A treatment, which prevent the membrane association of COPI coat. This eventually leads to the disappearance of the Golgi by quick resorption into the ER whenever a tubular connection is established between them, in a way suggesting a purely physical (tensiondriven) mechanism (Sciaky et al., 1997). According to our model, preventing COPI vesicle budding should lead to larger cisternae, as discussed above. The brefeldin A phenotype can be reconciled with our framework by invoking the fact that larger cisternae should have a higher surface to volume ratio and hence a lower membrane tension, decreasing the force required to nucleate membrane tubules and increasing the probability of direct tubular connection with the ER. Such tension effects are clearly beyond the scope of the present model, but this phenotype could be qualitatively reproduced by invoking an increase of the ER fusion parameter ${\alpha}_{\mathrm{ER}}$ with the size of the compartment due to membrane mechanics considerations.
The model does not account for the possibility of lateral partitioning between processing and exporting domains within compartments which could play an important role in regulating cisterna composition and cargo processing (Patterson et al., 2008). A detailed kinetic model of organisation based on this concept and supported by extensive quantitative live cell imaging experiments is proposed in Patterson et al., 2008 to account for the distinct composition of the different Golgi cisterna and the kinetics of cargo transport. This rapid partitioning model (RPM) is very different in scope with the model we propose here, since it focuses on the transport of lipids and cargo molecules in a preestablished Golgi apparatus with a prescribed structure, while the present model describes denovo Golgi formation and maintenance, through selforganised fusion and scission mechanisms between dynamic compartments. The two models share important similarities, such as the fact that no directionality is assumed a priori for intracisternal vesicular transport, and that export is allowed from every cisterna, which is crucial to reproduce the exponential kinetics observed in experiments. The RPM is much more precise in terms of the microscopic description of the biochemical processes, and consequently involves many parameters. Combining our simple model of Golgi selforganisation with the detailed description given by the RPM for the kinetics of lipid, enzymes and cargo within the Golgi is the logical next step to push forward our understanding of Golgi selforganisation.
In summary, we have analyzed a model of selforganisation and transport in cellular organelles based on a limited number of kinetic steps allowing for the generation and biochemical conversion of compartments, and vesicular exchange between them. Although we kept the complexity of individual steps to a minimum, our model gives rise to a rich diversity of phenotypes depending on the parameter values, and reproduces the effect of a number of specific protein mutations. We identify the concomitance of a structural transition (from large and mixed to small and pure subcompartments) and a dynamical transition (from retrograde to anterograde vesicular exchange). Within our model, anterograde vesicular transport is accompanied by many spurious events of vesicle backfusion. For very high budding rates, compartments are very pure and anterograde transport is dominated by vesicles undergoing biochemical conversion after budding from a compartment. In our model, such vesicle biochemical conversion events, which have been described in the secretory pathway as a way to direct vesicular traffic and to prevent vesicles backfusion (Lord et al., 2011), only occur when compartments are very small (high budding rate). Adding compositiondependent feedback involving specific molecular actors to tune budding and fusion rates  for instance to reduce or prevent the budding of the majority species  will displace the purity transition to lower budding rates and thus extend the regime of anterograde transport to larger subcompartment steadystate sizes. Cisternal stabilizers like GRASP or vesicles careers like Golgins, known to already be present in the ancestor of eukaryotes (Barlow et al., 2018), are obvious candidates. Although more complex models, and in particular the inclusion of spatial dependencies, are surely relevant to dynamics of the organelle, the fundamental relationship between kinetics, structure and transport highlighted by our model is a universal feature of the interplay between biochemical conversion and vesicular exchange in cellular organelles.
Appendix 1
Detailed model
General overview
The nomenclature we use in this section is presented in the Appendix 1—tables 1–3. We propose a coarsegrained, discrete model of the Golgi Apparatus. The smallest element is a vesicle, which defines the unit size of the system. Other compartments are assemblies of fused vesicles; their size is an integer corresponding to an equivalent number of vesicles that have fused to build this compartment. Each unit of membrane is of either cis, medial or transidentity. A vesicle has a unique identity, and a compartment of size $n$ is a assembly of $n$ discrete membrane patches of different identities.
We let the system selforganizes between two boundaries, the ER and the TGN, that are two biological boundaries of the Golgi (Presley et al., 1997; Klemm et al., 2009). This theoretical Golgi is selforganized as we do not impose any topological or functional constraints, such as the number of compartments or the directionality of vesicular trafficking. The evolution of the system is only dictated by three mechanisms that are budding, fusion and biochemical conversion (short description in the main text, complete description bellow). In order to keep the model as simple as possible, we neglect all mechanical and spatial dependencies. This means that we do not take into account motion of compartments in space, or mechanical properties of the membranes. Any two compartments have the possibility to fuse together, with a rate that depends on their relative composition (see below), but not on their size. An alternative model where compartments cannot fuse with one another is presented in Appendix 8.
Biochemical conversion
In the Golgi, it is known that cargo (Kelly, 1985) and cisternae (Day et al., 2013) undergo biochemical modifications over time. Changes in membrane identity are driven by maturation cascades of proteins such Rab GTPases (Suda et al., 2013). These biochemical conversions can be directly observed using live microscopy, and have been best characterized in the yeast Golgi (MatsuuraTokita et al., 2006). Biochemical conversion of membrane identity is a complex and possibly multistep process involving different kinds of enzymes. As we are mostly interested in the interplay between biochemical conversion (which mixes different identities within a single compartment) and sorting mechanisms, we coarsegrain the biochemical conversion events into simple, one step, stochastic processes: cis→medial and medial→trans. To simplify the model, we choose the same rate ${K}_{\mathrm{m}}$ for both conversions. Each patch of membrane (in a vesicle or a compartment, and whatever the surrounding patches are) has the same transition rate.
Fusion
Intracellular transport strongly relies on vesicular trafficking (Takamori et al., 2006). This involves at least two steps: the budding of a vesicle from a donor compartment (described in the next section), and its fusion with a receiving compartment. The fusion process itself requires close proximity between the receiving compartment and the vesicle, followed by the pairing of fusion actors resulting in membrane fusion. In the current model, the rate of encounter of two compartments is constant, and equal to ${K}_{\mathrm{f}}$, while the actual fusion event depends on the biochemical composition of the fusing compartments. In vivo, some of the key proteins involved in the fusion process, such as the SNAREs or tethering proteins are known to closely interact with membrane markers like Rab proteins (Cai et al., 2007). This is thought to accelerate fusion events between compartments of similar biochemical identities and decrease it between compartments of different identities, a process often called homotypic fusion (Marra et al., 2007). To take this into account, fusion is modulated by the probability that both compartments exhibit the same identity at the contact site, corresponding to Equation 1 in the main text.
The fusion rate is maximum (${K}_{\mathrm{f}}$) for two identical compartments, and vanishes between two compartments with no common identity. The exchange of vesicles between compartments of different identities (Glick and Luini, 2011), or the fusion of compartments of different identities (Marsh et al., 2004), are sometimes regarded as heterotypic fusion. Such processes are allowed within our homotypic fusion model, provided the different compartments share at least some patches of similar identities.
Budding
Budding is the formation of a new vesicles from a large compartment. We assume that each budding vesicle is composed of a single membrane identity, which is consistent with the high specificity of the in vivo budding machinery (Bonifacino and Glick, 2004). In our model, the flux of vesicles budding from a compartment depends on the size and composition of the compartment. Compartments smaller than 2 cannot bud a vesicle. Compartments of size 2 can split into two vesicles. For larger compartments, we consider a general budding flux for vesicles of identity $i$ (for $i=$cis, medial and trans) of the form:
with $n$ the size of the compartment’, ${K}_{\mathrm{b}}$ a constant and $f({\varphi}_{i})$ a function of the composition of the compartment.
One expects the budding flux to follow a MichaelisMenten kinetics; linear with the number of patches of a given identity ($f(\varphi )=\varphi $) at small concentration, and saturating to a constant at high concentration (Vagne and Sens, 2018b). This scenario is consistent with the fact that budding proteins interact very dynamically with the membrane, attaching and detaching multiple times before budding a vesicle (Hirschberg et al., 1998). Previous theoretical works suggest that vesicular sorting of different membrane species is more efficient in the saturated regime (Dmitrieff and Sens, 2011). Since one of the main features of the Golgi is to segregate different biochemical species into different compartments, we assume a very high affinity between membrane patches and their budding partners, setting $f({\varphi}_{\mathrm{i}})=1$ if ${\varphi}_{\mathrm{i}}>0$ and $f({\varphi}_{\mathrm{i}}=0)=0$ (corresponding to Equation 2 of the main text). The budding flux for a given identity thus depends on the total size of the compartment rather than on the amount of that particular identity. This leads to a total budding flux ${K}_{\mathrm{b}}\times n\times {n}_{\mathrm{id}}$ for a compartment carrying ${n}_{\mathrm{id}}$ different identities. The results of simulations with this linear budding scheme are discussed in Appendix 8.
Exits
In cells, the Golgi is placed between two intracellular structures that are the ER and the TGN (Presley et al., 1997; Klemm et al., 2009). Fluxes of material leaving the Golgi thus include retrograde fluxes toward the ER, carrying immature components such as cis and medialGolgi enzymes or recycling ER enzymes, and anterograde fluxes toward the TGN of mature components, such as processed cargo that properly underwent all their posttranslational maturation steps (Boncompain and Perez, 2013). The exiting fluxes are accounted for by allowing the different compartments in the system to fuse with the boundaries. All compartments, from vesicles to the largest ones, can fuse homotypically with the ER or the TGN to exit the system. Thus, these boundaries are modelled as stable compartments, containing a fraction ${\alpha}_{\mathrm{ER}}$ of cis components for the ER and ${\alpha}_{\mathrm{TGN}}$ of trans components for the TGN. These fraction are simply written $\alpha $ when they are assumed to be the same. This allows immature (cis) compartments to undergo retrograde exit, and mature (trans) components to undergo anterograde exit. In the main text we focus on the case $\alpha =1$. Lower values of $\alpha $ reduce fusion with the boundaries and increase the residence time of components in the system. This increases the average size of compartments, and increases fluctuations in the system, as large compartments exiting the system induce large fluctuations in the instantaneous size and composition of the system. The impact of $\alpha $ is studied in some details in Appendix 7.
Injection of components
Influx to the Golgi
The Golgi receives material from different compartments like the endosomal network, lyzosomes, etc… (Boncompain and Perez, 2013). As we are primarily interested in characterizing the relationship between the structure and dynamics of the Golgi, and in relating these to the rates of individual biochemical conversion and transport processes, we focus here on the secretory role of the Golgi, and only account for the incoming flux of immature components coming from the ER. We define a rate $J$ of injection of cisvesicles in the system. As the parameters of the system are varied, the injection rate is varied as well in order to keep the total system’s size to an almost constant value. This constraint is only approximately enforced using Equation 6, which is only strictly valid when all compartments are pure.
The impact of the system’s size on the structure of the Golgi is shown on Appendix 1—figure 1A. Increasing the total size increases the number of compartments and hence the total fusion flux between compartments. Consequently, compartments are larger in larger systems, and tend to be (slightly) less pure, as fusion increases mixing. In the main text, we restrict ourselves to a system size of $N=300$. One should remember that Equation 6, used to predict the system’s size $N$, assumes that compartments are pure (perfectly sorted). Appendix 1—figure 1B, shows that this assumption fails to predict $N$ for systems where the budding rate is low compared to the fusion rate (highly interacting compartments). In this regime, compartments have a great probability to fuse together, creating hybrid cis/medial/transstructures. This makes medialpatches sensitive to the interaction with the ER and the TGN, and creates an exit flux for these patches that is not observed in a pure regime. Such systems exhibit a larger exit flux and thus a lower $N$. Note this is only true for biochemical conversion rates lower than the fusion rate, as the system is saturated with transpatches for high ${k}_{\mathrm{m}}$ (see Equation 4).
Transport of passive cargoproteins
The direction of vesicular transport is assessed by following the transport of passive cargo injected from the ER as part of incoming cisvesicles. Cargo molecules are passive in the sense that their probability of joining a vesicle is insensitive to the vesicle membrane identity: when a compartment of size $n$ buds a vesicle, each cargo in this compartment has a probability $1/n$ to join the budding vesicle. In the simulation, the number of cargo molecules in the system is kept to a fixed number (typically 20) by injecting a new one each time one leaves the system trough fusion with the boundaries.
Gillespie algorithm and software implementation
All rates in the system are normalized by the fusion rate: ${k}_{\mathrm{b}}={K}_{\mathrm{b}}/{K}_{\mathrm{f}}$, ${k}_{\mathrm{m}}={K}_{\mathrm{m}}/{K}_{\mathrm{f}}$ and $j=J/{K}_{\mathrm{f}}$. The dynamics of the system can be exactly simulated (to yield the correct trajectories, stochastic noise, etc…) using a Gillespie algorithm (Gillespie, 1977), also known as stochastic simulation algorithm. After choosing the initial state of the system (see below), the algorithm is implemented as follows:
The rates of all possible events (be it the input of a new vesicle in the system, the fusion between two compartments or between a compartment and one boundary, the budding of a vesicle, or the biochemical conversion of a membrane patch) are calculated. The time interval $\delta t$ for the next reaction to occur is randomly picked from an exponential distribution with mean $1/{R}_{\mathrm{TOT}}$, where ${R}_{\mathrm{TOT}}$ is the sum of the reaction rates of all possible reactions. The reaction that actually occurs is randomly picked with a probability equal to the rate of this reaction, normalized by ${R}_{\mathrm{TOT}}$.
The state of the system is modified according to the picked reaction. The rates of all events in the modified system are calculated, and the current time of the simulation is incremented by $\delta t$.
The program comes back to step 1 IF the predefined maximum number of iterations has not been reached, AND there is at least one reaction that has a nonzero rate (${R}_{\mathrm{TOT}}\ne 0$).
With the Gillespie algorithm, there is no direct dependency between the number of simulation steps and the actual time that is simulated. The physical time between two consecutive steps decreases with increasing reaction rates. Since we are interested by steadystate quantities, we average over a fixed number of steps rather than a given physical time. As shown in Appendix C, the steadystate of the system is reached after a physical time of order $1/{K}_{\mathrm{f}}$. Depending on the parameters, the steadystate is reached after between 10^{3} and 10^{6}simulation steps. To characterise the steadystate, we disregards the transient regime and data are recorded after 10^{6} steps. The full simulation typically last at least 10^{7} steps. This arbitrary amount is a good compromise between computation time and the need to accumulate sufficient amount of data to obtain enough statistics on all the measured quantities. In practice, the time needed to reach steadystate can be shortened by starting the simulation from a vesicular system with the predicted amount of cis, medial and transspecies (with Equation 4).
Because simulation steps have different durations, one should be careful when computing time averages. Two different approaches can be used. Either we weight each configuration by the duration between two consecutive timesteps, or we resample data to get a fixed timestep between observations, $\mathrm{\Delta}t$. We checked that both procedures give the same results, and adopted the second approach for the analysis. As the steadystate of the system is reached after a physical time of order $1/{K}_{\mathrm{f}}$ , we take $\mathrm{\Delta}t=1/{K}_{\mathrm{f}}$.
Appendix 2
Quantification and statistical analysis
Computation of the purity
The purity of a compartment $P$ is defined such that its value 0 is for a perfectly mixed compartment containing the same amount of cis, medial, and transspecies, and it is equal to 1 for a pure compartment containing a single species. With ${\varphi}_{i}$ the fraction of the species $i$ in the compartment, the purity $P$ is defined as:
On a triangular composition space where each corner corresponds to a pure compartment, $P$ is the distance from the center of the triangle, see Appendix 2—figure 1. When we show snapshots of the system, each compartment is represented as a sphere whose area is proportional to the compartment size (defined as an equivalent number of vesicles). The composition of a compartment is represented as a color following the color code shown in Appendix 2—figure 1. The purity of the system is a global average (over time and over all compartments) of the purity of individual compartments, weighted by their sizes and ignoring vesicles (that are always pure as they are composed of one unique identity).
Computation of the typical size
Among other mechanisms, this theoretical organelle selforganizes by budding and fusion of components. In that sense, it is a scissionaggregation structure which should follow the laws dictating the behavior of this class of systems. One of these laws is the fact that the size distribution for small compartments should follow a powerlaw. Because of the scission, compartments cannot grow indefinitely and the powerlaw ends by an exponential cutoff (Turner et al., 2005; Vagne et al., 2015). Thus, the size distribution ${c}_{n}$ of compartments of size $n$ should, on a first approximation, follow this general formulation:
where ${n}_{c}$ is the cutoff size and $\beta $ is an exponent that has been calculated to be $\beta =3/2$ in a similar system (Turner et al., 2005). There are multiple ways to characterize the average size of the distribution, using ratios of moments $k+1$ over $k$ of the distribution:
It turns out that for ${n}_{c}\gg 1$, ${\sum}_{n=1}^{\mathrm{\infty}}{n}^{k\beta}{e}^{n/{n}_{c}}\sim 1$ if $1+k\beta \phantom{\rule{thinmathspace}{0ex}}<\phantom{\rule{thinmathspace}{0ex}}0$ and ${\sum}_{n=1}^{\mathrm{\infty}}{n}^{k\beta}{e}^{n/{n}_{c}}\sim {n}_{c}^{1+k\beta}$ if $1+k\beta \phantom{\rule{thinmathspace}{0ex}}>\phantom{\rule{thinmathspace}{0ex}}0$. In order to have ${\u27e8n\u27e9}_{k}\sim {n}_{c}$, we need to choose the exponent $k$ such that $k\phantom{\rule{thinmathspace}{0ex}}>\phantom{\rule{thinmathspace}{0ex}}\beta 1=1/2$ in the present case. In the main text, we adopt $k=1$ and define the characteristic size of the distribution as:
As shown in Figure 2A of the main text, the calculated average is in good agreement with simulations’ data.
Computation of the mean enrichment
To discriminate between the two models of Golgi’s dynamics that we can find in the literature, we need to quantify whether the vesicular transport is anterograde or retrograde. Indeed, the 'Vesicular transport’ model assumes that cargo is transported sequentially from cis to medial to transcompartments while resident enzymes remain in place, meaning that the vesicular flux is anterograde. On the other hand, 'Cisternal maturation’ models assume that cargo remain inside cisternae, and resident enzymes are recycled, thus requiring a retrograde vesicular flux. One way to measure this flux in our simulations is to follow passive cargo molecules and quantify whether they move to more or less mature compartments, as they get carried by vesicles. To do so, we record all events that affect cargo molecules. Every time a compartment buds a cargo, we store the composition of this compartment and compare it to the one in which the vesicle later fuses. Both compositions are a vector $\overrightarrow{C}$, with three components that are the fraction ${\varphi}_{i}$ ($i$ equals cis, medial and trans) of the donor and acceptor compartments. Defining $\overrightarrow{{C}_{d}}$ the composition of the donor compartment, and $\overrightarrow{{C}_{a}}$ the composition of the acceptor, we can compute the enrichment $\overrightarrow{E}$ as:
The sum of both $\overrightarrow{C}$ components equals 1 (as they are fractions of each identities), and the sum of $\overrightarrow{E}$ components equals $\overrightarrow{0}$, which simply means one cannot gain in fraction of any identity without loosing the same amount of the others. We can now compute $\u27e8\overrightarrow{E}\u27e9$ to calculate the mean enrichment in cis, medial and transspecies, between a budding event and the next fusion event.
However, and because of back fusion events (fusion into the same compartment that previously budded the vesicle), $\u27e8\overrightarrow{E}\u27e9$ components can be close to 0. This is particularly true for pure, sorted systems. In that case, a budded vesicle has the same identity as its donor compartment, and thus has a great probability to fuse back with the same compartment or one of very similar composition. That is why, nontreated data do not allow to discriminate well between an anterograde and a retrograde regime when the purity (and thus ${k}_{\mathrm{b}}$) is high (Appendix 6—figure 1). As we are primarily interested in the sign of these vectors’ components, we normalize $\overrightarrow{E}$ using the ${L}_{1}$ norm (${\sum}_{i}{E}_{i}=1$).
Computation of the vesicular flux vector field
The vesicular flux discussed in the previous section can be directly displayed on the triangle of composition we introduced in Appendix 2—figure 1. Instead of computing $\u27e8\overrightarrow{E}\u27e9$ for all $\overrightarrow{E}$, we can bin data with respect to the donor compartment composition (typically triangular bins of size $\mathrm{\Delta}\varphi =0.1$) and calculate the average enrichment for each bin of donor compartments. To get rid of back fusion effects in this quantification, we remove the vectors $\overrightarrow{E}$ for which all components are smaller (in absolute values) than the binning meshgrid. For each binned composition we can now compute the mean enrichment vector, and plot this vector on the 2D triangular composition space. To emphasize the dominant fluxes in this vector field, the opacity of vectors is set proportionally to the flux of transported cargo per unit time (normalized by the total amount of cargoproteins in the system). This can be seen on Figure 4  main text.
Appendix 3
Analytical model for the purity of compartments
We want to build a simple model to explain how the purity of the system is affected by the parameters. This is challenging as all events (budding, fusion and biochemical conversion) influence the purity of a compartment. We seek to determine the purity transition, namely the values of ${k}_{\mathrm{b}}$ and ${k}_{\mathrm{m}}$ for which a pure system becomes impure.
We make the following approximations, based on the assumption that the system is almost pure: (i) transcompartments are pure as they cannot be contaminated by biochemical conversion. (ii) cis and medialcompartments share the same purity. (iii) The total amount of the different species ${N}_{\text{\mathit{c}\mathit{i}\mathit{s}}}$, ${N}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}}$ and ${N}_{\text{\mathit{t}\mathit{r}\mathit{a}\mathit{n}\mathit{s}}}$ follows the steadystate repartition of cis, medial and transspecies in a pure limit, given by Equation 4. Most of the components of the system are within compartments.
The purity of a compartment  defined in the main text  is $P=\sqrt{\left(3/2{\sum}_{i}{({\varphi}_{i}1/3)}^{2}\right)}$. We consider that ciscompartments (medialcompartments) are contaminated by an average fraction $\u27e8\varphi \u27e9\ll 1$ of medialpatches (transpatches), while transcompartments are pure. The purity of cis and medialcompartments is thus ${P}_{\mathit{\text{cis, medial}}}=\sqrt{13\u27e8\varphi \u27e9(1\u27e8\varphi \u27e9)}$ and ${P}_{\text{\mathit{t}\mathit{r}\mathit{a}\mathit{n}\mathit{s}}}=1$. Weighting these with the fraction of the different species in the system (from Equation 4) gives
To compute the average contamination $\u27e8\varphi \u27e9$, we consider the different events that affect the purity of a ciscompartment of size $n$ contaminated by a fraction $\varphi $ of medialpatches, surrounded by a number ${n}_{\mathrm{v}}$ of cisvesicles, and a number ${n}_{\mathrm{v}}^{\prime}$ of medialvesicles. These events are sketched in Appendix 3—figure 1, and listed below:
Budding of cispatches (increases $\varphi $)
Budding of medialpatches (decreases $\varphi $)
Fusion of cisvesicles (decreases $\varphi $)
Fusion of medialvesicles (increases $\varphi $)
Conversion of cispatches (increases $\varphi $)
Conversion of medialpatches (decreases $\varphi $)
Fusion with a medialcompartment ($\varphi \to \frac{1}{2}$)
Fusion with the ER ($\varphi \to 0$)
The rate $p(i)$ (normalized by the fusion rate) for each mechanism, and the extent $\delta \varphi $ to which it affects $\varphi $, are:
We can now compute the temporal variation of $\varphi :\dot{\varphi}=\sum p(i)\delta \varphi (i)$. We concentrate on the limit of large compartments: $n\gg 1$. For simplicity, we also assume that $n\gg ({n}_{\mathrm{v}}{n}_{\mathrm{v}}^{\prime})$, so that the net contribution of vesicle fusion to the purity variation, which scales like $({n}_{\mathrm{v}}^{\prime}{n}_{\mathrm{v}})/n$, is negligible. Neglecting stochastic fluctuations, the meanfield evolution of the fraction of contaminating species is:
Below, we discuss separately the case where $\alpha =1$ (fast fusion with the ER boundary  discussed in the main text), and the case where $\alpha =0$ (no exit through the ER  discussed in Appendix 7).
$\alpha =1$. A pure compartment; $\varphi =0$, is stable ($\dot{\varphi}{}_{\varphi =0}\phantom{\rule{thinmathspace}{0ex}}<\phantom{\rule{thinmathspace}{0ex}}0$) when $k}_{\mathrm{m}}\phantom{\rule{thinmathspace}{0ex}}<\phantom{\rule{thinmathspace}{0ex}}{k}_{\mathrm{b}$. If $k}_{\mathrm{m}}\phantom{\rule{thinmathspace}{0ex}}>\phantom{\rule{thinmathspace}{0ex}}{k}_{\mathrm{b}$, the compartment is contaminated by a stable fraction of medialpatches given by the single stable root of Equation 17 (see Appendix 3—figure 1B). The purity transition lines shown in Figure 2C of the main text are obtained by inserting this root in Equation 16.
$\alpha =0$. If ${k}_{\mathrm{b}}<{k}_{\mathrm{m}}$, the single root of Equation 17 is for $\varphi =1/2$ and is stable. The system is always impure with mixed cismedial compartments. If ${k}_{\mathrm{b}}>{k}_{\mathrm{m}}$, pure compartments ($\varphi =0$) are stable at the meanfield level, but there exist an unstable fixed point for small values of $\varphi $. Stochastic fluctuations may bring the system passed the unstable fixed point toward the stable mixed solution $\varphi =1/2$. For ${k}_{\mathrm{b}}>{k}_{\mathrm{m}}+\frac{1}{8}$, $\varphi =1/2$ is an unstable fixed point and the only stable solution is for $\varphi =0$ (see Appendix 3—figure 1C). The purity transition is thus with the range ${k}_{\mathrm{m}}<{k}_{\mathrm{b}}<{k}_{\mathrm{m}}+\frac{1}{8}$.
A more precise characterization of the transition in the case ${k}_{\mathrm{b}}>{k}_{\mathrm{m}}$ requires to include stochastic effects. The different stochastic events are associated to very different fluctuations amplitudes. Budding and fusion of vesicles and biochemical conversion of membrane patches represent small fluctuations for a large compartment, while fusion between compartments and with the boundaries represent a very large alteration of the purity. In what follows, we propose a simplified treatment of this process, focusing on the case $\alpha =0$ for simplicity, and neglecting the contribution of vesicle fusion as discussed above. For a compartment of size $n$ contaminated by ${n}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}}$ patches, we identify 2 types of mechanisms (depending on their amplitude). First, a smooth drift in the contamination due to biochemical conversion (with a rate ${k}_{\mathrm{m}}*(n{n}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}})$) and budding (with a rate ${k}_{\mathrm{b}}n$). Second, a jump in the contamination due to fusion with the neighbor compartment (rate ${n}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}}(n{n}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}})/{n}^{2}$). We compute below the mean contamination, assuming that the ciscompartment spends a time ${\tau}_{\mathrm{conta}}$ undergoing small fluctuations $\delta \varphi \simeq 0$ close to the pure state, before fusing with a medialcompartment of composition $1\delta \varphi $. The resulting compartment with $\varphi \simeq 1/2$ then spends a time ${\tau}_{\mathrm{deconta}}$ undergoing a decontamination process, which depends on the balance between budding and biochemical conversion.
We first calculate the average value $\delta \varphi $ of the contamination around $\varphi =0$ due to biochemical conversion and budding, disregarding intercompartment fusion. The probability $\mathbb{P}({n}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}})$ of finding ${n}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}}$ in the cisterna satisfies:
In the limit $n\gg {n}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}}$, this becomes:
The steadystate solution is
and the average contamination $\delta \varphi $ due to budding and biochemical conversion is:
These fluctuations in composition allow fusion between compartments, at a rate $\delta \varphi (1\delta \varphi )\simeq \delta \varphi $ (for $\delta \varphi \ll 1$). The average waiting time for an intercompartment fusion event is thus ${\tau}_{\mathrm{conta}}=1/\delta \varphi $.
After intercompartment fusion, the compartment is at $\varphi =1/2$, which is an unstable fixed point according to Equation 17. Any small fluctuation in composition leads to a smooth decontamination that satisfies $\dot{\varphi}=({k}_{\mathrm{m}}{k}_{\mathrm{b}})(12\varphi )$, disregarding intercompartment fusion. The decontamination time ${\tau}_{\mathrm{deconta}}$ and the average contamination during the process ${\u27e8\varphi \u27e9}_{\mathrm{deconta}}$, can be estimated by integrating $\dot{\varphi}$, from $\varphi =1/21/n$ to $\varphi =0$:
We can now calculate $\u27e8\varphi \u27e9$ as a temporal average of $\varphi $:
When $\mathrm{log}(n/2)$ is close to 1 (typically true for the range of $n$ that is interesting here, $n\sim 100$), $\u27e8\varphi \u27e9$ can be simplified to:
Injecting this result into Equation 16 gives the approximate purity boundary for ${\alpha}_{\mathrm{ER}}=0$:
Despite numerous simplifications, the model gives good predictions on the purity, both for ${\alpha}_{\mathrm{ER}}={\alpha}_{\mathrm{TGN}}=1$ (Figure 2  main text) and for ${\alpha}_{\mathrm{ER}}=0$, ${\alpha}_{\mathrm{TGN}}=1$ (Appendix 7—figure 1). Depending on the normalized biochemical conversion rate ${k}_{\mathrm{m}}$, we can discriminate three regions in the purity phase diagram:
${k}_{\mathrm{m}}\gg 1$, the system is saturated in transspecies and is thus always pure.
${k}_{\mathrm{m}}\sim 1$, compartment’s contamination occurs by biochemical conversion, matured species are sorted by budding. The purity transition occurs for ${k}_{\mathrm{b}}\simeq {k}_{\mathrm{m}}$ (compartments are impure system if ${k}_{\mathrm{b}}<{k}_{\mathrm{m}}$).
${k}_{\mathrm{m}}\ll 1$, the system is saturated by cis and medialcompartments, and the transition of purity occurs for ${k}_{\mathrm{b}}\ll 1$: the fusion rate is dominant in this regime. Thus interaction with the boundaries or other compartments are crucial to understand the transition of purity. For ${\alpha}_{\mathrm{ER}}=1$, fusion with the ER decontaminates the system, thus the transition occurs for ${k}_{\mathrm{b}}<{k}_{\mathrm{m}}$. The transition becomes independent of ${k}_{\mathrm{b}}$ for very small biochemical conversion rates, where the system is saturated in cisspecies. For ${\alpha}_{\mathrm{ER}}=0$, purifying the system by recycling compartments is not possible, and intercompartment fusion decreases the purity, thus the transition of purity requires large values of the budding rate ${k}_{\mathrm{b}}>{k}_{\mathrm{m}}$, as given by Equation 25.
Appendix 4
de novo formation and cargo export kinetics
Transient evolution toward the steadystate: De novo formation
The evolution of the system’s size with time for different values of the dimensionless budding rate ${k}_{\mathrm{b}}$ in shown in Appendix 4—figure 1. The organelle forms denovo, and reaches the steadystate after a time of order a few $1/{K}_{\mathrm{f}}$. The steadystate is reached under any parameter values, but the fluctuations around the steadystate depend on the parameter, as discussed in Appendix 5. The analytic calculation presented in Equations 3 (red line in Figure 1 above) is very accurate for high budding rate, but gives a 10% error for low budding rate.
Transit kinetics of passive cargo
The cargo used as marker of vesicular transport directionality can be used to quantify the transit kinetics of a passive cargo across the selforganised organelle. Appendix 4—figure 1 shows the amount of cargo present in the system as a function of time during a typical ‘pulsechase’ experiment, in which a given number of cargo molecules is allowed to exit the E.R. at $t=0$, by joining the vesicles constituting the influx $j$. The amount of cargo in the system initially increases, then decreases in a closetoexponential fashion. The export kinetics resemble the one observed by iFRAP (Figure 2i of Patterson et al., 2008) or using the RUSH technics (Figure 3e of Boncompain et al., 2012).
Appendix 5
Detailed organisation
The fluctuations of the system around its steadystate are characterized by computing the temporal standard deviation of the total system’s size and purity, shown in Appendix 5—figure 1. Fluctuations decrease as the budding rate ${k}_{\mathrm{b}}$ increases, owing to the fact that compartments are on average smaller for larger budding rate, so that the removal of a compartment by fusion with the boundaries has a smaller impact on the state of the system. We note in Appendix 5—figure 1, that the total size of the system depends on the budding rate ${k}_{\mathrm{b}}$, despite the fact that the influx $j$ is varied according to Equation 6 to limit variations of the system’s size. This is due to the fact that Equation 6 is valid only in the limit of pure compartments (large ${k}_{\mathrm{b}}$). For smaller values of the budding rate, medialcompartments contaminated by cis and transspecies and may thus exit the system by fusing with the boundaries, which decreases the total size for low ${k}_{\mathrm{b}}$. For the regime we are interested in, namely ${k}_{\mathrm{m}}\sim \alpha \sim 1$, this phenomenon is practically negligible, but its impact on the average size is more severe for lower ${k}_{\mathrm{m}}$ (Appendix 1).
We also characterize the variation of the size of compartments and their abundance as a function of their composition. In the main text (Figure 4), variations of the total size (number of membrane patches) and the number of compartments is shown as a function of their composition, in the composition triangular space defined in Appendix 2. This is obtained by binning the composition space in such a way that the characterization is both relatively precise and statistically significant (binning $\varphi $ with a mesh of $1/30$). The size and number of compartments are averaged over the simulation time (once the steadystate is reached) for each composition bin. These results are reproduced in Appendix 5—figure 1, where the average size of compartment, computed as the ratio of average size over average number of compartments, is also shown.
For low values of the budding rate ${k}_{\mathrm{b}}$, both the size and the number of compartments follow the theoretical line obtained from the ‘cisternal maturation’ limit (Equation 9). Compartments also grow linearly in time, in agreement with the linear prediction of Equation 7. Their sizes are shown in Appendix 5—figure 1, as a function of their theoretical lifetime. The sizes are computed by measuring the average compartment size for each composition bin along the theoretical line. The compositions are then used to estimate the lifetime of compartments using Equation 9.
Appendix 6
Detailed characterization of vesicular transport
Quantification of the backfusion
In the main text, we discussed the importance of backfusion when the budding rate ${k}_{\mathrm{b}}$ increases. To determine the likelihood of a vesicle to fuse back with the compartment it budded from, or a compartment of similar composition, we analyzed the trajectories of passive cargoproteins.
Back fusion is defined as a vesicular transport event where the donor and acceptor compartments fall within the same composition bin. The fraction of the total vesicular flux leaving a compartment that undergoes back fusion is shown as a function of the composition of compartments in a triangular composition space in Appendix 6—figure 1. As expected, backfusion is dominant near the vertices of the composition triangle, where compartments are pure and mostly bud vesicles of identity similar to their own identity. These vesicles are very likely to fuse back with the donor compartment by homotypic fusion. The back fusion flux increases with the budding rate ${k}_{\mathrm{b}}$, going from $\sim [90\mathrm{\%},\text{}10\mathrm{\%},\text{}20\mathrm{\%}]$ for respectively [cis, medial and trans]compartments when ${k}_{\mathrm{b}}={10}^{2}$, to $\sim 100\%$ when ${k}_{\mathrm{b}}={10}^{2}$ for all compartments.
The mean enrichment (composition difference between donor and acceptor compartments) seen by a transport vesicle is plotted in the main text in a normalized fashion. The nonnormalized enrichment is plotted in Appendix 6—figure 1. The mean enrichment for all three identities vanishes when ${k}_{\mathrm{b}}$ is large. This is explained by the contribution of back fusion. When ${k}_{\mathrm{b}}$ is large, the compartments are very pure, and vesicles have a high probability to undergo back fusion, leading to a small value of the net enrichment.
Vesicular transport by composition
Figure 4 of the main text shows the average flux of vesicular transport as a vector field on the triangular composition phase space. The main conclusion is that the vesicular flux is retrograde (toward ciscompartments), when the compartments are mixed (low budding rate ${k}_{\mathrm{b}}$), anterograde when the compartments are pure (high ${k}_{\mathrm{b}}$), and centripetal (toward compartments of medialidentity) for intermediate budding rates. In the main text, only the total flux is described, which is the sum of the fluxes of vesicles of the three different identities. To better understand this net enrichment, the total flux can be dissected into the flux of of cis, medial and transvesicles. These results are shown in Appendix 6—figure 1. For low budding rate (${k}_{\mathrm{b}}\ll 1$), compartments are mixed with all three identities, with a composition distributed around a trajectory given by Equation 9 and shown in Appendix 6—figure 1. These mixed compartments emit vesicles of all identities. The cisvesicles predominantly fuse with pure ciscompartments, which always exist at steadystate, since pure cisvesicles are continuously being injected in the system. On the other hand, medial and transvesicles fuse with cis/medial and cis/transcompartments, as pure medialcompartments are absent, and pure transcompartments are few and transitory. In this limit, the net enrichment is completely dominated by the retrograde transport of cisvesicles, and is positive in cisidentity and negative in transidentity.
For intermediate budding rate (${k}_{\mathrm{b}}\sim 0.1\to 1$), and as ${k}_{\mathrm{b}}$ increases, components can be efficiently sorted as they mature, since the budding and biochemical conversion rates are of the same order. Medial compartments are rather unstable, as they can fuse with both cisrich and transrich compartments. Their reformation involves medial vesicles leaving cisrich and transrich compartments to fuse with each other or small medialrich compartments. Appendix 6—figure 1, shows that this is the main contribution to the overall vesicular flux, thereby defining the centripetal net vesicular flux shown in Figure 4  main text. Note that they primarily exit transcompartments for ${k}_{\mathrm{b}}\sim 0.1$, and ciscompartments for ${k}_{\mathrm{b}}\sim 1$, denoting a transition from retrograde to anterograde vesicular fluxes.
For high budding rate (${k}_{\mathrm{b}}\gg 1$), compartments are very pure, and most budding events lead to vesicle back fusion from the donor compartment (Appendix 6—figure 1). A new phenomenon can be observed, which is the anterograde transport of cis and medialvesicles fusing with medial and transcompartments, respectively. This is due to vesicle biochemical conversion after budding, as discussed below.
Role of vesicles biochemical conversion
The biochemical conversion of vesicles after their budding from an immature compartment has been presented as a mechanism to prevent back fusion with the donor compartment and to promote anterograde vesicular transport (see Discussion  main text). This mechanism is naturally included in our model, as vesicles have the same biochemical conversion rate than any membrane patch that belong to bigger compartments. One can notice on Appendix 6—figure 1, that there is a major increase of the anterograde cistomedial and medialtotrans vesicular fluxes between ${k}_{\mathrm{b}}=10$ and $k}_{\mathrm{b}}={10}^{2$. This is due to the biochemical conversion of the vesicles after their budding. The budding vesicle has the same identity than the donor compartment, but undergoes biochemical conversion before undergoing back fusion and fuses with a more mature compartment.
Appendix 6—figure 1, shows the net vesicular flux (for the same simulation than in the main text) ignoring the events where vesicles undergo biochemical conversion after their budding. The net vesicular flux vanishes in the vector field for high budding rates when these events are removed. We interpret this result considering dimers and trimers are dominant when $k}_{\mathrm{b}$ is large. Disregarding vesicle biochemical conversion, such compartments emit as many immature vesicles fusing with immature compartments as mature vesicles fusing with mature compartments, yielding a vanishing net vesicular flux. However, this is not true for the net enrichment which still displays an anterograde signature of the vesicular transport. Indeed, the explanation of the vesicular transport considering the minority identity (see Discussion section in the main text) is still relevant in this regime; compartments larger than three vesicles can generate an anterograde transport by budding patches of membrane that just underwent biochemical conversion. Consequently, the net enrichment displays the characteristics of an anterograde vesicular flux even in the absence of vesicle biochemical conversion.
Appendix 7
Boundaries composition
Our system contains two boundaries: the cis face (the ER) and the trans face (the TGN). cis vesicles are injected via the cis face, and compartments may exit the system by fusing with either boundary. The way the different compartments fuse with the boundaries has a strong impact on their lifetime and average size. In the model, this is characterized by the parameter $\alpha $, which represents the fraction of the boundary composed of species able to elicit homotypic fusion of Golgi’s components, that is, the fraction of cis species in the ER and of trans species in the TGN, and might be different for the two boundaries. In the main text, we focus on the case $\alpha =1$ (arguably the maximum possible value) for both cis and trans boundaries. Here, we discuss the impact of this parameter on the structure of the steadystate, and we compare the results of simulations with the analytical model derived in Appendix 3.
Decreasing this parameter has two major impacts on the steadystate organisation: it decreases the biochemical conversion rate for which the steadystate fraction in cis, medial and transspecies are equal and it increases the residence time of compartments in the system. This increases their size, as they have more time to aggregate, and increases the fluctuations in the system, as larger compartments have a stronger impact each time they exit or fuse together.
The impact of $\alpha $ on the fraction of cis, medial and transspecies in the system is rather straightforward if compartments are pure (Equation 4). In this case, the exit flux can be exactly computed, and Equation 4 shows that the different species are in equal amount at steadystate when ${k}_{\mathrm{m}}=\alpha $. If ${k}_{\mathrm{m}}\gg \alpha $, the membrane patches have time to undergo biochemical conversion before exiting the system, which is then dominated by transspecies. The purity of the compartment is high, but this is a rather uninteresting limit as the system is dominated by a single species. If ${k}_{\mathrm{m}}\ll \alpha $, compartments are recycled out of the system before transmembrane patches appear, and the system is dominated by cis and medialspecies, which are of equal amount at steady state if the system is perfectly sorted.
The same conclusion can be reached in the slow budding regime when compartments are not pure. Let’s consider a cisrich compartment of size $n$ contaminated by a fraction ${\varphi}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}}$ of medialspecies, and let’s disregard, for simplicity, budding and fusion with other compartments. The number of medialpatches increases by biochemical conversion with a flux ${k}_{\mathrm{m}}n(1{\varphi}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}})$, and decreases by compartment fusion with the boundary, which removes all medialpatches with an average flux $\alpha (1{\varphi}_{\mathit{\text{medial}}})\times n\text{}{\varphi}_{\mathit{\text{medial}}}$. At steadystate, one thus expects ${\varphi}_{\text{\mathit{m}\mathit{e}\mathit{d}\mathit{i}\mathit{a}\mathit{l}}}\simeq {k}_{\mathrm{m}}/\alpha $, suggesting that the lowest purity will be observed for ${k}_{\mathrm{m}}\simeq \alpha $. The results of simulations, shown in Appendix 7—figure 1, confirm this prediction.
The transition of purity occurs when sorting mechanisms become more efficient than mixing mechanisms. As a first approximation mixing occurs by biochemical conversion and sorting by budding, yielding the prediction that the purity transition occurs when ${k}_{\mathrm{b}}\simeq {k}_{\mathrm{m}}$. Additional mixing mechanisms include the fusion of two slightly impure compartments of different compositions. Such slow events become relevant if compartments remain in the system for a long time, that is, if fusion with the boundaries is slow (small values of $\alpha $). In this case, the transition of purity occurs for ${k}_{\mathrm{b}}>{k}_{\mathrm{m}}$. Phase diagram of the purity of the system for a system where $\alpha ={10}^{2}$ can be found in Appendix 7—figure 1. As expected the diagram is centered around ${k}_{\mathrm{m}}=\alpha $, with a transition of purity between ${k}_{\mathrm{b}}\sim {10}^{2}$ and ${k}_{\mathrm{b}}\sim {10}^{1}$ (and thus ${k}_{\mathrm{b}}>{k}_{\mathrm{m}}$).
As it tunes the residing time of compartments in the system, $\alpha $ also impact the average size of compartments. Decreasing the value of $\alpha $ increases the residence time and leads to larger compartments, as predicted by Equation 4. Note however that for small budding rates, when compartments are mixed and their exit through the boundary is difficult to estimate, the size distribution deviates from the singlecomponent ideal distribution.
The Golgi being a highly polarized organelle, it could be argued that the parameters controlling the fusion with the cis face (ER) and the trans face (TGN), ${\alpha}_{\mathrm{ER}}$ and ${\alpha}_{\mathrm{TGN}}$ could be different. An obvious way to reduce the retrograde exit of material from the ER is to reduce ${\alpha}_{\mathrm{ER}}$. To investigate the role of this asymmetry, we use an extreme regime of ${\alpha}_{\mathrm{ER}}=0$ and ${\alpha}_{\mathrm{TGN}}=1$. This prevents the recycling of impure compartments by fusion with the ER and broadens the low purity region of the parameter space, as shown in Appendix 7—figure 1, and explained by the analytical model presented in Appendix 3.
Appendix 8
Other budding/fusion schemes
Linear budding
To verify whether the assumption of saturated budding (Equation 2) has a strong impact on the results, we performed simulations with a linear budding mechanism ($f(\varphi )=\varphi $ in Equation 10), for which the budding flux for a given identity is linear with the number of patch with this identity (${J}_{\mathrm{b},i}={K}_{\mathrm{b}}n{\varphi}_{i}={K}_{\mathrm{b}}{n}_{i})$. As shown in Appendix 8—figure 1, the average compartment size and purity follow the trends discussed in the main text for the saturated budding regime. However, the purity transition occurs for larger ${k}_{\mathrm{b}}$, since the total budding flux is merely proportional to the size of the compartment for linear budding, and is thus smaller than for saturated budding, where it is proportional to the size times the number of species. Saturated budding thus promotes systems that are at the same time well sorted, and with large compartments.
Even though the purity transition is shifted toward higher values of ${k}_{\mathrm{b}}$ the link between purity and vesicular flux is qualitatively conserved (Appendix 8—figure 1). For low purity, the system is retrograde with an enrichment in cisspecies and a depletion of more mature species, and for high purity it is anterograde. However the retrograde flux is less marked than with a saturated budding rate (no clear depletion in transspecies). Indeed, efficient sorting relies on the capacity to export the minority component out of a compartment. Within the saturated budding scheme, the rate of export is proportional to the size of the compartment and does not depend on the number of minority components to export. Within the linear budding scheme however, it is proportional to the number of minority components. In the low ${k}_{\mathrm{b}}$ regime, transrich compartments are large. They emit a large vesicular flux of immature components within the saturated budding scheme, but this flux is much smaller within the linear budding scheme. This explains the qualitative difference between the enrichment curves in the low ${k}_{\mathrm{b}}$ regime for the saturated budding scheme (large depletion in trans identity  Figure 3  main text) and the linear budding scheme (almost no change in trans identity  Appendix 8—figure 1).
Nonfusing compartments
Fusion between compartments is an important ingredient of the model discussed in the main text. It could be argued that the fusion of large cisternae with one another could be a much slower process than fusion involving much smaller transport vesicles. We have performed simulations where fusion between compartments is prohibited if both compartments are larger than a vesicle. Note that in the model, such compartments are still able to fuse with the boundaries, so that the composition of the system can still be predicted by Equation 4 The results of this model are shown in Appendix 8—figure 1.
The size distribution of nonfusing compartment shows a broad peak, instead of the powerlaw seen in the case where compartments can fuse. As compartments can only grow by vesicle fusion and not by compartments fusion, their size is smaller than in the previous case. The compartments size is results from a balance between the rate at which compartments are recycled (fuse with one boundary) and the rate at which they grow by vesicle fusion. Consequently, the size distribution only weakly depends on ${k}_{\mathrm{b}}$, as compartments have very little time to fuse with vesicles before exiting the system (Appendix 8—figure 1).
As expected, and because we remove the possible interactions between compartments, the centripetal flux we observe between the previously described anterograde and retrograde is abolished. This is particularly visible if we plot the (normalized) mean enrichment in cis, medial and transspecies seen by cargo as they are transported via vesicular fluxes (Appendix 8—figure 1). As before, the system goes from a purely retrograde flux, characterized by a depletion in transspecies and an enrichment in cisspecies, to an anterograde flux, characterized by an enrichment in transspecies and a depletion in cisspecies. However, the transitional centripetal regime we previously described, characterized by an enrichment in medialspecies, vanishes. This suggests this particular regime is indeed resulting from the possible fusion between compartments, as it destabilizes medialcompartments that have a great probability to fuse with slightly impure cis or transcompartments. Thus the transport is anterograde as soon as the purity transition occurs (Appendix 8—figure 1).
Appendix 9
Experimental proposal
Data availability
All data generated from computer simulations and analysed during this study are included in the manuscript and supporting files. Source code has been provided.
References

The RSNARE endobrevin/VAMP8 mediates homotypic fusion of early endosomes and late endosomesMolecular Biology of the Cell 11:3289–3298.https://doi.org/10.1091/mbc.11.10.3289

Knockout of the Golgi stacking proteins GRASP55 and GRASP65 impairs Golgi structure and functionMolecular Biology of the Cell 28:2833–2842.https://doi.org/10.1091/mbc.e17020112

Golgi enlargement in Arfdepleted yeast cells is due to altered dynamics of cisternal maturationJournal of Cell Science 127:250–257.https://doi.org/10.1242/jcs.140996

The many routes of Golgidependent traffickingHistochemistry and Cell Biology 140:251–260.https://doi.org/10.1007/s0041801311247

A threestage model of Golgi structure and functionHistochemistry and Cell Biology 140:239–249.https://doi.org/10.1007/s0041801311283

Cooperative protein transport in cellular organellesPhysical Review E 83:041923.https://doi.org/10.1103/PhysRevE.83.041923

Exact stochastic simulation of coupled chemical reactionsThe Journal of Physical Chemistry 81:2340–2361.https://doi.org/10.1021/j100540a008

Models for Golgi traffic: a critical assessmentCold Spring Harbor Perspectives in Biology 3:a005215.https://doi.org/10.1101/cshperspect.a005215

Morphology and function of membranebound organellesCurrent Opinion in Cell Biology 26:79–86.https://doi.org/10.1016/j.ceb.2013.10.006

Golgi fragmentation is associated with ceramideinduced cellular effectsMolecular Biology of the Cell 16:1555–1567.https://doi.org/10.1091/mbc.e04070594

A model for the selforganization of vesicular flux and protein distributions in the golgi apparatusPLOS Computational Biology 9:e1003125.https://doi.org/10.1371/journal.pcbi.1003125

Vps74p controls Golgi size in an Arf1dependent mannerFEBS Letters 592:3720–3735.https://doi.org/10.1002/18733468.13266

Regulation of Golgi cisternal progression by ypt/Rab GTPasesDevelopmental Cell 36:440–452.https://doi.org/10.1016/j.devcel.2016.01.016

Segregation of sphingolipids and sterols during formation of secretory vesicles at the transGolgi networkThe Journal of Cell Biology 185:601–612.https://doi.org/10.1083/jcb.200901145

Evolution and diversity of the GolgiCold Spring Harbor Perspectives in Biology 3:a007849.https://doi.org/10.1101/cshperspect.a007849

A modeling approach to the selfassembly of the Golgi apparatusBiophysical Journal 98:2839–2847.https://doi.org/10.1016/j.bpj.2010.03.035

Secretory protein trafficking and organelle dynamics in living cellsAnnual Review of Cell and Developmental Biology 16:557–589.https://doi.org/10.1146/annurev.cellbio.16.1.557

Golgi maturation visualized in living yeastNature 441:1002–1006.https://doi.org/10.1038/nature04717

The biogenesis of the golgi ribbon: the roles of membrane input from the ER and of GM130Molecular Biology of the Cell 18:1595–1608.https://doi.org/10.1091/mbc.e06100886

Lipid transfer proteins and the tuning of compartmental identity in the Golgi apparatusChemistry and Physics of Lipids 200:42–61.https://doi.org/10.1016/j.chemphyslip.2016.06.005

The dynamics of engineered resident proteins in the mammalian Golgi complex relies on cisternal maturationThe Journal of Cell Biology 201:1027–1036.https://doi.org/10.1083/jcb.201211147

Roles for major histocompatibility complex glycosylation in immune functionSeminars in Immunopathology 34:425–441.https://doi.org/10.1007/s0028101203099

Golgi tubule traffic and the effects of brefeldin A visualized in living cellsJournal of Cell Biology 139:1137–1155.https://doi.org/10.1083/jcb.139.5.1137

(Re)modeling the GolgiMethods in Cell Biology 118:299–310.https://doi.org/10.1016/B9780124171640.000185

Golgi glycosylationCold Spring Harbor Perspectives in Biology 3:a005199.https://doi.org/10.1101/cshperspect.a005199

A missense mutation in the NSF gene causes abnormal golgi morphology in Arabidopsis thalianaCell Structure and Function 43:41–51.https://doi.org/10.1247/csf.17026

Nonequilibrium raftlike membrane domains under continuous recyclingPhysical Review Letters 95:168301.https://doi.org/10.1103/PhysRevLett.95.168301

Stochastic model of maturation and vesicular exchange in cellular organellesBiophysical Journal 114:947–957.https://doi.org/10.1016/j.bpj.2017.12.018

Stochastic model of vesicular sorting in cellular organellesPhysical Review Letters 120:058102.https://doi.org/10.1103/PhysRevLett.120.058102

Golgi ribbon disassembly during mitosis, differentiation and disease progressionCurrent Opinion in Cell Biology 47:43–51.https://doi.org/10.1016/j.ceb.2017.03.008
Decision letter

Raymond E GoldsteinReviewing Editor; University of Cambridge, United Kingdom

Suzanne R PfefferSenior Editor; Stanford University School of Medicine, United States

Jane KondevReviewer; Brandeis University, United States
In the interests of transparency, eLife publishes the most substantive revision requests and the accompanying author responses.
Acceptance summary:
Vagne et al. study a minimal model for the selforganisation of the Golgi, showing how the traditional models of "vesicular transport" and "cisternal maturation" arise as limiting cases of rich dynamical behaviour in a model of vesicle budding, maturation, fusion that leads to selforganised Golgi structure and cargo transport.
Decision letter after peer review:
[Editors’ note: the authors submitted for reconsideration following the decision after peer review. What follows is the decision letter after the first round of review.]
Thank you for submitting your work entitled "A minimal selforganization model of the golgi apparatus" for consideration by eLife. Your article has been reviewed by three peer reviewers, and the evaluation has been overseen by a Reviewing Editor and a Senior Editor. The following individual involved in review of your submission has agreed to reveal their identity: Jane Kondev (Reviewer #3).
You will see that the reviewers had differing views on the paper, with two highly critical and one supportive. Our decision has been reached after consultation between the reviewers. Based on these discussions and the individual reviews below, we regret to inform you that your work will not be considered further for publication in eLife.
The consensus among the reviewers is that there is no doubt in the potential utility of the kind of minimal model that you have developed, but there were serious questions about the connection between that model and existing data, about the presentation of the model itself, and the structure of the paper.
Reviewer #1:
The Sens group has been developing mathematical models of the Golgi apparatus since 2011 and this paper represents a continuation of this effort. Here, theoretical roles of vesicle budding, maturation and selforganization are used to describe the Golgi's various activities. While the mathematical computations and simulations are not in question, the paper provides no experiments on the Golgi to support or verify the proposed model, and the study fails to incorporate significant features of the Golgi (including the Golgi's unique lipid composition, its tubularvesicular character, and its ability to form de novo from ER). The paper reads as an exercise for modeling a minimal, selforganizing system that, while attempting to capture some features of organelles (i.e., fusion, budding and maturation), is far removed from the actual complex system comprising the Golgi apparatus. Because the proposed model does not significantly further our understanding of the Golgi apparatus, I cannot support its publication in eLife. Among the concerns of this study is that the model is built on faulty or unproven assumptions regarding the Golgi apparatus (some of these are listed below). It is also worrisome that several papers providing quantitative experimental data detailing Golgi dynamics are never cited, perhaps because these data cannot easily be explained by their model.
1) In their Introduction, the authors state that proteins passing through the Golgi interact with Golgi enzymes in a predetermined order dictated by the position of the enzyme in the Golgi stack (i.e., cis first and trans last). There is no evidence for this in the literature. Indeed, Farquhar et al., showed many years ago that some transacting enzymes are, in fact, in the cismost cisternae. We also know that carbohydrate processing by Golgi enzymes occurs in a sequential manner whether or not the enzymes are in the ER (i.e., in BFAtreated cells) or Golgi. In both cases, it's substrate availability that dictates the processing events.
2) The authors claim in their Introduction that there are only two major classes of models for explaining the Golgi's activities (i.e., the vesicle transport and maturation models). However, a third model involving tubular continuities and lipid partitioning has also been proposed by the groups of Luini and LippincottSchwartz. This third model can account for key experimentally observed features of the Golgi that neither the vesicle transport or maturation models explain. These features include: the finding that secretory cargo rapidly mixes throughout the Golgi stack upon arrival; that different types of cargo exhibit different export kinetics from the Golgi; and that cargo export from the Golgi follows mass action laws. Sens' model cannot account for any these observed features of the Golgi.
3) The authors claim that each membrane patch in their model undergoes maturation from a cis to medial to trans identity. This may have some type of support in yeast (in which membrane patches budding out from the ER clearly “mature” by exchanging cytosolic components over time), but it has no support in mammalian Golgi systems. There, cargo transport through the Golgi stack has never been visualized due to the extremely compact nature of the stack, whose individual cisterna are separated by <20 nm, well below the diffraction limit of imaging.
4) The authors' model predicts that there will be more cis/medial components associated with Golgi structures having low maturation rates, and more trans components for Golgi with high maturation rates. There is no evidence in the literature supporting this prediction. If anything, when one lowers the temperature to slow movement through the Golgi (which in the Sens' model would slow maturation), molecules tend to accumulate at the trans face of the Golgi, the opposite of their prediction.
5) The authors state that cargo flux through the Golgi in yeast is two orders of magnitude faster than in mammalian Golgi. This is unlikely. Bonfanti's paper examining flux of VSVG through the Golgi by EM, which they cite as supporting evidence, is inappropriate for this calculation at many levels from using repeated temperature shifts to counting a few immunogold particles per Golgi (a super crude quantification) to having to compare different fixed cells at different time points. The authors ignore the quantitative live cell imaging data of Golgi export of VSVGGFP provided by Hirschberg et al., 1998, which measured VSVG export out of the Golgi in real time, with thousands of data points per cell. There, VSVG was shown to be transported out the Golgi at an enormous rate (~7,000 molecules leaving the Golgi per sec) at peak flux after release from the ER. Importantly, the Hirschberg paper also revealed that the amount of VSVG cargo leaving the Golgi is directly correlated with the amount of cargo in the Golgi, with a single rate constant characterizing VSVG efflux kinetics. Sens' paper does not integrate any of these findings into their analysis.
6) The few predictions made by the author's model are not supported by findings in the literature. For example, the authors' state that their model predicts that deletion of Arf1 will decrease the number of Golgi compartments and increase their size. Experimentally, this has already been tested with brefeldin A (which inactivates Arf1) or with Arf1 mutants. Contrary to the Sens' model prediction, when Arf1 is inactivated in cells, the Golgi disappears, being resorbed quickly into the ER with no increase in any Golgi subcompartment. This occurs through retrograde transport back to the ER, a pathway not incorporated into the Sens' model.
I appreciate the efforts made by the Sens group in using modeling to describe selforganizing systems, but their attempt to describe the Golgi apparatus with this approach is missing key aspects of Golgi organization and dynamics. This has resulted in predictions by the model that are counter to experimental data.
Reviewer #2:
Poorly written, hard to understand.
I was excited to read this paper since the work of one of the authors is known to me and I expected an understandable and wellcrafted paper. I was rather crestfallen to find a poorly written paper, and one where I could not figure out what the authors were actually doing.
One major concern is that this paper is not at all easy to read. It is either i) not selfcontained, or ii) sloppily put together. The authors refer to a "steady state" but do not at any time show that such a state exists, as far as I can see. There is no plot with the time evolution of any quantities being measured. In fact, they admit that there is no steady state, since the size of the system is not constant. They state that "the maximum number of timesteps is typically set to 10^{6} or 10^{7} in order to reach steady state and accumulate enough statistics on all the measured quantities", but that in and of itself is meaningless if one doesn't correctly specify the conditions underlying the measurement.
The paper as a whole is an odd combination of being verbose AND telegraphic at the same time. Echoing my above comment, I do not feel that I am getting enough information, enough detail, despite the verbosity. An example: it takes a long time to get to the point of how heterogeneous compartments can form; they say "this homotypic fusion mechanism (relying on local interactions) allows vesicular transport between compartments of different identities – a process that may be regarded as heterotypic fusion. It merely requires that the receiving compartment contains some membrane patches of identity similar to that of the emitted vesicle". Why not just give the rule, rather than the metaphysical description? It would take up less space!
Does seemingly everything have to be in the Appendix? One can understand page limits, but this is ridiculous. More explanation of the model in the main text would have been very welcome, because in its present form, the paper is extremely difficult to understand.
The question of system size should appear in the main text. As mentioned above, they cannot really keep it constant, so they just adjust the injection rate so it is more or less constant when they change K. They say that they aim for N ~ 300, but the actual system size (sometimes called size, sometimes mass) in the "steady state" changes significantly with K (see Appendix 1—figure 1B and Appendix 4). Moreover, the fluctuations are very large (since no time evolution plot is ever shown, one suspects that they are not really fluctuations…). Of course, the results would be highly dependent on system size, and one gets no justification for their choice of N ~ 300 until the Appendix.
The plots that are shown in the paper are not very useful. They rather tend to obfuscate. Figure 2A intends to show a power law for the size of the compartments with an exponential cutoff, but by plotting a bar diagram (even though it is a loglog scale), it is very difficult to see how well this works. They should plot actual points and a fit, or a parallel line showing an actual power law. In a similar vein, they overuse density plots that are not really quantitative (even in the Appendices). In Figure 2BC, they could at the very least plot some contour lines (perhaps dotted) for the same values that are plotted with a continuous line for the analytical approximation – from the colour map, it is difficult to judge.
Nomenclature:
It is confusing to keep track for reasons that are fully avoidable. For instance, their parameters are (K_{i},K_{f},K_{b},K_{m}): the rates of injection, fusion, budding and maturation. All of them are normalised by Kf, which sets the time unit. But then they define k_{m} = K_{m}/K_{f} and, instead of maintaining the use of lowercase k for the others, they then define K=K_{b}/K_{f}, J=K_{i}/K_{f}. As if the other problems with the paper weren't enough, now the reader has to remember what each means.
Main result is obfuscated:
The main result, in my opinion, is the existence of an intermediate "sorted" regime. It would be trivial to have a regime with a large mixed compartment or with small pure ones. The potentially interesting aspect of the model is that it can create pure compartments of at least intermediate size, if the parameters are chosen correctly. This is highly dependent on the way budding is modeled  only explained in the Appendix, in an unclear way. In short, the budding flux for species i is defined as
J_{b,}i = K_{b} \times n \times f(\phii)
where n is the size of the compartment and phi_{i} the proportion of vesicles of type i (phi_{i} = n_{i} /n). This is already confusing nomenclature, since they use J_{b} for this flux, when J is also the normalised injection rate. They then say that if the budding rate for species i depends linearly on its concentration, so
J_{b,i} = K_{b} n \phi_{i} = K_{n} n_{i}
then you cannot get the intermediate regime. For that you need what they call "nonlinear budding". In their words: "In order to reproduce this feature we choose a highly nonlinear budding scheme f(phi_{species}) = 1 if phi_{species} > 0 and f(phi_{species}) = 0 otherwise". First of all, they are now using phi_{species} for what was before phi_{i}.
I can only imagine what biological readers might make of this labeled "highly nonlinear scheme". It is, after all, just a constant!
phi_{i} cannot be negative and if phi_{i} = 0, then there are no vesicles to bud. So the budding rate is just a constant, independent of everything and the same for all species in the same compartment.
Problematic statistical analysis:
The statistical analysis is problematic. First of all, because we don't know whether there is a steady state, can one even define averages? But even if there is a steady state, there are things like the figure in Appendix 4. There, they measure the "temporal standard deviation" of the mass and purity. How are they defining this standard deviation, when the individual measurements come from a correlated time series and are, therefore, not independent? Note that there is a very large dependence with K, which could be trivial for this reason (as the correlation time would obviously depend on K).
Conclusion:
I cannot recommend publication of the paper – it is too flawed at present – and importantly, it is not possible to evaluate whether the research presented is even correct. A revised version, should the editors wish to solicit one, at the very least, must show some evidence of a steady state and explain how averages are calculated. The authors should pay special attention and care to make the paper more readable. Much more readable. In particular, the model (the actual model) could be fit in the text and justification/discussion of terms could be left for the Appendix.
Reviewer #3:
This paper explores a simple model for the dynamics of vesicles as they progress through the Golgi apparatus and shows that a number of observed phenomena can be understood based on a small number of biologically, and physically, wellmotivated assumptions.
I thoroughly enjoyed reading this paper, even though I am not an expert on the Golgi. I think it's a great example of theory at its best, applied to an interesting biological problem. The authors start with some simple observations that have been reported in the literature, which they turn them into a mathematical model that upon analysis yields some interesting and experimentally testable predictions. I am optimistic that their results will lead to some new experiments.
The key result of the paper is the connection between the directionality of the vesicular transport and the chemical composition of the Golgi compartments, for which there seems to be experimental evidence. The connection is established within a very simple model whose assumptions seem reasonable. The paper also describes the emergence of a number of structural and dynamical properties of the Golgi, from the same small set of assumptions. A particularly interesting feature of the model is how different structures and dynamics can be accommodated within the model by changing one or more parameters (of which there are only four). For example, the model predicts that there is an inverse relationship between the speed of vesicle transport and fidelity of processing, which is consistent with published comparisons between the structure of Golgi in normal and starved yeast.
The paper is very clearly written, and there were only a few places where I was not quite sure about what the authors had in mind.
1) When describing the Budding mechanism, there is reference to a nonlinear budding scheme that depends on the compartment size and contamination state. These are described in detail later, but I would recommend that the authors say a few words here to give the reader an intuitive sense of the thing without having to skip ahead.
2) In describing maturation, the chemical transformations are described as going only one way. How is this justified? Is it known that the reverse reactions are very slow? Also, is the maturation rate independent of the composition of the compartment and if so, what is the justification for that?
3) A typical compartment size is defined as the ratio of the second and the first moment of the size distribution. I was struck by the fact that this would give a "typical size" even for an exponential distribution of sizes, which I would not say defines a typical size (at least not from the point of view that the mean and standard deviation of such a distribution are equal). Might be worth a short comment on why this is a good definition of "typical size".
4) In the "Steadystate organization" section there is reference to a cutoff size, which I don't think was defined.
5) In the Discussion there is a conclusion that the spatial information is less crucial than biochemical composition for the organization of the Golgi. I am not sure this is warranted since a comparison between the two cannot be made in a model lacking spatial information.
In general I was left wondering whether simple testable quantitative predictions can be made from the model. I very much liked the qualitative observations that the model neatly explains, but I would be even more excited if there was a prospect for putting it to a much more stringent test provided by quantitative experiments, which the authors might propose for an intrepid experimentalist to do.
[Editors’ note: further revisions were suggested prior to acceptance, as described below.]
Thank you for submitting your article "A minimal selforganization model of the golgi apparatus" for consideration by eLife. Your revised article has been reviewed by two peer reviewers, and the evaluation has been overseen by a Reviewing Editor and Suzanne Pfeffer as the Senior Editor. The reviewers have opted to remain anonymous.
The reviewers have discussed the reviews with one another and the Reviewing Editor has drafted this decision to help you prepare a revised submission.
We would like to draw your attention to changes in our revision policy that we have made in response to COVID19 (https://elifesciences.org/articles/57162). Specifically, we are asking editors to accept without delay manuscripts, like yours, that they judge can stand as eLife papers without additional data, even if they feel that they would make the manuscript stronger. Thus the revisions requested below only address clarity and presentation.
Following your submission of a revised manuscript we sought the opinion of the original reviewers and, more recently, that of an additional reviewer (a theorist chosen specifically to be outside the field of Golgi research). As you can appreciate, purely theoretical papers in eLife need to be accessible to the biological readership of the journal and, if on the abstract side, may require extra effort to be suitable for publication. Such is the case here. As detailed below, we fully appreciate the significance of your major conclusions regarding the emergence of two distinct models from a singe underlying dynamics. At the same time, the reviewers request changes in the structure of the paper to improve readability, and that you address some key questions about the role of the actual structure of the Golgi.
Summary:
Vagne et al. study a minimal model for the selforganisation of the Golgi, showing how the traditional models of "vesicular transport" and "cisternal maturation" arise as limiting cases of rich dynamical behaviour in a model of vesicle budding, maturation, fusion that leads to selforganised Golgi structure and cargo transport.
Revisions:
1) Regarding prior work, the reviewers note that as the work of Patterson et al., 2008, provided a detailed quantitative model of the Golgi using one set of parameters, and furthermore assessed and tested those parameters by doing extensive quantitative live cell imaging experiments, it would be appropriate to seek a comparison between your model (with its dimensionless parameters) and those results. This would likely enlarge the scope of potential predictions for direct experimental tests. Workers in the Golgi field would benefit from some suggested experiments expected to skew the pathway in one direction or another, with predicted outcomes on the size of the compartment or the speed of traversal etc. For example, the cisternal progression model predicted that larger cargo export rates would be more sensitive to nocodazole dissolution of microtubules and resulting creation of ministacks, a prediction that turned out to be correct: PMID 25103235.
2) To help the general reader understand the context of this work better, we would like to see an introductory figure illustrating the structure of the Golgi (and perhaps some of the models for transport in the Golgi).
3) The reviewers suggest that putting more of the model equations in the main text would streamline the argument; one way of making things more understandable might be to start by writing down simple meanfield equations (as in the early parts of Appendix 3), before delving into the analysis of the detailed model. The backandforth between the main text and the appendices makes for difficult reading.
4) An important concern is that the model shows how Golgi structure and vesicular transport can arise in a selforganised fashion, but it remains unclear to what extent the authors can claim that vesicular transport in actual Golgi is dominated by this transport arising from selforganisation. When reading the paper, we wondered about the contribution of spatial structure or biochemical processes (that might bias which vesicles cargo goes into, or which compartments a vesicle fuses into). It was only pages later, in the Discussion, that one finds a paragraph on the role of spatial structure: even if experiments show that "Golgi functionality is preserved under major perturbation of its spatial structure", how is transport affected? From a quick look at the reference (Dunlop et al., 2017), it seems that transport can be slowed down massively, hinting that spatial structure is just as important as selforganisation. Are there similar experimental results on biochemical effects?
5) One of the main results in the paper is the existence of different regimes of selforganisation, and in particular the existence of an intermediate, sorted regime with large pure compartments. However, the definition of these regimes remains rather qualitative. Perhaps some plots of mean compartment size against mean purity, to complement Figure 2, would put this part of the analysis on a clearer footing.
6) The model shows very clear anterograde transport of cargo in the limit k_{b}>>1, but the retrograde transport in the limit k_{b}<<1 is less clear; in particular, the model appears to show transport from the trans to wellmixed compartments, but not into the cis compartments: rather, there is weak flow from the cis compartments to wellmixed ones, i.e. anterograde flow in part of the system. Is it possible to show that flow cannot be purely retrograde? – This might be an important constraint on the "cisternal maturation" mechanism. (Perhaps part of this issue is addressed in Appendix 6.)
https://doi.org/10.7554/eLife.47318.sa1Author response
We thank you and the referees for taking the time and effort to review our paper. After carefully reading through the referee reports, we are convinced that we can provide satisfactory answers to all the criticisms offered by reviewers 1 and 2, and that we can use the comments from the three reviewers to improve our paper.
Reviewer #1:
The Sens group has been developing mathematical models of the Golgi apparatus since 2011 and this paper represents a continuation of this effort. Here, theoretical roles of vesicle budding, maturation and selforganization are used to describe the Golgi's various activities. While the mathematical computations and simulations are not in question, the paper provides no experiments on the Golgi to support or verify the proposed model, and the study fails to incorporate significant features of the Golgi (including the Golgi's unique lipid composition, its tubularvesicular character, and its ability to form de novo from ER). The paper reads as an exercise for modeling a minimal, selforganizing system that, while attempting to capture some features of organelles (i.e., fusion, budding and maturation), is far removed from the actual complex system comprising the Golgi apparatus.
Our model does not provide new experimental data, but discusses existing data is a new way, in particular regarding the correlation between size and purity of Golgi compartments. The reviewer is correct that our model is far from encompassing the complexity of the Golgi. This was stated explicitly in the original version, including in the title, which promises a “minimal” model. There will undoubtedly be experimental observations that our model will not be able to explain. In the revised version, we discuss the few observations mentioned by the reviewers that are indeed not reproduced by our model, we explain why our model fails to reproduce these data and how it can be improved to do so. But as we discuss below, many of the observations said to be at odds with our model are in fact entirely consistent with it, and none contradict the model. These observations are also discussed in the revised version of the manuscript.
Notwithstanding, we think that our model offers fundamental conceptual advances that are important to convey to the Golgi community, and more generally to cell biologists interested in intracellular transport and organelle dynamics. Most existing models discuss transport through a preexisting compartmentalised Golgi structure. The main goal of our model is to explain how this compartmentalised structure can spontaneously emerge through selforganisation. Our aim was and still is to investigate how the interplay between the biochemical maturation of membrane components and vesicular transport can lead, at the same time, to compartmentalisation and directed fluxes. The model involves three main parameters, so it definitely misses some of the Golgi’s complexity, but we believe that we have succeeded in answering that question from a theoretical point of view, and that our conclusions are very much relevant to Golgi dynamics.
Regarding lipid (and protein) composition of the Golgi, we only consider three types of membrane identity. It could be argued that this is not enough to account for the complex lipid composition of the Golgi, but we show that this is definitely enough to produce a rich and versatile dynamics. Our model Golgi is certainly able to form de novo and reaches a steadystate that does not depend on its initial state. Although this paper focusses on the study of the steadystate, the de novo formation is now explicitly shown in Appendix 4 (see also our response to reviewer 2 regarding the existence of a steadystate). Regarding the tubularvesicular nature, we agree that we have disregarded the existence of tubular connections between cisternae. We note that the relevance of tubules in Golgi dynamics is still debated, except in some extreme conditions such as Brefeldin A treatment (see our response to point (6) below). In the revised version, a paragraph has been added discussing in more details the limitation of our model, in particular regarding tubular transport and lateral segregation and partitioning within cisternae, and many references (provided below) have been added.
Because the proposed model does not significantly further our understanding of the Golgi apparatus, I cannot support its publication in eLife. Among the concerns of this study is that the model is built on faulty or unproven assumptions regarding the Golgi apparatus (some of these are listed below). It is also worrisome that several papers providing quantitative experimental data detailing Golgi dynamics are never cited, perhaps because these data cannot easily be explained by their model.
We strongly disagree that our model does not further the understanding of Golgi structure and dynamics. Whether Golgi cisternae are longlived structure exchanging material, or dynamical structures evolving over time is still fiercely debated after many decades. Our model is the first to show that, for the same set of individual mechanisms, there is a transition from one dynamics to the other upon varying the relative rates. Furthermore, we provide a complete analytical understanding of why this happens. We argue that the competition between sorting and maturation dictates the Golgi structure and dynamics. We have never seen this discussed before. This is an achievement, which allowed us to explain a number of experimental observations, in particular relating the size and the purity of Golgi compartments, which have been barely discussed and certainly not been explained before. This has been made clearer in the revised version of the manuscript.
We also disagree that our model is based on faulty assumptions. The reviewer criticised our treatment of (molecular) biochemical maturation (question 3). This is a fundamental assumption, which is indeed unproven for mammalian Golgi. This is explicitly mentioned as a hypothesis in the revised version. Our response to question 3 explains why we believe this is a reasonable assumption. We note that many influential models of the Golgi are based on it, so it is certainly worthwhile to quantitatively assess its potential role on Golgi dynamics (references are provided below). The other criticisms do not address model assumptions, but rather the discussion of model predictions. We hope that the reviewer will be convinced by our responses below.
Notwithstanding, we agree with the reviewer that important papers were not properly cited in the original version of the manuscript. This has been corrected in the revised version, where many references (listed in our response below) have been added. We now discuss how our model relates to the data provided in these papers, and how it should be improved (at the expense of increasing the number of parameters) when it fails to explain them.
1) In their Introduction, the authors state that proteins passing through the Golgi interact with Golgi enzymes in a predetermined order dictated by the position of the enzyme in the Golgi stack (i.e., cis first and trans last). There is no evidence for this in the literature. Indeed, Farquhar et al., showed many years ago that some transacting enzymes are, in fact, in the cismost cisternae. We also know that carbohydrate processing by Golgi enzymes occurs in a sequential manner whether or not the enzymes are in the ER (i.e., in BFAtreated cells) or Golgi. In both cases, it's substrate availability that dictates the processing events.
A more precise reference to the work of Farquhar et al. would be helpful to know which of her work the reviewer refers to. Nevertheless, many papers show that disrupting Golgi structure and dynamics affects the processing of proteins, and in particular glycosylation. Among those: (Hu et al., 2005), reporting that ceramide, which induces Golgi fragmentation, also induces glycosylation defects of integrin, disrupting cell adhesion. A similar phenotype with observed upon Brefeldin A treatment, strongly suggesting that BFA also induces glycosylation defects. (Shen et al., 2007) showed that siRNA suppression of BIG1 (a GEF also inhibited by BFA) disrupted Golgi structure and altered β1 integrin glycosylation, leading to aberrant cell adhesion and migration. (Xiang et al., 2013) and (Bekier et al., 2017) showed that knockout of the Golgi stacking proteins GRASP55 and GRASP65 impairs glycosylation and sorting of several proteins and lipids. In any case, this introductory sentence has no impact on the layout of our model what so ever, and can be removed if problematic. In the revised version, this sentence was modify to specifically mention glycosylation and references (Hu et al., 2005, Shen et al., 2007, Xiang et al., 2013, Bekier et al., 2017) are cited.
2) The authors claim in their Introduction that there are only two major classes of models for explaining the Golgi's activities (i.e., the vesicle transport and maturation models). However, a third model involving tubular continuities and lipid partitioning has also been proposed by the groups of Luini and LippincottSchwartz. This third model can account for key experimentally observed features of the Golgi that neither the vesicle transport or maturation models explain. These features include: the finding that secretory cargo rapidly mixes throughout the Golgi stack upon arrival; that different types of cargo exhibit different export kinetics from the Golgi; and that cargo export from the Golgi follows mass action laws. Sens' model cannot account for any these observed features of the Golgi.
There are indeed different qualitative models, including Luini’s cisternal maturation with tubular connections (Trucco et al., 2004, Rizzo et al., 2013), LippincottSchwartz’s rapid partitioning between processing and exporting subcompartments (Patterson et al., 2008), or Pfeffer’s cisternal progenitor model (Pfeffer, 2010). The strengths and weaknesses of these different models are reviewed by Glick and Luini in (Glick and Luini, 2011). Variations of these models, such as the “rim progression” model also exist (Lavieu et al., 2013). These papers should have been cited and are cited in the revised version.
Notwithstanding, none of these models provide a quantitative analysis of the generation and maintenance of Golgi compartments, which are generally assumed to be of fixed number and size. The strength of our model is that the structure spontaneously emerges from the dynamics. Our model is in that respect fairly close to the concepts qualitatively exposed in the cisternal progenitor model (Pfeffer, 2010), but it is quantitative and can explain denovo Golgi formation, as we now show in Appendix 4. One difference with the cisternal progenitor model is that we restrict intercisternae exchange to be vesiclebased, while Pfeffer suggests that larger cisterna fragments could also be exchanged. We have rephrase our description of existing models to explain that they mostly belong to two classes: one involving stable compartments exchanging components through fusionbased mechanisms, to which the vesicular transport, rim progression, rapid partitioning and cisternal progenitor models belong, and one involving transient compartment undergoing enblock maturation and in which fusion mechanisms are dispensable for cargo transport. The papers by Luini, LippincottSchwartz, Pfeffer and others (refs. (Trucco et al., 2004,Rizzo et al., 2013,Patterson et al., 2008,Pfeffer, 2010,Lavieu et al., 2013)) are now cited.
The features referred to the reviewer pertain to the dynamics of cargo transport through the Golgi, and our model can account for them. This is clearly a very important topic of much interests to us, but we emphasise that the present manuscript only discusses the steadystate organisation and vesicular exchange of Golgi components. The present article is already fairly dense, and we found it important to have a full description and quantification of the selforganisation model before tackling the issue of cargo transport. The dynamics of cargo sorting and transit through the Golgi is also very rich and deserve a focussed publication. This will be the subject of a followup “Research Advance” article which will be submitted to eLife later this year. Nevertheless, the simulations described in the present paper already includes cargo, which we use as marker of vesicular transport directionality. We emphasise that these cargo are “passive” in the sense that they do not feel the biochemical identity of the membrane patch they are bound to. Their dynamics thus correlates with membrane patches dynamics, averaged over all identities. Our simulations can already predict the Golgi export kinetics of such passive cargo. We have added a figure (Appendix 4—figure 2) showing that the export kinetics of cargo is indeed exponential. This dynamics is similar to the one observed by iFRAP (Figure 2 i of (Patterson et al., 2008)) or using the RUSH technics (Figure 3 e of (Boncompain et al., 2012)). The exponential exit kinetics was one of the main arguments in (Patterson et al., 2008) to dismiss a “strictly cisternal maturation” dynamics, a conclusion we fully support. The revised version now cites (Patterson et al., 2008, Boncompain et al., 2012).
The specific features mentioned by the reviewer are entirely consistent with our model.
1) Rapid mixing of cargo: The kinetics of cargo distribution amongst cisternae of different identities following their release from the ER is at present not available from our simulations, and we plan to do these simulations in the followup paper. Nevertheless, the fact that incoming cargo can distribute throughout the spectrum of compartment identity in a time scale not limited by biochemical maturation can already be inferred from Figure 4 (main text) and Appendix 6—figure 1B of our paper. Although these are steadystate distributions of compartment identity and intercompartment exchange (vesicular in our model), comparing the distribution of the system’s size and the distribution of vesicular exchange arrows show that compartments of all identity are dynamically connected through intercisternal exchange. Furthermore, compartments can fuse with one another in our model, which further speedsup the kinetics of cargo mixing. This leads to the presence of cargo in compartments of all identity within a time scale that is not limited by biochemical maturation. We realize that this remains qualitative at this point, but we hope that this will convince the reviewer that rapid mixing of cargo (which is not the topic of the present paper) is not at odd with our model.
2) Cargo export follow law of mass action: this is directly demonstrated by the exponential decay of cargo abundance in the Golgi following a pulsechaselike simulation. This is now shown in Appendix 4—figure 2 and discussed in the main text.
3) Different export kinetics for different cargo: our model naturally leads to different export kinetics for cargo that interact differently with the various membrane identities (which will be called “ colored cargo” in the followup article), related to the fact that these cargo would reside in different membrane domains (patches of different identities in our framework). Trans cargo (residing preferentially in Trans membrane patches), exit faster than Medial and Cis cargo because they can enter Trans vesicles that can fuse directly with the exit. Cis cargo are often recycled to the ER via their packaging into Cis vesicles and show the slowest exit rate.
We can make the data regarding colored cargo available to the reviewer, although additional work is needed to fully understand them theoretically. We stress that cargo dynamics is beyond the scope of the present work, which focusses on Golgi selforganisation. We strongly feel that adding results on cargo dynamics to the present paper could obfuscate our main message.
3) The authors claim that each membrane patch in their model undergoes maturation from a cis to medial to trans identity. This may have some type of support in yeast (in which membrane patches budding out from the ER clearly “mature” by exchanging cytosolic components over time), but it has no support in mammalian Golgi systems. There, cargo transport through the Golgi stack has never been visualized due to the extremely compact nature of the stack, whose individual cisterna are separated by <20 nm, well below the diffraction limit of imaging.
The maturation of membrane patches, which is one of the three fundamental mechanisms in our model, refers to the local biochemical maturation of membrane components. Importantly, it does not refer to the global maturation of entire Golgi cisternae. In fact, we have replaced the term maturation (which bears strong connotation in the Golgi field) by biochemical conversion in the revised version to avoid confusing the molecular mechanism with its putative largescale consequence.
Our inspiration for this model draws from the socalled Rab cascade which, as the reviewer correctly points out, has not yet been demonstrated in Mammalian Golgi. However, many of the proteins involved in the Yeast Rab cascade are conserved in mammalian cells (Klute et al., 2011), so it is important to explore the possibility that this mechanism is also relevant to mammalian Golgi. Furthermore, the Rab cascade concept is at the basis of the cisternal progenitor model (Pfeffer, 2010), and it is fundamentally important to test this model quantitatively, which is what we do here. One strong result of the present work is that this molecular biochemical conversion mechanism might lead to the overall maturation of cisternae (as in Yeast), but might also result in the export of the matured component from a cisternae whose identity remains mostly constant over time, which could correspond to mammalian Golgi. This is a conceptual advance, which is permitted by the quantitative model we propose. Again, our point of view is to adopt simple rules and to see the variety of structures that can issue as a function of the parameters. These rules are surely not sufficient to reproduce the entire diversity of mechanisms at play in the Golgi, but our results are robust and well controlled. Any model of the Golgi that includes some form of biochemical conversion and intercisternal exchange will inevitably encounter the universal phenomenon uncovered by our model.
We have emphasized that the biochemical conversion mechanism as a hypothesis in the revised version of the manuscript. We also have stressed that this does not have to rely on the Rab cascade, but could happen in many ways. It could also correspond to the modification of lipid components by enzymes in the Golgi, as nicely reviewed in (McDermott and Mousley, 2016).
4) The authors' model predicts that there will be more cis/medial components associated with Golgi structures having low maturation rates, and more trans components for Golgi with high maturation rates. There is no evidence in the literature supporting this prediction. If anything, when one lowers the temperature to slow movement through the Golgi (which in the Sens' model would slow maturation), molecules tend to accumulate at the trans face of the Golgi, the opposite of their prediction.
It is difficult to find clear experimental quantification of what would be the result of a decrease of maturation rate that would not affect other processes. So we agree with the reviewer that several predictions we make are not (yet) supported by experiments. We don’t feel that we should be criticised for making these predictions anyway. In the paper, we briefly discuss the result of starvation (glucose deprivation) experiments which slowed cisternal maturation by a factor of 2.5 and coincides with the Golgi becoming more polarized (Levi et al., 2010). This agrees with the general trends we put forward, but is in no way a falsifiable test, since one cannot be certain that starvation slowdown is due to a decrease of the maturation rate alone. We have since then become aware of more direct observations of the effect of varying the maturation rate (Kim et al., 2016). Interpreting this publication’s results is less straightforward as they can tune independently multiple steps of Golgi maturation, which is not the case in the present model implementation where all maturation steps occur with the same rate k_{m}. In this paper, they show that overexpression of Ypt1 (in Yeast) increases the transition rate from early to transitional Golgi, and increases the colocalization of early and late Golgi’s markers. We interpret this as a decrease of Golgi’s purity upon increasing the maturation rate by unbalancing the ratio of budding to maturation rates. This suggests that the wildtype Yeast Golgi is closed to the line k_{m} = k_{b} (see Figure 2C of the paper). In the same paper, they report that overexpression of Ypt31 has a more complex phenotype. It accelerates transition to Late Golgi identity, but also increases exit from the Golgi via clathrin vesicles. Within our model, this would correspond to a concomitant increase of k_{m} for the medialtrans conversion and increase the budding rate of trans species. We did not attempt to model this, which goes beyond our simplifying approximation that the different rates are the same for all identities.
We would not venture to predict the result of a temperature shift, which probably affects all the parameters of the model in a way hardly controllable. It is not clear to us why the reviewer says that decreasing temperature affect the maturation rate and not budding or fusion? Temperature decrease to 20^{◦}C slows down ER to Golgi export, and slows done Golgi to TGN export even more dramatically, leading to large Trans compartment (Ladinsky et al., 2002). This indicates that exchange is affected as much as, if not more than maturation. Our prediction are actually not based on the value of individual rates, but rather on the value of ratios of rates. Regarding the prediction alluded to by the reviewer, our claim is that increasing the ratio of maturation rate over budding rate should result in larger Trans cisternae. This fact agrees with results of Arf1 depletion, which likely decreases the budding rate (Bhave et al., 2014), and would also agree with the temperature shift data provided budding is slowed more than maturation, which seems reasonable given our current knowledge.
The revised version of the manuscript now proposes additional experimental strategies to test our predictions, detailed in our response to reviewer 3’s last comment.
5) The authors state that cargo flux through the Golgi in yeast is two orders of magnitude faster than in mammalian Golgi. This is unlikely. Bonfanti's paper examining flux of VSVG through the Golgi by EM, which they cite as supporting evidence, is inappropriate for this calculation at many levels from using repeated temperature shifts to counting a few immunogold particles per Golgi (a super crude quantification) to having to compare different fixed cells at different time points. The authors ignore the quantitative live cell imaging data of Golgi export of VSVGGFP provided by Hirschberg et al., 1998, which measured VSVG export out of the Golgi in real time, with thousands of data points per cell. There, VSVG was shown to be transported out the Golgi at an enormous rate (~7,000 molecules leaving the Golgi per sec) at peak flux after release from the ER. Importantly, the Hirschberg paper also revealed that the amount of VSVG cargo leaving the Golgi is directly correlated with the amount of cargo in the Golgi, with a single rate constant characterizing VSVG efflux kinetics. Sens' paper does not integrate any of these findings into their analysis.
We were well aware of the data reported in Hirschberg’s paper, which was cited in the original version and is not in contradiction with our statement. We did not state that the flux going through the Golgi is smaller in mammalian cells, but that the Golgi dynamics – which includes rates of budding, fusion and maturation – is slower. The number of 7000 molecules per sec. leaving the Golgi is a flux. This is controlled by the parameter J in our model. Doubling this flux would double the steadystate Golgi size without changing its dynamics, which is characterised by rates. Hirschberg reports Golgi to plasma membrane export rates for VSVG of order 3.0% per min, which corresponds to an exit rate constant of 1/30min. This paper reports Golgi residence time around 40 min. This is consistent with iFRAP (Patterson et al., 2008) and RUSH (Boncompain et al., 2012) experiments, which consistently report Golgi transit time of order 20 to 40 min for transit cargo (much longer times for Golgi resident proteins). On the other hand, the dynamics of Yeast Golgi, characterised by the maturation of cisternae, is of order 1 min (Losev et al., and MatsuuraTokita et al., Science 2006, both cited in the paper. We thus firmly stand by our claim that all available evidence we are aware of indicate a Golgi dynamics one to two orders of magnitude faster in Yeast than in mammalian cells. We have clarified our statement in the revised version, and included the references above (including Hirschberg et al.). As already discussed in our response to question 2, the cargo export kinetics is exponential, in agreement with Experimental observations. This is now discussed in the paper, with one added figure (Appendix 4—figure 2).
6) The few predictions made by the author's model are not supported by findings in the literature. For example, the authors' state that their model predicts that deletion of Arf1 will decrease the number of Golgi compartments and increase their size. Experimentally, this has already been tested with brefeldin A (which inactivates Arf1) or with Arf1 mutants. Contrary to the Sens' model prediction, when Arf1 is inactivated in cells, the Golgi disappears, being resorbed quickly into the ER with no increase in any Golgi subcompartment. This occurs through retrograde transport back to the ER, a pathway not incorporated into the Sens' model.
Retrograde transport back to the ER is definitely included in the model. The fusion with the entry (the ER) occurs according to the rules of homotypic fusion, and can be modulated through the parameter α_{ER} (see subsection “Model” and Appendix 7).
Nevertheless, this is an important question that is discussed in the revised version of the manuscript. BFA is fairly nonspecific, and as far as we know, there is no clear explanation of how it affects the Golgi structure. However, a similar phenotype is observed with Golgicide A (GCA), which inhibit GBF1, a GEF for Arf1. This leads to the dissociation of COPI coat and induces disassembly of the Golgi and the TGN (Saenz et al., 2009). There are some subtle differences, as BFA induces TGN tubulation, while GCA apparently causes TGN dispersal into small vesicles. This shows how difficult it would be for a model, even a highly complex multiparameter model, to capture these effects.
Nevertheless, we discuss this issue in a specific paragraph dedicated to the limitations of our model. One limitation is that we include a single exchange mechanism, based on vesicles. Our model does not include tubules. We argue in the revised version of the manuscript that transient tubular connections between compartments or tubules that fuse with one compartment after detachment (Sciaky et al., 1997,Trucco et al., 2004) can to some extent be modeled in similar ways than our vesicle exchange mechanism, with the need of defining tubule nucleation/fusion and scission rates. However the tubules induced by BFA are of a different nature. They are much more stable and do not detach. This results in the Golgi quickly redistributing in the ER as soon as a tubular connection is established between the Golgi and the ER, in a way suggesting a purely physical (tensiondriven) mechanism (Sciaky et al., 1997). Our model does not include this possibility, as is now clearly stated in the revised version. We actually believe that these observations could be made consistent with our predictions by adding a sizedependent tubulation mechanism based on membrane mechanics. The force required to form membrane tubules increases with membrane tension (Der´enyi et al., 2002). We predict that Arf1 inactivation increases the cisterna surface area S (which is what we call the size in the model). It should also increase the cisterna volume V , and should increase the socalled surfacetovolume ratio S^{3/2}/V ∼ N^{1/2} (where N is the number of vesicles that fused together to form the cisterna). The membrane tension of a vesicle decreases (in a fairly nonlinear way) upon increase of surfacetovolume ratio (Evans and Rawicz, 1990). One thus expects tubulation to increase with cisterna size, which (with the help of our model) could explain why Arf1 inactivation (indirectly) induces massive cisterna tubulation and strongly increases the fusion probability with the ER. Since large cisternae are under low tension (probably lower than that of the ER (Upadhyay and Sheetz, 2004)), the tension difference would drive their absorption into the ER upon tubular contact. We could thus reproduce the BFA phenotype by saying that the ER fusion parameter α_{ER} is size dependent. The argument above is not supposed to be quantitative, and we are well aware that the effect of BFA could be multifactorial. It could for instance prevent the scission of tubules that form continuously in normal conditions but never get long enough to reach the ER. The argument is merely intended to show that the BFA treatment is not in contradiction with our prediction, but requires improvement of the model (and the use of more parameters). We have added a shorter version of the comment above near the end of the Discussion section.
I appreciate the efforts made by the Sens group in using modeling to describe selforganizing systems, but their attempt to describe the Golgi apparatus with this approach is missing key aspects of Golgi organization and dynamics. This has resulted in predictions by the model that are counter to experimental data.
We hope we could convince the reviewer that most of the data she/he thought were conflicting with our predictions are actually in agreement with them, or at least not inconsistent with them, but beyond the scope of such a minimal model. We do agree that the limitations of the model need to be discussed in more detail, and we have done so in the revised version with the help of the reviewer’s comments. We want to stress again that this is a conceptual model, not designed to reproduce specific experiments, but rather to make a conceptual point that Golgi structure and dynamics are intimately linked and cannot be studied separately. We find remarkable that completely different phenotypes can be explained within such a simple model solely by changing rates. We feel that it is important to convey this message to the Golgi community.
Reviewer #2:
Poorly written, hard to understand.
I was excited to read this paper since the work of one of the authors is known to me and I expected an understandable and wellcrafted paper. I was rather crestfallen to find a poorly written paper, and one where I could not figure out what the authors were actually doing.
We are very disappointed that the reviewer found our paper difficult to read. We can assure her/him that this is not due to our sloppiness, as we actually worked hard, but apparently failed, to make our paper understandable without the need for an extensive theoretical background. We have endeavoured to improve the paper readability, with the help of the reviewers’ comments, in the revised version of the manuscript.
One major concern is that this paper is not at all easy to read. It is either i) not selfcontained, or ii) sloppily put together. The authors refer to a "steady state" but do not at any time show that such a state exists, as far as I can see. There is no plot with the time evolution of any quantities being measured. In fact, they admit that there is no steady state, since the size of the system is not constant.
Several of the reviewer’s comment relate to the existence of a steadystate. Let us be very clear that a steadystate does exist in our system. We have added Appendix 4 to the paper including a figure (Appendix 4—figure 1) showing the evolution of the system’s size with time. The organelle forms denovo, and reaches the steadystate – under any parameter values – after a time of order a few 1/K_{f}. There is apparently a confusion in the reviewer’s mind between the variation of the system’s size with time and its variation with the parameters. We aimed at fixing the steadysate systems size to a constant value (independent on budding, fusion and maturation rates) by adjusting the influx to the different parameter values, in order to precisely characterise compartmentalisation. Since this is done through an analytical model that is only rigorously valid under high budding rate, we observe a small deviation from the targeted system size (or order 10%) for low budding rates, as shown in Appendix 5—figure 1. Once the steadysate is reached, the system’s size fluctuates around the steadystate value (characterized by the standard deviation – Appendix 5—figure 1) in a parameterdependent way that we qualitatively understand, and explain in Appendix 5 (and also in the main text in the revised version). Fluctuations are large for small budding rate k_{b} because compartments are large and transient. This feature is realistic and relevant to Golgi organisation. We have expanded the Discussion of the steadystate fluctuations in the revised version of the paper.
They state that "the maximum number of timesteps is typically set to 10^{6} or 10^{7} in order to reach steady state and accumulate enough statistics on all the measured quantities", but that in and of itself is meaningless if one doesn't correctly specify the conditions underlying the measurement.
We agree with the reviewer that giving the typical number of simulation timesteps is not in itself sufficiently informative. We have modified Appendix 1 (and added a discussion of the transient regime in Appendix 4) to explain this better. The steadystate is reached after a time ∼ 1/K_{f}. Depending on the parameter, this corresponds to between 10^{3} and 10^{6} simulation timesteps. In practice, we start the simulations with a system already containing the steadysate amount of cis, medial and trans components to reduces the duration of the transient regime. Our simulations last at least 10^{7} timesteps, which is much longer than any correlation time that might exist in our system (see Appendix 4—figure 1 in the revised version). All quantities reported in the paper are averaged over the full simulation time, excluding the first 10^{6} timesteps to get rid fo the transient regime. We have checked that these quantities do not depend on the total simulation time.
The paper as a whole is an odd combination of being verbose AND telegraphic at the same time. Echoing my above comment, I do not feel that I am getting enough information, enough detail, despite the verbosity. An example: it takes a long time to get to the point of how heterogeneous compartments can form; they say "this homotypic fusion mechanism (relying on local interactions) allows vesicular transport between compartments of different identities – a process that may be regarded as heterotypic fusion. It merely requires that the receiving compartment contains some membrane patches of identity similar to that of the emitted vesicle". Why not just give the rule, rather than the metaphysical description? It would take up less space!
The Golgi field is a highly contentious one, in which there exist people that are quick to dismiss models that do not fit their preconceived ideas. We wrote this particular sentence related to homotypic fusion following comments made by Golgi experts that homotypic fusion was inconsistent with the fact that vesicles could be exchanged between compartments of different identities. While this is obviously false, it unfortunately can be sufficient for some people to reject the model. By trying to protect us against these types of criticisms, we apparently ran into the problem of verbosity. This sentence has been simplified, and the section describing the model has been modified in the revised version of the paper.
Does seemingly everything have to be in the Appendix? One can understand page limits, but this is ridiculous. More explanation of the model in the main text would have been very welcome, because in its present form, the paper is extremely difficult to understand.
We tried to give a qualitative feel for the model and to explain the main results in the main text, and to be rigorous in explaining the role of every parameters and of the model specific assumptions in the Appendices. In the revised version, we have made major changes to the “Model” section following the recommendations of the reviewers.
The question of system size should appear in the main text. As mentioned above, they cannot really keep it constant, so they just adjust the injection rate so it is more or less constant when they change K. They say that they aim for N ~ 300, but the actual system size (sometimes called size, sometimes mass) in the "steady state" changes significantly with K (see Appendix 1—figure 1B and Appendix 4). Moreover, the fluctuations are very large (since no time evolution plot is ever shown, one suspects that they are not really fluctuations…). Of course, the results would be highly dependent on system size, and one gets no justification for their choice of N ~ 300 until the Appendix.
There is apparently another confusion here regarding the system’s size (mentions of the system’s mass has been corrected in the revised version). Appendix 1—figure 1B shows how the size and purity phase diagrams varying the system’s size from N = 30 to N = 1000. The main message is that the structure of the phase diagrams does not depend on system size. The three regimes (mixed, sorted, vesicular) still exist and the transition between them happens for the similar values of parameters. Appendix 5 show a 10% variation of the system’s size with the parameter, which is small, and proves that we can successfully fix the system’s size while varying the parameters. The fluctuations characterised in Appendix 5 (previously Appendix 4) are real fluctuations (also observable in the timevariation of the system’s size shown in the new Appendix 4—figure 1), which depend on the parameter for reasons we understand and explain in Appendix 5 (and also in the main text in the revised version) (low budding rate, large compartments, hence large fluctuations related to compartment’s fusion with the exits). The standard deviation is indeed large in the small budding regime. It amounts to about a third of the system’s size, and corresponds to the typical size of compartments in this regime. This is explained better in the revised version, and is supported by the addition of Appendix 4—figure 1 showing the fluctuations in a time series. The revised version also includes a discussion of the effect of the system’s size in the main text.
The plots that are shown in the paper are not very useful. They rather tend to obfuscate. Figure 2A intends to show a power law for the size of the compartments with an exponential cutoff, but by plotting a bar diagram (even though it is a loglog scale), it is very difficult to see how well this works. They should plot actual points and a fit, or a parallel line showing an actual power law. In a similar vein, they overuse density plots that are not really quantitative (even in the Appendices). In Figure 2BC, they could at the very least plot some contour lines (perhaps dotted) for the same values that are plotted with a continuous line for the analytical approximation  from the colour map, it is difficult to judge.
The figures of the revised version has been modified with the help of these suggestions. The new version of Figure 2 shows the power law as a dashed line. The power law is very accurate up to the crossover size, shown (as before) as a vertical red bar. For Figure 2B and Figure 2C, discrete filled contour plots are used to ease the comparison with the analytical expression. We believe that the new figure clearly show that the continuous approximation qualitatively recapitulate the numerical behaviour quite well. The colors were chosen to be color blind friendly.
Nomenclature:
It is confusing to keep track for reasons that are fully avoidable. For instance, their parameters are (K_{i},K_{f},K_{b},K_{m}): the rates of injection, fusion, budding and maturation. All of them are normalised by Kf, which sets the time unit. But then they define k_{m} = K_{m}/K_{f} and, instead of maintaining the use of lowercase k for the others, they then define K=K_{b}/K_{f}, J=K_{i}/K_{f}. As if the other problems with the paper weren't enough, now the reader has to remember what each means.
This has been fixed in the revised version. We use systematic notations with dimensionless quantities designed by lower case: k_{m} = K_{m}/K_{f}, k_{b} = K_{b}/K_{f}, j = J/K_{f}…
Main result is obfuscated:
The main result, in my opinion, is the existence of an intermediate "sorted" regime. It would be trivial to have a regime with a large mixed compartment or with small pure ones. The potentially interesting aspect of the model is that it can create pure compartments of at least intermediate size, if the parameters are chosen correctly. This is highly dependent on the way budding is modeled  only explained in the Appendix, in an unclear way. In short, the budding flux for species i is defined as
J_{b,}i = K_{b} \times n \times f(\phii)
where n is the size of the compartment and phi_{i} the proportion of vesicles of type i (phi_{i} = n_{i} /n). This is already confusing nomenclature, since they use J_{b} for this flux, when J is also the normalised injection rate. They then say that if the budding rate for species i depends linearly on its concentration, so
J_{b,i} = K_{b} n \phi_{i} = K_{n} n_{i}
then you cannot get the intermediate regime. For that you need what they call "nonlinear budding". In their words: "In order to reproduce this feature we choose a highly nonlinear budding scheme f(phi_{species}) = 1 if phi_{species} > 0 and f(phi_{species}) = 0 otherwise". First of all, they are now using phi_{species} for what was before phi_{i}.
I can only imagine what biological readers might make of this labeled "highly nonlinear scheme". It is, after all, just a constant!
phi_{i} cannot be negative and if phi_{i} = 0, then there are no vesicles to bud. So the budding rate is just a constant, independent of everything and the same for all species in the same compartment.
The existence of a sorted regime with large and pure compartments is indeed an important result. The fact that structural transitions are directly related to transition in the direction of vesicular transport is another crucial result. Our original statement (that nonlinear budding was “essential” to the existence of a sorted regime) was too strong. The existence of a welldefined sorted regime with large compartments is indeed promoted by nonlinear budding. With a linear budding scheme, the purity transition appears for larger values of the budding rate (K_{b}/K_{f} ' 3 Appendix 8—figure 1B, instead of K_{b}/K_{f} ' 0.6 with the nonlinear scheme, Figure 3C). Consequently, compartments are smaller at the transition with a linear budding scheme (n ' 10 instead of n ' 50 with nonlinear budding when the system’s size is fixed at N = 300). The text has been revised to reflect this.
The rationale behind the nonlinear budding scheme is based on the MichaelisMenten kinetics, for which the rate of synthesis (the formation of a budded vesicle in our case) is linear with the substrate concentration at low concentration, and saturates at high concentration. This is a very common kinetics for enzymatic reactions in biological systems, and is relevant to the budding kinetics of coated vesicles. This can be shown with a model accounting for the kinetics of coat protein binding to the membrane, and the kinetics of selective recruitment of membrane components during coat formation (Vagne and Sens, 2018b). Nonlinear budding basically corresponds to a MichaelisMenten kinetics where the Michaelis constant is very low (high affinity between the coat and membrane components) so that saturation occurs for very small values of n_{i}. In the revised version, we have changed the name to “saturated budding” and φ_{species} has been changed into φ_{i}.
Problematic statistical analysis:
The statistical analysis is problematic. First of all, because we don't know whether there is a steady state, can one even define averages? But even if there is a steady state, there are things like the figure in Appendix 4. There, they measure the "temporal standard deviation" of the mass and purity. How are they defining this standard deviation, when the individual measurements come from a correlated time series and are, therefore, not independent? Note that there is a very large dependence with K, which could be trivial for this reason (as the correlation time would obviously depend on K).
As said above, there is always a steady state at which exiting fluxes balance the incoming flux. Fluctuations are correlated in the low budding regime, since they come from the fact that large compartments build up over time, and disappear by fusion with the boundary. This correlation occurs over time scales of order the inverse fusion rate 1/K_{f} and does not challenge the calculation of the temporal standard deviation because this measurement is made on time series over thousands of 1/K_{f}. The revised version contains a new Appendix 4 with a figure showing the fluctuations in a time series. We hope that this figure will make it clear that we are indeed looking at a steadystate, with real fluctuations which statistics is properly quantified. We have expanded the discussion on statistical analysis at the end of Appendix 1 in the revised version.
Conclusion:
I cannot recommend publication of the paper – it is too flawed at present – and importantly, it is not possible to evaluate whether the research presented is even correct. A revised version, should the editors wish to solicit one, at the very least, must show some evidence of a steady state and explain how averages are calculated. The authors should pay special attention and care to make the paper more readable. Much more readable. In particular, the model (the actual model) could be fit in the text and justification/discussion of terms could be left for the Appendix.
We have followed these recommendations to prepare the revised version of the manuscript. We hope that this revised version and the response above will convince the reviewer that our statistical analysis of the model is entirely correct.
Reviewer #3:
This paper explores a simple model for the dynamics of vesicles as they progress through the Golgi apparatus and shows that a number of observed phenomena can be understood based on a small number of biologically, and physically, wellmotivated assumptions.
I thoroughly enjoyed reading this paper, even though I am not an expert on the Golgi. I think it's a great example of theory at its best, applied to an interesting biological problem. The authors start with some simple observations that have been reported in the literature, which they turn them into a mathematical model that upon analysis yields some interesting and experimentally testable predictions. I am optimistic that their results will lead to some new experiments.
The key result of the paper is the connection between the directionality of the vesicular transport and the chemical composition of the Golgi compartments, for which there seems to be experimental evidence. The connection is established within a very simple model whose assumptions seem reasonable. The paper also describes the emergence of a number of structural and dynamical properties of the Golgi, from the same small set of assumptions. A particularly interesting feature of the model is how different structures and dynamics can be accommodated within the model by changing one or more parameters (of which there are only four). For example, the model predicts that there is an inverse relationship between the speed of vesicle transport and fidelity of processing, which is consistent with published comparisons between the structure of Golgi in normal and starved yeast.
The paper is very clearly written, and there were only a few places where I was not quite sure about what the authors had in mind.
We thank the reviewer for this positive comments. We are particularly pleased that the reviewer found our paper very clearly written considering the comments of reviewer 2. It is also comforting that he remarked on the small number of parameters involved in the model, which we hope will help temper the criticisms of reviewers 1 that some experimental observation are not directly explained by the model.
1) When describing the Budding mechanism, there is reference to a nonlinear budding scheme that depends on the compartment size and contamination state. These are described in detail later, but I would recommend that the authors say a few words here to give the reader an intuitive sense of the thing without having to skip ahead.
The budding mechanism is described in more detailed in the main text of the revised version of the paper. More generally, the quantitative description the model, which was initially relegated to the Appendices, now also appears in the main text.
2) In describing maturation, the chemical transformations are described as going only one way. How is this justified? Is it known that the reverse reactions are very slow? Also, is the maturation rate independent of the composition of the compartment and if so, what is the justification for that?
We have discussed our model of maturation in our response to point 3) of reviewer 1 above. The motivation for the maturation scheme is the observation of irreversible maturation of Golgi cisternae in Yeast, which has been explain by the Rab cascade (RiveraMolina and Novick, 2009). While there is no clear evidence that such mechanism is also at play in mammalian Golgi, most of the proteins involved in this mechanisms are conserved in mammals (Klute et al., 2011) and irreversible maturation (now called biochemical conversion) has been put forward as a possible organising principle there as well (Pfeffer, 2010). Our biochemical conversion scheme is designed to be generic, and could also represent the irreversible modification of lipid components, such as phosphoinositides, for which evidences exist in mammalian Golgi (McDermott and Mousley, 2016). We have expanded the paragraph justifying and describing our model of biochemical conversion in the revised version of the paper.
Whether the conversion rate could be cooperative and depend on the local composition is a very interesting question that we plan to examine in the followup paper. We have assumed a constant conversion rate, independent of the local composition. This is likely to be a simplification, as conversion could be a cooperative phenomenon. In a previous paper (Vagne and Sens, 2018a), we have shown that introducing cooperativity of biochemical conversion increases robustness. We choose not to include this possibility here because it requires the introduction of additional parameters. We discuss this (and cite this paper) in the revised version of the manuscript.
3) A typical compartment size is defined as the ratio of the second and the first moment of the size distribution. I was struck by the fact that this would give a "typical size" even for an exponential distribution of sizes, which I would not say defines a typical size (at least not from the point of view that the mean and standard deviation of such a distribution are equal). Might be worth a short comment on why this is a good definition of "typical size".
With our definition of the typical size as the second moment over the mean, there is indeed a typical size for an exponential distribution, which is the characteristic size over which this distribution decays. This is explained in Appendix 2 Eqs 57. This does not mean that most compartments will be of that size, but rather that it is unlikely to find compartments with a size larger that this typical size. The original version of the manuscript was imprecise in that respect, and has been corrected. In our case, the distribution is close to a powerlaw, with an exponential cutoff. As shown in by the red bars in Figure 2A, and mathematically in the Appendix 2, our way to define the typical size appropriately quantify the decay of the size distribution.
4) In the "Steadystate organization" section there is reference to a cutoff size, which I don't think was defined.
We have rephrased this. The cutoff size is intended to describe the size beyond which the distribution decays exponentially, which is also the characteristic size as defined by the ratio of second to first moment.
5) In the Discussion there is a conclusion that the spatial information is less crucial than biochemical composition for the organization of the Golgi. I am not sure this is warranted since a comparison between the two cannot be made in a model lacking spatial information.
This comment was intended to reflect on experimental observations, not on the model. This was indeed unclear and has been rephrased in the revised version. Data suggest that the transport dynamics is not strongly affected by major modifications of the spatial structure of the Golgi. Indeed, transport still occurs (although at a slightly lower rates) when cisternae are locked to mitochondria, more dispersed than in the usual Golgi stack, and unable to fuse with one another (Dunlop et al., 2017). This suggests that the spatial structure may not be crucial to the Golgi dynamics, and justify the relevance of our selforganisation model lacking spatial information. We have rephrased the text in the revised version.
In general I was left wondering whether simple testable quantitative predictions can be made from the model. I very much liked the qualitative observations that the model neatly explains, but I would be even more excited if there was a prospect for putting it to a much more stringent test provided by quantitative experiments, which the authors might propose for an intrepid experimentalist to do.
We thank the reviewer for this comment, which prompted us to think deeper about falsifiable tests of our theory by experimentally exploring the twodimensional phase diagrams we provide. Our main prediction regarding the Golgi structure is that compartment size and purity are linked and should be affected in a correlated fashion by physiological perturbations (as shown in the phase diagrams, Figure 2 BC of the main text). Importantly, these 2D phase diagrams can in principle be fully explored by independently varying the ratios of maturation over fusion rate k_{m} and of budding over fusion rate k_{b}. We present such possibility in Appendix 9 of the revised version of the manuscript.
Compartment purity can be altered without changing their size by acting on k_{m}, while it should be correlated with a change of size by acting on k_{b}. Experimentally, it was found that impairing COPI activity (decreasing k_{b}) decreases purity and increases size (Papanikou et al., 2015), as we predict (red arrow in Appendix 9—figure 1). We can predict further that this loss of purity phenotype could be reversed upon decreasing the maturation rate k_{m} (e.g. by impairing Ypt1 activity – as done in (Kim et al., 2016)) and that this would not change cisternae size (dashed red arrow in Appendix 9—figure 1). It was also found that Ypt1 overexpression (which increases the early to transitional maturation rate) decreases purity as well (Kim et al., 2016). We predict that this should not be associated to change of cisternae size (at steadystate) if altering Ypt1 does not modify the budding rate k_{b} (green arrow in Appendix 9—figure 1). We can predict further that the loss of purity phenotype could be reversed upon increasing the budding rate by overexpressing COPI, but that this would be accompanied by an decrease of cisternae size (dashed green arrow in Appendix 9—figure 1). We will add this proposal to the revised version of the manuscript. Of course, the arrows drawn in Appendix 9—figure 1 are somewhat arbitrary, but we hope that this proposal gives a good feel for the way our model can be used to make new and nontrivial predictions.
[Editors’ note: what follows is the authors’ response to the second round of review.]
Revisions:
1) Regarding prior work, the reviewers note that as the work of Patterson et al., 2008, provided a detailed quantitative model of the Golgi using one set of parameters, and furthermore assessed and tested those parameters by doing extensive quantitative live cell imaging experiments, it would be appropriate to seek a comparison between your model (with its dimensionless parameters) and those results. This would likely enlarge the scope of potential predictions for direct experimental tests. Workers in the Golgi field would benefit from some suggested experiments expected to skew the pathway in one direction or another, with predicted outcomes on the size of the compartment or the speed of traversal etc. For example, the cisternal progression model predicted that larger cargo export rates would be more sensitive to nocodazole dissolution of microtubules and resulting creation of ministacks, a prediction that turned out to be correct: PMID 25103235.
Regarding experimental suggestions, we already provide fairly detailed suggestions of experimental approaches to test some of our predictions, which are discussed in the Discussion section and in Appendix 9—figure 1. These experiments could test the correlation between compartment size and purity which we predict. The present paper focusses on Golgi selforganisation and does not investigate the different transport kinetics for different cargo. Such discussion requires additional assumptions regarding the coupling between the cargo and the local cisterna environment, and will be the subject of a subsequent paper (See also our response to Question 4). We strongly feel that adding results on cargo dynamics to the present paper could obfuscate our main message.
Regarding previous works, and in particular the rapid partitioning model (RPM). We would like to stress that while the RPM and our model share some similarities, they are fundamentally different in their scope. Our model describes denovo Golgi formation and maintenance, through selforganized fusion and scission mechanisms between dynamic compartments, while the RPM focuses on the transport of lipids and cargo molecules in a preestablished Golgi apparatus with a prescribed structure. While we assume that cisternae identity results from a balance between exchange and irreversible biochemical maturation, the RPM assumes that lipid composition (which can be a proxy for cisternal identity) is tuned by lipid transport mechanisms.
In spite of these fundamental differences, the two models converge on key aspects of cargo transport and export mechanisms. In both the RPM and our model, export is allowed from every cisterna, which is crucial to reproduce the exponential kinetics observed in experiments. In addition, no directionality is assumed a priori for intracisternal transport in the RPM, which is in accordance with the budding/fusion mechanism in our model. Finally, in the RPM, the transport rate towards other organelles is different from the intraGolgi transport rate. A similar effect is obtained in our model, through the parameters α_{ER} and α_{TGN} that tune the fusion rate with entry or exit compartments.
The RPM is much more precise than our model, in terms of the microscopic description of the biochemical processes. While we study here the effect of only 4 coarsegrained parameters, there are about 50 different parameters in the full version of the RPM. This makes further comparison difficult, but it shows how our simple theoretical model could evolve to become more realistic. Adding Golgi enzymes, cargo molecules, or lipids to our model, and describing their synthesis and/or sorting into membrane domains, as is done in the RPM, is the logical next step that will push forward our understanding of Golgi selforganisation.
A paragraph summarising the discussion above has been added to the Discussion section.
2) To help the general reader understand the context of this work better, we would like to see an introductory figure illustrating the structure of the Golgi (and perhaps some of the models for transport in the Golgi).
We followed the suggestion and we have extended Figure 1 that now presents an illustration of a classical representation of the Golgi structure, together with four sketches some of the main existing transport models.
3) The reviewers suggest that putting more of the model equations in the main text would streamline the argument; one way of making things more understandable might be to start by writing down simple meanfield equations (as in the early parts of Appendix 3), before delving into the analysis of the detailed model. The backandforth between the main text and the appendices makes for difficult reading.
The mean field model for the composition of the system in the wellsorted limit, and for the size and composition of compartments in the wellmixed limit are now presented at the beginning of the Results section.
4) An important concern is that the model shows how Golgi structure and vesicular transport can arise in a selforganised fashion, but it remains unclear to what extent the authors can claim that vesicular transport in actual Golgi is dominated by this transport arising from selforganisation. When reading the paper, we wondered about the contribution of spatial structure or biochemical processes (that might bias which vesicles cargo goes into, or which compartments a vesicle fuses into). It was only pages later, in the Discussion, that one finds a paragraph on the role of spatial structure: even if experiments show that "Golgi functionality is preserved under major perturbation of its spatial structure", how is transport affected? From a quick look at the reference (Dunlop et al., 2017), it seems that transport can be slowed down massively, hinting that spatial structure is just as important as selforganisation. Are there similar experimental results on biochemical effects?
There are two distinct questions above, one regarding the spatial structure, and one regarding transport biased through biochemical effect.
Regarding spatial structure, the fact that our model does not include space was not clearly discussed in the Introduction. This has been changed in the revised version of the manuscript. We have also made the statement in the Discussion section more precise. What we meant to say was that a Golgi with landlocked cisternae is still able to transport and process cargo. To quote Dunlop et al., 2017: “It is amazing that cargo processing and transport continue with normal efficiency and are merely slowed down when the normal topological relationships among Golgi cisternae are grossly altered […].” Our original statement lacked precision, and we now explicitly write that transport is only slowed down by a factor 2 (from 20 min; to 40 min. on average, see Table 3 in Dunlop et al., 2017), despite major spatial perturbation. In our opinion, this remarkable observation shows fairly convincingly that spatial proximity is not required for functional intraGolgi transport, although it does help increasing its rate. We have modified the discussion of Dunlop et al., 2017 in the Discussion section.
Regarding biochemical bias of intraGolgi transport, it is certainly included in the model as a key element for selforganisation through homotypic fusion of vesicles with compartment. The way we understand this question is whether Cargo transport through the Golgi is bound to follow the average vesicular transport which result from selforganisation. This is a very interesting question, which, as we said above, will be the subject of a followup paper. We have nevertheless added the following comment to the text: Note that vesicular transport of cargo through the Golgi is not bound to follow the net vesicular flux between compartments discussed above. Indeed, the net flux is the average of the flux of vesicles of the three identities. A given cargo follows this flux – on average – if it does not interact preferentially with membrane of particular identity. On the other hand, a cargo that is preferentially packaged into vesicles of trans identity (for example) will be transported toward transrich compartments even if the net vesicular flux is mostly retrograde.
5) One of the main results in the paper is the existence of different regimes of selforganisation, and in particular the existence of an intermediate, sorted regime with large pure compartments. However, the definition of these regimes remains rather qualitative. Perhaps some plots of mean compartment size against mean purity, to complement Figure 2, would put this part of the analysis on a clearer footing.
It is correct that the transition between the vesicular and mixed regimes is rather smooth, so that the definition of an intermediate, sorted regime is somewhat arbitrary. We now show a plot of the mean compartment purity vs. the mean size in Figure 2 E. This shows a clear inflexion point which corresponds to what we call the sorted regime.
6) The model shows very clear anterograde transport of cargo in the limit k_{b}>>1, but the retrograde transport in the limit k_{b}<<1 is less clear; in particular, the model appears to show transport from the trans to wellmixed compartments, but not into the cis compartments: rather, there is weak flow from the cis compartments to wellmixed ones, i.e. anterograde flow in part of the system. Is it possible to show that flow cannot be purely retrograde? – This might be an important constraint on the "cisternal maturation" mechanism. (Perhaps part of this issue is addressed in Appendix 6.)
It is correct that the net retrograde flux in the limit 1 in our model does not correspond to the stereotypical trans → medial → cis flux, but is dominated by fluxes “from transrich toward less mature compartments” as indicated in the original manuscript. The absence of net medial → cis flux in this limit can be understood by analysing the fluxes of vesicles of different identities, shown in Appendix 6. The compartments that contain the highest fraction of medial identity in this regime are well mixed (φ_{cis} ' φ_{medial} ' φ_{trans} ' 1/3). The vesicular flux leaving such compartment is equally split into cis, medial and trans vesicles, which fuse homotypically with compartments enriched in their identity, yielding a vanishing net flux averaged over the three identities. On the other hand, vesicles emitted from early compartments will either undergo back fusion with early compartments or anterograde fusion with more mixed compartments. Since backfusion is not included in the vector plot of Figure 4, the vesicular flux from cis compartment is necessarily anterograde. In our model, flow is never either purely retrograde or purely anterograde. The flow of cis vesicles is always retrograde (except if a vesicle undergoes maturation before fusion) and the flux of trans vesicles is always anterograde, as can be seen in Appendix 6—figure 1. What matters really is the intensity of this flux, which, as the reviewer points out, is very small (essentially negligible relatively to other fluxes) in the low budding regime. This directly translate into a clear retrograde enrichment in Figure 3B.
In the classical “cisternal maturation” mechanism, retrograde transport is needed to recycle cis Golgi components. Within our model, this can in principle be achieved by targeting these components to cis vesicles budding from more mature compartments.
These points are now made clear on the new sketches of vesicular transport in the revised version of Figure 3 and are now discussed in the text.
References
Bhave, M., Papanikou, E., Iyer, P., Pandya, K., Jain, B. K., Ganguly, A., Sharma, C., Pawar, K., Austin, J., Day, K. J., Rossanese, O. W., Glick, B. S., and Bhattacharyya, D. (2014). Golgi enlargement in Arfdepleted yeast cells is due to altered dynamics of cisternal maturation. Journal of Cell Science, 127(1):250–257.
Der´enyi, I., Ju¨licher, F., and Prost, J. (2002). Formation and interaction of membrane tubes. Phys. Rev. Let., 88:238101.
Dunlop, M. H., Ernst, A. M., Schroeder, L. K., Toomre, D. K., Lavieu, G., and Rothman, J. E. (2017). Landlocked mammalian Golgi reveals cargo transport between stable cisternae. Nature Communications, 8(1):432.
Evans, E. and Rawicz, W. (1990). Entropydriven tension and bending elasticity in condensedfluid membranes. Phys. Rev. Let., 64(17):2094–2097.
Glick, B. S. and Luini, A. (2011). Models for golgi traffic: A critical assessment. Cold Spring Harb Perspect Biol, 3:a005215.
Klute, M. J., Melan¸con, P., and Dacks, J. B. (2011). Evolution and diversity of the Golgi. Cold Spring Harbor Perspectives in Biology, 3(8):a007849.
Ladinsky, M., Wu, C., McIntosh, S., McIntosh, R., and Howell, K. (2002). Structure of the golgi and distribution of reporter molecules at 20^{◦}c reveals the complexity of the exit compartments. Mol. Biol. Cell., 13:2810–2825.
Levi, S. K., Bhattacharyya, D., Strack, R. L., Austin, J. R., and Glick, B. S. (2010). The yeast grasp grh1 colocalizes with copii and is dispensable for organizing the secretory pathway. Traffic, 11:1168–1179.
Papanikou, E., Day, K. J., Austin, J., and Glick, B. S. (2015). COPI selectively drives maturation of the early Golgi. eLife, 4.
RiveraMolina, F. and Novick, P. (2009). A rab gap cascade defines the boundary between two rab gtpases on the secretory pathway. P Natl Acad Sci Usa, 106:14408–14413.
Saenz, J., Sun, W., Chang, J. W., Li, J., Bursulaya, B., Gray, N., and Haslam, D. (2009). Golgicide A reveals essential roles for GBF1 in golgi assembly and function. Nat Chem Biol, 5:157–165.
Upadhyay, A. and Sheetz, M. (2004). Tension in tubulovesicular networks of golgi and endoplasmic reticulum membranes. Biophys. J., 86:2923–2928.
Vagne, Q. and Sens, P. (2018b). Stochastic model of vesicular sorting in cellular organelles. Phys. Rev. Lett., 120:058102.
https://doi.org/10.7554/eLife.47318.sa2Article and author information
Author details
Funding
Université Paris Descartes (Ecole Doctorale 474 FIRE  Programme Bettencourt)
 JeanPatrick Vrel
The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.
Acknowledgements
We acknowledge stimulating discussions with Bruno Goud, Cathy Jackson, Grégory Lavieu, Franck Perez and Stéphanie MisereyLenkei, and critical reading of the manuscript by Matthew S Turner.
Senior Editor
 Suzanne R Pfeffer, Stanford University School of Medicine, United States
Reviewing Editor
 Raymond E Goldstein, University of Cambridge, United Kingdom
Reviewer
 Jane Kondev, Brandeis University, United States
Publication history
 Received: April 1, 2019
 Accepted: July 7, 2020
 Version of Record published: August 5, 2020 (version 1)
Copyright
© 2020, Vagne et al.
This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.
Metrics

 1,739
 Page views

 272
 Downloads

 4
 Citations
Article citation count generated by polling the highest count across the following sources: Crossref, PubMed Central, Scopus.
Download links
Downloads (link to download the article as PDF)
Open citations (links to open the citations from this article in various online reference manager services)
Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)
Further reading

 Physics of Living Systems
Schooling in fish is linked to a number of factors such as increased foraging success, predator avoidance, and social interactions. In addition, a prevailing hypothesis is that swimming in groups provides energetic benefits through hydrodynamic interactions. Thrust wakes are frequently occurring flow structures in fish schools as they are shed behind swimming fish. Despite increased flow speeds in these wakes, recent modeling work has suggested that swimming directly inline behind an individual may lead to increased efficiency. However, only limited data are available on live fish interacting with thrust wakes. Here we designed a controlled experiment in which brook trout, Salvelinus fontinalis, interact with thrust wakes generated by a robotic mechanism that produces a fishlike wake. We show that trout swim in thrust wakes, reduce their tailbeat frequencies, and synchronize with the robotic flapping mechanism. Our flow and pressure field analysis revealed that the trout are interacting with oncoming vortices and that they exhibit reduced pressure drag at the head compared to swimming in isolation. Together, these experiments suggest that trout swim energetically more efficiently in thrust wakes and support the hypothesis that swimming in the wake of one another is an advantageous strategy to save energy in a school.

 Physics of Living Systems
When a fish beats its tail, it produces vortices in the water that other fish could take advantage of to save energy while swimming.