# Long-term coastal macrobenthic Community Trajectory Analysis reveals habitat-dependent stability patterns

Chirine Toumi, Miquel De Cáceres, Jacques Grall, Aurélien Boyé, Éric Thiébaut, Marion Maguer, Vincent Le Garrec, Caroline Broudin, Céline Houbin, Olivier Gauthier

<figure><img src="image_1.png" /></figure>

To cite this version:

Chirine Toumi, Miquel De Cáceres, Jacques Grall, Aurélien Boyé, Éric Thiébaut, et al.. Long-term coastal macrobenthic Community Trajectory Analysis reveals habitat-dependent stability patterns. *Ecography*, 2023, 10.1111/ecog.06489. hal-04095678

HAL Id: hal-04095678

https://hal.univ-brest.fr/hal-04095678

Submitted on 13 May 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
ECOGRAPHY

Research article

# Long-term coastal macrobenthic Community Trajectory Analysis reveals habitat-dependent stability patterns

Chirine Toumi<fnref n="1" />✉, Miquel De Cáceres<fnref n="1" />✉, Jacques Grall<fnref n="1" />,<fnref n="3" />, Aurélien Boyé<fnref n="1" />✉, Éric Thiébaut<fnref n="5" />, Marion Maguer<fnref n="3" />, Vincent Le Garrec<fnref n="3" />, Caroline Broudin<fnref n="6" />, Céline Houbin<fnref n="6" /> and Olivier Gauthier<fnref n="1" />✉,<fnref n="3" />

**Correspondence:** Chirine Toumi (chirine.toumi@univ-brest.fr)

## Ecography

2023: e06489

doi: 10.1111/ecog.06489

Subject Editor: Julia Baum

Editor-in-Chief: Miguel Araújo

Accepted 15 February 2023

<figure><img src="image_3.png" /></figure>

Long-term monitoring programs are fundamental to detect changes in ecosystem health and understand ecological processes. In the current context of increasing anthropogenic threats on marine ecosystems, understanding the dynamics and response of communities becomes essential. We used data collected over 14 years in the REBENT benthic coastal invertebrates monitoring program, at a regional scale in the North-East Atlantic, covering a total of 26 sites and 979 taxa. Four distinct habitats were studied: two biogenic habitats associated with foundation species in the intertidal and subtidal zones and two bare sedimentary habitats in the same respective tidal zones. We used community trajectory analysis (CTA), a statistical approach that allows for quantitative measures and comparisons of temporal trajectories of ecosystems. We compared observed community trajectories to trajectories simulated under a non-directional null model in order to better understand the dynamics of the communities, their potential drivers, and the role of the studied habitats in these dynamics. Despite strong differences in the community compositions between sites and habitats, the communities followed non-directional dynamics during the 14 years monitored, which suggested stability at the regional scale. However, the shape, size, and direction of the trajectories of benthic communities were more similar within than among habitats, also suggesting the influence of the nature of the habitat on community dynamics. Results showed a higher variability in community composition the first years of the monitoring in the intertidal bare habitat and confirmed the role of biogenic habitats in maintaining temporal stability. They also highlighted the need to apprehend the role of transient and rare species and the scale of observation in temporal beta diversity analyses. Finally, our study confirmed the usefulness of CTA to link observed trajectory patterns to fundamental ecological processes.

**Keywords:** benthic communities, biogenic habitats, long-term monitoring, null model, temporal β diversity, trajectory analysis

<figure><img src="image_4.png" /></figure>

www.ecography.org

© 2023 The Authors. Ecography published by John Wiley & Sons Ltd on behalf of Nordic Society Oikos

This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.

<footnote n="1">¹LEMAR, Univ. Brest, CNRS, IRD, Ifremer, Plouzané, France</footnote>

<footnote n="2">²CREAF, Campus de Bellaterra, Univ. Autònoma de Barcelona, Cerdanyola del Vallès, Barcelona, Spain</footnote>

<footnote n="3">³OSU IUEM, Univ. Brest, CNRS, IRD, Plouzané, France</footnote>

<footnote n="4">⁴Ifremer, Centre de Bretagne, DYNECO, Laboratory of Coastal Benthic Ecology, Plouzané, France</footnote>

<footnote n="5">⁵Sorbonne Univ., CNRS, Station Biologique de Roscoff, Roscoff, France</footnote>

<footnote n="6">⁶CNRS, Station Biologique de Roscoff, Sorbonne Univ., Roscoff, France</footnote>
Introduction

A key challenge in community ecology is to monitor, detect, quantify and predict temporal changes in biodiversity (Dornelas et al. 2013, Buckley et al. 2021a). Indeed, biodiversity is a key driver of ecosystem functions (Gamfeldt et al. 2015, Duffy et al. 2017) providing essential services to society (Kremen 2005). Moreover, understanding the drivers and consequences of changes in biodiversity is necessary to set up management strategies (Palumbi et al. 2008). Long-term monitoring programs allow for detecting changes in ecosystem health, and understanding fundamental ecological processes, community dynamics and their responses to environmental constraints (Giron-Nava et al. 2017, Kominoski et al. 2018). Species abundances and community composition are examples of essential biodiversity variables that arise from monitoring programs and serve as indicators of ecosystem change (Pereira et al. 2013). Indeed, species richness variations are insufficient to capture changes in biodiversity (Dornelas et al. 2014, Hillebrand et al. 2018, Blowes et al. 2019). Nonetheless, temporal alpha diversity has received more attention than temporal beta diversity or the shift in the identity and/or abundance of named taxa in communities over time (Magurran et al. 2019).

The methods applied in temporal community ecology have grown over the past decades, from descriptive (e.g. bar graphs, ordinations) to more computationally complex methods (e.g. Moran eigenvector maps, machine learning methods) (Buckley et al. 2021b). Community trajectory analysis (CTA) is a multivariate method specifically tailored to study temporal community dynamics (De Cáceres et al. 2019). Starting from a classical ordination, it performs a geometric analysis of temporal trajectories that allows for identifying temporal patterns and variations in community dynamics (De Cáceres et al. 2019, Sturbois et al. 2021a). CTA allows for describing single trajectories by quantifying the changes between consecutive observations, the direction of these changes or the overall dynamics of the community. Additionally, it allows for comparing trajectories and apprehending the spatial variability of community dynamics and their underlying drivers (Legendre and De Cáceres 2013, De Cáceres et al. 2019). This is fundamental as community responses are diverse and not consistent across locations and scales (Hewitt and Thrush 2009, Blowes et al. 2019). Matthews et al. (2013) followed by Lamothe et al. (2019) proposed a theoretical framework linking temporal trajectory patterns and ecological processes. Following Van Meerbeek et al. (2021), a system is considered stable if it retains its reference conditions (state or dynamic) under changing conditions. In CTA, subsequent temporal observations very close to each other would imply stable communities that follow non-directional and gradual changes. Directional (i.e. trajectory following a particular direction) or saltatory changes (sudden and abrupt change between consecutive observations) would imply succession or regime shifts after a disturbance (Matthews et al. 2013, Lamothe et al. 2019).

In marine ecosystems, macrobenthic species are useful to measure ecosystem changes as they are not very mobile, show various life span and a large range of sensitivity to disturbance (Bessa et al. 2014, Dauvin et al. 2017). Changes in structure and composition of these communities are mostly the result of interactions between drivers acting at different temporal and spatial scales (e.g. Schückel and Kröncke 2013, Kröncke et al. 2019, Thrush et al. 2021). Interestingly, long-term stability with little changes has been observed in macrobenthic communities (Hinz et al. 2011), even in areas under continuous anthropogenic pressures (Bacouillard et al. 2020). Habitat-dependent factors could also drive heterogeneous responses. For example, intertidal communities – subject to repetitive physical stresses – may be more often reset than subtidal ones (Defeo and McLachlan 2005, Gray and Elliott 2009, Quillien et al. 2018, Boyé et al. 2019). Previous research has focused on communities’ response after a particular disturbance to study their recovery potential (e.g. Fromentin et al. 1997, Gilkinson et al. 2005). These studies highlighted that the removal of habitat-forming species and ecosystem engineers lowers community resilience (Cimon and Cusson 2018). Indeed, biogenic habitats with high levels of recruitment or connectivity promote resilience (O’Leary et al. 2017). Moreover, they increase community stability through the multiplicity of niches they create, promoting species richness, populations stability and asynchronous fluctuations across species (Lamy et al. 2020), but also through the attenuation of physical disturbances (e.g. thermal buffering) (Jurgens and Gaylord 2018). Other factors might influence the temporal dynamics of the communities. For example, marine ecosystems show more singletons and more transient species than terrestrial ones because of more open and less isolated communities (Raffaelli et al. 2005, Snell Taylor et al. 2018). Low occurrences and low abundances species (hereafter referred to as LOLAS) are frequent in marine communities, but it is mostly impossible to distinguish rare (always present but not always sampled because of sampling effort) from transient species (observed occasionally as a result of dispersal from adjacent habitats). LOLAS do not interact with their biotic and abiotic environment as do core species (Snell Taylor et al. 2018). As such, LOLAS may have a real ecological impact on community trajectories or they may simply make the signal noisier and hinder our ability to analyze trajectories.

Temporal trajectories have already been used in studies on changes in macrobenthic communities, but the interpretation of the shapes of the trajectories was mostly subjective. Indeed, the multidimensional space used to display and interpret trajectories has most often been the output from a non-metric multidimensional scaling (nMDS) analysis which is not suited for geometric and quantitative comparisons (e.g. Fromentin et al. 1997, Warwick et al. 2002, Beuchel et al. 2006). Here, we used CTA to study the dynamics of macrobenthic communities in 26 sites monitored for 14 years. Sites were located in biogenic and bare benthic habitats in the intertidal and subtidal zones of the coast of Brittany (France). A null model was used as a reference for non-directional
dynamics and analyses were conducted with and without LOLAS to better understand their influence on the dynamics of ecosystems. Given the hypotheses that biogenic habitats enhance stability and resilience of their associated communities and that intertidal and subtidal communities would show different dynamics, we investigated the following questions: 1) How did macrobenthic communities of biogenic and bare habitats in the intertidal and subtidal zones change over 14 years? 2) Were there any similarities in community dynamics of the different habitats at the regional scale? 3) How did removing LOLAS influence the inferred temporal dynamics of the communities? We hypothesized that 1) temporal trajectories of communities associated with biogenic habitats would show gradual changes and short trajectories or rapid return towards a stable point after a potential perturbation whereas communities from bare habitats would have more directional dynamics and/or saltatory changes following a potential perturbation, 2) temporal trajectories of communities in subtidal areas would show less variability (e.g. shorter trajectory or segment lengths) compared to the intertidal, 3) temporal trajectories would show similar features (e.g. in the direction of changes) at the regional scale as large-scale environmental changes may outweigh local ones, 4) LOLAS would induce larger stochasticity in trajectories through higher species turnover.

# Material and methods

## Study area

At the northwestern tip of France, Brittany is a biogeographic transition zone between the English Channel and the Bay of Biscay. It is a hotspot for macrobenthic fauna richness, characterized by a high diversity of benthic habitats (Gallon et al. 2017). Brittany harbors habitats associated with foundation species, the most common being intertidal seagrass beds *Zostera marina* and *Zostera noltei* and subtidal maerl (or rhodolith) beds principally *Lithothamnion corallioides* and *Phymatholithon calcareum*.

## Sampling

Benthic communities have been monitored yearly since 2003 along the coast of Brittany (France) within the REBENT program (www.rebent.org). We focused on four habitats: intertidal seagrass beds (only *Zostera marina* beds are monitored within the REBENT), intertidal sandy beaches, subtidal maerl beds and subtidal soft sediments (respectively referred to as intertidal biogenic habitat (IBIO), intertidal bare habitat (IBAR), subtidal biogenic habitat (SBIO) and subtidal bare habitat (SBAR) from this point forward).

At each site three faunal samples were taken at each of three fixed sampling points distributed 200 m apart (0.03 m² cores in the intertidal and 0.1 m² Smith-McIntyre grabs in the subtidal; see Boyé et al. 2019), except for Pierre Noire (8) where up to 10 grabs were taken at the sampling site.

Sampling was performed between the end of February and the beginning of May, before recruitment of most species occurs in the region (Dauvin et al. 2007, Boyé et al. 2019). In the laboratory, specimens were sorted, counted and identified to the lowest possible taxonomic level (usually species). Since the acquisition and identification of specimens were not carried out by the same people over the years of the monitoring program, we proceeded to a taxonomic homogenization: each recorded taxon was scrutinized by experts in benthic taxonomy and their names were checked thanks to the World Register of Marine Species (WoRMS Editorial Board 2021) to ensure a consistent taxonomic resolution.

In order to minimize the prevalence – and potential effect – of missing data as much as possible, we only selected sites with complete time series and with at least three core or grab samples in any particular year. Samples were pooled to estimate abundances at the site level. In the end, this led to a selection of 26 sites monitored from 2005 to 2018 while keeping a spatial resolution covering the coasts of Brittany and encompassing most of the environmental settings found in this region (Boyé et al. 2017, 2019). Of these 26 sites, 5 were in IBAR, 7 in IBIO, 8 in SBAR and 6 in SBIO (Fig. 1). We conducted the analyses at the habitat level or within a same tidal level because it may not be relevant to run analyses including two different tidal levels since sampling gears differ between intertidal and subtidal sites.

Hereafter, the term ‘site’ refers to a given location in a given habitat and the term ‘observation’ refers to a sampling occasion at a given site in a given year (here there were 364 observations in total over the 26 sites and 14 years monitored).

# Numerical analyses

## Community trajectory analysis (CTA)

CTA is based on the dissimilarity between pairs of community observations in space and time (De Cáceres et al. 2019). We defined the multivariate space of resemblance between community observations using species abundances and the Hellinger distance coefficient, which is equivalent to the Euclidean distance computed on the square root of species relative abundances. It does not give an excessive weight to rare species and has the advantage of fulfilling the metric and Euclidean properties (Legendre and Gallagher 2001, Legendre and De Cáceres 2013). We used a principal coordinates analysis (PCoA) to display trajectories and compute CTA metrics without any distortions, for each habitat or each tidal zone separately. Each observation is represented in the multivariate space by its coordinates. Two consecutive observations of a same site are linked by a segment, the ensemble of the segments of a site representing its temporal trajectory. Geometrical metrics computed on site trajectories describe their ecological dynamics (see De Cáceres et al. 2019 and Sturbois et al. 2021a for the detailed formulas). Trajectories are usually represented on the first two axes of the ordination, but CTA can compute these metrics using all dimensions. Here we computed the following metrics for the whole multidimensional space:
<figure><img src="image_5.png" /><figcaption>Figure 1. Map of the monitored sites in the four distinct habitats along the coasts of Brittany (France) (sources: OpenStreetMap, European Environment Agency).</figcaption></figure>

1) *segment length* or the distance between two consecutive observations. With the Hellinger distance, the maximum value between two observations is $\sqrt{2} \approx 1.41$,

2) *trajectory length* or the total path length of the trajectory which is the sum of all the segment length of the trajectories. With the Hellinger distance, the maximum value of the trajectory length is $(n-1) \times \sqrt{2}$, where n is the number of observations in a site. Here 14 observations give trajectories of 13 segments and a maximum trajectory length of $13 \times \sqrt{2} \approx 18.38$,

3) *net change* or the distance between the first and the last observations. Here it is the distance between the observations made in 2005 and 2018 at a given site. With the Hellinger distance, the maximum value remains $\sqrt{2} \approx 1.41$,

4) *angle θ* between two consecutive segments takes values from 0° to 180°: 0° indicates three observations completely aligned with no change in direction, whereas 180° indicates two vectors with the same orientation but opposite directions,

5) *overall directionality* of the trajectory considers the angles between consecutive segments and their lengths. It reflects the degree to which the community is consistently moving in a particular direction. Directionality takes values

between 0 and 1, with 1 representing a trajectory following a completely directional pathway,

6) *resemblance* between a pair of trajectories is assessed using the symmetrized directed segment path dissimilarity ($D_{SDSP}$) (Besse et al. 2016, De Cáceres et al. 2019) which takes into account shape, size, direction and position of the trajectories. Because positions of the trajectories are highly influenced by species compositions at each site, we centered the trajectories prior to the calculation of trajectory distance to focus on compositional dynamics rather than spatial variation of species composition. After computing $D_{SDSP}$ between all centered trajectories of sites in each habitat, we assessed the overall variation of community dynamics in each habitat using the total dynamic beta diversity ($\text{dBD}_{\text{tot}}$) (Legendre and De Cáceres 2013, De Cáceres et al. 2019).

Analyses were repeated focusing on the core species of the communities of each of the four habitats by removing the LOLAS. Core species are defined as the most abundant and persistent ones whereas LOLAS occur more infrequently and often have low abundance (Magurran and Henderson 2003). To achieve a compromise between the maximum abundance
of each taxon and the number of observations in which they were detected, we set an arbitrary occurrence threshold (corresponding to 1/5 of the maximum occurrence in each habitat) under which species were classified as LOLAS and removed (Supporting information). This led to the removal of 50–75% of the species of each site in each habitat, as LOLAS were preponderant in the studied communities (Supporting information). However, removing LOLAS resulted in datasets containing from 79% to 96% of original total abundance.

We also tested the influence of habitat type (i.e. bare versus biogenic) on the dynamics of the communities by performing a PERMANOVA (Anderson 2017) on the symmetric matrix of distances between trajectories (i.e. with D<sub>SDSP</sub> values computed on centered trajectories of sites). For this, trajectories were first recomputed within tidal levels.

Finally, we were interested in comparing directionality values of sub-trajectories in IBAR especially because we identified two periods of time which seemed to show different dynamics. Sub-trajectories are obtained by splitting the trajectories into two or more sub-trajectories to compare time periods within the overall trajectory. The coordinates of the observations in the multidimensional space are the same as for the whole trajectories but the directionality metric is computed on each sub-trajectory independently.

## Simulations

We simulated communities under a pure non-directional dynamic to create a null model allowing for comparison with observed CTA metrics. Communities were simulated for each habitat separately while preserving the following habitat-specific properties:

1) Species pool composed of all species recorded in the habitat during the 14 years,

2) The empirical distribution of species abundances, based on species abundances recorded for each observation in the habitat,

3) The empirical distribution of species occurrences, based on the number of observations each species was recorded in the habitat,

4) The empirical distribution of total abundances, based on the total abundance recorded for each observation in the habitat,

5) A fixed mean extinction rate, i.e. the average proportion of species lost between consecutive years (0.46 for IBAR, 0.41 for IBIO, 0.38 for SBAR, 0.39 for SBIO).

The simulation procedure is detailed in the Supporting information and was designed to break the temporal dependency between consecutive observations except for the identity of the species kept. CTA metrics were computed on 100 trajectories of simulated communities in each habitat, with a 14 year dynamic as the observed communities. CTA metrics used to compare simulated and observed community dynamics were the directionality of trajectories, the net change, the total length of trajectories and the length of the consecutive segments. We used two-sided Kolmogorov–Smirnov tests to compare the simulated and observed distributions of the different metrics.

All analyses and simulations were conducted with the R programming language ver. 4.1.2 (www.r-project.org) and packages 'ecotraj' (De Cáceres et al. 2019, Sturbois et al. 2021a) and 'adespatial' (Dray et al. 2021).

## Results

### Regional scale non-directionality

Community trajectories on the first two axes of the PCoA represented from 35% (SBAR) to 49% (IBAR) of the total variance (Fig. 2). Site trajectories occupied clearly different positions, reflecting site-specific composition and structure within each habitat that were also reflected by high habitat-wise total beta diversity which represents the variance of the community matrix and takes a maximum value of 1 (Legendre and De Cáceres 2013): $BD_{tot}=0.70$ for IBAR, $BD_{tot}=0.55$ for IBIO, $BD_{tot}=0.60$ for SBAR and $BD_{tot}=0.59$ for SBIO. In addition to their different composition, the four habitats had different species richness: a total of 299 taxa was recorded over the 14 years in IBAR, 493 in IBIO, 527 in SBAR and 665 in SBIO. Overall, directionalities of sites trajectories were very similar even between habitats, and ranged from 0.34 to 0.39. These low values indicated a weak directionality in all monitored sites as did mean angles higher than 90° (Table 1). Trajectories rotated on themselves or oscillated around a point and can therefore be qualified as non-directional.

### Habitat-dependent dynamics

Mean trajectory length in IBIO was somewhat shorter than that of other habitats but all fell within the same range (Table 1). These shorter trajectories mean that community composition in IBIO was more similar and less variable from year to year. IBIO and SBIO had the smallest net changes across the survey period and thus a high similarity between their states in 2005 and 2018. On the contrary, IBAR trajectories exhibited the highest net changes, with high variability in community structure and composition between the first and last survey. Interestingly, IBAR also possessed the highest $dB_{tot}$ meaning that there was more variability between the dynamics of IBAR sites compared to the other habitats.

PERMANOVA computed on trajectories dissimilarity ($D_{SDSP}$) between centered trajectories for the subtidal and the intertidal zones separately showed that 10% of the variance of trajectories shapes can be attributed to habitat type (intertidal: $F=1.115, R^2=0.10, p=0.046$; subtidal: $F=1.228, R^2=0.093, p=0.001$), lending support to a habitat-dependent influence on the temporal dynamics of communities.

As expected, removing LOLAS and focusing on core species reduced $dB_{tot}$, net change and trajectory length (Table 1). However, directionality remained low (Supporting information) and angles were > 90° (Table 1), with sites still tending to return to their previous state. IBAR still presented the highest values of $dB_{tot}$, trajectory length and net change
<figure><img src="image_6.png" /><figcaption>Figure 2. Representation of community trajectories on the first two PCoA axes with their associated directionality values for each site in the four habitats monitored from 2005 to 2018. One point represents the community state of a site in a given year (one observation). Site specific consecutive community states are linked by a segment and taken together depicts the site trajectory. The arrow represents the final community state of a trajectory (i.e. the community state of a given site in 2018 here).</figcaption></figure>

while IBIO allowed for similar or higher levels of stability
than in the subtidal with or without LOLAS.

Stability despite variability

Distributions of directionality values of the simulated
and observed trajectories were not significantly differ-
ent (Kolmogorov–Smirnov test; D < 0.37, p > 0.05 in all
habitats), with directionality values oscillating around 0.36
in observed and simulated communities: both had a weak

directionality (Fig. 3). Overall, trajectory lengths, net changes (Fig. 3) and segment lengths (Fig. 4) of simulated communities were higher, with the exception of net change in IBAR (Fig. 3). Simulated distributions of net change and trajectory lengths were more symmetric and less spread out than observed ones (Fig. 3). In IBAR, simulated and observed distributions of CTA metrics had more similar distributions (Kolmogorov–Smirnov test; net change: D=0.48, p > 0.05; trajectory length: D=0.7, p < 0.01) than in the three other habitats (Kolmogorov–Smirnov test; D > 0.87, p < 0.001),
<table><caption>Table 1. Mean and standard error values of sites' trajectory metrics in the four habitats monitored from 2005 to 2018, considering the community with and without LOLAS. dBDtot = the total dynamic beta diversity computed on the dissimilarity D<sub>SDSP</sub> between centered trajectories, length = the total path length of the trajectory, net change = the distance between the starting point and the final point of the trajectory, angle = the mean angle between two consecutive segments of the trajectory.</caption><thead><tr><th rowspan="2">Habitat</th><th colspan="4">Global community</th><th colspan="4">Without LOLAS</th></tr><tr><th>dBD<sub>tot</sub></th><th>Mean length ± SE</th><th>Mean angle ± SE</th><th>Mean net ± SE</th><th>dBD<sub>tot</sub></th><th>Mean length ± SE</th><th>Mean angle ± SE</th><th>Mean net ± SE</th></tr></thead><tbody><tr><td>Intertidal bare</td><td>0.21</td><td>8.61 ± 0.42</td><td>114.83 ± 1.90</td><td>0.94 ± 0.06</td><td>0.18</td><td>7.84 ± 0.75</td><td>114.16 ± 2.30</td><td>0.84 ± 0.07</td></tr><tr><td>Intertidal biogenic</td><td>0.15</td><td>7.68 ± 0.51</td><td>115.88 ± 0.96</td><td>0.76 ± 0.03</td><td>0.09</td><td>5.63 ± 0.50</td><td>114.21 ± 1.63</td><td>0.58 ± 0.02</td></tr><tr><td>Subtidal bare</td><td>0.18</td><td>8.34 ± 0.41</td><td>117.73 ± 1.43</td><td>0.82 ± 0.04</td><td>0.11</td><td>6.72 ± 0.35</td><td>117.44 ± 1.92</td><td>0.56 ± 0.03</td></tr><tr><td>Subtidal biogenic</td><td>0.19</td><td>8.31 ± 0.91</td><td>117.23 ± 1.35</td><td>0.74 ± 0.07</td><td>0.10</td><td>6.22 ± 0.60</td><td>118.3 ± 1.88</td><td>0.58 ± 0.07</td></tr></tbody></table>

<figure><img src="image_7.png" /><figcaption>Figure 3. Comparison of three CTA metrics (directionality, trajectory length and net change) computed on the observed communities of the four habitats monitored from 2005 to 2018 and the communities simulated under the non-directional null model. Differences between observed and simulated distributions were tested by a two-sided Kolmogorov–Smirnov test (significance code: absence of code = p > 0.05, * = p < 0.05, ** = p < 0.01, *** = p < 0.001).</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>Figure 4. Comparison of trajectory segment lengths of the communities observed in the four habitats monitored from 2005 to 2018 and the simulated ones.</figcaption></figure>

confirming that IBAR trajectories were the most temporally variable among the studied habitats.

As expected, simulated segment lengths showed no partic-
ular trend and simulated segment length distributions were
very similar from year to year regardless of habitat (Fig. 4).
Simulated trajectories presented longer segments than obser-
vations, reflecting higher variability in community composi-
tions between consecutive years. Some differences appeared
between habitats in simulations: for example, the simulated
segment lengths in IBAR were more variable than in SBAR.
As simulations were based on habitat-wise observed distribu-
tions, this might reflect the observed heterogeneity of seg-
ment lengths in IBAR. The stability in the simulated segment
lengths diverged from the observed trajectories, especially

for IBAR. Indeed, in IBAR the first four observed segments
(2005–2009) were longer than the following 9 (2009–2018),
revealing a more stable dynamic of the communities these
last nine years (Fig. 2, 4). However, we did not find a cor-
responding change in directionality for these sub-trajectories
(2005–2009: 0.35; 2009–2018: 0.34).

Discussion

Non-directional regional dynamic

We used CTA to describe and understand 14 years of ecologi-
cal dynamics of macrobenthic communities located in four
distinct coastal habitats. Results suggested a global stability of the system at the regional scale with non-directional temporal trajectories. The low directionality observed is perhaps due to the nature of the considered systems. CTA was firstly applied on plant communities (De Cáceres et al. 2019), with more persistent and fewer species than expected in marine systems that have higher diversity and higher temporal turnover, and present less ordered ecological succession (Raffaelli et al. 2005, Henderson and Magurran 2014). The low directionality could also be due to the frequency of sampling in our study. It is likely that some ecological processes in the benthic macrofauna communities monitored are faster than our sampling frequency (i.e. seasonal variations). Furthermore, macrobenthic communities can exhibit different multi-year cycles (e.g. 3–9 years, 5–7 years, 10–12 years; Thrush et al. 2021) that could not be detected within a 14 year time series.

Long-term stability of macrobenthic communities has already been demonstrated in the English Channel (e.g. Fromentin et al. 1997, Bacouillard et al. 2020). Recently, CTA was used on intertidal communities monitored in a single bay, and results also demonstrated stability at the scale of the study area with changes mainly reflected by random population dynamics of structuring prevalent species under strong natural pressures (Sturbois et al. 2021b). CTA conducted on phytoplanktonic communities of the Eastern English Channel also demonstrated an overall stability in community composition (Lefran et al. 2021). Brittany is a biogeographic transition zone with a high diversity of benthic habitats (Gallon et al. 2017) and its geographical location makes it a very open system. This might enable large-scale transport of propagules from outside source patches (Ayata et al. 2010) and multiple sources of colonists thanks to the multiplicity of habitats. These factors could enhance recovery or persistence of communities in the region (Ellis et al. 2000). Moreover, Brittany is a hotspot of biodiversity (Gallon et al. 2017) and this could also have a stabilizing effect on communities (Downing et al. 2014, Craven et al. 2018), especially if species show asynchronous responses to environmental fluctuations and differences in the speed at which they respond to perturbations (Loreau and de Mazancourt 2013). However, the diversity-stability debate is still ongoing in ecology (Ives and Carpenter 2007, Kéfi et al. 2019).

Trajectory length and net change mostly had higher values in simulated communities compared to observed ones. Simulations smoothed out site-specific effects on trajectories: for instance, in observed communities, Trévignon (20) in SBIO expressed the lowest net change and trajectory length (Supporting information), since this site was largely dominated by a single and temporally persistent taxon *Porcellanidae* spp. (Supporting information). Moreover, temporal autocorrelation may be lower in simulated compared to observed communities: temporal turnover should be higher in simulated communities given the simulation procedure with abundances resampled at each step. Indeed, the within-site year-to-year BDtot was higher in simulated communities (Supporting information). However, the proportion of replacement and difference of richness (see Legendre and

De Cáceres 2013, Legendre 2014) were fairly homogeneous between simulated and observed communities in IBIO and SBAR but not in IBAR and SBIO (Supporting information). With habitat-scale constraints, the simulation did not only remove autocorrelation, but also broke co-occurrence patterns and site-specific constraints, making it a regional-scale neutral model. This also suggests that biotic or abiotic factors (which were not considered in the simulations), such as characteristics of the habitat itself, could act as filters preventing observed communities from adopting the dynamic produced in the simulation (especially the higher temporal turnover), and thus enhance stability at the regional scale.

## Habitat-dependent dynamics

### Biogenic versus bare

Despite non-directional regional scale dynamics, results revealed habitat-dependent community dynamics. First, within tidal zone differences exist between biogenic and bare habitats. In the intertidal, all CTA metrics were lower in IBIO compared to IBAR, translating a higher temporal stability in this former habitat. Hily and Bouteille (1999) demonstrated that the development of eelgrass meadows enhances abundances and biomasses and increases specific and functional diversity compared to IBAR. Indeed, the shelter created by seagrass patches can improve the recruitment, survival and diversity of species (Boström and Bonsdorff 2000). Moreover, their physical structures enhance habitat heterogeneity and complexity, increasing shelters and food resources availability and promoting biodiversity (Thomaz and Cunha 2010), hence stability given the mechanisms suggested by Loreau and de Mazancourt (2013). Biogenic habitats can also act as environmental buffers: in an environment facing harsh conditions, such as the intertidal zone, habitat-forming species can sustain biodiversity and ecosystem functioning through the reduction of physical stress (Bulleri et al. 2018), for example by altering local hydrodynamic conditions and sediment dynamics or reducing thermal stresses for inhabitant taxa thanks to their physical structures (Peterson et al. 2004, Bouma et al. 2009, Jurgens and Gaylord 2018).

In the subtidal zone, we expected more stable trajectory patterns in SBIO than SBAR. Indeed, the high architectural complexity of maerl beds may reduce predation stress in communities (Bouma et al. 2009), increase specific (Grall and Glémarec 1997, Grall and Hall-Spencer 2003) and functional diversity and redundancy (Boyé et al. 2019) which allow ecosystems to resist disturbances because multiple species could take on critical roles (Palumbi et al. 2008). SBIO communities had non-directional dynamics and the lowest average net change but not the lowest average trajectory length. Long trajectories returning toward their initial ecological state could make us lean toward a resilience hypothesis. However, here the low net change was strongly affected by an atypical site: Trévignon (20), the single site with the lowest net change and segment lengths, which presented a nearly persistent community (Gray and Elliott 2009, Thrush et al. 2021) (Supporting information). Contrary to our expectations, removing this
site from the analyses resulted in increased average net change and mean trajectory length in SBIO, the latter becoming the highest of all habitats (Supporting information). However, this did not radically change the observed patterns of segment lengths in SBIO (Supporting information). We hypothesized that the higher species richness, and hence the higher number of LOLAS, increased temporal turnover in SBIO because each LOLAS was only present over a small fraction of the time series (Magurran and Henderson 2010, Snell Taylor et al. 2018). Removing LOLAS decreased net change and trajectory length in all habitats. However, a high species diversity with numerous LOLAS could support the insurance hypothesis, where LOLAS can be more resistant or well suited to environmental change (Hewitt et al. 2016, Thrush et al. 2021) with asynchronous responses of species to fluctuations (Yachi and Loreau 1999). Core species could maintain local stability because they are suited to existing environmental conditions whereas LOLAS could maintain regional and long-term stability by replacing core species following an environmental change (Coyle et al. 2013, Henderson and Magurran 2014, Vermeij and Grosberg 2018).

## Intertidal versus subtidal

Biogenic habitats did not have a consistent effect on community dynamics in the two tidal zones, supporting the existence of different underlying mechanisms:

1) Foundation species and stress levels: the positive effect of foundation species on biodiversity can be emphasized in stressful environments (e.g. intertidal zone) and can be dampened in mild environments (e.g. subtidal zone) (Watt and Scrosati 2013).

2) Different foundation species: the high diversity of seagrass beds seem to be more dependent on transient species (Boyé et al. 2019) while we posit that it stems more from rare ones in maerl beds. Transient species can have a strong effect on the long-term return of communities through the insurance hypothesis (Arnoldi et al. 2018). However, transients rarely interact with other members of the community (Snell Taylor et al. 2018), thus weaker species interactions are expected in seagrass beds communities. Community composition might be temporally less variable with community dynamics governed by weak interactions (van Nes and Scheffer 2005, Magurran and Henderson 2010, Thrush et al. 2021) as illustrated in IBIO.

3) Richness artifact: shorter trajectories in IBIO compared to SBIO might be an artifact of their large difference in species richness (493 in IBIO versus 665 in SBIO). Higher species richness with a higher proportion of LOLAS in SBIO might induce a higher variability in community composition between years.

The two bare habitats also had distinct dynamics with SBAR being less variable in structure and composition of communities than IBAR (Supporting information). Subtidal macrobenthic communities have previously been shown to have heterogeneous dynamics. Some studies reported changes

and shifts in community composition mostly linked to environmental changes (Warwick et al. 2002, Frid 2011, Ghodrati Shojaei et al. 2016, Bonifácio et al. 2019), but results are inconsistent across spatial and temporal scales. Indeed, other studies demonstrated temporal stability of subtidal communities even under environmental or anthropogenic pressures (e.g. Hinz et al. 2011, Quillien et al. 2018, Bacouillard et al. 2020). In contrast to SBAR, sandy beaches (IBAR) are physically dynamic environment: instability is characteristic of this habitat where temporal environmental variability is a key driver compared to biotic interactions that are considered less important (Defeo and McLachlan 2005, Schlacher et al. 2008, Schlacher and Thompson 2013). Even if competition and predation in beach ecosystems are limited compared to all other littoral ones (Schlacher and Thompson 2013) which could have induced a higher stability (Loreau and de Mazancourt 2013), in this habitat the physical environment and habitat conditions are the main controlling factors (Defeo and McLachlan 2005, McLachlan and Dorvlo 2005). Our results showed that compositional dynamics of IBAR can be described as non-directional with saltatory changes in the first five years of the monitoring program. Directional saltatory changes usually translate into a shift of the community from one stable state to a different one (Scheffer and Carpenter 2003, deYoung et al. 2008, Lamothe et al. 2019). However, directionality of the first four segments was indistinguishable from that of the last nine segments, suggesting that the system did not shift from one regime to another. Non-directional saltatory changes could occur when systems are facing multiple short term and delineated disturbances (Lamothe et al. 2019). Gray (1977) suggested that soft-sediments systems have multiple stable-states and tend to return to an equilibrium point after any perturbation of limited extent. Gray and Elliott (2009) also argued that marine benthic communities can exhibit poly-climax and neighborhood stability with several alternate dominant species. Interestingly, communities in IBAR presented greater changes in dominant taxa during the monitored period than any other habitats (Supporting information). This could be an artifact of the low number of IBAR sites. Indeed, global stability is mostly observed when measures are conducted over large spatial scales while neighborhood stability is a more appropriate model over relatively small spatial scales (Gray and Elliott 2009). Notwithstanding, other habitats had similar number of sites (only one more for SBIO compared to IBAR) and we observed more of a global stability pattern, hence we assume the neighborhood pattern observed is not an artifact of the number of sites.

Because in IBAR an ecological change is mainly the result of an environmental change (Schoeman et al. 2014), we hypothesized that, along the first five years of the survey, communities had to face multiple pulse disturbances that led to changes in taxa dominance reflected in the observed trajectory patterns (saltatory changes between consecutive years). Moreover, CTA metrics from observed and simulated communities were the most similar in IBAR. The effect of pulse perturbations could break temporal autocorrelation in the observed communities as in the null model, because
communities with low species richness can often be reset when facing harsh environmental conditions (Boyé et al. 2019). dBD\_tot was the highest in IBAR suggesting that temporal dynamics differed more between sites: either the same type and frequency of disturbance on different sites produced different responses at the community level (because of the different species composition or history between sites) or the sources and timing of disturbances may be different (which is often the case in intertidal habitats where natural and anthropogenic disturbances are confounded) (Whomersley et al. 2010). Indeed, there seemed to be site-specific effects on community dynamics in each habitat and maybe even more in IBAR: distributions of CTA metrics were always more heterogeneous compared to simulations where the sites characteristics and differences are not taken into account. To verify whether the observed variability could be linked to greater environmental variability during the beginning of the monitoring (e.g. alternating colder winters and warmer summers; Beukema and Dekker 2020) we looked at temperature, wind speed and precipitation anomalies (Supporting information) for each year of the time series and every site in IBAR. We could not detect any trend changes between 2005–2009 and 2009–2018, nor any pulse perturbations that could explain the observed trajectories. However, partitioning the components of beta diversity (Legendre and De Cáceres 2013, Legendre 2014) revealed a higher proportion of species replacement between 2006 and 2007 than ever observed over the whole time series (Supporting information). This peak in replacement is followed by consecutive years of higher richness differences from 2007 to 2009, and a lower mean species richness after 2009. Interestingly, temperature and wind speed anomalies started to oscillate between positive and negative values from 2006. One hypothesis could be that the system changed from a system with relatively stable environmental conditions to a system with a varying environment. This might have led to a loss of intolerant species and colonization of species benefiting the empty space (high replacement), followed by a stabilization phase with competitive exclusion and disappearance of opportunistic species (high richness differences), leading to a new non-directional stable state constituted of species more tolerant to a varying environment. This hypothesis is coherent with the neighborhood stability hypothesis with alternate dominant species. However, we would have expected the succession of replacement and richness differences would lead to a higher directionality of the trajectory for the first four segments which was not the case. It raises the perspective to look at the different beta diversity components to better apprehend the changes in trajectory metrics in further work, and also to look at community changes preceding the observed trajectory changes because of a potential lag in species response and dynamics.

## Linking trajectory shapes to ecological dynamics

To summarize, community state can be represented as the position of a ball in a cup-shaped landscape (Lamothe et al. 2019, Dakos and Kéfi 2022). In such a landscape, the ball

represents the current state of the system, the cup represents its current domain of attraction. Communities showing high stability are portrayed as balls within deep cups with steep walls. In this case, important perturbations are needed to change the community state, and the community will quickly revert to its previous state. In our study:

1) Observed variations in IBAR could correspond to a neighborhood stability (Gray and Elliott 2009): within a wide cup, community state alternates between shallow neighboring cups – each corresponding to a different dominant taxon – during the first five years of the survey and stays in another cup for the remaining nine years (Fig. 5). Such a landscape can be pictured by multiple narrow and shallow attraction basins in the ball-and-cup analogy (Fig. 5a).

2) SBAR showed more of a general stability pattern with moderate variability between consecutive years compared to SBIO or IBIO. Such a moderate variability stability landscape can be pictured by a deep attraction basin in the ball-and-cup analogy (Fig. 5b).

3) IBIO showed a general stability pattern with the lowest variability of communities between consecutive years among all the studied habitats. This may be due to the influence of transient species, which weaken the role of species interactions in community dynamics, thus narrowing the attraction basin. Such a low variability stability landscape can be pictured by a narrow and deep attraction basin in the ball-and-cup analogy (Fig. 5c).

4) SBIO showed a general stability pattern with a high variability between consecutive years compared to SBAR and IBIO. We hypothesize that the high diversity and number of rare species promote long term stability in this habitat but that they are also factors that can induce a higher temporal turnover, thus a wider attraction basin. Such a high variability stability landscape can be pictured by a wide and deep attraction basin in the ball-and-cup analogy (Fig. 5d).

## Conclusions

Our study demonstrated that habitat-dependent dynamics exist despite a general stability of macrobenthic communities at the regional scale. Indeed, the four monitored habitats could be matched to four different patterns of stability and we hypothesize 1) that these are driven by different mechanisms, related to biotic and abiotic factors involved such as habitat characteristics and 2) that the presence of numerous LOLAS can maintain the long-term stability of the systems even if it increases the temporal turnover of communities compared to the dynamics of core species only. Our study corroborates that community dynamics are not consistent across habitats and scales, which could have important consequences in the context of global change. Thus, it highlights the importance of taking the scale of observation into account in temporal studies because an observed local heterogeneity can be part of a system that is fairly homogeneous (Gray and Elliott 2009).
<figure><img src="image_9.png" /><figcaption>Figure 5. Synthetic figure of the global pattern of temporal trajectories (on the left in each box) translated into a ball and cup representation (on the right in each box) and the hypothetic drivers of the patterns observed for each of the four monitored habitats. The red ball represents the community state of a given observation and the arrows indicate the directionality of movement across the cup-shaped landscape (i.e. all possible states the community could take in this landscape). Tr. spp. stands for transient species and rare spp. for rare species.</figcaption></figure>

It would be interesting in future works to investigate the threshold of variability in communities beyond which change would be of concern and fall outside the basin of attraction. Finally, our study confirmed the usefulness of CTA, coupled to null models, in order to describe, understand and draw hypotheses or conclusions on ecological dynamics of long-term monitored macrobenthic communities. Still, further work needs to be done in order to better understand and discriminate the biotic and abiotic drivers of the different temporal dynamics observed as well as the dynamic of communities when looking at other dimensions of diversity (e.g. functional or phylogenetic diversities) which are fundamental to better assess temporal beta diversity in ecological communities (Magurran et al. 2019).

**Acknowledgements** – We thank the crew of the research vessels ‘Thalia’, ‘Côtes de la Manche’, ‘Neomysis’ and ‘Albert Lucas’ for their valuable assistance during the different field surveys. We also gratefully acknowledge all students, technicians and colleagues who contribute to the sampling and laboratory works since 2003. We thank Mathias Rouan and Michel Aidonidis for providing us with meteorological data from the Météo France models. We also thank Mathieu Helias for his help with the figures’ arrangement. Finally, we are grateful for the valuable and insightful comments from the reviewers of this manuscript that helped us to improve its quality.

**Funding** – This study was supported by the REBENT program coordinated by Sandrine Derrien (MNHN) and its funding partners (Agence de l'eau Loire-Bretagne, Région Bretagne, DREAL Bretagne). This work was supported by ISblue project, Interdisciplinary graduate school for the blue planet (ANR-17-EURE-0015) and co-funded by a grant from the French government under the program ‘Investissements d’Avenir’ embedded in France 2030. This work was made possible thanks to the campaigns conducted by the Flotte Océanographique Française (https://doi.org/10.18142/160).

## Author contributions

**Chirine Toumi:** Conceptualization (equal); Data curation (lead); Formal analysis (lead); Methodology (equal); Visualization (lead); Writing – original draft (lead); Writing – review and editing (lead). **Miquel De Cáceres:** Conceptualization (equal); Methodology (equal); Writing – review and editing (equal). **Jacques Grall:** Conceptualization (equal); Investigation (equal); Supervision (equal); Writing – review and editing (equal). **Aurélien Boyé:** Writing – review and editing (equal). **Éric Thiébaut:** Writing – review and editing (equal). **Marion Maguer:** Investigation (equal). **Vincent Le Garrec:** Investigation (equal). **Caroline Broudin:** Investigation (equal). **Céline Houbin:** Investigation (equal). **Olivier Gauthier:** Conceptualization (equal); Supervision (equal); Writing – review and editing (equal).

## Transparent peer review

The peer review history for this article is available at https://publons.com/publon/10.1111/ecog.06489.

## Data availability statement

Data are available from the Dryad Digital Repository: https://doi.org/doi:10.5061/dryad.8gth76rx (Toumi et al. 2023).

## Supporting information

The Supporting information associated with this article is available with the online version.
## References

Anderson, M. J. 2017. Permutational multivariate analysis of variance (PERMANOVA). — In: Balakrishnan, N., Colton, T., Everitt, B., Piegorsch, W., Ruggeri, F. and Teugels, J. L. (eds), Wiley StatsRef: statistics reference online, pp. 1–15.

Arnoldi, J.-F., Bideault, A., Loreau, M. and Haegeman, B. 2018. How ecosystems recover from pulse perturbations: a theory of short- to long-term responses. — J. Theor. Biol. 436: 79–92.

Ayata, S.-D., Lazure, P. and Thiébaut, É. 2010. How does the connectivity between populations mediate range limits of marine invertebrates? A case study of larval dispersal between the Bay of Biscay and the English Channel (North-East Atlantic). — Prog. Oceanogr. 87: 18–36.

Bacouillard, L., Baux, N., Dauvin, J.-C., Desroy, N., Geiger, K. J., Gentil, F. and Thiébaut, É. 2020. Long-term spatio-temporal changes of the muddy fine sand benthic community of the Bay of Seine (eastern English Channel). — Mar. Environ. Res. 161: 105062.

Bessa, F., Gonçalves, S. C., Franco, J. N., André, J. N., Cunha, P. P. and Marques, J. C. 2014. Temporal changes in macrofauna as response indicator to potential human pressures on sandy beaches. — Ecol. Indic. 41: 49–57.

Besse, P. C., Guillouet, B., Loubes, J.-M. and Royer, F. 2016. Review and perspective for distance-based clustering of vehicle trajectories. — IEEE Trans. Intell. Transp. Syst. 17: 3306–3317.

Beuchel, F., Gulliksen, B. and Carroll, M. L. 2006. Long-term patterns of rocky bottom macrobenthic community structure in an Arctic fjord (Kongsfjorden, Svalbard) in relation to climate variability (1980–2003). — J. Mar. Syst. 63: 35–48.

Beukema, J. and Dekker, R. 2020. Half a century of monitoring macrobenthic animals on tidal flats in the Dutch Wadden Sea. — Mar. Ecol. Prog. Ser. 656: 1–18.

Blowes, S. A. et al. 2019. The geography of biodiversity change in marine and terrestrial assemblages. — Science 366: 339–345.

Bonifácio, P., Grémare, A., Amouroux, J.-M. and Labrune, C. 2019. Climate-driven changes in macrobenthic communities in the Mediterranean Sea: a 10-year study in the Bay of Banyuls-sur-Mer. — Ecol. Evol. 9: 10483–10498.

Boström, C. and Bonsdorff, E. 2000. Zoobenthic community establishment and habitat complexity—the importance of seagrass shoot-density, morphology and physical disturbance for faunal recruitment. — Mar. Ecol. Prog. Ser. 205: 123–138.

Bouma, T. J., Olenin, S., Reise, K. and Ysebaert, T. 2009. Ecosystem engineering and biodiversity in coastal sediments: posing hypotheses. — Helgol. Mar. Res. 63: 95–106.

Boyé, A., Legendre, P., Grall, J. and Gauthier, O. 2017. Constancy despite variability: local and regional macrofaunal diversity in intertidal seagrass beds. — J. Sea Res. 130: 107–122.

Boyé, A., Thiébaut, É., Grall, J., Legendre, P., Broudin, C., Houbin, C., Le Garrec, V., Maguer, M., Droual, G. and Gauthier, O. 2019. Trait-based approach to monitoring marine benthic data along 500 km of coastline. — Divers. Distrib. 25: 1879–1896.

Buckley, H. L., Day, N. J., Case, B. S. and Lear, G. 2021a. Measuring change in biological communities: multivariate analysis approaches for temporal datasets with low sample size. — PeerJ 9: e11096.

Buckley, H. L., Day, N. J., Lear, G. and Case, B. S. 2021b. Changes in the analysis of temporal community dynamics data: a 29-year literature review. — PeerJ 9: e11250.

Bulleri, F., Eriksson, B. K., Queirós, A., Airoldi, L., Arenas, F., Arvanitidis, C., Bouma, T. J., Crowe, T. P., Davoult, D., Gui-

zien, K., Iveša, L., Jenkins, S. R., Michalet, R., Olabarria, C., Procaccini, G., Serrão, E. A., Wahl, M. and Benedetti-Cecchi, L. 2018. Harnessing positive species interactions as a tool against climate-driven loss of coastal biodiversity. — PLoS Biol. 16: e2006852.

Cimon, S. and Cusson, M. 2018. Impact of multiple disturbances and stress on the temporal trajectories and resilience of benthic intertidal communities. — Ecosphere 9: e02467.

Coyle, J. R., Hurlbert, A. H. and White, E. P. 2013. Opposing mechanisms drive richness patterns of core and transient bird species. — Am. Nat. 181: E83–E90.

Craven, D. et al. 2018. Multiple facets of biodiversity drive the diversity-stability relationship. — Nat. Ecol. Evol. 2: 1579–1587.

Dakos, V. and Kéfi, S. 2022. Ecological resilience: what to measure and how. — Environ. Res. Lett. 17: 043003.

Dauvin, J.-C., Ruellet, T., Desroy, N. and Janson, A.-L. 2007. The ecological quality status of the Bay of Seine and the Seine estuary: use of biotic indices. — Mar. Pollut. Bull. 55: 241–257.

Dauvin, J.-C., Lucas, S., Navon, M., Lesourd, S., Mear, Y., Poizot, E. and Alizier, S. 2017. Does the hydrodynamic, morphometric and sedimentary environment explain the structure of soft-bottom benthic assemblages in the eastern Bay of Seine (English Channel)? — Estuar. Coast. Shelf Sci. 189: 156–172.

De Cáceres, M., Coll, L., Legendre, P., Allen, R. B., Wiser, S. K., Fortin, M., Condit, R. and Hubbell, S. 2019. Trajectory analysis in community ecology. — Ecol. Monogr. 89: e01350.

Defeo, O. and McLachlan, A. 2005. Patterns, processes and regulatory mechanisms in sandy beach macrofauna: a multi-scale analysis. — Mar. Ecol. Prog. Ser. 295: 1–20.

deYoung, B., Barange, M., Beaugrand, G., Harris, R., Perry, R. I., Scheffer, M. and Werner, F. 2008. Regime shifts in marine ecosystems: detection, prediction and management. — Trends Ecol. Evol. 23: 402–409.

Dornelas, M., Magurran, A. E., Buckland, S. T., Chao, A., Chazdon, R. L., Colwell, R. K., Curtis, T., Gaston, K. J., Gotelli, N. J., Kosnik, M. A., McGill, B., McCune, J. L., Morlon, H., Mumby, P. J., Øvreås, L., Studeny, A. and Vellend, M. 2013. Quantifying temporal change in biodiversity: challenges and opportunities. — Proc. R. Soc. B 280: 20121931.

Dornelas, M., Gotelli, N. J., McGill, B., Shimadzu, H., Moyes, F., Sievers, C. and Magurran, A. E. 2014. Assemblage time series reveal biodiversity change but not systematic loss. — Science 344: 296–299.

Downing, A. L., Brown, B. L. and Leibold, M. A. 2014. Multiple diversity-stability mechanisms enhance population and community stability in aquatic food webs. — Ecology 95: 173–184.

Dray, S., Bauman, D., Blanchet, G., Borcard, D., Clappe, S., Guenard, G., Jombart, T., Larocque, G., Legendre, P., Madi, N. and Wagner, H. H. 2021. adespatial: multivariate multiscale spatial analysis. — R package ver. 0.3-14, https://cran.r-project.org/package=adespatial.

Duffy, J. E., Godwin, C. M. and Cardinale, B. J. 2017. Biodiversity effects in the wild are common and as strong as key drivers of productivity. — Nature 549: 261–264.

Ellis, J. I., Norkko, A., Thrush, S. F. 2000. Broad-scale disturbance of intertidal and shallow sublittoral soft-sediment habitats; effects on the benthic macrofauna. — J. Aquat. Ecosyst. Stress Recovery 7: 57–74.

Frid, C. L. J. 2011. Temporal variability in the benthos: does the sea floor function differently over time? — J. Exp. Mar. Biol. Ecol. 400: 99–107.

Fromentin, J. M., Ibanez, F., Dauvin, J. C., Dewarumez, J. M. and Elkaim, B. 1997. Long-term changes of four macrobenthic
assemblages from 1978 to 1992. – J. Mar. Biol. Assoc. 77: 287–310.

Gallon, R. K. 2017. Regional and latitudinal patterns of soft-
bottom macrobenthic invertebrates along French coasts:
results from the RESOMAR database. – J. Sea Res. 130:
96–106.

Gamfeldt, L., Lefcheck, J. S., Byrnes, J. E. K., Cardinale, B. J., Duffy, J. E. and Griffin, J. N. 2015. Marine biodiversity and ecosystem functioning: what's known and what's next? – Oikos 124: 252–265.

Ghodrati Shojaei, M., Gutow, L., Dannheim, J., Rachor, E., Schröder, A. and Brey, T. 2016. Common trends in German Bight benthic macrofaunal communities: assessing temporal variability and the relative importance of environmental varia- bles. – J. Sea Res. 107: 25–33.

Gilkinson, K. D., Gordon, D. C., MacIsaac, K. G., McKeown, D. L., Kenchington, E. L. R., Bourbonnais, C. and Vass, W. P. 2005. Immediate impacts and recovery trajectories of macro-faunal communities following hydraulic clam dredging on Ban-quereau, eastern Canada. – ICES J. Mar. Sci. 62: 925–947.

Giron-Nava, A., James, C., Johnson, A., Dannecker, D., Kolody, B., Lee, A., Nagarkar, M., Pao, G., Ye, H., Johns, D. and Sug-ihara, G. 2017. Quantitative argument for long-term ecological monitoring. – Mar. Ecol. Prog. Ser. 572: 269–274.

Grall, J. and Glémarec, M. 1997. Biodiversité des fonds de maerl en Bretagne: approche fonctionnelle et impacts anthropiques. – Vie et Milieu / Life and Environ. 47: 339–349.

Grall, J. and Hall-Spencer, J. M. 2003. Problems facing maerl con-
servation in Brittany. – Aquat. Conserv. 13: S55–S64.

Gray, J. S. 1977. The stability of benthic ecosystems. Helgolander Wiss. Meeresunters 30: 427–444.

Gray, J. S. and Elliott, M. 2009. Ecology of marine sediments: from science to management, 2nd ed. – Oxford Univ. Press.

Henderson, P. A. and Magurran, A. E. 2014. Direct evidence that density-dependent regulation underpins the temporal stability of abundant species in a diverse animal community. – Proc. R. Soc. B 281: 20141336.

Hewitt, J. E. and Thrush, S. F. 2009. Reconciling the influence of global climate phenomena on macrofaunal temporal dynamics at a variety of spatial scales. – Global Change Biol. 15: 1911–1929.

Hewitt, J. E., Thrush, S. F. and Ellingsen, K. E. 2016. The role of time and species identities in spatial patterns of species richness and conservation. – Conserv. Biol. 30: 1080–1088.

Hillebrand, H., Blasius, B., Borer, E. T., Chase, J. M., Downing, J. A., Eriksson, B. K., Filstrup, C. T., Harpole, W. S., Hodapp, D., Larsen, S., Lewandowska, A. M., Seabloom, E. W., Van de Waal, D. B. and Ryabov, A. B. 2018. Biodiversity change is uncoupled from species richness trends: consequences for con- servation and monitoring. – J. Appl. Ecol. 55: 169–184.

Hily, C. and Bouteille, M. 1999. Modifications of the specific diver-
sity and feeding guilds in an intertidal sediment colonized by
an eelgrass meadow (Zostera marina) (Brittany, France). – C. R.
Acad. Sci., Ser. III 322: 1121–1131.

Hinz, H., Capasso, E., Lilley, M., Frost, M. and Jenkins, S. 2011. Temporal differences across a bio-geographical boundary reveal slow response of sub-littoral benthos to climate change. – Mar. Ecol. Prog. Ser. 423: 69–82.

Ives, A. R. and Carpenter, S. R. 2007. Stability and diversity of ecosystems. – Science 317: 58–62.

Jurgens, L. J. and Gaylord, B. 2018. Physical effects of habitat-
forming species override latitudinal trends in temperature. –
Ecol. Lett. 21: 190–196.

Kéfi, S., Domínguez-García, V., Donohue, I., Fontaine, C., Thébault, E. and Dakos, V. 2019. Advancing our understanding of ecological stability. – Ecol. Lett. 22: 1349–1356.

Kominoski, J. S., Gaiser, E. E. and Baer, S. G. 2018. Advancing theories of ecosystem development through long-term ecological research. – BioScience 68: 554–562.

Kremen, C. 2005. Managing ecosystem services: what do we need to know about their ecology? – Ecol. Lett. 8: 468–479.

Kröncke, I., Neumann, H., Dippner, J. W., Holbrook, S., Lamy, T., Miller, R., Padedda, B. M., Pulina, S., Reed, D. C., Reini-kainen, M., Satta, C. T., Sechi, N., Soltwedel, T., Suikkanen, S. and Lugliè, A. 2019. Comparison of biological and eco-logical long-term trends related to northern hemisphere cli-mate in different marine ecosystems. – Nat. Conserv. 34: 311–341.

Lamothe, K. A., Somers, K. M. and Jackson, D. A. 2019. Linking the ball-and-cup analogy and ordination trajectories to describe ecosystem stability, resistance, and resilience. – Ecosphere 10: e02629.

Lamy, T., Koenigs, C., Holbrook, S. J., Miller, R. J., Stier, A. C. and Reed, D. C. 2020. Foundation species promote community stability by increasing diversity in a giant kelp forest. – Ecology 101: e02987.

Lefran, A., Hernández-Fariñas, T., Gohin, F. and Claquin, P. 2021. Decadal trajectories of phytoplankton communities in contrasted estuarine systems in an epicontinental sea. – Estuar. Coast. Shelf Sci. 258: 107409.

Legendre, P. 2014. Interpreting the replacement and richness difference components of beta diversity: replacement and richness difference components. – Global Ecol. Biogeogr. 23: 1324–1334.

Legendre, P. and De Cáceres, M. 2013. Beta diversity as the variance of community data: dissimilarity coefficients and partitioning. – Ecol. Lett. 16: 951–963.

Legendre, P. and Gallagher, E. D. 2001. Ecologically meaningful transformations for ordination of species data. – Oecologia 129: 271–280.

Loreau, M. and de Mazancourt, C. 2013. Biodiversity and ecosystem stability: a synthesis of underlying mechanisms. – Ecol. Lett. 16: 106–115.

Magurran, A. E., Dornelas, M., Moyes, F. and Henderson, P. A. 2019. Temporal β diversity—A macroecological perspective. – Global Ecol. Biogeogr. 28: 1949–1960.

Magurran, A. E. and Henderson, P. A. 2003. Explaining the excess of rare species in natural species abundance distributions. – Nature 422: 714–716.

Magurran, A. E. and Henderson, P. A. 2010. Temporal turnover and the maintenance of diversity in ecological assemblages. – Phil. Trans. R. Soc. B 365: 3611–3620.

Matthews, W. J., Marsh-Matthews, E., Cashner, R. C. and Gel-
wick, F. 2013. Disturbance and trajectory of change in a stream
fish community over four decades. – Oecologia 173: 955–969.

McLachlan, A. and Dorvlo, A. 2005. Global patterns in sandy beach macrobenthic communities. – J. Coast. Res. 214: 674–687.

O'Leary, J. K., Micheli, F., Airoldi, L., Boch, C., De Leo, G., Elahi, R., Ferretti, F., Graham, N. A. J., Litvin, S. Y., Low, N. H., Lummis, S., Nickols, K. J. and Wong, J. 2017. The resilience of marine ecosystems to climatic disturbances. – BioScience 67: 208–220.

Palumbi, S. R., McLeod, K. L. and Grünbaum, D. 2008. Ecosystems in action: lessons from marine ecology about recovery, resistance, and reversibility. – BioScience 58: 33–42.
Pereira, H. M. 2013. Essential biodiversity variables. – Science 339: 277–278.

Peterson, C., Luettich, R., Micheli, F. and Skilleter, G. 2004. Attenuation of water flow inside seagrass canopies of differing structure. – Mar. Ecol. Prog. Ser. 268: 81–92.

Quillien, N., Nordström, M. C., Le Bris, H., Bonsdorff, E. and Grall, J. 2018. Green tides on inter- and subtidal sandy shores: differential impacts on infauna and flatfish. – J. Mar. Biol. Ass. 98: 699–712.

Raffaelli, D., Solan, M. and Webb, T. J. 2005. Do marine and terrestrial ecologists do it differently? – Mar. Ecol. Prog. Ser. 304: 283–289.

Scheffer, M. and Carpenter, S. R. 2003. Catastrophic regime shifts in ecosystems: linking theory to observation. – Trends Ecol. Evol. 18: 648–656.

Schlacher, T. A., Schoeman, D. S., Dugan, J., Lastra, M., Jones, A., Scapini, F. and McLachlan, A. 2008. Sandy beach ecosystems: key features, sampling issues, management challenges and climate change impacts. – Mar. Ecol. 29: 70–90.

Schlacher, T. A. and Thompson, L. 2013. Environmental control of community organisation on ocean-exposed sandy beaches. – Mar. Freshwater Res. 64: 119.

Schoeman, D. S., Schlacher, T. A. and Defeo, O. 2014. Climate-change impacts on sandy-beach biota: crossing a line in the sand. – Global Change Biol. 20: 2383–2392.

Schückel, U. and Kröncke, I. 2013. Temporal changes in intertidal macrofauna communities over eight decades: a result of eutrophication and climate change. – Estuar. Coast. Shelf Sci. 117: 210–218.

Snell Taylor, S. J., Evans, B. S., White, E. P. and Hurlbert, A. H. 2018. The prevalence and impact of transient species in ecological communities. – Ecology 99: 1825–1835.

Sturbois, A., De Cáceres, M., Sánchez-Pinillos, M., Schaal, G., Gauthier, O., Mao, P. L., Ponsero, A. and Desroy, N. 2021a. Extending community trajectory analysis: new metrics and representation. – Ecol. Modell. 440: 109400.

Sturbois, A., Cormy, G., Schaal, G., Gauthier, O., Ponsero, A., Le Mao, P., Riera, P. and Desroy, N. 2021b. Characterizing spatiotemporal changes in benthic communities: taxonomic and functional trajectories of intertidal assemblages in the bay of Saint-Brieuc (English Channel). – Estuar. Coast. Shelf Sci. 107603.

Thomaz, S. M. and da Cunha, E. R. 2010. The role of macrophytes in habitat structuring in aquatic ecosystems: methods of measurement, causes and consequences on animal assemblages' composition and biodiversity. – Acta 22: 218–236.

Thrush, S., Hewitt, J., Pilditch, C. and Norkko, A. 2021. Ecology of coastal marine sediments: form, function, and change in the Anthropocene, 1st ed. – Oxford Univ. Press.

Toumi, C., De Cáceres, M., Grall, J., Boyé, A., Thiébaut, E., Maguer, M., Le Garrec, V., Broudin, C., Houbin, C. and Gauthier, O. 2023. Data from: Long-term coastal macrobenthic community trajectory analysis reveals habitat-dependent stability patterns – Dryad Digital Repository. – Dryad Digital Repository, https://doi.org/doi:10.5061/dryad.8gt76rx.

Van Meerbeek, K., Jucker, T. and Svenning, J.-C. 2021. Unifying the concepts of stability and resilience in ecology. – J. Ecol. 109: 3114–3132.

van Nes, E. H. and Scheffer, M. 2005. Implications of spatial heterogeneity for catastrophic regime shifts in ecosystems. – Ecology 86: 1797–1807.

Vermeij, G. J. and Grosberg, R. K. 2018. Rarity and persistence. – Ecol. Lett. 21: 3–8.

Warwick, R., Ashman, C., Brown, A., Clarke, K., Dowell, B., Hart, B., Lewis Re Shillabeer, N., Somerfield, P. and Tapp, J. 2002. Inter-annual changes in the biodiversity and community structure of the macrobenthos in Tees Bay and the Tees estuary, UK, associated with local and regional environmental events. – Mar. Ecol. Prog. Ser. 234: 1–13.

Watt, C. A. and Scrosati, R. A. 2013. Bioengineer effects on understory species richness, diversity, and composition change along an environmental stress gradient: experimental and mensurative evidence. – Estuar. Coast. Shelf Sci. 123: 10–18.

Whomersley, P., Huxham, M., Bolam, S., Schratzberger, M., Augley, J. and Ridland, D. 2010. Response of intertidal macrofauna to multiple disturbance types and intensities – An experimental approach. – Mar. Environ. Res. 69: 297–308.

WoRMS Editorial Board 2021. World register of marine species. – www.marinespecies.org at VLIZ.

Yachi, S. and Loreau, M. 1999. Biodiversity and ecosystem productivity in a fluctuating environment: the insurance hypothesis. – Proc. Natl Acad. Sci. USA 96: 1463–1468.