Modularity and stability in ecological communities

Networks composed of distinct, densely connected subsystems are called modular. In ecology, it has been posited that a modular organization of species interactions would benefit the dynamical stability of communities, even though evidence supporting this hypothesis is mixed. Here we study the effect of modularity on the local stability of ecological dynamical systems, by presenting new results in random matrix theory, which are obtained using a quaternionic parameterization of the cavity method. Results show that modularity can have moderate stabilizing effects for particular parameter choices, while anti-modularity can greatly destabilize ecological networks.

cological communities are structured both in space, as in a stratified lake and in time, as for migratory birds. This temporal and spatial organization is reflected in the strength of species interactions, where we expect populations dwelling in the same location, or being active in the same season, to interact more frequently than those that are not.
When drawing ecological interaction networks-where species are the nodes and the edges connecting the species stand for interactions (for example, consumption, pollination and competition)-we therefore expect to find that species can be partitioned into distinct 'groups', such that the frequency of interaction largely depends on group-membership 1,2 . Networks with this property are said to be 'block-structured': a modular structure is a particular block structure, in which a network is divided into subsystems, and within-subsystem interactions are much more frequent than those between subsystems 3,4 . If, on the other hand, interactions occur exclusively between groups, we obtain a bipartite network-another type of block structure with many ecological applications 5 .
The appealingly simple idea of a block structure has been formalized in the measure of modularity, Q 3,4 , calculated as the difference between observed and expected within-group interactions, divided by the total number of interactions. Positive values indicate that interactions occur predominantly within-groups, while negative values that interactions are more frequent between-than within-group. Modularity has become one of the most investigated network metrics, with applications spanning biological, social and technological networks [6][7][8][9][10][11] .
Unveiling the relationship between the network structure of biological systems and their dynamical properties has been a long-sought goal of the discipline, and many authors have hypothesized that biological networks must be shaped by evolution [12][13][14][15] and co-evolution 16 , favouring configurations yielding controlled dynamics [17][18][19] .
In ecology, the idea that a modular organization would be beneficial for the local stability of ecological communities (that is, their ability to recover from small perturbations) dates back to work on complexity and stability by May 20 , where he suggested 'that our model multi-species communities [...] will do better if the interactions tend to be arranged in blocks'. This hypothesis was challenged by a number of authors 21,22 , who produced simulations showing the opposite effect. However, recent studies found that modularity can indeed enhance species persistence 23,24 , further complicating the picture. Though the relationship between block structure and dynamics on networks has since been investigated in many fields, including epidemiology 25 , neuroscience 26 and complex systems in general 27 , a systematic classification is still lacking.
Here we study how the local stability of ecological systems is influenced by modularity, providing new results on the theory of random matrices that can be used to draw a direct relationship between modularity Q and local stability. The mathematical results are briefly stated in the Methods section, while a calculation based on the cavity method [28][29][30][31] with quaternions is carried out in the Supplementary Notes.
We show that, with respect to the corresponding unstructured case, modularity can have moderate stabilizing effects for particular parameter choices, while anti-modularity can greatly destabilize networks. The rich range of possible effects associated with the same level of modularity stresses the fact that a given network structure is not 'stabilizing' or 'destabilizing' per se, but only for particular regimes.

Results
Building community matrices. We study the stability of a community matrix M, modelling a continuous-time, dynamical ecological system composed of S populations, resting at a feasible equilibrium point. We remove self-interactions from the matrix (setting M ii ¼ 0), so that we can concentrate on inter-specific effects (adding intra-specific effects would not qualitatively alter the results 32 ).
M is obtained by multiplying element-by-element two matrices, a matrix of interaction strengths, W, where W ij expresses the effect of species j on species i around equilibrium, and the adjacency matrix of an undirected graph, K. The community matrix is therefore M ¼ W J K (Fig. 1).
Initially, we independently sample the coefficients in W in pairs, drawing (W ij , W ji ) from a bivariate distribution with identical marginals, defined by the mean m¼E W ij Â Ã ¼E W ji Â Ã , the variance s 2 ¼E½ W ij À m À Á 2 ¼E½ W ji À m À Á 2 and the correlation r¼ E W ij W ji Â Ã À m 2 À Á =s 2 . By varying these parameters, we can model different types of interactions between the species from preponderantly predator-prey to dominated by competition or mutualism 33 .
The binary matrix K dictates 'who interacts with whom', and is symmetric because we assume pairwise interactions. Thus, K determines which interactions in W are activated, and which are suppressed. Here we study the case of a block-structured K: we assume that the community is composed of two subsystems, of sizes aS and (1 À a)S, respectively (with ar1/2), and that species in the same subsystem interact with The matrix of interaction strengths W, whose diagonal elements are all zero and the off-diagonal elements are independently sampled in pairs from a bivariate distribution (in this case, a bivariate normal with identical marginals defined by m ¼ 0, s ¼ 1/2 and correlation r ¼ À 3/4; red stands for negative coefficients, blue for positive one and the intensity of the colour is proportional to the coefficient value). (b) The adjacency matrix K, which has a block structure; in this case, we have that a ¼ 1/4 of the species belong to the first subsystem, and that the modularity Q ¼ 0.3, meaning that interactions tend to occur within-subsystem. (c) The community matrix M is obtained by multiplying W and K element-by-element.
probability C w (within-subsystem connectance), while species in different subsystems with probability C b (between-subsystem connectance).
This parameterization is especially intuitive, because for C b ¼ C w , we recover the well-studied case of a random ecological community 20,33 . Hence, by varying C b and C w , we can isolate the effect of having a modular or anti-modular structure (Fig. 2), with the two extreme cases being networks composed of two separate subsystems (perfectly modular), and those in which interactions occur exclusively between subsystems (perfectly anti-modular or bipartite). For simplicity, we speak of a 'modular' structure whenever C w 4C b , and of an 'anti-modular' structure when C w oC b . The case of C w ¼ C b represents 'unstructured' systems, such as those studied by May and other authors 20,33 . Equivalently (Methods), we can express C w and C b in terms of the overall connectance C (that is, the overall density of interactions in K) and Q, the modularity of the network 3,4 , defined as: where L w is the observed number of interactions within the subsystems, L b the number of between-subsystem interactions and E L w ½ is the number of within-subsystem interactions we would expect by chance. Values of Q40 (Qo0) mean that we observe within-subsystem interactions more (less) frequently than expected by chance. To calculate E L w ½ , we need to choose a reference model for network structure, and here we use the Erd+ os-Rényi random graph 34 , as we want to contrast our results with those found for an unstructured, random network. Though Q is bounded by À 1 and 1, only a smaller interval of values might be achievable for a given choice of a reference model, a and C (Methods; for example, in our case, when a ¼ 1/2, then QA[ À 1/2, 1/2], provided that Co1/2).
In summary, the parameterization M ¼ W J K allows us to separate the effect on stability of a block structure (modelled by K through a, C and Q) from those due to the distribution of interaction strengths (modelled by W). Given that the case of Q ¼ 0 (unstructured network) has been studied intensively, and that the calculation of the stability of these matrices can be achieved analytically 20,32,33,35 , we use it as a reference point to determine the effect of Q on stability.
Effect of modularity on stability. We want to study the effect of modularity on stability. Therefore, we contrast the value for the real part of the 'rightmost' eigenvalue of M, Re(l M,1 ), with ReðlM ;1 Þ, the value found for e M, a matrix with exactly the same coefficients, but re-arranged according to a random network structure (Q ¼ 0). Re(l M,1 ), is a measure of stability, as it expresses the amount of self-regulation we would need to stabilize the equilibrium 20,32 .
Our analysis ( Fig. 3; Methods) highlights that there are three main parameterizations we need to consider: (a) mean interaction strength close to zero (mE0); (b) strongly negative mean interaction strength; and (c) strongly positive mean interaction strength. Without loss of generality, we can set s 2 ¼ 1, and study the effect of the modularity Q on the stability of the community, measured as the ratio G ¼ Reðl M;1 Þ=ReðlM ;1 Þ, for a given choice of a (controlling the size of the smaller subsystem), r (correlation of interaction strengths) and C (overall connectance of the system). Values Go1 are found when imposing the block structure helps stabilizing the community (for example, in Fig. 3, for Qo0 and m ¼ 0), while ratios G 41 stand for destabilizing effects (for example, any Qa0 for positive mean). For an unstructured matrix (Q ¼ 0), the ratio is exactly 1.
In Fig. 4, we show the effect of modularity on stability in a community composed of 1,000 species, when we set C ¼ 0.2. Take the case of two equally sized subsystems (a ¼ 1/2), for which we derive new results allowing us to express the ratio analytically (Methods; Supplementary Information): when mZ0, we have no effect of modularity on stability; when mo0, on the other hand, a bipartite structure is highly destabilizing, while a modular structure is moderately stabilizing. Both effects are stronger in the case of a negative correlation.
When the two subsystems have different sizes (ao1/2), the stabilizing effect of modularity found for the case of mo0 is greatly diminished, and eventually also modularity can become destabilizing (especially for positive r and a ( 1=2). For m40, any Qa0 is destabilizing, while for mE0, a modular structure is always destabilizing, and a anti-modular structure can be stabilizing, provided that r is sufficiently negative. These effects hold qualitatively for different levels of C (Methods; Supplementary Figs 2-3), with higher connectances leading to more marked effects. Though we cannot predict the ratio in full generality for the case of ao1/2, we can treat the extreme cases of a perfectly modular and perfectly bipartite structure, and this is sufficient to understand the qualitative behaviour of all cases (Methods).
The picture emerging from these results is much more nuanced and complex than what was previously hypothesized [20][21][22] : modularity can have a moderate stabilizing effect when the two subsystem have about the same size, and the mean m is negative or destabilizing, when mZ0 and the subsystems have different sizes. Similarly, anti-modularity is highly destabilizing for ma0, but can be stabilizing for m ¼ 0.
The qualitative behaviour of these systems can be understood quite simply, when considering the distribution of the eigenvalues of the block-structured matrices in the complex plane. As shown in Fig. 5, when there are two subsystems the spectrum of M is composed of a 'bulk' of eigenvalues, and up to two 'outlier', real eigenvalues. When mE0, there are no outliers, and thus stability is determined in all cases by the rightmost eigenvalue(s) in the bulk. When ma0, we have only one outlier in the case of unstructured networks: if mo0, then the outlier lies to the left of the bulk and thus has limited effects on stability; if m40, on the other hand, the outlier lies to the right of the bulk and therefore solely determines stability. The modular case is similar to the unstructured one, though we now have two outilers, in that both lie either to the right (m40) or the left (mo0) of the bulk. In the bipartite case, however, for any ma0 the spectrum presents an outlier to the right (determining stability) and one to the left of the bulk.
These simple observations are sufficient to understand the very strong destabilizing effect of a bipartite structure when mo0: in this case, the stability of the unstructured network, ReðlM ;1 Þ, is determined by the bulk of the eigenvalues, while that of the blockstructured network, Re(l M,1 ), by the outlier to the right of the bulk (Fig. 5, red). When both Re(l M,1 ) and ReðlM ;1 Þ are determined by the bulk (for example, modular case with mo0 or any structure with mE0), the either stabilizing or destabilizing effect is going to be moderate. Moderate effects are also observed when both Re(l M,1 ) and ReðlM ;1 Þ are associated with an outlier lying to the right of the bulk (m40). When both Re(l M,1 ) and ReðlM ;1 Þ are determined by the same type (bulk, outlier) of eigenvalue, the precise stabilizing or destabilizing effect depends nonlinearly on the parameters a, C, Q, m, s and r (Methods).
To summarize, a block structure for an otherwise random ecological system can help stabilization in only two cases: (a) when the structure is modular, and mo0 (though small a or a positive r could reverse this effect); and (b) when the structure is bipartite, mE0, and the correlation is negative. For all the other cases, the effect of a block structure ranges from neutral to highly destabilizing.
Food-web structure. Clearly, ecological systems do not follow a random graph structure, for example, displaying a directionality in the flow of energy from producers to consumers. This directionality proved important in our previous study 36 , where we showed that when the mean of the negative effects dominates that of the positive effects, systems built according to the cascade 37 model (in which 'larger' species consume 'smaller' ones) are more likely to be stable than their random counterparts. We therefore analysed matrices constructed using a variation of the cascade model, where we assign a 'size' to each species and each species can only consume smaller species, and has a preference for those in the same subsystem (Q40), or for those in the other subsystem (Qo0). Note that in this case, we need to set a mean for the positive effects and another one for the negative effects (Methods). Figure 6 shows that the stabilizing effect of modularity found before for the case of mo0 is practically negligible, while the other results are qualitatively the same.

Discussion
We have studied the effect of a modular or anti-modular network structure on the stability of an otherwise random ecological system. Our parameterization makes it easy to compare the effect of the network structure with that one would obtain for unstructured, random systems such as those studied in the past 20,33 : a ratio Go1 stand for stabilizing and G41 for destabilizing effect of network structure. À Á =ReðlM ;1 Þ (y axis), measuring the stabilizing/destabilizing effect of the modularity Q. The ratio is obtained dividing the real part of the leading eigenvalue of the block-structured matrix, Re(l M,1 ), by that of the corresponding unstructured matrix, ReðlM ;1 Þ. Ratios greater than one indicate destabilization, lower than one stabilization and whenever Q ¼ 0 (unstructured matrices) we should find G ¼ 1. For a given connectance (C ¼ 0.1 in this case), correlation (r ¼ À 3/4) and size of the smaller subsystem (a ¼ 1/4), we fix the variance s 2 ¼ 1 and study three different cases, in which the mean m is negative (red), positive (blue) or zero (green). The figure highlights that the same modularity Q, all other parameters being fixed, can have completely different effects on stability, depending on the value of m. Each point is obtained averaging the ratio over 50 replicates of a community matrix containing 1,000 species (Methods).
For block-structured matrices, we showed that modularity can have a positive effect on stability only when (a) the system is composed of two subsystems of about the same size and (b) the overall mean interaction strength is negative. The stabilizing effect is stronger for negative correlations. Anti-modularity, on the other hand, is typically strongly destabilizing, besides the case in which the average interaction strength is close to 0 (a well-studied case 20,33 , though of limited biological realism). When the mean interaction strength is positive, both modularity and anti-modularity are destabilizing.
Through numerical simulations, we have investigated the more complex case of an interaction between modularity and food-web structure, and found that the results are qualitatively unchanged.
The picture emerging from both simulations and mathematical analysis is much more complex than previously hypothesized. Block structure can have an effect on the local asymptotic stability of the underlying system. However, unless we are in particular areas of the parameter space, the effect tends to be destabilizing. Our results stress the fact that, when discussing the relationship between network structure and local stability, we need to qualify our statements, as a given structure is not stabilizing or destabilizing per se, but is only so under certain specific conditions.
Though we have illustrated this point by studying the modular structures, we believe this phenomenon to hold generally: any network structure can have different effects on stability, depending on the choice of parameters. To reinforce this message, in Fig. 7, we show three cases in which an empirical network structure makes the system more or less stable than its random counterpart, depending on the parameterization of the coefficients.
Practically, this means that the challenge of proving that biological network structure emerges because of a selective process, removing configurations yielding unfavourable dynamics [17][18][19] is much harder than expected: network structure, without estimates of the distribution of the coefficients, cannot be used to determine the effect on dynamical properties.

Methods
Constructing the community matrix. M is the S Â S community matrix, representing the population dynamics of an unknown dynamical system around a feasible equilibrium point. We consider two cases: (a) random ecological networks with block structure and (b) food webs with block structure.
For case (a), we sampled the pairs (W ij , W ji ) independently from a bivariate normal distribution with means (m, m) T , and covariance matrix For case (b), we first assigned a 'size' to each species (randomly sampling it from a uniform distribution between 0 and 1), and then sampled the pairs (W ij , W ji ) independently from a normal bivariate distribution with means ((1 þ x)m, (1 À x)m) T (ensuring that E W ij Â Ã ¼m, for Fig. 6, x ¼ 3), and covariance matrix S, whenever i was larger than j. This means that species can only consume 'smaller' prey, such as in the cascade model. In particular, in this case we could order the rows and columns of matrix W so that all the positive effects would be confined to the upper-triangular part and the negative effects to the lower-triangular part. For each value of a (determining the size of the smaller subsystem; rows) and r (expressing the correlation between interaction strengths; columns), we vary Q and record the ratio G. Given that the maximum and minimum Q that can be attained depend on C and a (Methods), we take 20 equally spaced points between the minimum and maximum Q for each configuration. We set C ¼ 0.2, S ¼ 1,000 and s 2 ¼ 1, and track the case of m ¼ 0 (green), m ¼ À 1 (red) and m ¼ 1 (blue). The dots represent numerical simulations, obtained by averaging over 50 replicates. The open circles are the corresponding analytical predictions (Methods). Because for negative means the destabilizing effect of bipartite structures is so strong that plotting it would make the other effects difficult to see, we only plot the region of G bounded by 1/2 and 3/2. See Supplementary Fig. 1 for the complete graph.
In both cases, the matrix M is obtained via the Hadamard (element-by-element) product of W and the adjacency matrix K. The matrix K is characterized by four parameters: S, a, C and Q. S is the size, a is the proportion of species belonging to the first subsystem, Q is the modularity and C the overall connectance of K (density of the nonzero elements). The first aS species are assigned to the first subsystem, and the remaining (1 À a)S species to the second subsystem. The vectorg encodes the subsystem-membership of each species. Then, we set (K ij , K ji ) to (1, 1) with probability C w , when g i ¼ g j and with probability C b when g i ag j . The 'within-subsystem connectance' C w is: and the 'between-subsystem connectance' C b is: Note that, given the Erd+ os-Rényi reference model, the values of Q that are attainable depend on both a and C: Numerical simulations. For each a and r, we set S ¼ 1,000, C ¼ 0.2, s ¼ 1 and m ¼ 0 (green), m ¼ À 1 (red) or m ¼ 1 (blue), and varied Q from its minimum to its maximum possible value in twenty equally sized steps. For each parameter set, we produced 50 block-structured matrices M, and 50 unstructured matrices e M, obtained by setting Q ¼ 0. The ratio G ¼ Re l M;1 À Á =ReðlM ;1 Þ was computed by averaging over the replicates. In many cases, one can obtain the expectation for the ratio analytically (below; Supplementary Information). The simulations were repeated for both random ecological networks (Fig. 4) and cascade-based food webs (Fig. 6).
The spectrum of block-structured matrices. For our derivations, we adopt a slightly more general notation, which includes that discussed above as a special case. We consider the matrix M, with M ii ¼ 0 and the off-diagonal coefficients independently sampled in pairs from either of two distributions: Hence, the pairs come from a certain bivariate distribution Z w , when i and j belong to the same subsystem, while from a different distribution Z b , when i and j belong to different subsystems. We do not need to specify the exact form of the distributions Z w and Z b , given that, as for many results in random matrix theory, our findings are consistent with the 'universality' property 38,39 : once fixed the mean and covariance matrices of Z w and Z b , and provided that the fourth moment of each is bounded, any choice of distributions yields the same result, for S-N.
When considering the case examined in the main text, where the elements M ij ¼ W ij K ij , the universality property helps us in two ways. First, consider that the pairs (M ij , M ji ) are zero with probability 1 À C w (case g i ¼ g j ) or probability 1 À C b (case g i ag j ), and that the nonzero pairs are sampled from a bivariate distribution (for the elements of W), defined by the parameters m, s 2 and r. This is sufficient to Q > 0 (modular) Figure 5 | Effects of modularity on the shape of the spectrum. The spectrum of block-structured matrices is given by a 'bulk', and-when there are two subsystems-up to two 'outlier' eigenvalues (highlighted points). The stability of the system is determined by the rightmost eigenvalue(s) (crosses), which can either be one of the outliers (cases with m40 and when Qo0, mo0), or the rightmost eigenvalue(s) of the bulk (cases with m ¼ 0, and mo0 with Qr0). Large effects of Q on stability are found when varying Q changes the type of eigenvalue determining stability (mo0), while, when the type of eigenvalue determining stability does not depend on Q, effects will be moderate. The panels show the eigenvalues of a single matrix with S ¼ 1,000,  Fig. 4, but with the matrix W built such that when two species i and j interact, and species i is 'smaller' than j, then the effect of j on i, W ij is negative (on average), while that of i on j positive. This means that the food web is structured, as in a cascade model with block structure (Methods). While the results are generally similar to those in Fig. 4, food-web structure greatly decreases the stabilizing effect of modularity found when mo0 and aE1/2 in the random case. The same qualitative results are found when varying the connectance C ( Supplementary Figs 4-6).  We show the results of numerical simulations where we parameterize empirical networks, and we study the stabilizing/destabilizing effect of network structure (G, y axis) when varying a 'critical' parameter y (x axis). When parameterizing a high-school social network 42 (left), y is the average of the nonzero coefficients: for negative y, we observe that matrices with empirical structure are easier to stabilize than their random counterparts (Go1); the reverse is found for positive means. For the parameterization of the Little Rock food web 43 (centre), y expresses the relative magnitude of positive effects (effects of prey on predators) with respect to the negative ones (effects of predators on prey), so that when yE1 positive effects are about as strong as the negative ones, and network structure is destabilizing; for y much lower (higher) than one, on the other hand, network structure is stabilizing. Finally, for a pollination network 44 This means that the effect of the connectances is somewhat trivial: we can 'absorb' the connectances in the parameters m w , m b , s w , s b , r w and r b , which are our 'effective' parameters, dictating the shape of the distribution of the eigenvalues of M.
The second advantage of having universal results is that we are free to choose any distribution for the pairs (W ij , W ji ). In all our examples, these are sampled from a bivariate normal distribution.
Decomposition. We want to study the limiting (that is, when S is large) distribution of the eigenvalues of M for the case of a random community (for the food webs following the cascade model, we rely exclusively on simulations). Following the approach by Allesina et al. 36 , we write the matrix M as the sum of two matrices, where A is a matrix with block structure whose elements are and B is obtained by difference: B ¼ M À A. Thus, the diagonal elements of . This parameterization is very convenient, as the spectrum of matrix B describes the bulk of the eigenvalues of M, while the outlier eigenvalues of M are given by the nonzero eigenvalues of A, modified by a small correction 40 .
The eigenvalues of A. The eigenvalues of A are easy to obtain for any choice of S, a, m w and m b , with all being zero besides which can be both zero as well (m w ¼ m b ¼ 0), both different from zero (aa1/2, m w am b ), or one zero and one nonzero (m w ¼ m b a0). Thus, there are going to be up to two outlier eigenvalues. These are the approximate locations of the two outlier eigenvalues of M (only one outlier when m b ¼ m w a0, as found for example in the 'unstructured' case). The exact location of the outliers depends also on B, as explained below.
and a variance that is the arithmetic mean of the variances in B, In this case, it is convenient to express these values in terms of m, s 2 , r, C and Q, because this makes the role of modularity in modulating the stability much clearer: e s 2 ¼ Cðs 2 þ ð1 À C À 4CQ 2 Þm 2 Þ ð 16Þ From the two equations above, it is clear that (a) for m ¼ 0, modularity has no effect on the spectrum, while for ma0 the sign of m does not affect the spectrum; (b) modular and bipartite structures have the same effect: the eigenvalues of B will be approximately the same, when we have Q ¼ q or Q ¼ À q; (c) the effect of modularity is going to be more marked for large C or |m|; and (d) the radius of B is always lower or equal than that we would find by setting Q ¼ 0.
Summarizing, for a ¼ 1/2, the eigenvalues of B are contained in an ellipse whose horizontal semi-axis is always smaller or equal than that found for the corresponding unstructured matrix. This explains the stabilizing effect of a modular structure (Qo0) we observed for mo0, as in that case the rightmost eigenvalue of M is the rightmost eigenvalue of B.
New case: r w ¼ 0 (perfectly bipartite). When s w ¼ 0 (that is, C w ¼ 0), the nonzero coefficients of B are exclusively contained in the two blocks, representing the interactions between subsystems, as in a bipartite network. Hence, we can write the matrix B in block form: where X is a aS Â (1 À a)S matrix and Y is a (1 À a)S Â aS matrix. The two matrices on the diagonal contain all zeros, and have size aS Â aS and (1 À a)S Â (1 À a)S, respectively. The eigenvalues of B 2 , are the eigenvalues of B, squared: if l i is an eigenvalue of B, then l 2 i is an eigenvalue of B 2 . Squaring B, we obtain: XY has aS eigenvalues, while YX has (1 À a)S eigenvalues. The eigenvalues of YX are the same as those of XY, with the exception of (1 À a)S À aS eigenvalues which are exactly 0 (take v to be an eigenvector of XY. Then XYv ¼ lv. Let w ¼ Yv, hence, XYv ¼ Xw ¼ lv. Finally, consider YXw ¼ Y(lv) ¼ lYv ¼ lw. Thus, if la0, it is an eigenvalue of both XY and YX). Hence, we can study the eigenvalues of XY (the smaller matrix) without loss of generality.
In the Supplementary Information, we show that the eigenvalues of XY are contained in an ellipse in the complex plane, with centre in x c ¼Sr b s 2 b ; 0 À Á , horizontal semi-axis r x ¼Ss 2 b ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi a 1 À a ð Þ p 1 þ r 2 b À Á and vertical semi-axis r y ¼Ss 2 b ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi a 1 À a ð Þ p 1 À r 2 b À Á ( Supplementary Fig. 7, top). If the support of the distribution of the eigenvalues of XY is an ellipse in the complex plane, then the support of the eigenvalues of B for perfectly bipartite interaction matrices is obtained via a square-root transformation of the ellipse in the complex plane ( Supplementary Fig. 7, bottom), with the addition of the point (0, 0), stemming from the extra eigenvalues of YX.
To find the real part of the rightmost eigenvalue of B, we thus need to consider the square-root transformation of the ellipse found for XY. The eigenvalues of XY, which are the squared eigenvalues of B, are contained in the ellipse: x À x c ð Þ 2 r 2 x þ y 2 r 2 where x is the real part and y the imaginary part of the point z ¼ x þ iy (we can consider the case y40 without loss of generality, given that the spectrum of B is symmetric about the real and the imaginary axis). Then, the square root of z has