1 Introduction
Ecological theory has long sought general principles linking the architecture of interaction networks to ecosystem stability. Early theoretical work framed this relationship in terms of system complexity, most famously through the analysis of random community matrices by Robert May, who showed that increasing species richness, interaction density, and interaction strength variance should reduce local dynamical stability in large ecological systems (May 1972, 1974). Subsequent work refined this prediction by demonstrating that specific forms of network organisation such as interaction structure, trophic hierarchy, and compartmentalisation can substantially alter stability outcomes (McCann 2000; Allesina and Tang 2012; Stouffer and Bascompte 2010). More recently, structural theories such as trophic coherence have proposed that food webs occupy constrained regions of structural space in which particular architectural features strongly influence system stability (Johnson et al. 2014). At the same time, ecological stability itself is now widely recognised as a multidimensional concept encompassing distinct dynamical properties such as persistence, resistance to disturbance, and recovery dynamics (Domínguez-García et al. 2019; Ives and Carpenter 2007; Loreau and de Mazancourt 2013; Donohue et al. 2016; Chen et al. 2024). Together these developments suggest that the relationship between network structure and ecological stability may be more complex than originally assumed, specifically that both structure and stability are multidimensional properties of ecological systems. If so, apparent inconsistencies in the structure–stability literature may arise not because theoretical predictions are incorrect, but because different studies measure different aspects of network architecture and different components of ecological stability.
Ecological networks are commonly characterised using a large and diverse set of structural metrics, yet there remains little consensus on how these metrics should be interpreted or compared across studies (Vermaat et al. 2009; Lau et al. 2017). Measures of connectance, trophic level, modularity, omnivory, motif structure, and many others are routinely invoked as indicators of ecosystem organisation or stability. However, these descriptors often covary strongly, capture overlapping information, or reflect different scales of network organisation (Thompson et al. 2012). As a result, studies using different subsets of metrics frequently arrive at conflicting conclusions about the relationship between network structure and ecological stability. Rather than reflecting inconsistency or redundancy, these differences arise because different metrics implicitly target different components of stability (Lau et al. 2017).
A key, but often implicit, assumption in this literature is that stability is defined primarily through the temporal behaviour of populations (Pimm 1984; Ives and Carpenter 2007). Here, we extend this view by considering stability not only as whether populations return to equilibrium, but whether the network that constrains those dynamics continues to exist and function. This perspective aligns with a broader view emerging across complex systems, in which stability is associated with the persistence of network structure and function under perturbation, including in financial systems, infrastructure networks, and social systems (Albert et al. 2000; Burleson-Lesser et al. 2020). This perspective is particularly relevant for empirical food webs, which are typically observed as static interaction networks and therefore encode structural constraints on dynamics rather than dynamics directly.
From this perspective, both the structure and stability of a network can be understood as multiscale properties. Structural metrics capture different aspects of how interactions are organised, from species-level roles to global network architecture, and should not be expected to provide universal predictors of stability (Newman 2010). Likewise, different stability components emerge from processes operating at different levels of organisation. Node-level properties relate to species persistence and redundancy (Dunne et al. 2002), path-level structure governs the propagation of perturbations and thus resistance (Albert et al. 2000), and global network properties constrain system-wide behaviour and recovery potential (May 1972; Allesina and Tang 2012).
Within this framework, commonly used metrics can be reinterpreted as complementary descriptors of distinct stability mechanisms rather than competing proxies for a single concept of stability. For example, controllability identifies influential species at the node level (Liu et al. 2011; Cowan et al. 2012), robustness to species loss captures the maintenance of connectivity under perturbation (Dunne et al. 2002), subgraph structure reveals cohesive subcomponents of the network (Kitsak et al. 2010), and emergent properties reflect emergent constraints on system-level dynamics (Allesina and Tang 2012). Many of these metrics span multiple scales, linking local structure to global behaviour.
Adopting a multiscale perspective may provide a mechanistic explanation for why different studies report contrasting structure–stability relationships. Rather than reflecting inconsistency, such differences arise because different metrics implicitly target different stability components operating at different organisational scales. Identifying whether commonly used metrics organise into a smaller set of coherent structural dimensions therefore offers a way to reduce redundancy while preserving ecological meaning.
Here we investigate whether the large number of structural metrics used to characterise food webs can be reduced to a smaller set of coherent structural dimensions, and whether these dimensions relate to different components of ecological stability. Using a compilation of empirical food webs, we quantified a broad suite of commonly used network descriptors spanning species roles, interaction pathways, and global network organisation. We first examined whether these metrics self-organise into statistically robust modules based on their covariation across ecosystems. We then evaluated how these modules align with dominant multivariate axes of structural variation. Finally, we compared alternative representations of network structure to assess how different structural dimensions relate to three stability properties: robustness to species loss, resilience, spectral radius, and controllability. By explicitly mapping the multivariate structure of network architecture to stability outcomes, this analysis provides a principled framework for interpreting network metrics and for selecting structural descriptors that match specific ecological questions.
2 Materials & Methods
2.1 Data Compilation
We compiled food web data from Mangal (Poisot et al. 2016), Web of Life (Fortuna et al. 2014), and the canonical networks used by Vermaat et al. (2009), resulting in a total of XX networks. All networks were treated as binary and characterised using a suite of XX structural metrics (see Table 1).
| Label | Definition | Structural interpretation | Reference |
|---|---|---|---|
| Basal | Proportion of taxa with zero vulnerability (no consumers). | Quantifies the proportion of species representing basal energy inputs to the network. | |
| Top | Proportion of taxa with zero generality (no resources). | Describes the relative prevalence of terminal consumers in the network. | |
| Intermediate | Proportion of taxa with both consumers and resources. | Captures the proportion of species participating in both upward and downward energy transfer. | |
| Richness (S) | Number of taxa (nodes) in the network. | Describes network size. | |
| Links (L) | Total number of trophic interactions (edges). | Describes interaction density independent of network size. | |
| Connectance | \(L/S^2\), where \(S\) is the number of species and \(L\) the number of links | Measures the proportion of realised interactions relative to all possible interactions. | Dunne et al. 2002 |
| L/S | Mean number of links per species. | Captures average interaction density per taxon. | |
| Cannibal | Proportion of taxa with self-loops. | Quantifies the prevalence of cannibalistic interactions. | |
| Herbivore | Proportion of taxa feeding exclusively on basal species. | Describes the representation of primary consumers. | |
| Trophic level (TL) | Prey-weighted trophic level averaged across taxa. | Captures the vertical organisation of energy transfer. | Williams and Martinez (2004) |
| MaxSim | Mean maximum trophic similarity of each taxon to all others. | Quantifies functional similarity based on shared predators and prey. | Yodzis and Winemiller (1999) |
| Centrality | Node centrality averaged across taxa (definition-dependent). | Captures the distribution of influence or connectivity among species. | Estrada and Bodin (2008) |
| ChLen | Mean length of all food chains from basal to top taxa. | Describes the average number of steps in energy-transfer pathways. | |
| ChSD | Standard deviation of food chain length. | Captures variability in pathway lengths. | |
| ChNum | Log-transformed number of distinct food chains. | Quantifies the multiplicity of alternative energy pathways. | |
| Path | Mean shortest path length between all species pairs. | Describes the average distance between taxa within the network. | |
| Diameter | Maximum shortest path length between any two taxa. | Captures the largest network distance between species. | |
| Omnivory | Proportion of taxa feeding on resources at multiple trophic levels. | Describes vertical coupling of energy channels. | McCann (2000) |
| Loop | Proportion of taxa involved in trophic loops. | Quantifies the prevalence of cyclic interaction pathways. | |
| Prey:Predator | Ratio of prey taxa (basal + intermediate) to predator taxa (intermediate + top). | Describes the overall shape of the trophic structure. | |
| Clust | Mean clustering coefficient. | Measures the tendency for taxa sharing interaction partners to also interact with each other. | Watts and Strogatz (1998) |
| GenSD | Normalised standard deviation of generality. | Captures heterogeneity in the number of resources per taxon. | Williams and Martinez (2004) |
| VulSD | Normalised standard deviation of vulnerability. | Captures heterogeneity in the number of consumers per taxon. | Williams and Martinez (2004) |
| LinkSD | Normalised standard deviation of total links per taxon. | Quantifies variation in species connectivity. | |
| Intervality | Degree to which taxa can be ordered along a single niche dimension. | Measures the extent of niche ordering in trophic interactions. | Stouffer et al. (2006) |
| S1 (Linear chain) | Frequency of three-node linear chains (A → B → C) with no additional links. | Captures the prevalence of simple, unbranched energy-transfer pathways. | Stouffer et al. (2007) Milo et al. (2002) |
| S2 (Omnivory) | Frequency of three-node motifs forming a feed-forward loop (A → B → C, A → C). | Describes vertical coupling of trophic levels within small subnetworks. | Stouffer et al. (2007) Milo et al. (2002) |
| S4 (Apparent competition) | Frequency of motifs where one consumer feeds on two resources (A → B ← C). | Captures the prevalence of shared-predator structures among resources. | Stouffer et al. (2007) Milo et al. (2002) |
| S5 (Direct competition) | Frequency of motifs where two consumers share a single resource (A ← B → C). | Describes the occurrence of shared-resource structures among consumers. | Stouffer et al. (2007) Milo et al. (2002) |
| trophicVar | Measure of how much the trophic positions of species deviate from the mean trophic level. | Variance is linked to chain length, with low variance indicating few long chains and high variance indicates long chains and varied interactions | Pimm (1982) |
| TrophicCoherence | Measure how well the species in a food web fit into discrete trophic levels. | Highly coherent, neatly layered food webs are more stable to perturbations | Johnson et al. (2014) |
In addition to network structural descriptors, we quantified multiple complementary measures of ecological stability that capture distinct aspects of system response to perturbation. First, robustness to species loss (\(R_{50}\)) quantifies the proportion of primary extinctions required to reduce network richness to 50% of its initial value, providing a threshold-based measure of the system’s tolerance to species removal and its ability to maintain connectivity under perturbation (Dunne et al. 2002). Second, we quantify resilience as the area under the extinction curve, relating primary species loss to secondary extinctions, which captures how rapidly structural collapse propagates through the network as perturbations accumulate and thus reflects the rate of system degradation under sequential disturbance. Third, spectral radius (\(\rho\)), defined as the largest real part of the eigenvalues of the adjacency matrix(Staniczenko et al. 2013) making it a proxy for network persistence (Bascompte et al. 2003). This provides a global constraint on system behaviour and can be interpreted as a proxy for the system’s dynamical stability potential, governing the amplification or decay of perturbations. Finally, structural controllability quantifies the minimum number of driver species required to control the system’s dynamics, capturing the capacity for directed reorganisation following perturbation and identifying species that exert disproportionate influence on system behaviour (Liu et al. 2011). Together, these metrics capture complementary aspects of ecological stability: (i) resistance to perturbation, quantified by robustness as tolerance to species loss; (ii) recovery/persistence after perturbation, measured as the propagation of cascading extinctions following sequential species removal; (iii) stability potential, captured by the spectral radius as a global constraint governing the amplification or decay of perturbations; and (iv) controllability, reflecting the minimum set of driver species required to steer system dynamics and reconfigure network structure.
Finally, we incorporated two complementary measures that capture higher-order aspects of network organisation. First, we quantified complexity based on the singular value spectrum of the adjacency matrix. Specifically, we computed the entropy of the singular values obtained via singular value decomposition (SVD), which captures how structural variance is distributed across orthogonal modes. Networks with more evenly distributed singular values exhibit higher complexity, reflecting greater heterogeneity in interaction pathways and reduced dominance of any single structural mode (Strydom et al. 2021). This provides a spectral analogue to traditional structural descriptors, linking network organisation to the distribution of interaction strength across latent dimensions.
Second, we included a composite scaling term of the form \(-\frac{1}{2}(\log{(Richness)} + \log{(Connectance)})\) This term captures the joint scaling of network size and density and can be interpreted as a log-transformed geometric mean of these quantities. This formulation is motivated by classical stability results, where system behaviour depends on the combined scaling of species richness and interaction density, rather than either quantity in isolation. As such, this term provides a compact representation of size–complexity tradeoffs and complements both local (e.g., degree distributions) and global (e.g., spectral) descriptors of network structure.
Together with SVD-based complexity, this additional term allows us to (respectively) capture the distribution of structural information across latent spectral modes and the scaling relationships governing network size and interaction density..
2.2 Identification of Structural Modules
We tested the hypothesis that structural descriptors of food webs are organised into statistically coherent modules reflecting shared ecological function or scaling relationships, rather than forming arbitrary clusters driven by sampling noise. If such modular organisation exists, then metrics within a module should exhibit strong internal correlation relative to metrics assigned to different modules. The resulting clusters should be robust to resampling of the data. The identified modular structure should exceed expectations under null models of random association among metrics.
We quantified pairwise associations among the XX structural metrics using Pearson correlations computed across food webs. Because ecological metrics may covary either positively or negatively depending on scaling relationships, we constructed two alternative distance matrices: \(1−r\), which preserves the sign of correlations and distinguishes positive from negative association and \(1−∣r∣\), which groups variables based on the magnitude of their association regardless of sign. These two definitions allow us to test whether modular structure depends on directional relationships or simply on the strength of coupling among metrics. Hierarchical clustering was performed using average linkage on each distance matrix to identify candidate structural modules.
The optimal number of clusters was evaluated across a range of partition sizes (\(k = 2–10\)) using average silhouette width to assess within-cluster cohesion and between-cluster separation. To further evaluate cluster robustness, we implemented bootstrap resampling with 1,000 bootstrap replicates. This procedure estimates approximately unbiased (AU) p-values for each cluster, quantifying the probability that a cluster is supported under repeated resampling of the data. Clusters were considered statistically robust when they exhibited high silhouette support, and approximately unbiased bootstrap support ≥ 0.95. This dual criterion ensures that identified modules are both structurally coherent and stable to sampling variation.
2.3 Multivariate Structure and Module–Axis Alignment
To evaluate whether structural metrics of food webs organise into coherent multivariate modules and whether those modules define principal axes of variation, we conducted a principal component analysis (PCA) followed by a permutation-based test of module–axis alignment. Under this framework, principal components define dominant structural gradients in the metric space, whereas modules represent hypothesised mechanistic groupings of structurally related metrics. Demonstrating alignment between these two representations would suggest that modular decomposition captures fundamental axes of ecological variation.
Skewed count-based metrics (e.g., link counts, interval counts, and related quantities) were log-transformed using \(log(x + 1)\) to reduce right skew. All remaining variables were standardised to zero mean and unit variance prior to analysis to ensure that metrics measured on different scales contributed equally to the ordination.
To quantify how structural modules align with principal axes, we decomposed variance in each principal component according to module membership. For each module \(m\) and principal component \(k\), we computed:
\[ A_{mk} = \sum_{i \in m} l^{2}_{ik} \]
where \(l_{ik}\) is the loading of metric \(i\) on \(PC_k\) and \(A_{mk}\) represents the fraction of variance in \(PC_k\) attributable to metrics in module \(m\). Because squared loadings sum to one within each principal component, this provides a direct partitioning of PC variance across modules. This produces a module × PC matrix describing geometric alignment between modular structure and multivariate axes.
To evaluate whether observed module–PC alignment exceeded expectations under random module structure, we implemented a permutation test. We randomly permuted module labels among metrics 1 000 times while preserving the number and size distribution of modules, and the PCA loadings. For each permutation \(r\), we recomputed \(A^{(r)}_{mk}\) so as to form a null distribution for each module–PC pair. From this we computed p-values and corresponding z-scores, which was then used to infer significant alignment when \(p_{mk} < 0.05 \land |Z_{mk}| > 1.96\)
where:
\[ p_{mk} = P(A^{(r)}_{mk} > A_{mk}) \]
and:
\[ Z_{mk} = \frac{A_{mk} - \mu_{mk}}{\sigma_{mk}} \]
Additionally we also evaluated if modules collectively aligned with PCA structure beyond random expectation. Overall concentration of variance within module–PC space was calculated as \(T = \sum_{m, k} A^{2}_{mk}\) and significance was determined by comparing the observed (\(T\)) the permutation distribution (\(T^{(r)}\)). Where \(p_{global} = P(T^{(r)} \geq T)\). A significant result indicates that module structure explains multivariate variance better than random groupings.
2.4 Structural Metric Reduction and Representation of Network Space
To evaluate how network structure predicts ecological stability, we first reduced a high-dimensional set of topological metrics into alternative, biologically interpretable predictor sets. Because dimensionality reduction can fundamentally shape inference, we explicitly compared two conceptually distinct representations of network structure: (1) structural domain representatives and (2) latent axes of network variation.
2.4.1 Identification of Structural Domains
We previously grouped network metrics into clusters based on their pairwise correlations, such that each cluster represented a structurally coherent domain (e.g., trophic composition, centralisation, path structure). Clustering was performed using correlation-based distance, ensuring that metrics grouped together reflected shared structural information rather than raw scale. To construct a reduced predictor set from these domains, we selected a single representative metric per cluster using a medoid approach. Within each cluster, we computed the absolute correlation matrix among member metrics and defined distance as \(1−∣r∣\). The medoid was identified as the metric minimising mean pairwise distance to other metrics in the cluster. This approach preserves structural diversity while minimising redundancy and does not rely on principal component loadings, which can bias selection toward dominant variance axes. The resulting ‘cluster medoid’ set retained one interpretable metric per structural domain. Importantly, this procedure prioritises preservation of structural domains rather than overall variance magnitude, allowing low-variance but potentially mechanistically relevant features of network topology to be retained.
2.4.2 Extraction of Latent Structural Axes
As a complementary representation of network structure, we performed principal component analysis (PCA) on the scaled full metric matrix. Principal components were retained until cumulative explained variance exceeded 80%, yielding a set of orthogonal axes describing the dominant gradients of variation in network topology. These retained PC scores constitute a low-dimensional latent representation of network space. Unlike cluster medoids, which preserve structural categories, PC scores preserve maximal variance and ensure orthogonality among predictors. We also explored a third reduction approach—selecting the single highest-loading metric per retained PC. However, because principal components are linear combinations of multiple metrics, this PC-dominant metric strategy substantially reduces the variance represented by each axis. Consequently, this approach captures less of the total structural variance.
2.4.3 Conceptual Contrast Between Representations
The cluster-medoid and PCA-score representations reflect different theoretical assumptions about how structure relates to stability. The cluster-medoid approach assumes that ecological stability responds to discrete structural domains that may vary independently and need not align with the dominant gradients of network variation. In contrast, the PCA-score approach assumes that stability responds primarily to the major axes of variation in network topology, regardless of their interpretability in terms of structural domains. By explicitly comparing predictive performance across these representations, we test whether ecological stability is better explained by domain-specific structural mechanisms or by global geometric gradients of network organisation.
2.5 Structural Control of Stability
To evaluate how hierarchical representations of network structure explain variation in stability, we modelled four complementary stability outcomes: Robustness (\(R_{50}\)), Resilience, Spectral Radius (\(\rho\)), and Structural Controllability. These measures capture distinct dynamical processes governing ecosystem persistence, resistance, and controllability. Elastic net regression models were fit separately for each stability metric using three structural representations (cluster medoids, PC-dominant metrics, and full PCA scores) and XXX univariate, canonical proxies of stability (complexity, connectance, trophic coherence), allowing us to compare the structural determinants of different stability processes. All predictors and response variables were standardised prior to analysis to allow direct comparison of coefficient magnitudes.
Elastic net regularisation was chosen because it permits inference under correlated predictors by interpolating between ridge (distributed shrinkage) and lasso (sparse selection) penalties. The mixing parameter \(\alpha\) (0–1) determines the degree of sparsity, with lower \(\alpha\) values indicating distributed influence across many correlated structural descriptors (ridge-like behaviour), whereas higher \(\alpha\) values indicate that stability is driven by a smaller subset of structural predictors (lasso-like sparsity).
2.5.1 Model Evaluation
Predictive performance was assessed using repeated v-fold cross-validation (5 folds × 10 repeats) this improves the robustness of hyperparameter selection, particularly under correlated predictors and variable sample sizes. Repeating the partitioning reduces sensitivity to any single fold split and stabilises the selection of both \(\alpha\) and \(\lambda\), ensuring that the reported coefficients and variance contributions reflect consistent structural effects rather than idiosyncrasies of a particular data split. This approach enhances confidence that observed switching patterns in predictor importance across stability metrics represent genuine structural relationships rather than artefacts of sampling variability.. Within each training partition, models were fit across a grid of \(\alpha\) values (0–1 in increments of 0.25), and for each \(\alpha\), the optimal penalty strength (\(\lambda\)) was selected via internal cross-validation on the training data. Specifically, \(\lambda\) values were chosen from a dense path generated by the elastic net algorithm, and the \(\lambda\) that maximised cross-validated \(R^2\) within the training fold was recorded. This process ensured that \(\lambda\) selection was independent of the held-out test data. Model performance was quantified as the mean cross-validated \(R^2\) computed on held-out test partitions across all repeats and folds.
After cross-validation, final models were refit to the full dataset using the mean selected \(\alpha\), and the λ value closest to the cross-validated optimum was used to extract standardised regression coefficients. These coefficients characterise the direction and relative magnitude of structural effects on each stability metric.
The primary objective was interpretation rather than predictive optimisation, with the aim to quantify how structural organisation governs different stability processes. The elastic net mixing parameter \(\alpha\) provides insight into the architecture of structural control, indicating whether stability is influenced by many predictors (ridge-like) or a sparse subset (lasso-like). To quantify the relative contribution of different structural modules, standardised regression coefficients from the final models were squared and aggregated within modules. The proportion of variance explained by each module was then scaled by the model’s cross-validated \(R^2\), yielding a measure of the absolute variance in stability explained by each structural component.
3 Results
3.1 Structural Metrics Organise into Robust Modules
Hierarchical clustering of structural descriptors revealed clear modular organisation within food-web architecture Figure 2. Silhouette analysis indicated an optimal partition of the signed correlation matrix at k=7 modules, with a maximum average silhouette width of [insert value]. Clustering based on absolute correlations produced a slightly coarser but comparable solution (k=5), demonstrating that the modular structure is robust to the treatment of correlation sign and reflects strong underlying association patterns among metrics.
Bootstrap resampling further supported the stability of the major clusters. Most principal modules exhibited high approximately unbiased (AU) support values (AU ≥ [insert value]), indicating that the identified clusters are not artefacts of sampling variability but represent statistically robust groupings of structural descriptors.
The resulting modular partition identified coherent groups of metrics describing distinct aspects of network architecture Table 2. These included a large-scale architectural module integrating size and motif-based descriptors, modules capturing trophic organisation and energy routing, and two modules consisting of isolated structural descriptors that function independently from broader metric groupings. The six-module solution Table 2 reveals a hierarchical decomposition of food-web structural space into interpretable ecological dimensions. Collectively, these modules illustrate that network organization is shaped by multiple statistically supported axes of structural variation, reflecting size scaling, integration, energy flow, and control asymmetry.
| Module Name | Metrics Included | Ecological Interpretation |
|---|---|---|
| Macro Complexity | richness, links, l_S, S4, S5, intervals | Captures large-scale expansion of network structure and combinatorial complexity associated with increasing system size |
| Trophic Integration | connectance, diameter, intermediate, omnivory, cannibal, TL, ChLen, ChSD, S1, S2, loops, Clust, trophicCoherence, trophicVar | Describes how ordered the trophic structure is and how tightly species are embedded across trophic levels |
| Energy Transport | mean distance, path length | Quantifies geometric structure of energy routing independent of network size or trophic structure |
| Trophic Role Asymmetry | basal, predpreyRatio, GenSD, LinkSD, MaxSim | Captures structural imbalance across trophic roles, but it also incorporates heterogeneity in consumer behaviour and redundancy |
| Control Heterogeneity | top, VulSD, ChNum | Reflects heterogeneity in top-down regulation and predation pressure distribution |
| Flow Control | herbivory, centrality | Captures energy entry points and dominant structural hubs |
3.2 Structural Modules Align with Dominant Axes of Network Variation
Principal component analysis of standardised structural metrics revealed strong dimensional compression in food-web structure Figure 4. The first three principal components accounted for 68.5% of total variance (PC1 = 31.2%, PC2 = 21.5%, PC3 = 15.6%), with 80.7% explained by the first five components, indicating that a small number of dominant gradients capture most structural variation. Variance decomposition of principal components by module revealed a structured alignment between the modular organisation of metrics and the dominant axes of variation Figure 5. The first two principal components were jointly structured by Macro Complexity and Trophic Integration, but in contrasting ways indicating that variation primarily separates highly integrated trophic structure from size-driven architecture. PC3 represented a distinct structural dimension dominated by Control Heterogeneity, with additional contributions from Trophic Role Asymmetry, reflecting variation in the distribution of predation pressure and imbalances among trophic roles. This axis captures differences in how strongly interactions are concentrated among species and how unevenly trophic roles are distributed across the network.
Beyond the first three axes, additional modules aligned with more specialised components of structural variation. In particular, Energy Transport exhibited strong alignment with higher-order components (notably PC6), indicating that the geometry of trophic pathways varies largely independently of both network size and trophic integration. Basal Control showed weaker and more diffuse alignment across higher-order components, suggesting that variation in energy entry points and centrality represents a secondary structural dimension.
A global permutation test confirmed that the observed alignment between modules and principal components was significantly stronger than expected under random assignment (\(p = 0.006\)), demonstrating that the modular structure captures meaningful geometric organisation in the multivariate space of network metrics Figure 8. Alignment was concentrated in a subset of principal components, indicating that modules map onto specific structural gradients rather than uniformly spanning the full space.
Collectively, these results demonstrate that food-web structure is organised along a small number of dominant, ecologically interpretable axes. A primary structural gradient reflects the interplay between network size and trophic integration, while secondary axes capture variation in control heterogeneity, trophic asymmetry, and energy transport geometry. This hierarchical organisation indicates that food-web architecture is governed by multiple, partially independent dimensions of structural variation.
3.3 Predictors of stability
3.3.1 Dimensional Reduction Characterisation
| Representation | Dimensionality | Variance Preserved | Structure Type |
|---|---|---|---|
| Medoids | 7 | 23% | Domain sampling |
| PC-dominant metrics | 4 | 13% | Axis proxy sampling |
| PC scores | 5 | 80% | True latent space |
To evaluate how alternative structural representations differed in their information content and redundancy, we quantified variance retention, internal correlation structure, effective dimensionality, and geometric similarity among representations.
Variance Retention: The three representations differed substantially in the proportion of total structural variance preserved (Fig. S1A). The cluster-medoid representation retained ~23% of the total variance in the full metric space despite reducing dimensionality to seven predictors. In contrast, the PC-dominant metric set retained only ~13% of total variance, reflecting the fact that individual metrics capture only a fraction of each principal component’s multivariate structure. As expected by construction, the retained PC-score representation preserved approximately 80% of total variance. Thus, domain-based reduction preserved more distributed structural information than selecting one metric per PC axis, whereas the PC-score representation maximally preserved dominant gradients of network variation.
Internal Redundancy: Representations also differed in their internal correlation structure (Fig. S1B). PC scores were orthogonal by definition (mean |r| ≈ 0), indicating complete statistical independence among predictors. In contrast, cluster medoids exhibited moderate residual correlation (mean |r| = X), suggesting partial overlap among structural domains. PC-dominant metrics showed comparable (or higher/lower — insert result) redundancy relative to medoids. These differences indicate that the three approaches vary not only in information retention but also in predictor independence.
Effective Dimensionality: We next quantified effective dimensionality as the number of axes required to explain 80% of variance within each reduced predictor set (Fig. S1C). The PC-score representation required X axes, reflecting its design to capture dominant structural gradients. The cluster-medoid representation required X axes to reach the same threshold, indicating that despite containing seven predictors, structural variation was concentrated along fewer effective dimensions. The PC-dominant set exhibited the lowest effective dimensionality (X axes), consistent with its reduced variance retention. Together, these results show that dimensional compression differed across approaches not only in magnitude but in the distribution of variance across axes.
Overall, the three structural representations differed substantially in variance retention, redundancy, and effective dimensionality. Cluster medoids preserved moderate variance while maintaining domain interpretability, PC-dominant metrics retained minimal total variance, and PC scores preserved dominant structural gradients while ensuring predictor orthogonality. These differences establish that the representations encode distinct aspects of network topology, justifying their empirical comparison in predictive analyses of stability.
3.3.2 Predictive Performance Across Structural Representations
3.3.2.1 Structural Modules and Differential Explanatory Power
Across stability components, variance partitioning revealed that different structural modules contribute unevenly to explaining ecological stability (Figure 10). This indicates that stability is not governed by a single dominant structural regime, but instead emerges from module-specific control structures that map onto distinct dynamical processes. For Stability Potential, explanatory power was concentrated in modules associated with global structural organisation, consistent with its strong alignment with spectral structure. In contrast, Controllability drew heavily from modules capturing trophic heterogeneity (e.g., GenSD), indicating that control capacity is shaped by distributed variation in species-level connectivity. Recovery/Persistence showed strong contributions from modules associated with intermediate trophic structure and interaction pathways, suggesting that extinction cascades are governed by structural features that regulate how perturbations propagate through the network. By contrast, Resistance exhibited relatively weak and diffuse module contributions, reinforcing the finding that tolerance to species loss is not strongly encoded in any single structural domain. Together, these patterns demonstrate that different stability components are explained by distinct structural modules, highlighting a modular organisation of ecological stability in which each dynamical process is governed by a different subset of network features.
3.3.2.2 Role Reversals Across Stability Metrics
Beyond differences in explanatory power, we observed a striking pattern of sign switching in predictor effects across stability components, indicating that the same structural feature can have opposing influences depending on the stability process considered Figure 11. The most notable example is network complexity, which exhibits a clear shift in its relationship with stability. In the context of Resistance, complexity is negatively associated with stability (estimate = −0.410), suggesting that increased structural heterogeneity reduces tolerance to species loss. However, for Stability Potential, complexity becomes strongly negative (estimate = −0.598), indicating that greater spectral heterogeneity constrains global stability properties. Conversely, in Controllability, complexity shifts to a positive effect (estimate = 0.342), suggesting that increased structural heterogeneity enhances the system’s ability to be externally controlled and reconfigured. This pattern of sign reversal extends to other structural features, including trophic composition and connectivity measures. For example, variables such as intermediate species and top taxa proportions switch direction across stability metrics, reflecting their context-dependent roles in either stabilising or destabilising different dynamical processes. These sign changes indicate that structural features are not universally stabilising or destabilising, but instead play process-dependent roles. A given network property may enhance one aspect of stability while simultaneously undermining another, revealing inherent trade-offs in ecological organisation.
Taken together, these results support a view of ecological stability as a multi-process system structured along two axes: Modular organisation: Stability components are governed by different structural modules, reflecting distinct mechanisms of control, propagation, and resistance. Sign-dependent roles of structure: Structural features exhibit role reversals across stability metrics, indicating that their influence is conditional on the specific dynamical process being considered. This dual structure suggests that ecological networks do not optimise for a single notion of stability. Instead, they balance competing stability demands, where the same structural properties may contribute positively to one aspect of stability while constraining another. This pattern of sign switching suggests that ecological stability is governed by trade-offs rather than universal stabilising principles, with structural features occupying different roles across distinct dynamical regimes.
3.3.2.3 Revisiting Classical Complexity and the May Stability Framework
A central implication of these results is that neither classical complexity-based theory nor general spectral measures fully capture the structure–stability relationships observed in empirical food webs. In the framework of Robert May, stability is expected to decline with increasing species richness, interaction density, and interaction strength variance, leading to the well-known prediction that large, complex systems should be less locally stable. This prediction has provided a foundational baseline for ecological stability theory, but it implicitly assumes that complexity acts as a single aggregated driver of stability.
Our results challenge this assumption in two important ways. First, while our general complexity metric captures aspects of global structural heterogeneity, its explanatory power is uneven across stability components, and it fails to account for processes such as recovery and controllability. Second, and more importantly, we observe that the relationship between complexity and stability is not fixed, but instead changes sign across stability metrics. Complexity is negatively associated with resistance and stability potential, consistent with May’s prediction that increased structural complexity can destabilise local dynamics. However, it becomes positively associated with controllability, indicating that greater structural heterogeneity can enhance the system’s capacity for directed reorganisation and control.
These findings suggest that the classical complexity–stability relationship identified by May emerges as a projection of a higher-dimensional relationship between structure and stability. When stability is decomposed into multiple dynamical components, complexity no longer acts as a universal constraint, but instead participates in a set of process-specific trade-offs that govern different aspects of ecological dynamics.
Collectively, these results demonstrate that stability components occupy distinct positions within structural space Figure 12. ρ emerges from coordinated multivariate gradients, and robustness remains weakly constrained by static architecture. Rather than reflecting a single structural driver, ecosystem stability is organised along layered, facet-specific dimensions of network architecture.
4 Discussion
4.1 Network structure is multidimensional
For several decades, ecological theory has sought to understand how the structure of food webs influences ecological stability. Yet despite extensive empirical and theoretical work, consensus has remained elusive. Different studies frequently identify different structural predictors of stability, and results often appear contradictory. Connectance has been linked both positively and negatively to stability [REF], omnivory has been argued to stabilise or destabilise communities [REF], and measures of modularity, trophic organisation, or motif structure have each been proposed as key determinants of ecosystem dynamics [REF].
Our results suggest that much of this inconsistency arises from a conceptual simplification: network structure is rarely a single property. Instead, food webs exhibit strong modular organisation among commonly used network metrics Figure 2, with descriptors grouped into statistically supported modules spanning macro-architectural complexity, trophic integration, routing geometry, and trophic imbalance. Principal component analysis further shows that these modules align with a small number of dominant structural gradients Figure 4, Figure 7, Figure 5. Rather than forming a single continuum, food-web structure is organised along multiple partially independent axes.
This modular organisation implies that commonly used metrics are not interchangeable proxies for a single structural attribute, but rather partial projections of a higher-dimensional structural space. Consequently, many apparent inconsistencies in the structure–stability literature arise because different studies implicitly sample different structural dimensions.
4.2 Structural modules reveal dimensions of network architecture
The clustering analysis identified several coherent structural domains Table 2. A dominant module captured macro-architectural complexity (richness, link density, motif structure), reflecting how structural complexity scales with system size. A second module captured trophic integration, describing how tightly interactions are organised across trophic levels. Additional modules described energy transport geometry, trophic role asymmetry, and control heterogeneity.
These modules reveal that food-web architecture is not flat, but hierarchically structured. Some properties co-vary as networks expand, while others vary independently, reflecting distinct ecological mechanisms. Importantly, these domains provide a basis for understanding why different structural metrics are associated with different ecological outcomes.
4.3 Structural modules and mechanisms of stability
Rather than reflecting a single property, ecological stability emerges from multiple dynamical processes operating across scales. Our results show that different stability components are explained by different structural modules Figure 10, supporting a modular view of stability.
Robustness to species loss (resistance to extinction cascades) is only weakly constrained by any single structural domain, suggesting that tolerance to perturbation depends on distributed redundancy and context-dependent interactions. Stability potential, in contrast, is strongly governed by modules capturing global structural organisation, consistent with the role of large-scale architecture in constraining system-wide dynamics. Controllability is dominated by modules associated with trophic role asymmetry and heterogeneity, indicating that system steerability depends on uneven distributions of influence across species. Recovery and persistence are shaped by intermediate trophic structure and pathways, highlighting the importance of how perturbations propagate through interaction chains. Together, these results demonstrate that ecological stability is not a single property, but a collection of process-specific responses to perturbation, each governed by different structural dimensions.
4.4 Role reversals and sign switching across stability metrics
A central result of this study is that the same structural feature can have opposing effects depending on the stability metric considered Figure 11. This sign switching reveals that structural descriptors are not inherently stabilising or destabilising, but instead participate in process-dependent trade-offs.
The clearest example is network complexity. In line with classical expectations, increasing complexity is associated with reduced resistance and stability potential, consistent with the theoretical framework developed by Robert May, which predicts that increasing species richness and interaction density can destabilise systems (May 1972, 1974). However, complexity becomes positively associated with controllability, indicating that structurally diverse systems can be more amenable to directed intervention and reorganisation.
This pattern indicates that complexity is not a universal driver of instability, but a multifunctional structural property. Its effects depend on which aspect of stability is being considered. This resolves a long-standing tension between theoretical predictions and empirical observations: both are correct, but they refer to different stability components.
4.5 Revisiting May and classical complexity theory
These results suggest that classical complexity–stability theory represents only one projection of a more general structural–dynamical relationship. May’s framework captures local dynamical stability, where increasing complexity tends to reduce equilibrium stability. However, ecological stability is multidimensional, and once additional components are considered—such as robustness, recovery, and controllability—the relationship between structure and stability becomes more nuanced.
In this broader context, our general complexity metric does not act as a universal predictor of stability. Instead, its effects vary across stability components, and it participates in sign reversals similar to those observed for other structural descriptors. This indicates that the classical complexity–stability relationship emerges from aggregating multiple, potentially opposing, effects across different dimensions of stability.
More generally, these findings suggest that complexity is not a single axis of stability, but part of a higher-dimensional structural space in which different dimensions of complexity govern different dynamical outcomes. In this view, May’s result is recovered as a special case within a broader, multivariate framework.
4.6 Predictive structure: modules vs axes
Our predictive models further demonstrate that the relationship between structure and stability depends on how network structure is represented. Across stability components, different structural modules exhibit different explanatory power, indicating that stability is governed by module-specific mechanisms rather than a single dominant predictor.
At the same time, the choice of representation strongly affects which structural features appear important. Fine-scale (cluster-based) representations distribute explanatory power across many correlated features, while PCA-based representations concentrate influence along a few dominant axes. Importantly, these representations emphasise different aspects of structure: domain-based representations highlight mechanistic components, whereas PCA representations highlight global gradients.
This distinction is critical for interpreting structure–stability relationships. Apparent contradictions in the literature may arise not only from differences in ecological systems, but also from differences in how structure is represented and reduced.
4.7 Toward a modular, multidimensional theory of stability
Together, these results support a shift from univariate to multivariate thinking about ecological stability. Rather than asking whether complexity or connectance stabilises ecosystems, a more informative question is which structural dimension governs which component of stability. Here we show that different modules explain different stability processes, and that structural features change sign across stability metrics, revealing trade-offs rather than universal effects. This dual structure suggests that ecological networks do not optimise for a single notion of stability. Instead, they balance competing demands: resistance, recovery, persistence, and controllability. Structural features may enhance one aspect of stability while constraining another.
4.8 Implications for ecological inference
Our findings highlight the importance of how network structure is represented. Dimensional reduction choices—whether via modules, medoids, or PCA—shape which structural features are emphasised and therefore which ecological relationships are recovered. Recognising this allows for more principled comparisons across studies and reduces ambiguity in interpreting structure–stability relationships. More broadly, these results suggest that ecological networks should be understood as multivariate systems structured along a small number of coherent dimensions. Future work should focus on aligning structural representations with specific ecological processes, rather than relying on individual metrics as universal predictors.
5 Conclusion
The structure–stability relationship in ecology is not governed by a single rule, but by a set of process-specific interactions between structural dimensions and stability components. Classical results such as those of May emerge as special cases within this broader framework, which reveals that ecological stability is inherently multidimensional and governed by trade-offs across structural space.