# Cross-scale environmental impacts across persistent and dynamic aggregations within a complex population: implications for fisheries management

Georgios Kerametsidis, James Thorson, Vincent Rossi, Diego

Álvarez-Berastegui, Cheryl Barnes, Gregoire Certain, Antonio Esteban, Encarnación García, Angélique Jadaud, Safo Piñeiro, et al.

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

## ► To cite this version:

Georgios Kerametsidis, James Thorson, Vincent Rossi, Diego Álvarez-Berastegui, Cheryl Barnes, et al.. Cross-scale environmental impacts across persistent and dynamic aggregations within a complex population: implications for fisheries management. Canadian Journal of Fisheries and Aquatic Sciences, 2023, 00, pp.1-17. 10.1139/cjfas-2023-0120. http://04480896

HAL Id: hal-04480896

https://hal.science/hal-04480896v1

Submitted on 2 Mar 2024

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>
# Cross-scale environmental impacts across persistent and dynamic aggregations within a complex population: implications for fisheries management

Georgios Kerametsidis<sup>a,b</sup>, James Thorson<sup>c</sup>, Vincent Rossi<sup>d</sup>, Diego Álvarez-Berastegui<sup>a</sup>, Cheryl Barnes<sup>e</sup>, Gregoire Certain<sup>f</sup>, Antonio Esteban<sup>g</sup>, Encarnación García<sup>g</sup>, Angélique Jadaud<sup>f</sup>, Safo Piñeiro<sup>d,a</sup>, Miguel Vivas<sup>g</sup>, and Manuel Hidalgo<sup>d,a</sup>

<sup>a</sup>Centro Oceanográfico de Baleares (IEO, CSIC), Moll de Ponent s/n, 07015, Palma, Balearic Islands, Spain; <sup>b</sup>University of the Balearic Islands, Carretera de Valldemossa, km 7.5, 07122, Palma, Balearic Islands, Spain; <sup>c</sup>Habitat and Ecological Processes Research Program, Alaska Fisheries Science Center, NOAA, Seattle, WA, USA; <sup>d</sup>Mediterranean Institute of Oceanography (UM110, UMR 7294), CNRS, Aix Marseille Univ., Univ. Toulon, IRD, Marseille, France; <sup>e</sup>Department of Fisheries, Wildlife, and Conservation Sciences, Oregon State University, Newport, OR, USA; <sup>f</sup>MARBEC, Univ. Montpellier, CNRS, Ifremer, IRD, Sète, France; <sup>g</sup>Centro Oceanográfico de Murcia (IEO-CSIC), C/ Varadero, 1, 30740 Lo Pagan, Murcia, Spain

Corresponding author: Georgios Kerametsidis (email: georgios.kerametsidis@ieo.csic.es; georgios.kerametsidis@uib.eu)

## Abstract

Accounting for marine stocks spatiotemporal complexity has become one of the most pressing improvements that should be added to the new generation of stock assessment. Disentangling persistent and dynamic population subcomponents and understanding their main drivers of variation are still stock-specific challenges. Here, we hypothesized that the spatiotemporal variability of density in two adjacent fish stocks is associated with spatially structured environmental processes across multiple spatiotemporal scales. To test this, we applied a generalized empirical orthogonal function and dynamic factor analysis to fishery-independent and -dependent data of red mullet, a highly commercial species, in the Western Mediterranean Sea. Areas with persistent and dynamic high aggregations were detected for both stock units. A large-scale climatic index and local open-ocean convection were associated with both stocks, while other variables exhibited stock-specific effects. We also revealed spatially structured density dynamics within the examined management units. This suggests a metapopulation structure and supports the future implementation of a spatial stock assessment. Considering the common—generally unrealistic—assumptions of panmictic structure and absence of connectivity with neighbouring stock units, our methodology can be applied to other species and systems with putative spatial complexity to inform more accurate population dynamics and structure.

**Key words:** Mullus barbatus, empirical orthogonal function, dynamic factor analysis, Mediterranean Sea, open-ocean convection, stock identification

## 1. Introduction

Understanding species spatiotemporal distributions and densities is critical to linking ecological mechanisms and conservation measures. For commercially exploited species, obtaining information about the spatiotemporal dynamics of their distribution and density is a prerequisite to providing robust management advice (Rufener et al. 2021), while it is also central to the Ecosystem-Based Fisheries Management (EBFM). EBFM is widely viewed as a set of tools for successful management in different contexts (Trochta et al. 2018) and constitutes a fundamental component of sustainability and balance between environmental, social, and economic objectives (Marshall et al. 2018). Current stock assessment frameworks often do not include key features such as ecosystem components and ecological complexity of popu-

lations, although it has already been shown that environmental and climatic variation are closely associated with the abundance variation of commercially important species (e.g., Schlenker et al. 2023). However, the ongoing development of next-generation stock assessments highly prioritizes the incorporation of missing critical information (Punt et al. 2020). This inherently requires the identification of spatial structure to model more realistic population dynamics and develop reliable assessment and management frameworks (Punt 2019; Punt et al. 2020). Despite the recognition of complex spatial structures in many stocks (Reiss et al. 2009; Matić-Skoko et al. 2018; Hidalgo et al. 2019a), management unit boundaries are not yet generally reconciled, and stock assessments rely on the assumption of single panmictic stocks for the majority of resources. A plethora of studies have shown that not
<figure><img src="image_3.png" /></figure>

accounting for spatial complexity can unintentionally lead to overfishing, with often devastating consequences for the stocks (Goethel and Berger 2016; Kerr et al. 2017; Cadrin 2020; Punt 2023). Conversely, modelling studies suggest that incorporating spatial complexity improves the performance of population models used in stock assessment (Goethel et al. 2011, 2021; Punt 2019).

A stock unit is generally defined as a subset of the whole population of a species that shares the same growth, recruitment, spawning, and mortality characteristics and inhabits the same geographic area (Sparre and Venema 1998). However, all these biological parameters and processes can exhibit non-stationary relationships with environmental variables. Not accounting for the latter may result in unrealistic outputs of the fisheries assessment models, which can hamper efficient fisheries management. This risk is even more pronounced in stocks with complex spatial structures where local environmental drivers may differ across different stock subunits. As an essential component of local habitats, dynamic environmental conditions can result in spatial restructuring of the stocks with further demographic implications (Szuwalski and Hollowed 2016; Kerr et al. 2017). The environment influences critical ecological processes such as recruitment (Houde 2016), spawning (Di Stefano et al. 2023) or mortality (Kerametsidis et al. 2023), as well as several life history parameters of fish, such as growth or reproductive scheduling (Perry et al. 2005; Clark et al. 2020). In fact, in some cases, no clear stock-recruitment patterns can be evidenced, which leaves the environment as the main driver of recruitment with minimal influence of the spawning stock (Szuwalski et al. 2015; Hidalgo et al. 2019b). Therefore, considering the spatial heterogeneity that exploited stocks are subject to, the hypothesis that spatially structured environmental processes (i.e., environmental processes varying across multiple spatial scales) explain critical ecological and demographic processes seems plausible and needs to be explored.

Various life stages and traits are linked to cross-scale environmental processes (Twiname et al. 2020), and thus both fine- and broad-scale environmental processes must be considered to properly comprehend the full life cycle of fish. Complex spatiotemporal variation among local habitats is often aggregated into broader-scale measurements such as climate indices (Stenseth and Mysterud 2005; Thorson et al. 2020a). Besides this cross-scale nature of environmentally driven processes, there are other important factors to consider when investigating the dynamics of animal populations. These include movements towards favourable habitats and persistent climatically and environmentally driven population hotspots. Hence, providing a mechanistic understanding of the links between environmental and density variability in space and time is essential, particularly for commercial stocks. This is especially critical considering the current climate change scenarios (e.g., Lotze et al. 2019) and environment-driven shifts in density and distribution that have been noted for a variety of organisms within the marine realm (Bowler et al. 2017; Champion et al. 2021; Thorson et al. 2021).

More than 30% of fish stocks are exploited beyond biologically sustainable limits globally (FAO 2022a). In certain basins, this percentage might be considerably higher. For instance, 73% of all fish stocks are overexploited in the Mediterranean Sea (FAO 2022b). Despite the recent implementation of highly promising management plans, such as the Multinational Plan for Demersal Fish Stocks in the Western Mediterranean Sea (Sánchez Lizaso et al. 2000), there are still issues inhibiting the proper management of commercial stocks (Cardinale et al. 2021). While recent studies have recommended larger species-specific units for the better assessment of some Mediterranean stocks (e.g., Fiorentino et al. 2015; Spedicato et al. 2021; STECF 2021), there is also cumulative evidence of the importance of local and regional dynamics (Hidalgo et al. 2019b; Paradinas et al. 2022), challenging these recommendations and encouraging the implementation of spatially structured stock assessment frameworks (Cadrin 2020; Punt 2023). Obtaining high-resolution information on population spatial structure and connectivity is a prerequisite for spatial stock assessments (Goethel et al. 2023) and is logistically feasible for data-rich species and systems.

For stocks with long, adequate and spatially resolved time series of fishery-dependent and -independent data, scientists and stakeholders might benefit from analyzing spatiotemporal patterns and/or connectivity between different stock units. In this study, we used red mullet (Mullus barbatus) in the north-western (NW) Mediterranean Sea as a case study, where its populations have been historically managed in two separate units, the Northern Spanish Coast and the Gulf of Lions. We hypothesized that the spatiotemporal variability among and within the two management units would be associated with specific spatially structured environmental processes across multiple scales. To test this hypothesis, we first described persistent and dynamic density hotspots, and then we identified the environmental variables that drive the spatiotemporal variability of red mullet density across different areas. We employed a two-fold spatiotemporal approach: (i) we used abundance-sampling data from scientific trawling to explore the principal modes of density variability, the persistent density hotspots, as well as the putative environmental drivers of the spatiotemporal variability in the two management units, and (ii) we analyzed monthly landings per unit effort (LPUE) data from commercial fishing operations to explore whether the seasonal and long-term trends are spatially structured.

# 2. Materials and methods

## 2.1. Target area and species

Our study area encompasses two management units that are used by the General Fisheries Commission for the Mediterranean—the Northern Spanish Coast (Geographic Sub-Area (GSA) 06) and the Gulf of Lions (Geographic Sub-Area (GSA) 07) in the NW Mediterranean—and covers the trawlable waters of these areas (Fig. 1).
**Fig. 1.** Study area encompassing Geographic Subareas 6 (GSA06, in dark blue) and 7 (GSA07, in red). Dark blue and red dots correspond to the MEDITS hauls in GSA06 and GSA07, respectively. The three subregions that GSA06 was subdivided into in this study are also shown: south of the Valencia Channel, north of the Valencia Channel, and Catalan Coast. The 20 ports (dark red triangles) that were used for the dynamic factor analysis (DFA) and the 200 and 1000 bathymetric contours (light grey lines) are also depicted. The two polygons around the MEDOC point that were used to calculate the mixed-layer depth (MLD) are present too. Map projection is WGS84 UTM.

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

*Mullus barbatus* (hereafter referred to as red mullet) is a demersal fish distributed in the eastern Atlantic Ocean, from the North Sea to Senegal, and throughout the Mediterranean and Black Seas (Fischer et al. 1987). Red mullet commonly resides over soft sandy and muddy substrates and is mainly distributed along the continental shelf in depths of up to 300–400 m, with significant declining trends in waters deeper than 200 m (Tserpes et al. 2019). It is a high-value commercial species and one of the most important and overexploited resources in Mediterranean demersal fisheries, predominantly targeted by bottom-trawl fisheries. It is one of the two demersal fish species that have long been assessed separately on an annual basis in the GSA06 and GSA07, while a complex demographic structure partially shaped by the elevated dispersal abilities in early life stages has been evidenced in other regions (Gargano et al. 2017). In the Western Mediterranean Sea, three spatially segregated and persistent density hotspots have been detected along the Northern Spanish Coast in the Mediterranean Sea (i.e., within the management unit GSA06) (Paradinas et al. 2020). This might indicate a more complex population structure resembling that of a metapopulation system, as suggested for other harvested species in the Mediterranean Sea (Hidalgo et al. 2019a; Gargano et al. 2022).

## 2.2. Data sources

### 2.2.1. Biological data

Data on red mullet were obtained from two different sources to apply two complementary modelling approaches: (1) standardized density data from scientific trawling and (2) LPUE from commercial fisheries. Annual standardized abundance data for GSA06 and GSA07 were collected during the EU-funded International Mediterranean Bottom Trawl Survey (MEDITS) that was carried out between spring and early summer (April–June) between 1994 and 2019 (Anonymous 2017; Spedicato et al. 2019). The MEDITS project follows a stratified sampling design based on the coverage of five bathymetric strata (10–50, 51–100, 101–200, 201–500, and 501–800 m) in each Mediterranean GSA and uses a GOC-73 net with a 20 mm mesh size. Sampling stations were randomly placed within each stratum at the beginning of the project, and in all subsequent years, sampling was carried out in similar randomized locations. The duration of the hauls has been specified at 30 min in depths of <200 m and 60 min in >200 m waters. Since the abundance of red mullet declines significantly over 200 m (Tserpes et al. 2019), data from hauls deeper than 200 m were excluded from this study. Preliminary
analysis confirmed this declining trend in our dataset. A total of 1499 and 1419 trawl operations were observed over a period of 26 years, averaging 58 and 55 fishing hauls per year for the GSA06 and the GSA07, respectively. Mean densities were approximately 651 and 494 ind. km⁻² haul⁻¹ in the GSA06 and GSA07, respectively, and red mullet was encountered in 47%–93% of hauls in each year in both areas.

Since MEDITS only takes place during a short prespecified period in late spring/early summer, monthly LPUE (during 2004–2019) from commercial fishing operations was used as a complementary data source to capture the intra- and inter-annual density variability. Monthly LPUE were standardized considering the number of vessels and the number of fishing days per month to make data from different ports comparable as landings in kg per fishing trip (Puerta et al. 2016). Fishing operations in the area are always developed on a daily basis within short distances from the port. LPUE data were used to explore the dynamics across the GSA06, as previous studies have already proposed three persistent hotspots for the species in the GSA06 (Paradinas et al. 2020). LPUE data for the GSA07 were not available. There are 52 ports with commercial fisheries operations in the GSA06; however, only those ports for which landing data were available for at least 10 months, of which 6 consecutives, every year during 2004–2019, were selected to be included in our analyses. The selected ports were grouped based on the three known aggregations: the area south of the Valencia Channel, the area north of the Valencia Channel, and the Catalan Coast (Fig. 1). In addition, to closely capture the dynamics in the persistent and dynamic density hotspots of the species (see Section 2.3.1), only ports near the Ebro Delta area (north of Valencia Channel) and in the very south of the area south of Valencia Channel were retained. After applying the above criteria, a total of 20 ports were included in the LPUE analyses (Fig. 1).

### 2.2.2. Environmental data

Eight environmental indices with known influence in the same region and other species were utilized in the present study (Table 1 and Supplementary materials S1). To capture the effects of broad-scale coupled atmosphere-ocean climatic processes, three indices were used. Data on the Atlantic Multidecadal Oscillation (AMO) and the North Atlantic Oscillation (NAO) were downloaded from the National Oceanic and Atmospheric Administration (NOAA, Physical Sciences Laboratory: https://www.psl.noaa.gov/data/correlation/amon.us.lo ng.data and https://www.cpc.ncep.noaa.gov/products/precip/ CWlink/pna/norm.nao.monthly.b5001.current.utf8.table; accessed on February 2022). AMO represents changes in the sea surface temperature (SST) in the North Atlantic Ocean from 0° to 70° N that are characterized by multidecadal variability (Dijkstra et al. 2006; Knight et al. 2006). A negative AMO phase is associated with an anomalously low North Atlantic SST. Several impacts of AMO on climate (Knight et al. 2006) and the fish communities of the NE Atlantic (Zimmermann et al. 2019) and the Mediterranean Sea (Tsikliras et al. 2019) have been documented, including influence on recruitment success and community ratios. Second, the NAO was employed.

NAO represents an alternation in pressure between the sub-tropical atmospheric high-pressure zone centered over the Azores and the atmospheric low-pressure zone over Iceland. A positive phase of NAO results in a higher frequency and stronger winter storms crossing the Atlantic Ocean in a more northerly track. Inversely, a negative phase of NAO is associated with fewer and weaker winter storms crossing on a more west-east pathway (Ottersen et al. 2001). Finally, we used the Western Mediterranean Oscillation index (WeMOi), which is defined as the difference in the standardized values in sea level pressures between Cádiz-San Fernando (Spain) and Padua (Italy) (Martin-Vide and Lopez-Bustins 2006) and has been shown to influence fish species in the Western Mediterranean (Martín et al. 2012). Data on WeMOi were obtained from the Climatology Group of the University of Barcelona (http://www.ub.edu/gc/documents/We b_WeMOi-2019.txt, accessed in February 2022).

Surface chlorophyll-a (chl-a) concentration was used as a proxy of primary productivity, and SST was used to inform the regional thermal conditions in each management unit over two contrasting seasons. For both variables, seasonal spatial averages corresponding to winter (December–February) and late spring/summer (May–July) were computed. To assess the influence of these two variables on a local scale and specifically on the aggregations in GSA06, besides an average for the whole GSA06, seasonal spatial averages were also calculated for the subareas of the Catalan Coast, the Ebro Delta, and the Valencia Channel (Fig. 1). For the GSA07, only the unique spatial average for the whole region was calculated and used. Data on chl-a and SST were downloaded from the E.U. Copernicus Marine Service (chl-a: https://doi.org/10.48670/moi-00300; SST: https://doi.org/10.25423/CMCC/ MEDSEA_MULTIYEAR PHY_006_004_E3R1). Further, a Local Climatic Index (LCI) was also included. LCI quantifies integrated hydroclimatic variability at the regional scale, synthesizing the air temperature, SST, atmospheric sea level pressure, 500 hPa geopotential height, and precipitation rates at a monthly scale by means of the first axis of a principal component analysis (Molinero et al. 2005). Positive LCI values are associated with higher regional atmospheric sea level pressure and 500 hPa geopotential height, whereas negative values are associated with a high precipitation rate. The Ebro River runoff (m³ s⁻¹), calculated at its mouth, was also included in the study as a proxy of the nutrient discharge from the river into the sea (Ebro Hydrographic Confederation, https://www.chebro.es/). Sea bottom temperature (SBT) was not included in this study as it is generally more conservative than temperatures in the rest of the water column (Hiddink and ter Hofstede 2008). Specifically, for GSA06 and GSA07 combined, we examined the monthly average SBT from 1992 to 2019, and we found very little monthly variation (i.e., mean SBT ranging from 12.8 to 13.5 °C) and we therefore decided to exclude it from the rest of the analyses.

Finally, we also assessed the potential influence of open-ocean convection in the NW Mediterranean, as observed in other species (Martin et al. 2016; Hidalgo et al. 2019a). This phenomenon, whose magnitudes and extents vary over time (Houpert et al. 2016), occurs during winter, typically between December and March, mainly in the Gulf of Lions and off the
<table><caption>**Table 1.** Summary of the environmental variables that were used in this study.</caption><thead><tr><th>Env. variable</th><th>Years</th><th>Model</th><th>Reference</th></tr></thead><tbody><tr><td colspan="4"><em>Large-scale climatic indices</em></td></tr><tr><td>AMO</td><td>1994–2019</td><td>EOF and DFA</td><td>NOAA—Physical Sciences Laboratory</td></tr><tr><td>NAO</td><td>1994–2019</td><td>EOF and DFA</td><td>NOAA—National Centers for Environmental Prediction</td></tr><tr><td>WeMOi</td><td>1994–2019</td><td>EOF and DFA</td><td>Martin-Vide and Lopez-Bustins (2006)</td></tr><tr><td colspan="4"><em>Regional and local indices</em></td></tr><tr><td>Chl-a</td><td>1998–2019</td><td>EOF and DFA</td><td>E.U. Copernicus Marine Service Information</td></tr><tr><td>SST</td><td>1994–2019</td><td>EOF and DFA</td><td>E.U. Copernicus Marine Service Information</td></tr><tr><td>LCI</td><td>1994–2019</td><td>EOF and DFA</td><td>Molinero et al. (2005)</td></tr><tr><td>Ebro runoff</td><td>1994–2019</td><td>EOF and DFA</td><td>Ebro Hydrographic Confederation</td></tr><tr><td>MLD</td><td>1994–2019</td><td>EOF</td><td>E.U. Copernicus Marine Service Information</td></tr></tbody></table>

Note: AMO, Atlantic Multidecadal Oscillation; NAO, North Atlantic Oscillation; WeMOi, Western Mediterranean Oscillation index; Chl-a, chlorophyll-a; SST, sea surface temperature; LCI, Local Climatic Index; MLD, mixed-layer depth; EOF, generalization of empirical orthogonal function (with a spatiotemporal model); DFA, dynamic factor analysis; NOAA, National Oceanic and Atmospheric Administration. The period for which the information was available, the models that the variable was included, and the references are given in the table.

northern Catalan coast (Fig. 1). It is associated with the regional incidence of northerly cold and dry winds and local oceanographic features, which favour winter vertical mixing and deepening of the mixed layer. During some years, atmospheric forcing can be particularly strong, and stratification erosion can reach great depths or even the seafloor (Mertens and Schott 1998). Open-ocean deep convection is a major regional oceanographic process that greatly contributes to the primary productivity and nutrient exchange fluctuations in the area (Lavigne et al. 2015; Macias et al. 2018). To model this process, the spatial mean mixed-layer depth (MLD) in the MEDOC area and in a greater polygon around the MEDOC area (MEDOC group 1970) were used as proxies of the vertical extent of winter mixing resulting from open-ocean deep convection in the region (Fig. 1). The MEDOC point (42°N, 5°E) is considered to be the center of convection (Martin et al. 2016).

In accordance with the standardized density-based and LPUE-based approaches, annual and monthly estimates of the above environmental variables were obtained and used, respectively (Table 1). For either approach, these variables were examined on 0-, 1-, and 2-year lags. The eight environmental variables that were utilized in this study were independent and non-collinear, as indicated by the variance inflation factor score (VIF < 5), and thus they were examined separately. An influence of AMO on NAO has been documented (Börgel et al. 2020), but here we examined the two indices independently, prompted by the variance inflation factor score that we obtained. Regarding chl-a and SST, monthly spatial averages were obtained on a regional scale (i.e., a mean for the whole GSA06 and GSA07) as well as on a local scale for the subareas of the Catalan Coast, the Ebro Delta, and the Valencia Channel, as described above.

## 2.3. Modelling approaches

### 2.3.1. Generalization of empirical orthogonal function

We seek to estimate one or more modes of spatiotemporal variability in red mullet density and investigate the latent local and nonlocal environmental conditions that may

drive this variability. To do so, we employ a statistical generalization of EOF using a spatiotemporal model (Thorson et al. 2020b). Red mullet spatiotemporal density patterns include three components: (i) temporal variation ($\beta$), where the intercept is specified to follow a random walk and represents fluctuations at all locations from year to year; (ii) spatial variation ($\omega$), which represents a static long-term spatial pattern; and (iii) spatiotemporal variation ($\epsilon$), estimated as unmeasured variation that is expressed with one or more dominant modes of variability as well as a map representing the spatial response of density to these estimated modes of variability. The third component is analogous to “conventional” empirical orthogonal function (EOF) analysis and generates time series representing modes (i.e., indices) of variability as well as the associated maps of response. However, in the present study, the spatiotemporal generalization to EOF is estimated after strictly accounting for expected spatial and temporal components (Thorson et al. 2021).

The applied model is a generalized linear mixed model that approximates the dependent variable of interest using a log link and linear predictors, including Gaussian Markov random fields representing spatial ($\omega$) and spatiotemporal ($\epsilon$) variation:

$$ \text{Log (Density}(s_i, t_i)) = \beta(t_i) + \omega(s_i) \sum_{f=1}^{N_f} \lambda(t_i, f) \varepsilon(s_i, f) $$

where $s_i$ and $t_i$ are, respectively, the location and year associated with sample i, whereas $N_f$ is the number of estimated modes of spatiotemporal variability. Finally, the $\lambda(t_i, f)$ estimate indicates whether a given year $t$ has a positive phase $\lambda(t_i, f) > 0$ or a negative phase ($\lambda(t_i, f) < 0$) during the positive phase of mode f. The $\varepsilon(s_i, f)$ estimate provides the map associated with the time series $\lambda(t_i, f)$ and represents whether a given location $s$ has a positive or negative value. The number of modes of spatiotemporal variability was prespecified at 2 ($N_f = 2$) in this study, as preliminary sensitivity analyses showed that by increasing their number $f$, the explanatory power of the model only increased slightly (Supplementary materials S2).
Red mullet density data were fit with a Poisson-link delta-
gamma distribution (Thorson 2018), which is appropriate
for zero-inflated biological data where two linear predic-
tors are estimated, the product of which gives the depen-
dent variable of interest (for further details, see Grüss et
al. (2021)). This Poisson-link delta-gamma model specifies a
Bernoulli distribution using a complementary log-log link
for encounter/non-encounters (i.e., the probability mass at
zero for a delta-model), and simultaneously specifies that pos-
itive catches follow a gamma distribution (the probability dis-
tribution for nonzeros in a delta model) (Thorson et al. 2021).

Following Thorson et al. (2020a), our modelling approach
fixes λf(t) = 0 for all f > t, with the additional constraint ap-
plied in Grüss et al. (2021) that all indices have a sum of zero,
i.e., $\sum_{t=0}^{N_t} \lambda(t, f) \epsilon = 0$, where $N_t$ is the number of years in the
time series of interest, such that an “average” year (defined as
a hypothetical year $t*$ when $\sum_{t*=0}^{N_t} \lambda(t*, f) = 0$, for all factors)
has spatial distribution $\omega$. Then, we rotate the results so that
they are directly interpretable, as is commonly done for prin-
cipal components analysis (PCA) (Zuur et al. 2003a; Thorson
et al. 2015, 2020a). A rotation matrix is defined such that $\Lambda$
has columns identical to the eigenvectors of $\Lambda^t\Lambda$, and then
define $\Lambda$ as the factor index and $\Lambda\epsilon$ as the map associated
with it (Thorson et al. 2020a). This “PCA rotation” maximizes
the variance for each axis in sequential order (Thorson 2019a;
Thorson et al. 2020a).

The model parameters were estimated using R package
"VAST" (Thorson 2019b) release number 3.10.0, with code
and associated materials publicly available online at https://gi
thub.com/James-Thorson-NOAA/VAST. VAST estimates spatial
variables for $N_x$ "knots" for computational efficiency (Shelton
et al. 2014) and employs bilinear interpolation to obtain
model predictions between knot locations. For both models
performed in this study and after carrying out some prelim-
inary sensitivity tests (Supplementary Materials S3), we used
$N_x$ = 100 knots that were uniformly distributed over the two
management units and predicted response variables across
2500 and 1080 grid cells for GSA06 and GSA07, respectively.

After having calculated the annual predicted densities and
the static long-term spatial pattern (ω), we computed the 95%
quantiles of density in GSA06 and GSA07 and plotted them
onto the areas that were consistently characterized by high
densities. This was done to further investigate the annual
fluctuations of red mullet density, better locate and visualize
the areas of high aggregations, and, ultimately, examine the
environmental variables that may drive these fluctuations
(i.e., contraction or expansion of persistent hotspot areas).

For each fitted model and to detect environmental pro-
cesses that may be associated with the spatiotemporal vari-
ability in the two areas, we fit an autoregressive model (AR)
to the first mode of variability in the GSA06 (that explained
nearly all variation for GSA06; see Section 3.1) and the first
and second modes of variability in the GSA07. Specifically, a
lag-1, AR(1), model was fitted in each case since the modes
of variability presented high first-order autocorrelation (not
shown). To examine the linkage between the environment
and the modes of variability (Grüss et al. 2021), all envi-
ronmental covariates were incorporated separately in each
model. The best model among all models with one environ-

mental covariate was selected based on Akaike's informa-
tion criterion (AIC). This process was repeated for lag-1 and
lag-2 years as well as for those environmental variables that
were only available for a subset of the time series examined
(Table 1). To test if the annual fluctuations in consistently
high-density areas were environmentally driven, we also per-
formed Pearson's correlation tests between environmental
variables and the number of cells with an annual density
higher than the 95% quantile. If for a given environmental
variable, the Pearson's *r* was high in absolute value (>0.4),
it was concluded that this variable has an influence on the
mode of variability (Grüss et al. 2021). As the data on chl-a
were not available during 1994–1997 (Table 1), to test for sig-
nificant differences between two *r* that derived from these
different-length datasets, we followed the methodology of
(Diedenhofen and Musch 2015) using the “cocor” package
they developed. All Pearson’s correlations were carried out
in lag-0, lag-1, and lag-2 years.

2.3.2. Dynamic factor analysis

To identify underlying common trends among LPUE time
series and to explore the effects of environmental variables
on these long-term trends, a dynamic factor analysis (DFA)
was implemented (Zuur et al. 2003b). This is a dimension
reduction technique where a set of *n* observed time series
is modelled as a linear combination of *m* common trends
(where *m* ≠ *n*), factor loadings, and error terms to explain tem-
poral variability. Factor loadings are used to detect the associ-
ation between the time series (ports) and the obtained trends.
Covariates can optionally be included in the model. DFA is
suitable for non-stationary datasets with missing values (Zuur
et al. 2003b). The correlation of observation errors can be
modelled using different error matrices: (i) same variance and
no covariance (diagonal-equal); (ii) different variances and no
covariance (diagonal-unequal); (iii) same variance and covari-
ance (equalvarcov); and (iv) different variances and covari-
ances (unconstrained) (Keller et al. 2017; Holmes et al. 2021).
The correlations of observation errors were fitted to all pos-
sible model structures in the time series, including 1, 2, and
3 common trends.

Model selection was based on the standard correction
to Akaike's information criterion (AICc) as a measure of
goodness-of-fit where the model with the lowest AICc was
considered to be the best (Zuur et al. 2003b). Standard errors
and confidence intervals of the regression were calculated us-
ing the Hessian matrix approach, whereas significance was
assessed using a t-test derived from standard errors. DFA was
realized in the multivariate autoregressive state-space mod-
elling (MARSS) package developed for R software (Holmes et
al. 2012).

In this study, we used monthly time-series data that can
be decomposed into seasonal, long-term, and residual varia-
tions. The decomposition was realized in base-R with a time-
series analysis using the functions "ts" and "decompose". It
is thus important to clarify that the seasonal component
was removed from the time series since preliminary analyses
showed no differentiation in the patterns of LPUE seasonal
trends in the GSA06 (Supplementary materials S4). Data were
standardized to the mean zero and variance 1, while DFA was
applied to an additive combination of long-term and residual
components (Holmes et al. 2021).

3. Results

**3.1. Estimating modes of variability**

For the Spanish Coast, the first mode of the spatiotemporal
variability (ε) explained 93.9% of the total variance and repre-
sented a multi-decadal trend (i.e., the values of the mode were
positive from 1994 to 2004 and negative thereafter), whereas
the second mode explained 6.1% and was characterized by
higher interannual variability (Fig. 2). We also generated spa-
tial maps that are associated with each of the two modes of
variability in dynamic density hotspots. These areas, which
were closely associated with the positive phase of the mode,
included the waters off the Catalan Coast and the coastal wa-
ters at the south of the Valencia Channel. Conversely, two ar-
eas in the southern part of the Spanish Coast were highly as-
sociated with the negative phase of the mode. Finally, persis-
tent hotspots based on the long-term spatial prediction map
(i.e., spatial variation, ω) were also detected, with the most
prominent found off the Ebro Delta area, where the conti-
nental shelf is wider.

For the Gulf of Lions, the two modes of spatiotemporal vari-
ability explained 67.5% and 32.5% of the total variance (Fig. 3).
The first mode of variability for the Gulf of Lions—as it was
the case for the first mode for the Spanish Coast—represents
a multi-decadal trend (Fig. 3). The areas associated with the
positive phase of the first mode were the south-eastern part
of the inner Gulf, while a south-eastern extension of the same
area was also found to be associated with the positive phase
of the second mode of spatiotemporal variability. The coastal
waters at the western limit of the Gulf of Lions were associ-
ated with the negative phase of this second mode. The long-
term spatial average predicted by the model revealed two per-
sistent hotspots of abundance in the Gulf of Lions, one at the
offshore area of the westernmost edge and one at the north-
western limit of the wide coastal plain (Fig. 3).

The two models predicted the spatial density on an an-
nual basis for the two examined regions, where for the Span-
ish Coast, a high-density area off the Ebro Delta prevails ev-
ery year, with smaller and more dynamic high-density ar-
eas off the Catalan Coast and at the southern tip of the
GSA06 emerging in certain years (Supplementary materials
S5). Regarding the Gulf of Lions, the two persistent hotspots
identified by the long-term spatial average appear every
year throughout the time series examined, with a lower
contribution of the dynamic components in the prediction
maps.

To assess temporal variation of contraction/expansion over
time of persistent high/aggregation areas, we retrieved from
the long-term spatial prediction map (i.e., spatial variation, ω)
the cells with a density higher than the 95% quantile of total
density in the Spanish Coast and the Gulf of Lions and then
summarized the number of these cells on the grid as an indi-
cator of the extent of these areas. We show that the red mul-

let density on the Spanish Coast has been steadily high since
2012, while the years 2008 and 2011 exhibit minimal relative
density (Figs. 4a and 4b). Likewise, in the Gulf of Lions, we
observe two periods of high density, during 2006–2010 and
2014–2019 (Figs. 4c and 4d).

**3.2. Links to environmental variables**

The first mode of variability for the Spanish Coast was
highly related to the AMO on lag-1, the MLD in the greater
area around MEDOC, and the Ebro runoff point on lag-2 years.
When considering data only from 1998, no significant rela-
tionship with chl-a was identified. For the Gulf of Lions, the
first mode was highly associated with the winter SST and
the AMO index on lag-0 years. The second mode was asso-
ciated with the MLD in the MEDOC area on lag-1 years and
the Western Mediterranean Oscillation index on lag-2 years
(Table 2). No significant effect of chl-a data was detected in
models fit in the shorter time series (from 1998 to 2019)
for which chl-a data were available. Similar results were ob-
tained with Pearson's correlation tests when not account-
ing for first-order autocorrelation (Table S6). It is noteworthy,
though, that more environmental variables were found to be
significantly correlated with the modes of variability when
applying Pearson's correlations (complementarily and not ac-
counting for autocorrelation) rather than when we fit AR
models.

With regard to the drivers of contraction/expansion
hotspots, Pearson's correlation tests between the size of the
persistent areas and environmental variables revealed mod-
erate but significant links (*p* < 0.05) between the density on
the Spanish Coast and the WeMOi and NAO lag-1 indices. In-
terestingly, the highest correlation for the Spanish Coast was
detected to be with the lag-0 winter SST of the Gulf of Li-
ons (*r* = 0.57, *p* < 0.05). Finally, the only variable that was
found to be associated with the annual variations of den-
sity in the Gulf of Lions was the lag-1 AMO index (*r* = 0.43,
*p* < 0.05). Besides the four aforementioned drivers, all re-
maining environmental variables presented no significant
correlations.

**3.3. DFA on LPUE data**

After we had examined all possible model structures re-
garding the number of trends, the variance, and the error
matrices, we found that the best fit for the long-term trends
of LPUE included models with three trends (m = 3) (Fig. 5a)
and an unconstrained matrix (different variance and covari-
ance) (Table 3). The incorporation of any of the environmental
variables did not improve the model (not shown). The first
trend exhibits a steep increase in LPUE from 2006, peaking
in 2008, followed by a 10-year period of relatively low fluctu-
ations that ended with a decline in 2018. Ports around the
Ebro Delta were associated with the first trend, for which
they presented negative factor loadings (Fig. 5b). The second
trend reveals a steady upward tendency in LPUE starting in
2008, while the third trend shows a decline from 2005 un-
til 2010, followed by a relatively constant period of low LPUE
that ended in 2015 with a slow increase (Fig. 5b). Ports along
the Catalan Coast were associated with the second and third
**Fig. 2.** Portrayal of (a, c) the two dominant modes (factors) of density spatiotemporal variability for the Spanish Coast (GSA06), (b, d) the spatial maps associated with the two modes, and (e) the spatial variation ω, which depicts the long-term average spatial pattern of red mullet density. For each axis (a, c), the proportion of explained variance is indicated. In panels b and d, yellow and blue tones represent areas that are associated with the positive and negative phases of the mode of variability, respectively. In panel e, yellow and blue tones represent high- and low-density areas. Map projection is WGS84 UTM.

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

trends and presented substantial variability with negative
and positive factor loadings for either of the trends. Finally,
the four southernmost landing ports in our study area, south
of the Valencia Channel, were associated with the first and

second trends, with negative and positive factor loadings, re-
spectively (Fig. 5b). A preliminary analysis of the seasonal and
residual variation showed no differentiation in the seasonal
trends of LPUE across the Spanish Coast.
**Fig. 3.** Portrayal of (a, c) the two dominant modes (factors) of density spatiotemporal variability for the Gulf of Lions (GSA07), (b, d) the spatial maps associated with the positive phase of the two modes, and (e) the spatial variation ω, which depicts the long-term average spatial pattern of red mullet density. For each axis (a, c), the proportion of explained variance is indicated. In panels b and d, yellow and blue tones represent areas that are associated with the positive and negative phases of the mode of variability, respectively. In panel e, yellow and blue tones represent high- and low-density areas. Map projection is WGS84 UTM.

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

**4. Discussion**

In this study, we examined spatiotemporal variation in the distributions and densities of a commercially important species to identify environmental drivers of population structure and assess the potential utility of implementing a spatially explicit stock assessment. We detected persistent and dynamic aggregations of a red mullet population from fisheries-independent data and explored whether these aggregations are linked to spatial patterns derived from fisheries-dependent data of landings in the NW Mediterranean Sea. To do so, we applied a combination of a generalized EOF analysis (Thorson et al. 2020a; Grüss et al. 2021) and a DFA (Zuur et al. 2003b) to complement fine spatiotemporal variation at the annual scale with seasonal dynamics along the geographic gradient. Our study also informs us of the effect of latent environmental processes on the spatial biocomplexity of the red mullet population system in the NW Mediterranean Sea. The NW Mediterranean open-ocean convection was shown to influence the dynamics of both management units considered in this study, GSA06 and GSA07. Model predictions indicated that red mullet populations along the Spanish coast were equally characterized by dynamic and persistent aggregations, whereas those in the Gulf of Lions primarily consisted of persistent aggregations.

Regarding landings data, three trends were detected in the long-term monthly LPUE of the Spanish Coast, with a primary segregation of the ports that was closely related to the persistent population unit off the Ebro Delta.

Identifying persistent areas of high density for marine organisms is of paramount importance for conservation and management, as these could be candidate areas for the implementation of specific measures such as spatiotemporal closures through the detection of essential fish habitats (STECF 2022) or the designation of marine protected areas (Agardy 2000; Erisman et al. 2017). On the Spanish Coast, we found that the continental shelf off the Ebro Delta is the main high-density and high-persistence area for the red mullet, which is consistent with previous studies (e.g., Paradinas et al. 2020). Discharges from the Ebro River and the prevailing Northern current jointly increase the primary productivity around the Delta area (Estrada 1996), which could then favour the productivity at higher trophic levels (Donoso et al. 2017). Since Coll et al. (2016) pointed out that the greater Ebro Delta is an important biomass area for different species, recent studies have evidenced that this area is of central importance for a variety of demersal resources in the NW Mediterranean (Vilas et al. 2020; Paradinas et al. 2022), including the red mullet (Paradinas et al. 2020). As such, the Ebro Delta may be
<figure><img src="image_7.png" /><figcaption>**Table 2.** Summary of the best models with a significant environmental effect for each mode of spatiotemporal variability, geographic sub-area (GSA), and lag-k years, where k = 0, 1, 2.</figcaption></figure>

<table><thead><tr><th>Mode of spatiotemporal variability</th><th>Environmental variables of the best model</th><th>Lag</th><th>AIC</th><th>Number of years (total number of fits)</th></tr></thead><tbody><tr><td rowspan="3">GSA06</td><td>1</td><td>AMO (+)</td><td>1</td><td>57.41</td><td>25 (16)</td></tr><tr><td>1</td><td>Ebro runoff (-)</td><td>2</td><td>56.56</td><td>24 (16)</td></tr><tr><td>1</td><td>MLD (Dec-Feb) (+)</td><td>2</td><td>57.99</td><td>24 (16)</td></tr><tr><td rowspan="4">GSA07</td><td>1</td><td>SST_GSA07 (Dec-Feb) (-)</td><td>0</td><td>48.61</td><td>26 (12)</td></tr><tr><td>1</td><td>AMO (-)</td><td>0</td><td>50.04</td><td>26 (12)</td></tr><tr><td>2</td><td>MLD (Dec-Feb) (+)</td><td>1</td><td>41.41</td><td>25 (12)</td></tr><tr><td>2</td><td>WeMOi (-)</td><td>2</td><td>40.08</td><td>24 (12)</td></tr></tbody></table>

Note: GSA06, geographic sub-area 06; GSA07, geographic sub-area 07; AIC, Akaike information criterion; AMO, Atlantic Multi-decadal Oscillation; MLD, mixed-layer depth; SST, sea surface temperature; WeMOi, Western Mediterranean Oscillation index; (+), positive effect; (-), negative effect. In the case where the difference in AIC between the best two models was below 2 (ΔAIC < 2), the second-best model is also reported in the table. The direction of the effect is also provided. The number of years of the data that was used in each model as well as the total number of fitted models depending on the number of environmental covariates for each GSA are shown in the right column.
<figure><img src="image_8.png" /><figcaption>**Fig. 5.** Results of the best model of dynamic factor analysis (DFA) for landings per unit effort (LPUE) time series in the Spanish Coast. (a) Common trends, (b) Factor loadings.</figcaption></figure>

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

viewed as a demographic engine for regional persistence and population dynamics in the Western Mediterranean (Kerametsidis et al. 2023). In the Gulf of Lions, by contrast, two areas with persistent aggregations of red mullet were identified. The Gulf as a whole is indeed one of the most productive areas of the entire Mediterranean Sea due to complex coupled ocean-atmosphere processes (Millot 1990; Bosc et al. 2004). One of the two persistent high-density areas was located on the continental shelf off the westernmost part of the gulf (Fig. 3), which is consistent with the very high densities of red mullet reported by (Morfin et al. 2016).

The second persistent area is located in the eastern part of the Inner Gulf and could also be related to the Rhône River discharges (Vandenbulcke and Barth 2019), a prominent hydrographical feature that greatly affects productivity in the area (Ulses et al. 2008) similar to the Ebro River on the Spanish Coast.

We identified dynamic, high-density areas in both management units that are equally important as the persistent ones that were described above. These areas can also be used to explore the link between environmental and biological components of the ecosystem. A notable result of our study
<table><caption>**Table 3.** Summary of the models with different combinations of the covariance matrix structure of observation errors (R) and the number of trends (m).</caption><thead><tr><th>R</th><th>m</th><th>AICc</th></tr></thead><tbody><tr><td>Unconstrained</td><td>3</td><td>8760</td></tr><tr><td>Unconstrained</td><td>2</td><td>8841</td></tr><tr><td>Unconstrained</td><td>1</td><td>8934</td></tr><tr><td>Diagonal and unequal</td><td>3</td><td>9044</td></tr><tr><td>Diagonal and unequal</td><td>2</td><td>9179</td></tr><tr><td>Equalvarcov</td><td>3</td><td>9184</td></tr><tr><td>Diagonal and equal</td><td>3</td><td>9229</td></tr><tr><td>Equalvarcov</td><td>2</td><td>9296</td></tr><tr><td>Diagonal and equal</td><td>2</td><td>9349</td></tr><tr><td>Diagonal and unequal</td><td>1</td><td>9415</td></tr><tr><td>Equalvarcov</td><td>1</td><td>9461</td></tr><tr><td>Diagonal and equal</td><td>1</td><td>9513</td></tr></tbody></table>

Note: The models are presented in ascending order, based on the correction to the Akaike information criterion (AICc).

that holds true for both examined management units was the high correlation detected between the spatiotemporal variability of red mullet density and the NW Mediterranean open-ocean convection. Despite the latter being a very influential oceanographic phenomenon in the Western Mediterranean, greatly affecting primary productivity (Heimbürger-Boavida et al. 2013), zooplankton communities (Donoso et al. 2017), and even small-pelagic stocks (Feuilloley et al. 2020), its relation to demersal fish stock dynamics has to date only been demonstrated for two species in the area, hake and blue whiting (Massutí et al. 2008; Martín et al. 2012; Hidalgo et al. 2019a). Since deep convection strongly drives zooplankton dynamics, which represent a large portion of the red mullet diet over different ontogenetic stages (Bautista-Vega et al. 2008; Esposito et al. 2014), it is not surprising that we detected significant correlations between convective events and red mullet density in both management units. In addition, NW Mediterranean open-ocean convection events dramatically affect primary productivity dynamics in the region (Estrada 1996; Sabatés et al. 2015). Sabatés et al. (2015) maintained that chl-a concentrations are indirectly associated with the distribution and feeding dynamics of red mullet larvae, mediated via the fluctuations of one of its main preys, the cladocera. This is in line with our findings regarding the critical effect of open-ocean convection.

The AMO has been shown to play an important role in demersal fish and invertebrate communities (Nye et al. 2014; Zimmermann et al. 2019). In the Mediterranean Sea, small pelagic fish communities (Tsikliras et al. 2019) have also been associated with AMO, with the western basin being particularly affected. Our results suggest that AMO played a significant role in shaping red mullet dynamics across a wide geographical gradient in the NW Mediterranean Sea. Notably, AMO emerged as the predominant process affecting the extent of persistently high-density areas in the Gulf of Lions. Furthermore, another regional climatic index, the WeMO index, was associated with the spatiotemporal variability of red mullet in the Gulf of Lions. To our knowledge, this is

the first study to document the effect of this process on the red mullet, although the effect of WeMO on another important demersal species, the hake, has been previously documented on the Spanish Coast (Martín et al. 2012; Ordines et al. 2019). The relationship between red mullet density variations and two broad-scale coupled atmosphere-ocean climatic indices documented here underlines the utility of such synthesized indices to assess population dynamics on local scales as, in certain cases, they can be better tools to predict ecological processes than are local environmental variables (Hallett et al. 2004; Stenseth and Mysterud 2005). However, we must acknowledge that the first modes of variability for both areas display long-term trends and low-frequency cycles that could also be due to the recruitment spread among age classes, generating autocovariance between cohorts (Fromentin and Fonteneau 2001) and inter- or intra-cohort interactions (Bjørnstad et al. 1999). Also, the periods of climatic oscillations (e.g., >50 years for AMO) are often larger than the information available. These factors often limit the mechanistic understanding between climatic indices and ecological responses.

Our study also revealed the influence of other local environmental drivers. SST is a variable that is routinely used as an environmental driver for stock dynamics as it affects fish directly (e.g., via positive relationships with physiological processes such as metabolism) and indirectly (e.g., via fluctuations in prey availability) (Lloret et al. 2014; Brosset et al. 2015). Here, we found that winter SST is the primary driver of red mullet dynamics in the Gulf of Lions. This is in accordance with Kerametsidis et al. (2023), who documented that higher winter SST favours the recruitment success of red mullet in the Western Mediterranean. Furthermore, Levi et al. (2003) showed that SST anomalies significantly influence red mullet in the Central Mediterranean Sea. The density and biomass of demersal communities in the NW Mediterranean have been previously shown to be positively impacted by higher SST (Vilas et al. 2020).

The monthly LPUE data from commercial fisheries allowed for the conclusion that the intra-annual dynamics of the species and the supporting fishery are consistent with the spatial structure and dynamics of the sub-units identified here. The DFA performed on LPUE from ports along the Northern Spanish Coast revealed that the long-term variability can be explained by three underlying trends. No environmental influence on either of the trends was detected, as none of the examined environmental processes improved the performance of the model. Since DFA was performed on fishery-dependent data, it is very likely that other non-environmental (e.g., socio-economic factors), such as market-related variables, or intrinsic stock characteristics, may determine such trends (Bjørnstad et al. 1999; Aguilera et al. 2015). These trends show that the ports around the Ebro Delta were clustered together, separately from the remainders on the Catalan coast and the south of the Valencia Channel. This once again underlines the importance of the Ebro Delta for red mullet and its partially independent dynamics. This is also in accordance with EOFs that identified the area off the Ebro Delta as a persistent hotspot for red mullet density as well as for a variety of other species, as previously explained.
Non-environmental variables were not included in this study, as our primary objective was to identify the environmental components that drive the density variability of red mullet in the NW Mediterranean Sea across different scales. Although red mullet has constantly been harvested unsustainably (Colloca et al. 2017), obtaining reliable information on fishing mortality on a local level (i.e., inferior to the management unit level) proved challenging, so we opted not to incorporate this information into our study. Additionally, we refrained from including habitat degradation, as this has been shown to be relevant to the organisms of meiobenthos and superbenthos (Coll et al. 2010). These, in addition to other variables such as pollution, trophic interactions, etc. could be explored in future cross-scale studies and/or under various climate change scenarios. This would be particularly important for the trends in LPUE data, which seemed not to be related to any of the candidate environmental drivers that we explored.

**5. Conclusion, broad implications, and future research**

Here we demonstrate the existence of both persistent and dynamic high-density areas for a well-studied high-value species across two adjacent management units of the western Mediterranean Sea. Persistent dynamic aggregations were associated with different environmental drivers that varied by spatial scale and temporal dynamics. This shows that the population of the species in the examined region is spatially segregated among and within the two management units. This suggests a metapopulation structure for this species in the area (Cadrin and Dickey-Collas 2014), which should reasonably be considered in its assessment. Additional studies applying other stock identification techniques could provide conclusive evidence on the dynamics of such a metapopulation structure and on the (bi)directional connections among the sub-units. Thus, the combination of various tools (i.e., generalized EOF and DFA) applied to complementary datasets (i.e., CPUE and LPUE) to explore and quantify demographic connectivity at both local and regional scales is strongly encouraged for future studies. Some of these techniques might include tagging, analyses of otolith shape and microchemistry, white muscle stable isotope analysis, morphometrics, or population dynamics simulations assuming complex population structures (Goethel et al. 2011, 2021; ICES 2022). This study is the first on non-gadoid species, supporting a possible generalization of the impacts of open-ocean convective events on the local and/or regional population dynamics in the Mediterranean Sea. Taking into account that convective events outside the Mediterranean Sea have been associated directly and indirectly with salinity and primary productivity variations (e.g., Ferrari et al. 2015; Chan et al. 2017; Lowry et al. 2018), it could be of elevated importance to investigate the link between open-ocean convection and stock dynamics in areas with considerable convective events, such as the Greenland and Labrador Seas (Marshall and Schott 1999).

As we revealed persistent, high-density areas, our findings could be used as a basis for the development of area-

based fisheries management measures. Area-based management schemes can help achieve maximum sustainable yield (MSY) and maximize socio-economic benefits over the long term for the fishers (STECF 2022), and they are regarded as an effective framework to achieve the sustainability of marine resources. In our case study, for instance, one of the persistent areas that was identified for red mullet is the Ebro Delta, an area that has already been highlighted as important for other exploited demersal species (Vilas et al. 2020). It is noteworthy that time closures, as an area-based measure, have recently been enforced within the Ebro Delta area (STECF 2021). In addition, understanding the spatial structure of marine stocks is vital for effective management. The segregation of different subpopulations within the same management unit indicates that the implementation of spatial stock assessment frameworks (Goethel et al. 2011, 2023; Punt 2019) might be more suitable to inform the management of red mullet stock in the region. Finally, this study shows that besides spatial biocomplexity in marine stocks, it is equally crucial that other ecosystem components, such as environmental processes, be integrated into assessment frameworks in the context of EBFM.

**Acknowledgements**

We thank all the dedicated people of the Instituto Español de Oceanografía (IEO, Spain), who have worked on the research vessels and MEDITS surveys. Additionally, we would like to thank Dr. Patricia Puerta for her insightful comments during the implementation of DFA. Finally, we thank the anonymous reviewers for their careful reading of our manuscript and their highly useful remarks and suggestions.

**Article information**

**History dates**

Received: 2 May 2023

Accepted: 25 September 2023

Accepted manuscript online: 1 December 2023

Version of record online: 12 February 2024

**Copyright**

©2024 Authors Kerametsidis, Rossi, Álvarez-Berastegui, Barnes, Certain, Esteban, García, Piñeiro, Vivas, Hidalgo, and MARBEC, Univ Montpellier, CNRS, Ifremer, IRD, Sète, France. This work is licensed under a [Creative Commons Attribution 4.0 International License](https://creativecommons.org/licenses/by/4.0/) (CC BY 4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author(s) and source are credited.

**Data availability**

The data analyzed during this study are available from the corresponding author upon reasonable request.

**Author information**

**Author ORCIDs**

Georgios Kerametsidis https://orcid.org/0000-0002-6502-8122

James Thorson https://orcid.org/0000-0001-7415-1010
<figure><img src="image_10.png" /></figure>

Vincent Rossi https://orcid.org/0000-0001-7291-0415

Diego Álvarez-Berastegui https://orcid.org/
0000-0003-1304-7185

Cheryl Barnes https://orcid.org/0000-0002-2650-5299

Gregoire Certain https://orcid.org/0000-0002-5242-5268

Antonio Esteban https://orcid.org/0000-0002-2896-7972

Encarnación García https://orcid.org/0000-0001-5179-0098

Angélique Jadaud https://orcid.org/0000-0001-6858-3570

Safo Piñeiro https://orcid.org/0000-0003-0142-2398

Miguel Vivas https://orcid.org/0000-0003-4684-8202

Manuel Hidalgo https://orcid.org/0000-0002-3494-9658

## Author contributions

**Conceptualization:** MH

**Data curation:** GK, JT, VR, DAB, CB, GC, AE, EG, AJ, SP, MV, MH

**Formal analysis:** GK, JT, CB, MH

**Funding acquisition:** GK, AE, EG, MV, MH

**Investigation:** GK, MH

**Methodology:** GK, JT, CB, MH

**Project administration:** GK, MH

**Resources:** GK, JT, VR, DAB, CB, GC, AE, EG, AJ, SP, MV, MH

**Supervision:** MH

**Validation:** GK, JT, CB, MH

**Visualization:** GK, JT, GC, CB

**Writing – original draft preparation:** GK, MH

**Writing – review and editing:** GK, JT, VR, DAB, CB, GC, AE, EG, AJ, SP, MV, MH

All authors have read and agreed to the submitted version of the manuscript.

## Competing interests

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Funding information

MEDITS surveys are funded by IEO and the European Union (EU) through the European Maritime and Fisheries Fund (EMFF) within the National Program of collection, management, and use of data in the fisheries sector and support for scientific advice regarding the Common Fisheries Policy. GK, AE, EG, MV, and MH acknowledge funding from the COCOCHA project (Grant PID2019-110282RA-100 funded by MCIN/AEI/10.13039/501100011033). GK acknowledges funding from the Grant for PhD candidates (FPI) PRE2020-094321 funded by MCIN/AEI/10.13039/501100011033 by the "European Union NextGenerationEU/PRTR" and co-funded by the European Social Fund.

## Supplementary material

Supplementary data are available with the article at https://doi.org/10.1139/cjfas-2023-0120.

## References

* Agardy, T. 2000. Information needs for marine protected areas: scientific and societal. Bull. Mar. Sci. **66**: 875–888

* Aguilera, S.E., Cole, J., Finkbeiner, E.M., Le Cornu, E., Ban, N.C., Carr, M.H., et al. 2015. Managing small-scale commercial fisheries for adaptive capacity: insights from dynamic social-ecological drivers of change in monterey bay. PLoS One, **10**(3). doi:10.1371/journal.pone.0118992.

* Anonymous. 2017. MEDITS-handbook. Version n. 9, 2017. MEDITS Working Group. 106p.

* Bautista-Vega, A.A., Letourneur, Y., Harmelin-Vivien, M., and Salen-Picard, C. 2008. Difference in diet and size-related trophic level in two sympatric fish species, the red mullets *Mullus barbatus* and *Mullus surmuletus*, in the Gulf of Lions (north-west Mediterranean Sea). J. Fish Biol. **73**(10): 2402–2420. doi:10.1111/j.1095-8649.2008.02093.x.

* Bjørnstad, O.N., Ims, R.A., and Lambin, X. 1999. Spatial population dynamics: analyzing patterns and processes of population synchrony. Trends Ecol. Evol. **14**(11): 427–432.

* Börgel, F., Frauen, C., Neumann, T., and Meier, H.E.M. 2020. The Atlantic Multidecadal Oscillation controls the impact of the North Atlantic Oscillation on North European climate. Environ. Res. Lett. **15**(10). doi:10.1088/1748-9326/aba925.

* Bosc, E., Bricaud, A., and Antoine, D. 2004. Seasonal and interannual variability in algal biomass and primary production in the Mediterranean Sea, as derived from 4 years of SeaWiFS observations. Global Biogeochem. Cycles, **18**(1). doi:10.1029/2003gb002034.

* Bowler, D.E., Hof, C., Haase, P., Kröncke, I., Schweiger, O., Adrian, R., et al. 2017. Cross-realm assessment of climate change impacts on species' abundance trends. Nat. Ecol. Evol. **1**(3). doi:10.1038/s41559-016-0067-7.

* Brosset, P., Ménard, F., Fromentin, J.M., Bonhommeau, S., Ulses, C., Bourdeix, J.H., et al. 2015. Influence of environmental variability and age on the body condition of small pelagic fish in the Gulf of Lions. Mar. Ecol. Prog. Ser. **529**: 219–231. doi:10.3354/meps11275.

* Cadrin, S.X. 2020. Defining spatial structure for fishery stock assessment. Fish Res. **221**: doi:10.1016/j.fishres.2019.105397.

* Cadrin, S.X., and Dickey-Collas, M. 2014. Stock assessment methods for sustainable fisheries. ICES J. Mar. Sci. 1–6. doi:10.1093/icesjms/fsu228.

* Cardinale, M., Colloca, F., Bonanno, A., Scarcella, G., Arneri, E., Jadaud, A., et al. 2021. The Mediterranean fishery management: a call for shifting the current paradigm from duplication to synergy. Mar. Policy, **131**: doi:10.1016/j.marpol.2021.104612.

* Champion, C., Brodie, S., and Coleman, M.A. 2021. Climate-driven range shifts are rapid yet variable among recreationally important coastal-pelagic fishes. Front. Mar. Sci. **8**: doi:10.3389/fmars.2021.622299.

* Chan, P., Halfar, J., Adey, W., Hetzinger, S., Zack, T., Moore, G.W.K., et al. 2017. Multicentennial record of Labrador Sea primary productivity and sea-ice variability archived in coralline algal barium. Nat. Commun. **8**: doi:10.1038/ncomms15543.

* Clark, M., Macdiarmid, J., Jones, A.D., Ranganathan, J., Herrero, M., and Fanzo, J. 2020. The role of healthy diets in environmentally sustainable food systems. Food Nutr. Bull. **41**(Suppl. 2): 31S–58S. doi:10.1177/0379572120953734.

* Coll, M., Piroddi, C., Steenbeek, J., Kaschner, K., Lasram, F.B.R., Aguzzi, J., et al. 2010. The biodiversity of the Mediterranean Sea: estimates, patterns, and threats. doi:10.1371/journal.pone.0011842.

* Coll, M., Steenbeek, J., Sole, J., Palomera, I., and Christensen, V. 2016. Modelling the cumulative spatial-temporal effects of environmental drivers and fishing in a NW Mediterranean marine ecosystem. Ecol. Modell. **331**: 100–114. doi:10.1016/j.ecolmodel.2016.03.020.

* Colloca, F., Scarcella, G., and Libralato, S. 2017. Recent trends and impacts of fisheries exploitation on Mediterranean stocks and ecosystems. Front. Media S. A. doi:10.3389/fmars.2017.00244.

* Di Stefano, M., Legrand, T., Di Franco, A., Nerini, D., and Rossi, V. 2023. Insights into the spatio-temporal variability of spawning in a territo-rial coastal fish by combining observations, modelling and literature review. Fish. Oceanogr. **32**(1): 70–90. doi:10.1111/fog.12609.

* Diedenhofen, B., and Musch, J. 2015. Cocor: a comprehensive solution for the statistical comparison of correlations. PLoS One, **10**(4). doi:10.1371/journal.pone.0121945.

* Dijkstra, H.A., Te Raa, L., Schmeits, M., and Gerrits, J. 2006. On the physics of the Atlantic Multidecadal Oscillation. Ocean Dyn. **56**(1): 36–50. doi:10.1007/s10236-005-0043-0.

* Donoso, K., Carlotti, F., Pagano, M., Hunt, B.P.V., Escribano, R., and Berline, L. 2017. Zooplankton community response to the winter 2013 deep convection process in the NW Mediterranean.
ranean Sea. J. Geophys. Res. Oceans, **122**(3): 2319-2338. doi:10.1002/2016JC012176.

Erisman, B., Heyman, W., Kobara, S., Ezer, T., Pittman, S., Aburto-Oropeza, O., and Nemeth, R.S. 2017. Fish spawning aggregations: where well-placed management actions can yield big benefits for fisheries and conservation. Fish Fish. **18**(1): 128-144. doi:10.1111/faf. 12132.

Esposito, V., Andaloro, F., Bianca, D., Natalotto, A., Romeo, T., Scotti, G., and Castriota, L. 2014. Diet and prey selectivity of the red mullet, Mullus barbatus (Pisces: Mullidae), from the southern Tyrrhenian Sea: the role of the surf zone as a feeding ground. Mar. Biol. Res. **10**(2): 167-178. doi:10.1080/17451000.2013.797585.

Estrada, M. 1996. Primary production in the northwestern Mediterranean. In *The European anchovy and its environment*.

FAO. 2022a. The state of world fisheries and Aquaculture 2022. In *Towards blue transformation*. Rome, FAO. doi:10.4060/cc0461en.

FAO. 2022b. The state of Mediterranean and Black Sea Fisheries 2022. General Fisheries Commission for the Mediterranean, Rome. doi:10.4060/cc3370en.

Ferrari, R., Merrifield, S.T., and Taylor, J.R. 2015. Shutdown of convection triggers increase of surface chlorophyll. J. Mar. Syst. **147**: 116-122. doi:10.1016/j.jmarsys.2014.02.009.

Feuilloley, G., Fromentin, J.-M., Stemmann, L., Demarcq, H., Estournel, C., and Saraux, C. 2020. Concomitant changes in the Environment and small pelagic fish community of the Gulf of Lions. Prog. Oceanogr. **186**: 102375. doi:10.1016/j.pocean.2020.102375.

Fiorentino, F., Massutè, E., Tinti, F., Somarakis, S., Garofalo, G., Russo, T., et al. 2015. Stock units: identification of distinct biological units (stock units) for different fish and shellfish species and among different GFCM-GSA. STOCKMED Deliverable 03: Final report. January 2015, 310 p.

Fischer, W., Bauchot, M.L., and Schneider, M. 1987. Fiches FAO d'identification des espèces pour les besoins de la peche. (Révision 1). Mediterranée et Mer Noire. Zone de peche 37. 2. Vertébrés. Publication préparée par la FAO (Project GCP/INT/422/EEC). FAO: Rome. pp. 761-1530.

Fromentin, J.M., and Fonteneau, A. 2001. Fishing effects and life history traits: a case study comparing tropical versus temperate tunas. Fish. Res. **53**: 133-150. doi:10.1016/S0165-7836(00)00299-X.

Gargano, F., Garofalo, G., and Fiorentino, F. 2017. Exploring connectivity between spawning and nursery areas of Mullus barbatus (L., 1758) in the Mediterranean through a dispersal model. Fish Oceanogr. **26**(4): 476-497. doi:10.1111/fog.12210.

Gargano, F., Garofalo, G., Quattrocchi, F., and Fiorentino, F. 2022. Where do recruits come from? Backward lagrangian simulation for the deep water rose shrimps in the Central Mediterranean Sea. Fish Oceanogr. **31**(4): 369-383. doi:10.1111/fog.12582.

Goethel, D.R., and Berger, A.M. 2016. Accounting for spatial complexities in the calculation of biological reference points: effects of misdiagnosing population structure for stock status indicators. Can. J. Fish. Aquat. Sci. **74**(11): 1878-1894. doi:10.1139/cjfas-2016-0290.

Goethel, D.R., Berger, A.M., and Cadrin, S.X. 2023. Spatial awareness: good practices and pragmatic recommendations for developing spatially structured stock assessments. Fish Res. **264**. doi:10.1016/j.fishres.2023.106703.

Goethel, D.R., Bosley, K.M., Langseth, B.J., Deroba, J.J., Berger, A.M., Hanselman, D.H., and Schueller, A.M. 2021. Where do you think you're going? Accounting for ontogenetic and climate-induced movement in spatially stratified integrated population assessment models. Fish Fish. **22**(1): 141-160. doi:10.1111/faf.12510.

Goethel, D.R., Quinn, T.J., and Cadrin, S.X. 2011. Incorporating spatial structure in stock assessment: movement modeling in marine fish population dynamics. Rev. Fish. Sci. **19**(2): 119-136. doi:10.1080/10641262.2011.557451.

Grüss, A., Thorson, J.T., Stawitz, C.C., Reum, J.C.P., Rohan, S.K., and Barnes, C.L. 2021. Synthesis of interannual variability in spatial demographic processes supports the strong influence of cold-pool extent on eastern Bering Sea walleye pollock (*Gadus chalcogrammus*). Prog. Oceanogr. **194**. doi:10.1016/j.pocean.2021.102569.

Hallett, T. B., Coulson, T., Pilkington, J. G., Clutton-Brock, T. H., Pemberton, J. M., and Grenfell, B. T. 2004. Why large-scale climate indices seem to predict ecological processes better than local weather. Nature, **430**(6995): 71-75. doi:10.1038/nature02708.

Heimbürger-Boavida, L.-E., Lavigne, H., Migon, C., Estournel, C., Coppola, L., Miquel, J.-C., et al. 2013. Temporal variability of vertical export flux at the DYFAMED time-series station (Northwestern Mediterranean Sea). doi:10.1016/j.pocean.2013.08.005i.

Hidalgo, M., Ligas, A., Bellido, J.M., Bitetto, I., Carbonara, P., Carlucci, R., et al. 2019b. Size-dependent survival of European hake juveniles in the Mediterranean sea. Sci. Mar. **83**(S1): 207-221. doi:10.3989/scimar. 04857.16A.

Hidalgo, M., Rossi, V., Monroy, P., Ser-Giacomi, E., Hernández-García, E., Guijarro, B., et al. 2019a. Accounting for ocean connectivity and hydroclimate in fish recruitment fluctuations within transboundary metapopulations. Ecol. Appl. **29**(5). doi:10.1002/eap.1913.

Hiddink, J.G., and ter Hofstede, R. 2008. Climate induced increases in species richness of marine fishes. Global Change Biol. **14**(3): 453-460. doi:10.1111/j.1365-2486.2007.01518.x.

Holmes, E., Villanueva, N.M., Sest elo, M., and Zhang, W. 2012. MARSS: multivariate autoregressive state-space models for analyzing time-series data. R J. **4**(1): p.11. doi:10.32614/RJ-2012-002.

Holmes, E.E., Ward, E.J., Scheuerell, M.D., and Wills, K. 2021. MARSS: multivariate autoregressive state-space modeling. R package version 3.11.4. Available from https://CRAN.R-project.org/package=MARSS.

Houde, E.D. 2016. Recruitment variability. In Fish reproductive biology: implications for assessment and management. pp. 98-187. doi:10.1002/9781118752739.ch3.

Houpert, L., Durrieu de Madron, X., Testor, P., Bosse, A., D'Ortenzio, F., Bouin, M.N., et al. 2016. Observations of open-ocean deep convection in the northwestern Mediterranean Sea: seasonal and inter-annual variability of mixing and deep water masses for the 2007-2013 period. J. Geophys. Res. Oceans, **121**(11): 8139-8171. doi:10.1002/2016JC011857.

ICES. 2022. Stock identification methods Working group (SIMWG). ICES Sci. Rep. 4: 72. 66 p. doi:10.17895/ices.pub.20937001.

Keller, S., Quetglas, A., Puerta, P., Bitetto, I., Casciaro, L., Cuccu, D., et al. 2017. Environmentally driven synchronies of Mediterranean cephalopod populations. Prog. Oceanogr. **152**: 1-14. doi:10.1016/j.pocean. 2016.12.010.

Kerametsidis, G., Rueda, L., Bellido, J.M., Esteban, A., García, E., Gil de Sola, L., et al. 2023. The trade-off between condition and growth shapes juveniles' survival of harvested demersal fish of the Mediterranean sea. Mar. Environ. Res. **184**. doi:10.1016/j.marenvres.2022. 105844.

Kerr, L.A., Hintzen, N.T., Cadrin, S.X., Clausen, L.W., Dickey-Collas, M., Goethel, D.R., et al. 2017, Lessons learned from practical approaches to reconcile mismatches between biological population structure and stock units of marine fish. Oxford University Press. doi:10.1093/icesjms/fsw188.

Knight, J.R., Folland, C.K., and Scaife, A.A. 2006. Climate impacts of the Atlantic multic�edal oscillation. Geophys. Res. Lett. **33**(17). doi:10.1029/2006GL026242.

Lavigne, H., D'Ortenzio, F., Ribera D'Alcalà, M., Claustre, H., Sauzède, R., and Gacic, M. 2015. On the vertical distribution of the chlorophyll a concentration in the Mediterranean Sea: a basin-scale and seasonal approach. Biogeosciences, **12**(16): 5021-5039. doi:10.5194/ bg-12-5021-2015.

Levi, D., Andreoli, M.G., Bonanno, A., Fiorentino, F., Garofalo, G., Mazzola, S., et al. 2003. Embedding sea surface temperature anomalies into the stock recruitment relationship of red mullet (Mullus barbatus L. 1758) in the Strait of Sicily. Sci. Mar. **67**: 259-268. doi:10.3989/scimar.2003.67s1259.

Lloret, J., Shulman, G., and Love, R.M. 2014. Condition and health indicators of exploited marine fishes. Vol. **263**.

Lotze, H.K., Tittensor, D.P., Bryndum-Buchholz, A., Eddy, T.D., Cheung, W.W.L., Galbraith, E.D., et al. 2019. Global ensemble projections reveal trophic amplification of ocean biomass declines with climate change. Proc. Natl. Acad. Sci. U.S.A. **116**(26): 12907-12912. doi:10.1073/pnas.1900194116.

Lowry, K.E., Pickart, R.S., Selz, V., Mills, M.M., Pacini, A., Lewis, K.M., et al. 2018. Under-ice phytoplankton blooms inhibited by spring convective mixing in refreezing leads. J. Geophys. Res. Oceans, **123**(1): 90-109. doi:10.1002/2016JC012575.

Macias, D., Garcia-Gorriz, E., and Stips, A. 2018. Deep winter convection and phytoplankton dynamics in the NW Mediterranean Sea under
<figure><img src="image_11.png" /></figure>

present climate and future (horizon 2030) scenarios. Sci. Rep. 8(1). doi:10.1038/s41598-018-24965-0.

Marshall, J., and Schott, F. 1999. Open-ocean convection: observations, theory, and models. Rev. Geophys. 37(1): 1-64. doi:10.1029/98RG02739.

Marshall, K.N., Levin, P.S., Essington, T.E., Koehn, L.E., Anderson, L.G., Bundy, A., et al. 2018. Ecosystem-based fisheries management for social-ecological systems: renewing the focus in the United States with next generation fishery ecosystem plans. Conserv. Lett. 11(1). doi:10.1111/conl.12367.

Martin-Vide, J., and Lopez-Bustins, J.A. 2006. The Western Mediterranean oscillation and rainfall in the Iberian Peninsula. Int. J. Climatol. 26(11): 1455-1475. doi:10.1002/joc.1388.

Martin, P., Maynou, F., Recasens, L., and Sabatés, A. 2016. Cyclic fluctuations of blue whiting (Micromesistius poutassou) linked to open-sea convection processes in the northwestern Mediterranean. Fish. Oceanogr. 25(3): 229-240. doi:10.1111/fog.12147.

Martín, P., Sabatés, A., Lloret, J., and Martin-Vide, J. 2012. Climate modulation of fish populations: the role of the Western Mediterranean Oscillation (WeMO) in sardine (Sardina pilchardus) and anchovy (Engraulis encrasicolus) production in the north-western Mediterranean. Clim. Change, 110(3-4): 925-939. doi:10.1007/s10584-011-0091-z.

Massutí, E., Monserrat, S., Oliver, P., Moranta, J., López-Jurado, J.L., Marcos, M., et al. 2008. The influence of oceanographic scenarios on the population dynamics of demersal resources in the western Mediterranean: hypothesis for hake and red shrimp off Balearic Islands. J. Mar. Syst. 71(3-4): 421-438. doi:10.1016/j.jmarsys.2007.01.009.

Matić-Skoko, S., Šegvić-Bubić, T., Mandić, I., Izquierdo-Gomez, D., Arneri, E., Carbonara, P., et al. 2018. Evidence of subtle genetic structure in the sympatric species *Mullus barbatus* and *Mullus surmuletus* (Linnaeus, 1758) in the Mediterranean Sea. Sci. Rep. 8(1). doi:10.1038/s41598-017-18503-7.

MEDOC Group 1970. Observation of formation of deep water in the Mediterranean Sea, 1969. Nature, 227: 1037-1040. doi:10.1038/2271037a0.

Mertens, C., and Schott, F. 1998. Interannual variability of deep-water formation in the Northwestern Mediterranean. J. Phys. Oceanogr. 28(7): 1410-1424. doi:10.1175/1520-0485(1998)028%3c1410:IVODWF%3e2.0;CO;2.

Millot, C. 1990. The gulf of Lions' hydrodynamics. Continental Shelf Res. 10(9-11): 885-894. doi:10.1016/0278-4343(90)90065-T.

Molinero, J.C., Ibanez, F., Nival, P., Buecher, E., and Souissi, S. 2005. Molinero, Juan Carlos, Frédéric Ibanez, Paul Nival, Emmanuelle Buecher, and Sami Souissi. The North Atlantic climate and the northwestern Mediterranean plankton variability. Limnol. Oceanogr. 50(4): 2005, 1213-1220. Available from www.pml.ac.uk/. doi:10.4319/lo.2005.50.4.1213.

Morfin, M., Bez, N., and Fromentin, J.M. 2016. Habitats of ten demersal species in the Gulf of Lions and potential implications for spatial management. Mar. Ecol. Prog. Ser. 547: 219-232. doi:10.3354/meps11603.

Nye, J.A., Baker, M.R., Bell, R., Kenny, A., Kilbourne, K.H., Friedland, K.D., et al. 2014. Ecosystem effects of the Atlantic Multidecadal Oscillation. J. Mar. Syst. 133: 103-116. doi:10.1016/j.jmarsys.2013.02.006.

Ordines, F., Lloret, J., Tugores, P., Manfredi, C., Guijarro, B., Jadaud, A., et al. 2019. A new approach to recruitment overfishing diagnosis based on fish condition from survey data. Sci. Mar. 83(S1): 223-233. doi:10.3989/scimar.04950.03A.

Ottersen, G., Planque, B., Belgrano, A., Post, E., Reid, P.C., and Stenseth, N.C. 2001. Ecological effects of the North Atlantic Oscillation. Oecologia, 128: 1-14. doi:10.1007/s00442-0100655.

Paradinas, I., Conesa, D., López-Quílez, A., Esteban, A., López, L.M.M., Bellido, J.M., and Pennino, M.G. 2020. Assessing the spatiotemporal persistence of fish distributions: a case study on two red mullet species (*Mullus surmuletus* and *M. barbatus*) in the western mediterranean. Mar. Ecol. Prog. Ser. 644: 173-185. doi:10.3354/meps13366.

Paradinas, I., Giménez, J., Conesa, D., López-Quílez, A., and Pennino, M.G. 2022. Evidence for spatiotemporal shift in demersal fishery management priority areas in the western Mediterranean. Can. J. Fish. Aquat. Sci. doi:10.1139/cjfas-2021-0327.

Perry, A.L., Low, P.J., Ellis, J.R., and Reynolds, J.D. 2005. Climate change and distribution shifts in marine fishes. Available from https://www.science.org.

Puerta, P., Quetglas, A., and Hidalgo, M. 2016. Seasonal variability of cephalopod populations: a spatio-temporal approach in the Western Mediterranean Sea. Fish. Oceanogr. 25(4): 373-389. doi:10.1111/fog.12159.

Punt, A.E. 2019. Spatial stock assessment methods: a viewpoint on current issues and assumptions. Fish. Res. 213: 132-143. doi:10.1016/j.fishres.2019.01.014.

Punt, A.E. 2023. Those who fail to learn from history are condemned to repeat it: a perspective on current stock assessment good practices and the consequences of not following them. Fish. Res. 261. doi:10.1016/j.fishres.2023.106642.

Punt, A.E., Dunn, A., Elvarsson, B.P., Hampton, J., Hoyle, S., Maunder, M.N., et al. 2020. Essential features of the next-generation integrated fisheries Stock Assessment Package: a perspective. Fish. Res. 229, p. 105617. doi:10.1016/j.fishres.2020.105617.

Reiss, H., Hoarau, G., Dickey-Collas, M., and Wolff, W.J. 2009. Genetic population structure of marine fish: mismatch between biological and fisheries management units. Fish. Fish. 10(4): 361-395. doi:10.1111/j.1467-2979.2008.00324.x.

Rufener, M.C., Kristensen, K., Nielsen, J.R., and Bastardie, F. 2021. Bridging the gap between commercial fisheries and survey data to model the spatiotemporal dynamics of marine species. Ecol. Appl. 31(8). doi:10.1002/eap.2453.

Sabatés, A., Zaragoza, N., and Raya, V. 2015. Distribution and feeding dynamics of larval red mullet (*Mullus barbatus*) in the NW Mediterranean: the important role of cladocera. J. Plankton Res. 37(4): 820-833. doi:10.1093/planktf/fbv040.

Sánchez Lizaso, J.L., Goñi, R., Reñones, O., García Charton, J.A., Galzin, R., Bayle, J.T., et al. 2000. Density dependence in marine protected populations: a review. doi:10.1017/S0376892900000187.

Schlenker, L.S., Stewart, C., Rock, J., Heck, N., and Morley, J.W. 2023. Environmental and climate variability drive population size of annual penaeid shrimp in a large lagoonal estuary. PLoS One, 18(5 May). doi:10.1371/journal.pone.0285498.

Scientific, Technical and Economic Committee for Fisheries (STECF). 2021. Stock assessments: demersal stocks in the western Mediterranean Sea (STECF-21-11). Publications Office of the European Union, Luxembourg. EUR 28359 EN, ISBN 978-92-76-46116-6. JRC127744. doi:10.2760/046729.

Scientific, Technical and Economic Committee for Fisheries (STECF). 2022. Evaluation of maximum catch limits and closure areas in the Western Mediterranean (STECF-22-01). Publications Office of the European Union, Luxembourg. EUR 28359 EN. ISBN 978-92-76-51980-5. JRC129243. doi:10.2760/657891.

Shelton, A.O., Thorson, J.T., Ward, E.J., and Feist, B.E. 2014. Spatial semi-parametric models improve estimates of species abundance and distribution. Can. J. Fish. Aquat. Sci. 71(11): 1655-1666. doi:10.1139/cjfas-2013-0508.

Sparre, P., and Venema, S.C. 1998. Introduction to fish stock assessment. Part 1: manual. FAO Fish. Tech. Paper, 306(1).

Spedicato, M.T., Cannas, R., Mahe, K., Morales, B., Tsigenopoulos, C., Zane, L., et al. 2021. Study on advancing fisheries assessment and management advice in the Mediterranean by aligning biological and management units of priority species. MED_UNITS.

Spedicato, M.T., Massutí, E., Mérigot, B., Tserpes, G., Jadaud, A., and Reini, G. 2019. The medits trawl survey specifications in an ecosystem approach to fishery management. Sci. Mar. 83(S1): 9-20. doi:10.3989/scimar.04915.11X.

Stenseth, N.C., and Mysterud, A. 2005. Weather packages: finding the right scale and composition of climate in ecology. doi:10.1111/j.1365-2656.2005.01005.x.

Szuwalski, C.S., and Hollowed, A.B. 2016. Climate change and non-stationary population processes in fisheries management. ICES J. Mar. Sci. 1297-1305. doi:10.1093/icesjms/fsv229.

Szuwalski, C.S., Vert-Pre, K.A., Punt, A.E., Branch, T.A., and Hilborn, R. 2015. Examining common assumptions about recruitment: a meta-analysis of recruitment dynamics for worldwide marine fisheries. Fish Fish. 16(4): 633-648. doi:10.1111/faf.12083.
<figure><img src="image_12.png" /></figure>

Thorson, J.T. 2018. Three problems with the conventional delta-model for biomass sampling data, and a computationally efficient alternative. Can. J. Fish. Aquat. Sci. **75**(9): 1369-1382. doi:10.1139/cjfas-2017-0266.

Thorson, J.T. 2019a. Guidance for decisions using the Vector Autoregressive Spatio-Temporal (VAST) package in stock, ecosystem, habitat and climate assessments. Fish. Res. **210**: 143-161. doi:10.1016/j.fishres.2018.10.013.

Thorson, J.T. 2019b. Measuring the impact of oceanographic indices on species distribution shifts: the spatially varying effect of cold-pool extent in the eastern Bering Sea. Limnol. Oceanogr. **64**(6): 2632-2645. doi:10.1002/lno.11238.

Thorson, J.T., Arimitsu, M.L., Barnett, L.A.K., Cheng, W., Eisner, L.B., Haynie, A.C., et al. 2021. Forecasting community reassembly using climate-linked spatio-temporal ecosystem models. Ecography, **44**(4): 612-625. doi:10.1111/ecog.05471.

Thorson, J.T., Cheng, W., Hermann, A.J., Ianelli, J.N., Litzow, M.A., O'Leary, C.A., and Thompson, G.G. 2020b. Empirical orthogonal function regression: linking population biology to spatial varying environmental conditions using climate projections. Global Change Biol. **26**(8): 4638-4649. doi:10.1111/gcb.15149.

Thorson, J.T., Ciannelli, L., and Litzow, M.A. 2020a. Defining indices of ecosystem variability using biological samples of fish communities: a generalization of empirical orthogonal functions. Prog. Oceanogr. **181**. doi:10.1016/j-pocean.2019.102244.

Thorson, J.T., Scheuerell, M.D., Shelton, A.O., See, K.E., Skaug, H.J., and Kristensen, K. 2015. Spatial factor analysis: a new tool for estimating joint species distributions and correlations in species range. Methods Ecol. Evol. **6**(6): 627-637. doi:10.1111/2041-210X.12359.

Trochta, J.T., Pons, M., Rudd, M.B., Krigbaum, M., Tanz, A., and Hilborn, R. 2018. Ecosystem-based fisheries management: perception on definitions, implementations, and aspirations. PLoS One, **13**(1). doi:10.1371/journal.pone.0190467.

Tserpes, G., Massutí, E., Fiorentino, F., Facchini, M.T., Viva, C., Jadaud, A., et al. 2019. Distribution and spatio-temporal biomass trends of red mullets across the mediterranean. Sci. Mar. **83**(S1): 43-55. doi:10.3989/scimar.04888.21A.

Tsikliras, A.C., Licandro, P., Pardalou, A., McQuinn, I.H., Gröger, J.P., and Alheit, J. 2019. Synchronization of Mediterranean pelagic fish populations with the North Atlantic climate variability. Deep Sea Res. 2 Top. Stud. Oceanogr. **159**: 143-151. doi:10.1016/j.dsr2.2018.07.005.

Twiname, S., Audzijonyte, A., Blanchard, J.L., Champion, C., de la Chesnais, T., Fitzgibbon, Q.P., et al. 2020. A cross-scale framework to support a mechanistic understanding and modelling of marine climate-driven species redistribution, from individuals to communities. Ecography, **43**(12): 1764-1778. doi:10.1111/ecog.04996.

Ulses, C., Estournel, C., De Madron, X. D., and Palanques, A. 2008. Suspended sediment transport in the Gulf of Lions (NW Mediterranean): impact of extreme storms and floods. Continental Shelf Res. **28**(15): 2048-2070. doi:10.1016/j.csr.2008.01.015.

Vandenbulcke, L., and Barth, A. 2019. Upscaling of a local model into a larger-scale model. Ocean Sci. **15**(2): 291-305. doi:10.5194/os-15-291-2019.

Vilas, D., Pennino, M.G., Bellido, J.M., Navarro, J., Palomera, I., and Coll, M. 2020. Seasonality of spatial patterns of abundance, biomass, and biodiversity in a demersal community of the NW Mediterranean Sea. ICES J. Mar. Sci. **77**(2): 567-580. doi:10.1093/icesjms/fsz197.

Zimmermann, F., Claireaux, M., and Enberg, K. 2019. Common trends in recruitment dynamics of north-east Atlantic fish stocks and their links to environment, ecology and management. Fish Fish. **20**(3): 518-536. doi:10.1111/faf.12360.

Zuur, A.F., Fryer, R.J., Jolliffe, I.T., Dekker, R., and Beukema, J.J. 2003b. Estimating common trends in multivariate time series using dynamic factor analysis. Environmetrics, **14**(7): 665-685. doi:10.1002/env.611.

Zuur, A.F., Tuck, I.D., and Bailey, N. 2003a. Dynamic factor analysis to estimate common trends in fisheries time series. Can. J. Fish. Aquat. Sci. **60**(5): 542-552. doi:10.1139/F03-030.