Mitochondrial network complexity emerges from fission/fusion dynamics

Mitochondrial networks exhibit a variety of complex behaviors, including coordinated cell-wide oscillations of energy states as well as a phase transition (depolarization) in response to oxidative stress. Since functional and structural properties are often interwinded, here we characterized the structure of mitochondrial networks in mouse embryonic fibroblasts using network tools and percolation theory. Subsequently we perturbed the system either by promoting the fusion of mitochondrial segments or by inducing mitochondrial fission. Quantitative analysis of mitochondrial clusters revealed that structural parameters of healthy mitochondria laid in between the extremes of highly fragmented and completely fusioned networks. We confirmed our results by contrasting our empirical findings with the predictions of a recently described computational model of mitochondrial network emergence based on fission-fusion kinetics. Altogether these results offer not only an objective methodology to parametrize the complexity of this organelle but also support the idea that mitochondrial networks behave as critical systems and undergo structural phase transitions.

where X i (i = 1, 2, 3) are the number of nodes with degree i. Tip to tip reactions happen with association (dissociation) rate a 1 (b 1 ) between a random chosen pair of nodes with degree one (association) or for a random chosen site with degree 2 (dissociation). Tip to side reactions happen with association (dissociation) rate a 2 (b 2 ) between a random chosen pair of nodes with degrees one and two (association) or for a random chosen site with degree 3 (dissociation). Following Sukhorukov et al. [1], we assumed b 2 = (3/2)b 1 and varied the relative rates c i = a i /b i .
Monte Carlo simulations were performed using Gillespie algorithm [3]. Most of the simulations were performed for L = 15000, the estimated average number of edges in the control cells images. Preliminary tests showed that the the number of nodes with degree k become stationary when the number of iterations is ≈ 2L. So we run every simulation 3L iterations after which we measured different quantities. This procedure was repeated 100 times for different sequences of random numbers and the different quantities were averaged over this sample. The different quantities measured were: the average degree k (where the average is taken both over all the nodes in the network and over different runs), the average fraction of nodes in the largest cluster N g /N (order parameter of the percolation transition) and the average cluster size excluding the largest cluster s (again the average is taken both over the network and over runs). s was calculated using the expression from classical percolation theory [2], namely if N s is the number of clusters of size s and n s = N s /N , then s = s s 2 n s s sn s where the sums exclude the largest cluster in the network. In Fig.1 we illustrate the typical behavior of the different quantities as a function of c 2 for a fixed value of c 1 . As described in Ref. [1], for a given value of c 1 , c 2 controls the percolation transition, characterized as a peak in s , which happens when the order parameter (probability of a "giant component or cluster") is ≈ 0.3. Notice that the average degree is smaller than 2 around the critical region (onset of percolation). Such behavior was observed for all the range of values of c 1 of interest (i.e., those that correspond to values of the measured quantities observed in the experiments).
In Fig.2 we show a parametric plot of the order parameter vs. the average degree for a wide range of values of c 1 . We see that all the curves intersect approximately at the same point, corresponding roughly to N g /N = 0.8 and k = 2. Below those values, every point have associated well defined values of the parameters (c 1 , c 2 ).

Experimental supplementary information
This section contains additional evidence supporting the main hypothesis as well as details about the network extraction procedure.
• The performance of the branching detection method is presented in Table 1. To evaluate the ability of our algorithm to detect branching points we selected multiple ROIs from confocal images corresponding to ctl, pqt and mfn treated cells (n=10). The number of true brances (TP) was obtained by visual inspection and constituted our gold standard set (GS). We then ran the algorithm on each ROI using three different thresholds (th = 0.1, 0.125 and 0.15) and extracted the predicted branching points in each case. Predicted branching points not present in the GS were counted false positives (FP). Conversely, branches present in the GS but absent in the predicted set were considered false negatives (FN). Using the described parameters, we computed the sensitivity of the procedure, defined as s = T P T P +F N , and its precision, defined as p = T P T P +F P . Values are presented in in Table 1 • It is possible that a 2D projection of the 3D mitochondrial structure might introduce an artificial superposition of network branches. The most probable artifact shall be the presence of nodes with degree four, byproduct of the intersection (in 2D) of linear segments. This issue was investigated by computing the number of degree four nodes in all conditions. Figure 3 contains the results of the calculations done to demonstrate that the frequency of the artifactual branching points is negligible.
• Figures 4 and 5 present data from additional experimental manipulations. BAPTA-AM, a cell-permeant Ca 2+ chelator that boost ER-mediated mitochondrial fission, or FCCP, an uncoupling agent that cause mitochondrial membrane potential loss, were used to promote mitochondrial fission [4,5]. Also, based on the fact mitochondrial structure rely on the cytoeskeleton, pharmacological modulation of actin filaments and microtubules was performed  Table 1: Performance of the branching detection procedure. Sample ROIs from ctl, pqt and mfn images were analysed. Three different thresholds were used to convert images into skeletons. Sensitivity and precision were computed using T P T P +F N and T P T P +F P , respectively. by using Jasplakinolide or Taxol, respectively [6,7,8]. Finally, we over-expressed dynamin related protein 1 (DRP1), a pro-fission protein, or mitofusin 2 (MFN2), a pro-fusion protein. Figure 4 presents examples of phenotypes obtained using each of the treatments and Figure  5 the structural analysis of the mitochondrial networks derived from those images. Notice that the results in Figure 5 are fully consistent with those presented in the main manuscript, confirming the generality of the conclusions obtained using paraquat treatment and mitofusin 1 over-expression.
• Each CCDF shown in Figure 6 was computed by selecting random sets of 5 networks of the corresponding condition. The result demonstrates that differences depicted in Figure 4 of the main manuscript are statistically significant.  (d) (c) Figure 5: Changes in mitochondrial network structural properties upon fission/fusion balance alterations. Black circles correspond to control condition, while red squares and blue triangles correspond to drugs that promote mitochondrial fission and fusion, respectively. (a) Complementary cumulative distributions of cluster sizes. (b) Average fractal dimensions of the networks across the threshold range. Asterisks indicate statistically significant differences between ctl and treated networks (p < 0.05). (c) The deviation from a random structure was calculated using Kolmogorov complexity. Asterisks indicate statistically significant differences between ctl and treated networks (p < 0.05). Inset panel shows the deviation from the control structure of fissioned and fused networks. (d) Shanon entropy of cluster mass distributions. Shifts in the critical threshold th * upon treatment are depicted (red and blue arrows).  A random permutation procedure was implemented to calculate the average cumulative distribution of cluster sizes for sets of 5 networks. Within each group, every combination of 5 out of 10 networks was used to build the average distributions.