Resumen
The Biodiversity-Ecosystem Functioning (BEF) and the Food Web Structure (FWS) theories are cornerstones of contemporary ecology. However, while several theoretical hypotheses predict a link between BEF and FWS, the integration of both frameworks has only recently been considered. In this study, we applied structural equation models to evaluate 73 sink food webs of predatory fish from the Paraná River, which encompass a wide gradient of community richness. Our analysis revealed a well-supported causal model where species richness drives food web structure, increasing link density, modularity, intermodular connections, and weak interactions, while decreasing nestedness. Both link density and modularity were positively associated with standing biomass, suggesting that communities with multiple energy pathways and strong complementarity effects tend to support higher biomass. Surprisingly, while species richness had the largest overall effect on biomass, this effect was indirect, mediated through three positive and one negative pathways. These findings highlight the complex associations between BEF and FWS, suggesting that anthropogenic impacts on modularity, such as community functional homogenization, could shift positive BEF effects to negative, with cascading consequences for entire ecosystems.
Métodos
The Paraná River belongs to the Rio de la Plata basin, which has the sixth largest discharge in the world and boasts around a thousand species of fish (Albert & Reis et al., 2011). Within its middle reach, the Paraná River features an extensive floodplain characterized by high diversity of microhabitats and elevated fish species richness (Marchetti et al. 2013; Scarabotti et al., 2017). Simultaneously, its pronounced thermal seasonality and hydrological variability deeply modify environmental heterogeneity, riverscape connectivity, community composition, and interaction networks (Scarabotti et al., 2017; Borzone Mas et al., 2022). This combination of biodiversity, spatial heterogeneity, and temporal variability positions the Middle Paraná River as an interesting model to test hypotheses about the relationship between diversity, food web structure, and ecosystem functioning (Figure 1). Fishes were sampled from 27 water bodies grouped into four habitat types according to their flow, size and connectivity: major rivers; secondary channels, connected lakes and isolated lakes (see details in Borzone Mas et al. 2022). Sampling surveys were carried out at one of four different hydro-climatic conditions that resulted from the combination of hydrometric level and temperature seasonality (high waters cold season, high waters warm season, low waters cold season and low waters warm season). Piscivorous fish were identified, measured, and weighed. The stomachs were preserved in 10% formalin and examined in the laboratory using a binocular microscope, with prey items identified at the species level. Community biomass was estimated as the total weight of captured individuals standardized by sampling effort (see supplemental table 1). Estimation of trophic interactions and construction of food webs Food webs were constructed based on the cooccurrence of predator and fish prey species, predator consumption patterns, predator biomass and prey biomass, metabolic losses, and energy efficiencies. These unipartite food webs represent the sink webs of the 15 most abundant piscivorous fishes in the middle Paraná River. Communities with fewer than three predators were excluded from the analysis. The predator consumption profiles were derived from the analysis of the stomach content of previous studies (Borzone Mas et al., 2022), with each profile represented as a vector indicating the biomass contribution of each species of prey to the predator’s diet. In each community, we built the topological food web using these consumption profiles, standardizing the preferences for available prey to a value of 1. Co-occurrence matrices were developed for each combination of site and sampling survey (27 sites in four surveys). Each fish prey species was considered a potential resource for consumers if at least one trophic interaction with a predator was observed in any survey. Biomass was standardized by sampling effort at each site (Supplementary Table 1). Interaction strength was estimated using an energetic food web approach following Barnes et al. (2018). Energy fluxes (F) between nodes were calculated by balancing the energetic demands of biomass stocks with energy outflow assuming system equilibrium, (Barnes et al. 2018): F= 1/e ×(X+L) where is the diet-specific assimilation efficiency, is the metabolic demands of the predator node, and is the energy loss. Since we only analyzed fish consumption, we assumed an energy efficiency of 0.906 (Gauzens et al. 2018). The mass specific energy loss of each individual was calculated with a mass specific metabolic scaling as L= 〖wt[kg]〗^(-0.25). Therefore, the energy loss of each node was the sum of the energy losses of the total biomass of the node (Barnes et al. 2018). The flow between predators and their prey was calculated using the fluxing function of the fluxweb package (Gauzens et al. 2018), using the standardized consumption profile of each community as the preference of the predators. The resulting estimated flow matrices were logarithmically transformed (log(F+1)) for subsequent analysis, hereafter referred to as food webs. Structural properties of food webs Species richness, quantitative weighted link density, nestedness, interaction strength, modularity, and the strength of intermodular connection were calculated for each food web (Banasek -Ritcher et al. 2009; Borzone Mas et al. 2022). Quantitative weighted link density was estimated as the “effective” number of prey and predators. This was computed by first applying Shannon´s index and then calculating the exponential of this value (following Banasek-Ritcher et al. 2009; see Supplementary table 1 for more details). Nestedness was calculated as Weighted Nestedness metric based on Overlap and Decreasing Fill (WNODF, Almeida-Neto et al. 2011), where values range from 0 (no nesting) to 100 (total nesting). The skewness in the interaction strength (hereafter weak interactions) examines the bias in the distribution, where positive values indicate a higher proportion of weak interactions while negative values correspond to a higher proportion of strong interactions (Hair et al. 2019, see Figure 2). The modularity coefficient (Qw) was calculated using the leading non-negative eigenvector of the modularity matrix of the graph with Newmann's algorithm (Newman 2006). The intermodular connector role was estimated from c-values (Guimera et al. 2010, Borzone Mas et al. 2022), which indicate how the species connections are distributed throughout the network, where the values range from 0 (all connections within their module) and 1 (perfect distribution among all modules in the network). Here, we use the 90th percentile of c-values as the maximum strength of intermodular connection, since it is less sensitive to extreme values than the maximum c-value. Since the WNODF and Qw values are sensitive to the size of the network, values were corrected by null model expectations through standardization with z-scores. In both cases, 2000 null communities were generated with fixed sums for columns and rows (quasi-swap count). See Supplementary Table 1 for the functions and packages used to compute the structural variables.