Spatial autocorrelation of local patch extinctions drives recovery dynamics in metacommunities
Abstract
Human activities put ecosystems under increasing pressure, often resulting in local extinctions. However, it is unclear how local extinctions affect regional processes, such as the distribution of diversity in space, especially if extinctions show spatial patterns, such as being clustered. Therefore, it is crucial to investigate extinctions and their consequences in a spatially explicit framework. Using highly controlled microcosm experiments and theoretical models, we ask here how the number and spatial autocorrelation of extinctions interactively affect metacommunity dynamics. We found that local patch extinctions increased local diversity (α-diversity) and inter-patch diversity (β-diversity) by delaying the exclusion of inferior competitors. Importantly, recolonization dynamics depended more strongly on the spatial distribution than on the number of patch extinctions: clustered local patch extinctions resulted in slower recovery, lower α-diversity and higher β-diversity. Our results highlight that the spatial distribution of perturbations should be taken into account when studying and managing spatially structured communities.
1. Introduction
Understanding the causes and consequences of local extinctions and how they affect biological systems at larger spatial scales lies at the heart of metapopulation and metacommunity ecology. Natural metapopulations and metacommunities—sets of local populations and communities linked by dispersal [1]—naturally experience local extinctions [2–4], for instance, due to demographic stochasticity, natural disasters or disease outbreaks. In addition, global changes—including climate change, habitat loss and fragmentation due to land-use changes, deforestation and urbanization—put increasing stress on ecological communities [5,6] which contributes to local patch extinctions.
Local patch extinctions, which we here define as the disappearance of all species from a patch, can have various consequences. In trophic systems, sustained local patch extinctions can induce regional species extinctions [7,8] and thus reduce regional diversity. Top predators are more likely to go extinct than intermediate and basal species. As a consequence, prey species can benefit at the regional scale from local patch extinctions due to release from predation. In competitive communities with a competition–colonization trade-off occasional local patch extinctions can even prevent regional extinctions and increase regional diversity by allowing less competitive species to persist [9].
Clearly, the relationship between local processes, such as extinctions, and regional patterns, such as distributions of biodiversity, is non-trivial. This is because metacommunities consist of more or less independent units, patches, harbouring local communities, which are linked in space by dispersal events. The coupling of spatially distinct communities can reduce the effect of local extinctions if these are asynchronous: extinct patches can be recolonized from occupied ones. However, dispersal between local communities can also have detrimental effects by synchronizing populations and thereby decreasing spatial insurance effects [10]. Under strong dispersal, the effects of local extinctions can even spread throughout a metacommunity such that local events have a regional effect [11,12]. For example, biomass could decrease in unperturbed patches as they receive reduced biomass fluxes from perturbed patches, or their species composition could change if they exchange individuals with perturbed patches whose species composition differs during the recolonization process. At the metacommunity level, strong dispersal might homogenize the composition of perturbed and unperturbed patches, thus reducing inter-patch diversity.
One probably important factor that modulates the effects discussed above is the spatial distribution of local patch extinctions, that is, whether they are clustered in space or not. An increase in the spatial autocorrelation of local extinction events could have a destabilizing effect at the metacommunity scale by coupling local dynamics and thus increasing global extinction risk [13,14]. Indeed, climate models have predicted an increase in the spatial and temporal autocorrelation of temperature [15], implying an increase in the environmental similarity between communities in space and time. This is expected to result in more climate extremes, such as heatwaves, droughts or frosts, affecting increasingly larger areas and for a longer time. Such climatic extremes can lead to local extinctions of populations of organisms sensitive to temperature changes, as seen in episodes of coral bleaching [16] or forest die-offs [17].
Despite these trends that foreshadow greater numbers and especially stronger spatial autocorrelation of climate-induced local extinctions, few studies have taken an appropriate, spatially explicit view of disturbances and their effects on metacommunity dynamics. This leaves a gap in our understanding of how spatially clustered extinctions may affect the dynamics of ecological systems.
Here, we investigate how the number and spatial distribution of local patch extinctions affect recolonization dynamics in metacommunities. We are particularly interested in the effects on local diversity (α-diversity) and inter-patch diversity (β-diversity), which act together to determine regional diversity. Previous theory on landscape moderation of biodiversity [18] led us to speculate that local extinctions could increase β-diversity if different species thrive in perturbed and unperturbed patches, and α-diversity if mass effect dynamics take place between perturbed and unperturbed patches with different species. Using a full factorial design crossing three levels of extinction numbers and two levels of spatial autocorrelation, we forced local patch extinctions in experimental and simulated metacommunities and followed metacommunity dynamics. We focused on the dynamics of the recolonization process (i.e. shortly after the extinctions) to capture the transient effects of extinctions. We were able to show that local patch extinctions can increase both α- and β-diversity, and that this effect depends strongly on the spatial autocorrelation of extinctions: dispersed extinctions increase α-diversity more, while clustered extinctions increase mainly β-diversity.
2. Methods
We used a combination of laboratory experiments with metacommunities of three freshwater ciliates (Tetrahymena thermophila, Colpidium sp. and Blepharisma sp.) in microcosm landscapes and mathematical modelling of metacommunities to address our main research question. To do so, we forced local patch extinctions (not sustained in time, i.e.‘pulse’ perturbations [19]) in experimental microcosm landscapes [20] and followed metacommunity recovery in terms of species diversity and biomass as a function of the intensity (amount of extinctions) and spatial distribution (clustered versus dispersed) of the extinctions. Experiments and simulations followed the dynamics of metacommunities in landscapes made of 16 patches arranged in a square lattice and connected by tubes allowing active dispersal. We present here a summary description of our methods. A more detailed account is provided in electronic supplementary material, S3.
(a) Experiments
Experimental landscapes were made of 16 vials (volume: 20 ml) arranged on a grid and connected to their four nearest neighbours, allowing individuals to disperse from one patch to another. Local patch extinctions consisted in removing all individuals of all species in a given patch. Each patch was initially inoculated with one of the three species at half its carrying capacity. Extinctions were implemented once, two weeks after inoculation to allow community assembly to have taken place. Subsequently, we observed the recovery of the landscapes (figure 1). Since we expected the extinctions to have only a transient effect before the metacommunity reached an equilibrium dominated by the best competitor (Blepharisma sp.), we followed the recovery dynamics just after the extinctions for a duration of two weeks (which is the time it takes for Blepharisma sp. to exclude the other species in a single patch co-culture; electronic supplementary material, figure S1.12h–j). In order to explore the effects of the number of local patch extinctions and their spatial autocorrelation on the dynamics of metacommunities, we used a full factorial design crossing three levels of local patch extinctions (0, 4 or 8 extinctions out of 16 patches) with two levels of spatial autocorrelation (clustered: electronic supplementary material, figure S1.8, landscapes 7–9 and 13-15; dispersed: electronic supplementary material, figure S1.8, landscapes 4–6 and 10–12). This design yielded a total of five treatments (no extinction, four clustered extinctions, four dispersed extinctions, eight clustered extinctions, eight dispersed extinctions) that were each replicated in three landscapes, for a total of 15 landscapes and 240 patches.
Figure 1. Overview of the experimental design: Species density over time of Blepharisma sp. (Ble, red), Colpidium sp. (Col, green) and Tetrahymena thermophila (Tet, blue), in experiments (a,c,e) and simulations under the ‘competition–colonization trade-off’ scenario (b,d,f). (a,b) Dynamics in patches from control landscapes. (c,d) Dynamics in perturbed patches from landscapes with four clustered extinctions. (e,f) Dynamics in perturbed patches from landscapes with four dispersed extinctions. Black arrows represent the extinctions. Note that two treatments (eight clustered extinctions and eight dispersed extinctions) are not shown here. (Online version in colour.)
We followed metacommunity dynamics through time by measuring the density of each species in each patch using automated video analysis. Three times per week, 2 ml of medium were sampled from all microcosms and replaced with fresh medium. For each microcosm, a sub-sample of 250 μl was placed between two microscope slides (height: 500 μm) and filmed using an optical stereo-microscope (Perfex Pro 10) coupled with a camera (Perfex SC38800) for 10 s (150 frames). We used the Bemovi R package (v. 1.0) [21] to extract individuals characteristics (shape, speed, size etc.) from the videos. We identified individuals from their characteristics using a random forest algorithm (R package randomForest v. 4.6-14) trained on videos of the monocultures filmed on the same day [22]. We rejected all individuals with an identification confidence (proportion of trees leading to that identification) lower than 0.8 as a compromise between the number of observations discarded and the confidence of identification (electronic supplementary material, figure S1.11).
α-diversity was measured as the inverse of Simpson’s index, which represents an effective number of species [23], and takes the relative abundance of species into account. We used the function beta.div.comp (R package adespatial v. 0.3-8, Ruzicka-based index) to compute the total β-diversity among the patches of a landscape [24].
All statistical analyses were conducted in R (v. 4.0.2). To test the relative effects of spatial autocorrelation and number of local extinctions on metacommunitiy properties, we studied four metrics (biomass, α-diversity, β-diversity and biomass recovery time) using mixed-effects models (R package lme4 v. 1.1-23) with the measurement point and landscape ID (for patch level metrics) as random effects to account for non-independence of measures taken the same day and measures taken within one landscape. Fixed effects were the autocorrelation of extinctions, the number of extinctions, as well as their interaction. We normalized the response variables using the R package bestNormalize (v. 1.6.1): we used the function bestNormalize (which finds the best transformation to render some data Gaussian while losing the fewest degrees of freedom using a Pearson P statistic) on each response variable (β-diversity: no normalization needed; α-diversity, biomass and recovery time normalized using the ordered quantiles technique, function orderNorm). The biomass in each patch was estimated using the bioarea per volume, a measure of the total surface of organisms visible in a video divided by the volume of medium in the camera field. The biomass recovery from extinction was estimated as the time needed to reach a biomass higher that the quantile of pre-extinction biomass in a given patch. This time span is hereafter referred to as recovery time.
For each statistical model, we performed AICc-based model selection on all models from the intercept to the full model. We used the weighted average of the model predictions for visualization.
Because we use measurements taken at different times, the temporal autocorrelation of the data acquired in a patch could lead us to artificially increase our statistical power and report non-significant results as significant. To assess the robustness of our analysis, we reproduced it using a total pooling [25] at the patch level by analysing the average (electronic supplementary material, figure S1.4) or the median (electronic supplementary material, figure S1.5) of post-extinction data from a given patch. Both approaches are very conservative and free from issues of autocorrelation, yet they closely reproduce the patterns observed in the main text (figure 2).
Figure 2. Observed response variables in the experiments (dots) and averaged mixed model predictions (medians and 95% confidence intervals; electronic supplementary material, table S2.3) from the extinction events to the end of the experiments. (a) Biomass in perturbed patches (blue: dispersed extinctions, red: clustered extinctions) and patches from landscapes with no extinctions (green), (b) biomass recovery time in perturbed patches, (c) α-diversity (measured as Simpson’s index) in perturbed patches and patches from landscapes with no extinctions and (d) β-diversity in all landscapes. (Online version in colour.)
(b) Metacommunity model
We also developed a metacommunity model to replicate and generalize the experiment in silico. We used a set of ordinary differential equations to describe the dynamics of metacommunities, where the terms describe the local dynamics (f), the emigration (g) and the immigration (h) of species i in patch k, with Ni,k as the density of species i in patch k:
The local dynamics are described by a competitive Lotka–Volterra equation where Ni,k grows logistically (ri: growth rate, αi,i: intraspecific competition) and is additionally impacted by interspecific competition (αi,j):
The number of individuals emigrating from a patch k is defined by a constant dispersal rate mi:
In analogy, we obtain the number of individuals immigrating into patch k as follows:
We used four different parameterizations (see electronic supplementary material, S3) to investigate which biological processes may explain the patterns observed experimentally, hereafter described as ‘scenarios of species interactions’ (electronic supplementary material, table S2.7). The scenarios empirical interactions and competition–colonization trade-off use growth rate and interaction parameters fitted from experimental data (electronic supplementary material, figures S1.12 and S1.27) with equal dispersal rates for all species in the former and dispersal rate inversely proportional to the competitive ability in the latter. The scenario randomized interactions used the same parameters as the ‘empirical interactions’ scenario but with randomized interspecific interaction rates in order to test whether our results held for other community structures. Lastly, we used a scenario without competition (no interspecific interactions) as a null model. We also used this model to conduct sensitivity analysis on the landscape size and dispersal rate in order to test the generality of our results.
The simulations were run in R (v. 4.0.2), using the ode45 solver from the library deSolve (v. 1.30). We simulated dynamics using the same extinction plans as in the microcosm experiments with 100 replicates for each treatment. While the simulations are deterministic, the initial distribution of species was drawn randomly in each replicate of each treatment, leading to variability between replicates at a given timestep. To check that our integrator choice did not result in numerical errors, we also reproduced 10% or the simulations with a second integrator (lsoda) with very low error tolerance parameters (rtol = 10−9, atol = 10−9). These simulations (electronic supplementary material, figures S1.6 and S1.7) match the main simulations (figures 2 and 4), ruling out integrator-related numerical errors.
3. Results
(a) The effect of the spatial distribution of extinctions
In the experiments, both local and regional effects of local patch extinctions were mainly determined by the spatial autocorrelation of extinctions. Except for β-diversity, the number of extinctions alone only had a marginal effect on the outcome of the experiment as indicated by model selection (figure 2; electronic supplementary material, tables S2.3 and S2.5). For the local variables studied (α-diversity, biomass and biomass recovery time), the autocorrelation of extinctions was found to be more important than the number of extinctions (electronic supplementary material, tables S2.3 and S2.5). Both α-diversity in unperturbed patches (electronic supplementary material, tables S2.4b and S2.6b) and β-diversity (electronic supplementary material, tables S2.3b and S2.5b) were mostly explained by the interaction between autocorrelation and number of extinctions (statistical models without the interactions had either a null (for β-diversity) or low (for α-diversity) weight).
Numerical simulations of our metacommunity model with the same spatial configuration and extinctions patterns confirmed this important effect of the spatial arrangement of extinctions compared to that of their number for all competition scenarios (figures 3 and 5).
Figure 3. Observed response variables in numerical simulations of the metacommunity model displaying different metrics after the extinction events (all biomass and diversity values of all perturbed patches between Textinction + 50 and Textinction + 150). (a) Biomass in perturbed patches (blue: dispersed extinctions, red: clustered extinctions) and patches from landscapes with no extinctions (green), (b) biomass recovery time in perturbed patches, (c) α-diversity (measured as Simpson’s index) in perturbed patches and patches from landscapes with no extinctions and (d) β-diversity in all landscapes. The top labels denote the scenarios of species interactions: ‘emp.’ for ‘empirical interactions’, ‘comp.-col.’ for ‘competition–colonization trade-off’, ‘rand.’ for ‘randomized interactions’ and ‘no int.’ for ‘no interspecific interactions’. Note that the effective number of replicates in the ‘no interspecific interactions’ is only one, as all species quickly reach their equilibrium density in all patches in the absence of competition, erasing any initial variability between patches and landscapes. (Online version in colour.)
(b) Direct effects: recolonization dynamics in perturbed patches
(i) Biomass
The biomass in a given patch after local patch extinctions was slightly higher in perturbed patches from landscapes with dispersed extinctions than in perturbed patches from landscapes with clustered extinctions (figure 2a; median predictions: 6076 μm2 ml−1 versus 4855 μm2 ml−1, 5–95% quantiles: 4274 − 7436 μm2 ml−1 versus 2224 − 6372 μm2 ml−1). Note that this effect is weak as indicated by model selection which ranks the intercept model second with an AICc weight of 0.27 (electronic supplementary material, table S2.3). The recovery time needed to reach a biomass higher than the 2.5% quantile of the pre-extinction biomass was shorter in case of dispersed extinctions compared to clustered extinctions, and it slightly increased with the number of extinctions (electronic supplementary material, tables S2.3 and S2.5, figure 2b; electronic supplementary material, S1.2; median (5–95% quantiles) mixed model predictions: four dispersed: 122 h (85–152), 8 dispersed: 129 h (94–164), four clustered: 139 h (100–172) and eight clustered: 152 h (122–185)).
In simulations of the metacommunity model, recovery times (figure 3b) qualitatively matched the experimental patterns in all scenarios. Quantitatively, the recovery times were much shorter (less than 100 time units) than what we found experimentally, probably because dispersal in the experiments happened over discrete time intervals (4 h periods, three times per week) resulting in a lag in recolonization dynamics.
In simulations with fitted interaction terms (‘empirical interactions’ and ‘competition–colonization trade-off’), the biomass of perturbed patches during the recolonization process was on average higher than the biomass of patches from control landscapes (figure 3a) because of the fast recolonization and higher carrying capacity of the less competitive species (T. thermophila and Copidium sp.) compared to the most competitive species (Blepharisma sp.). In the other simulations, the biomass during recolonization did not differ much between the perturbed patches and the patches from control landscapes (figure 3a) and was only slightly lower in perturbed patches.
(ii) α-diversity
In patches from control landscapes, α-diversity increased at first as species dispersed between patches but quickly fell to 1 as Blepharisma sp. finally excluded the two other species and dominated the community (electronic supplementary material, figure S1.1). In perturbed patches of the landscapes with extinction treatments, α-diversity was higher during the recolonization process in comparison to patches from control landscapes since the species were present in more even densities in the former (figure 2c; electronic supplementary material, S1.1). This effect was stronger for dispersed extinctions than for clustered extinctions (figure 2c).
In simulations from the metacommunity model, the empirical α-diversity pattern was best recaptured by the ‘empirical’ and ‘randomized’ scenarios (figure 3c), as well as transiently in the ‘competition–colonization trade-off’ scenario (electronic supplementary material, figure S1.23). In the absence of interspecific interactions, all species coexisted locally and the α-diversity was high in all patches.
(c) Indirect effects: spread of extinctions effects to unperturbed patches and at the regional scale
(i) Biomass
In both experiments and simulations, we observed no strong difference in biomass between treatments (figures 4a and 5a; electronic supplementary material, S1.28a–c).
Figure 4. Observed response variables in the experiments (dots) and averaged mixed model predictions (medians and 95% confidence intervals; electronic supplementary material, table S2.4) in unperturbed patches adjacent to at least one perturbed patch (blue: dispersed extinctions, red: clustered extinctions) and in control landscapes (green). (a) Biomass in unperturbed patches (for the two measurements following the extinctions). (b) α-diversity (measured as Simpson’s index) in unperturbed patches at the last two measurements. (Online version in colour.)
(ii) α-diversity
Experimentally, α-diversity was higher in unperturbed patches than in patches from control landscapes, particularly for dispersed extinctions (figure 4b; electronic supplementary material, S1.28d–f). Most of the variation between treatments was explained by the spatial autocorrelation of extinctions rather than the number of extinctions (electronic supplementary material, tables S2.4b and S2.6b). Interestingly, the effect of the number of extinctions depended on their spatial organization: under clustered extinctions, the α-diversity in unperturbed patches decreased with the number of extinctions but it increased under dispersed extinctions (figure 4b; electronic supplementary material, S1.28d–f). Note that this was not observed in simulations. This discrepancy could be due to either condition-dependent dispersal not accounted for in simulations, or to the low statistical power when it comes to indirect effect.
In simulations lacking interspecific competition, α-diversity levels were similar in unperturbed patches (across all treatments) and patches from control landscapes. In all simulations that included interspecific competition, α-diversity increased with both the number of extinctions and their spatial autocorrelation (figure 5b). Nevertheless, the effect sizes were variable: empirical interactions yielded effect sizes consistent with the experimental results (according to qualitative visual inspection), while randomized interactions yielded smaller effects and the ‘competition–colonization trade-off’ scenario yielded stronger effects.
Figure 5. Observed response variables in numerical simulations of the metacommunity model showing (a) biomass and (b) α-diversity (measured as Simpson’s index) in unperturbed patches adjacent to at least one perturbed patches (blue: Dispersed extinctions, red: clustered extinctions) and in control landscapes (green) after extinction events (all biomass and diversity values of all unperturbed patches adjacent to a perturbed patch between Textinction + 50 and Textinction + 150). The top labels denote the scenarios of species interactions: ‘emp.’ for ‘empirical interactions’, ‘comp.-col.’ for ‘competition–colonization trade-off’, ‘rand.’ for ‘randomized interactions’ and ‘no int.’ for ‘no interspecific interactions’. Note that the effective number of replicates in the ‘no interspecific interactions’ is only one, as all species quickly reach their equilibrium density in all patches in the absence of competition, erasing any initial variability between patches and landscapes. (Online version in colour.)
(iii) β-diversity
In control landscapes, β-diversity was fairly low because the patches ended up being homogeneous and dominated by Blepharisma sp. (electronic supplementary material, figure S1.1). β-diversity was higher in landscapes with extinctions than in control landscapes because of differences in species composition and density between perturbed and unperturbed patches (electronic supplementary material, figure S1.1). This effect was stronger for eight extinctions than for four extinctions, particularly for clustered extinctions (figure 2d).
In simulations of the metacommunity model, these results held qualitatively for all scenarios (figure 3d). These effects were strong and on par with experimental effect sizes for realistic interaction matrices (scenarios ‘empirical interactions’ and ‘competition–colonization trade-off’). They were weaker for randomized interaction matrices and negligible in the absence of interspecific interactions.
(d) Sensitivity to landscape size and dispersal rates
The simulations on larger landscapes (16 × 16 patches) yielded results (electronic supplementary material, figures S1.13 and S1.14) remarkably consistent with those discussed above. Our results were more sensitive to dispersal rates, but most qualitative patterns described for the ‘empirical interactions’ and ‘competition–colonization trade-off’ scenarios (e.g. stronger influence of the spatial autocorrelation than the number of extinctions, higher β-diversity for clustered extinctions, higher α-diversity spillover and faster biomass recovery for dispersed extinctions) were coherent for dispersal rates up to two times stronger/weaker than our standard simulations (electronic supplementary material, figures S1.15 to S1.22).
4. Discussion
(a) The role of the spatial distribution of the extinctions
We found that the spatial autocorrelation of extinctions had a stronger effect than the number of extinctions per se on all metrics measured, both in experiments and in simulations. Since our simulations suggest that this effect is independent of community structure, this result must be explained by the connectivity and distance between perturbed and unperturbed patches: if extinctions are dispersed, perturbed patches are closer and better connected to unperturbed patches than when extinctions are clustered (electronic supplementary material, table S2.2; figures S1.2 and S1.3).
The analysis of the simulations of large landscapes indicates that patch-level metrics in perturbed patches (recovery time and α-diversity; electronic supplementary material, figures S1.24 and S1.25) depend only on the distance to the closest unperturbed patch and not on the connectivity to unperturbed patches. This is coherent with our experimental results, where perturbed patches that were two links away from unperturbed patches had a longer recovery time (electronic supplementary material, figure S1.2, in red) and a lower α-diversity (electronic supplementary material, figure S1.3, in red) than the patches adjacent to at least one unperturbed patch. These findings also explain why the number of extinctions had a marginal effect in dispersed treatments compared to clustered treatments (figures 2 and 3): increasing the number of extinctions did not increase the distance from perturbed to unperturbed patches for dispersed extinctions (electronic supplementary material, table S2.2). On the contrary, more clustered extinctions resulted in larger clusters and thus in a greater distance from perturbed to unperturbed patches (electronic supplementary material, table S2.2).
(b) Direct effects of extinctions
(i) Biomass recovery
Experimental data and simulations support the conclusion that simultaneously increasing the number and autocorrelation of extinctions increases the time needed for a metacommunity to recover its pre-extinction biomass (figures 2b and 3b). These results were surprisingly consistent between the experiments and the various simulations scenarios, highlighting that this pattern does not depend on species interactions but rather on the geometry of the patches to be recolonized. A high number of spatially clustered extinctions increases the recovery time by creating large areas of perturbed patches, thus increasing the average distance and reducing the average connectivity between perturbed and unperturbed patches (electronic supplementary material, table S2.2). Clustered extinctions therefore result in what Zelnik et al. [12] have termed ‘rescue recovery regime’ where biomass recovery relies mainly on local population growth and is thus slower.
Additionally, both experimentally and in model simulations, perturbed patches had a slightly higher biomass after recovery than patches from unperturbed landscapes (figures 2a and 3a). This is because unperturbed patches mainly had the better competitor left (Blepharisma sp., electronic supplementary material, figure S1.1), while all three species persisted in perturbed patches. Since poorly competitive species (especially Colpidium sp.) reached a higher biomass than Blepharisma sp., perturbed patches had a higher biomass. This result should hold for communities dominated by highly competitive but slowly reproducing species that do not reach high densities (e.g. if there is a trade-off between population growth rate and competitive ability [26]) or when populations are able to overshoot their equilibrium density. This should however not be the case for communities where the dominant species happens to reach higher equilibrium densities, as it is the case in forests, for instance, where transiently recolonizing species (e.g. grasses or shrubs) do not accumulate biomass and are slowly replaced by dominant species that do (trees).
(ii) α-diversity
Local patch extinctions generally increased α-diversity as delayed competitive exclusion of inferior competitors. The persistence of less competitive species in perturbed patches during the recolonization process can be explained both by the decrease in population density and by a competition–colonization trade-off across the three species: the low population density after extinction events decreases the intensity of competition, while the competition–colonization trade-off delays the recolonization by Blepharisma sp., both processes resulting in the delay of competitive exclusion. These results are similar to the effect described in the intermediate disturbance hypothesis which predicts that some degree of perturbation should result in a higher local and regional biodiversity by reducing the abundance of competitively dominant species and allowing the persistence of early succesional species [27,28]. However, previous experiments on similar systems found that local patch extinctions decreased local diversity [9]. This can be explained by differences in metacommunity composition: metacommunities skewed towards early-succesional species should exhibit the α-diversity increase observed here, while metacommunities skewed towards late-succesional species (as in [9]) should see α-diversity decrease with local patch extinctions.
Clearly, these effects may be relevant in the context of ecosystem management: while local perturbations decrease biomass, they can also allow the persistence of species that would otherwise be excluded and lead to an increased local diversity.
(c) Indirect effects
Besides the direct effects discussed above, local patch extinctions may also have indirect effects at the regional scale by altering species densities and composition in unperturbed patches [11,12].
(i) α-diversity
Unperturbed patches in landscapes with extinctions had a higher α-diversity than unperturbed patches from control landscapes (figure 4b). This is because dispersal of less competitive species (T. thermophila and Colpidium sp.) from perturbed patches, where they were present in high density during the recolonization process, allowed persistence in both patches (electronic supplementary material, figure S1.1) with perturbed patches acting as sources and unperturbed patches as sinks. These source-sink dynamics correspond to the cross-habitat spillover hypothesis discussed by Tscharntke et al. [18]. The increase of α-diversity was stronger in unperturbed patches from dispersed extinction treatments, as these patches were connected to more perturbed patches and thus received an increased number of less competitive dispersers than unperturbed patches from clustered extinction treatments.
The patterns observed experimentally were recovered in all simulations that included interspecific competition (figure 5b), showing that local diversity maintenance by local extinctions is not restricted to our particular experimental community but can occur as long a some species excludes others.
It is worth noting that the increase in α-diversity was only observed in patches adjacent to perturbed patches, which could be described as an edge effect. This means that isolated extinction events do not have large scale effects in our setting, as perturbed patches only have an effect on their local neighbourhood. Indirect effects, however, can affect large proportions of the landscape if extinctions are numerous and spatially dispersed. Dispersed extinctions thus have both a stronger effect on unperturbed patches and affect a greater number of unperturbed patches.
(ii) β-diversity
β-diversity was higher in landscapes that experienced local patch extinctions in comparison to control landscapes, both in experiments and in simulations including interspecific competition (figures 2d and 3d). This can be explained by the fact that perturbed patches had a different species composition than unperturbed patches. In unperturbed patches, communities were mainly composed of Blepharisma sp., while perturbed patches allowed less competitive species to persist during the recolonization process. While we find a strictly increasing relationship between the number of extinctions and β-diversity (figures 2d and 3d), [9] found a unimodal relationship between β-diversity and local patch extinction number. While this seems contradictory, it is also possible that we did not use enough extinctions to reveal a unimodal relationship, as β-diversity could decrease when extinctions affect more patches.
(d) Perspectives
The strong effect of the spatial distribution of extinctions we report can be interpreted as differences in recovery regimes across spatial treatments: clustered extinctions, characterized by a weak connectivity between perturbed and unperturbed patches, result in what Zelnik et al. [12] described as a ‘rescue recovery regime’, while dispersed extinctions, characterized by a strong connectivity between perturbed and unperturbed patches, result in a ‘mixing recovery regime’. Under the ‘rescue’ regime, dispersal between perturbed and unperturbed patches is marginal compared to local dynamics. Perturbed and unperturbed patches are strongly differentiated, and the recovery dynamics mainly rely on local growth. Because of this strong differentiation, β-diversity was higher than in the ‘clustered extinctions’ treatment, but the high α-diversity of perturbed patches did not spill over much to unperturbed patches. Under the ‘mixing’ regime, dispersal between perturbed and unperturbed patches is on a par with local dynamics. Perturbed and unperturbed patches are well mixed, and both local growth and dispersal from perturbed patches participate substantially to the recovery. Because of the mixing between perturbed and unperturbed patches, α-diversity in the ‘dispersed extinctions’ treatment in unperturbed patches increased greatly (due to dispersal from perturbed patches), but β-diversity was lower than in the ‘clustered extinctions’ treatment.
Strictly speaking, our experimental settings, with discrete patches, homogeneous conditions and only three non-redundant species, may be thought to conform best to the patch dynamics paradigm [29], making extrapolations potentially difficult. However, as Thompson et al. [30] point out, metacommunity dynamics are more complex that what is captured by the four archetypes described by Leibold et al. [29]. Here, by looking at the transient recolonization dynamics, we were able to observe patterns consistent with both species sorting (good competitors are found mainly in unperturbed patches, good colonizers in perturbed patches), and mass effects (perturbed patches act as a source of less competitive species), highlighting that these mechanisms may often act simultaneously [31]. Our work also showcases the importance of transient dynamics in shaping biodiversity patterns, especially when we consider that local patch extinctions in nature should be recurring and asynchronous, leaving patches at different stages of recolonization and potentially enhancing metacommunity stability [32] and β-diversity. Moreover, the spatial treatment strongly influenced which patterns we observed during the recolonization: landscapes with clustered extinctions verged more on species sorting while landscapes with dispersed extinctions were more in line with the mass effects paradigm because the spatial autocorrelation of extinctions decreased the overall dispersal between perturbed and unperturbed patches. The spatial patterns of local perturbations can thus deeply alter the functioning of a metacommunity, here driving it from one metacommunity paradigm to another. This is particularly concerning when we consider that climate change could increase the spatial and temporal autocorrelation of climatic events [15], as observed in the metapopulation of Melitaea cinxia [13]. Our study thus warrants the inclusion of finer processes into existing metacommunity theory in order to better understand how the spatial structure of perturbations and the following transient dynamics affect the functioning of metacommunities.
5. Conclusion
Overall, our study shows that the effects of local patch extinctions in metacommunities strongly depend on the spatial distributions of extinctions. Local patch extinctions can increase both α-diversity and β-diversity by allowing weak competitors to persist in the metacommunity and by forcing a differentiation between perturbed and unperturbed patches.
Dispersal and connectivity between patches are central to recovery as they allow the recolonization of perturbed patches but also a mixing between perturbed and unperturbed patches, which can result in the spread of local extinction effects to unperturbed patches. In our setting, this spread was characterized by an increase in α-diversity in unperturbed patches through dispersal from species-rich, previously perturbed patches to species poor, unperturbed patches.
By determining the connectivity between perturbed and unperturbed patches, the spatial autocorrelation of extinctions modulates the dynamics after extinction events: when extinctions are clustered, perturbed and unperturbed patches are weakly connected. This results in a slower biomass recovery, a weak spread of α-diversity and high β-diversity as perturbed and unperturbed patches are differentiated. On the contrary, dispersed extinctions imply higher connectivity between perturbed and unperturbed patches which translates into a faster biomass recovery, a stronger spread of α-diversity and a lower β-diversity as perturbed and unperturbed patches are better mixed.
Our highly controlled experiment in combination with the theoretical model provide a proof-of-concept for the importance of taking into account the spatial distribution of disturbances in biodiversity research. Of course, applying our findings to specific, real-world ecosystems will require a combination of field data and system-specific models to obtain better estimates of the effects of local extinctions in more realistic settings. Nevertheless, our work highlights the relevance of the spatial distribution of local extinctions when doing so.
Data accessibility
Data and code are available on GitHub via Zenodo: https://doi.org/10.5281/zenodo.6364903. Electronic supplementary material is available online at rs.figshare.com [33].
Authors' contributions
C.S.: conceptualization, formal analysis, investigation, methodology, resources, software, writing—original draft; S.K.: conceptualization, methodology, supervision, writing—original draft; C.G.-B.: investigation, writing—review and editing; B.R.: resources, software, writing review and editing; E.A.F.: conceptualization, methodology, supervision, resources, software, writing—original draft.
All authors gave final approval for publication and agreed to be held accountable for the work performed therein.
Conflict of interest declaration
We declare we have no competing interests.
Funding
The study was funded by a grant of the ENS to C.S. and CNRS funds to E.A.F. B.R. acknowledges the support of iDiv funded by the German Research Foundation (DFG–FZT 118, 202548816). This is publication ISEM-2022-073 of the Institut des Sciences de l’Evolution–Montpellier.