Mediterranean demersal resources and ecosystems: 25 years of MEDITS trawl surveys
M.T. Spedicato, G. Tserpes, B. Mérigot and E. Massutí (eds)

Stability of the relationships among demersal fish assemblages and environmental-trawling drivers at large spatio-temporal scales in the northern Mediterranean Sea

Bastien Mérigot 1, Jean-Claude Gaertner 2, Anik Brind’Amour 3, Pierluigi Carbonara 4, Antonio Esteban 5, Cristina Garcia-Ruiz 6, Michele Gristina 7, Taha Imzilen 1, Angélique Jadaud 1, Aleksandar Joksimovic 8, Stefanos Kavadas 9, Jerina Kolitari 10, Irida Maina 9, Porzia Maiorano 11, Chiara Manfredi 12, Reno Micallef 13, Panagiota Peristeraki 14, 15, Giulio Relini 16, Mario Sbrana 17, Maria Teresa Spedicato 4, Ioannis Thasitis 18, Stefania Vittori 19, Nedo Vrgoc 20

1 MARBEC - Université de Montpellier, CNRS, Ifremer, IRD, Sète, France.
(BM) (Corresponding author) E-mail: bastien.merigot@umontpellier.fr. ORCID iD: https://orcid.org/0000-0001-5264-4324
(TI) E-mail: taha.imzilen@ird.fr. ORCID iD: https://orcid.org/0000-0001-8257-9113
(AJ) E.mail: Angelique.Jadaud@ifremer.fr. ORCID iD: http://orcid.org/0000-0001-6858-3570
2 Institut de Recherche pour le Développement (IRD), UMR 241 Ecosystèmes Insulaire Océanique (EIO), Université de Polynésie française, Punaauia, BP 6570, 98702 Faa’a, Tahiti, Polynésie française.
(JCG) E-mail: jean-claude.gaertner@ird.fr. ORCID-iD: https://orcid.org/0000-0001-9248-6260
3 Institut Français de Recherche pour l’Exploitation de la Mer (Ifremer), Département Ecologie et Modèles pour l’Halieutique (EMH), 44311 Nantes, France.
(ABA) E-mail: Anik.Brindamour@ifremer.fr. ORCID-iD: https://orcid.org/0000-0003-4072-2588
4 COISPA Tecnologia e Ricerca, Stazione Sperimentale per lo Studio delle Risorse del Mare, Bari, Italy.
(PC) E-mail: carbonara@coispa.it. ORCID iD: https://orcid.org/0000-0002-2529-2535
(MTS) E-mail: spedicato@coispo.it. ORCID iD: https://orcid.org/0000-0001-9939-9426
5 Instituto Español de Oceanografía (IEO), Centro Oceanográfico de Murcia, Murcia, Spain.
(AE) E-mail: antonio.esteban@ieo.es. ORCID-iD: https://orcid.org/0000-0002-2896-7972
6 Instituto Español de Oceanografía (IEO), Centro Oceanográfico de Málaga , Fuengirola, Málaga, Spain.
(CGR) E-mail: cristina.garcia@ieo.es. ORCID-iD: https://orcid.org/0000-0003-2767-4200
7 Institute of Anthropic Impacts and Sustainability in marine environment (CNR), Castellammare del Golfo, Italy.
(MG) E-mail: michele.gristina@cnr.it. ORCID-iD: https://orcid.org/0000-0003-3639-7655
8 Institute of Marine Biology, University of Montenegro, 85330 Kotor, Montenegro.
(AJ) E-mail: acojo@ac.me. ORCID iD: https://orcid.org/0000-0002-0796-4641
9 Institute of Marine Biological Resources and Inland Waters, Hellenic Centre for Marine Research (IMBRIW/HCMR, 19013 Anavyssos, Attiki, Greece.
(SK) E-mail: stefanos@hcmr.gr. ORCID iD: https://orcid.org/0000-0003-3473-9084
(IM) E-mail: imaina@hcmr.gr. ORCID iD: https://orcid.org/0000-0001-5244-2722
10 Aquaculture and Fishery Laboratory Durres, Agricultural University Tirana, Durres, Albania.
(JK) E-mail: jerina_juka@yahoo.com. ORCID iD: https://orcid.org/0000-0002-6116-8726
11 Department of Biology, University of Bari Aldo Moro, Bari, Italy.
(PM) E-mail: porzia.maiorano@uniba.it. ORCID iD: https://orcid.org/0000-0001-5737-3025
12 Laboratorio di Biologia Marina e Pesca di Fano, Dipto. B.E.S., Università di Bologna, Fano, Italy.
(CM) E-mail: chiara.manfredi3@unibo.it. ORCID iD: https://orcid.org/0000-0002-2852-4856
13 Department of Fisheries and Aquaculture, Ministry for Sustainable Development, the Environment and Climate Change (MSDEC), Marsa, Malta.
(RM) E-mail: reno.micallef@gov.mt. ORCID iD: https://orcid.org/0000-0003-1921-508X
14 Institute of Marine Biological Resources and Inland Waters, Hellenic Center for Marine Research, 71003 Heraklion, Crete, Greece.
15 University of Crete, Biology Department, Stavrakia, Heraklion, Crete, Greece.
(PP) E-mail: notap@her.hcmr.gr. ORCID iD: https://orcid.org/0000-0002-8608-078X
16 Dipartimento di Scienze della Terra, dell’Ambiente e della Vita (DISTAV), Università 16132 Genova, Italy.
(GR) E-mail: biolmar@unige.it. ORCID iD: https://orcid.org/0000-0003-2502-5581
17 Consorzio per il Centro Interuniversitario di Biologia Marina ed Ecologia Applicata, Livorno, Italy.
(MS) E-mail: msbrana@cibm.it. ORCID iD: https://orcid.org/0000-0002-9562-8372
18 Department of Fisheries and Marine Research, Ministry of Agriculture, Natural Resources and Environment, Nicosia, Cyprus.
(IT) E-mail: ithasitis@dfmr.moa.gov.cy. ORCID iD: https://orcid.org/0000-0002-0940-2212
19 Dipartimento di Scienze della vita e dell’ambiente, Università di Cagliari, Cagliari, Italy.
(SV) E-mail: ste.vittori@gmail.com. ORCID iD: https://orcid.org/0000-0001-6857-0470
20 Institute of Oceanography and Fisheries, Split, Croatia.
(NV) E-mail: vrgoc@izor.hr. ORCID iD: https://orcid.org/0000-0002-5208-4512

Summary: Trawling pressure and environmental changes may affect the composition of fish assemblages. Our knowledge on large spatio-temporal patterns of demersal fish composition remains incomplete for the Mediterranean Sea. We investigated (1) the spatio-temporal stability of demersal assemblages, (2) the relationships between these assemblages and potential structuring factors (trawling pressure and environmental conditions) in order to assess the dynamic of the assemblage structure at the scale of the northern Mediterranean Sea. We analysed a dataset of 18062 hauls from 10 to 800 m depth performed annually during the last two decades across 17 Geographical Sub-Areas (GSAs) (MEDITS program). A multi-table analysis (STATICO-CoA) evidenced a strong inter-GSAs stability in the organization of assemblages, with specificities for some GSAs. The most stable structuring factors were linked to combined gradients of chlorophyll a, phytoplancton carbon biomass and temperature, inversely correlated with depth, salinity and nutrient gradients (axis 1 of the STATICO-CoA compromise, 93.74% of the total variability). A common pattern linking the distribution of species to these environmental gradients was evidenced for most of the 17 GSAs. Estimate of trawling pressure showed a minor role in the organization of the assemblages for the spatial scale and years investigated (axis 2, 4.67%).

Keywords: species composition; stability; demersal assemblages; environment; fishing pressure; large scale; co-inertia analysis; STATICO-CoA.

Estabilidad de las relaciones entre las agrupaciones de peces demersales y los factores ambientales-actividad pesquera a gran escala espacio-temporal en el norte del mar Mediterráneo

Resumen: La presión pesquera y los cambios ambientales pueden afectar a la composición de las agrupaciones de peces. En el Mediterráneo, nuestro conocimiento a gran escala sobre los patrones espacio-temporales de la composición de especies de peces demersales sigue siendo incompleto. Investigamos (1) la estabilidad espacio-temporal de las agrupaciones demersales (2) las relaciones entre dichas agrupaciones y los posibles factores estructurantes (presión pesquera de arrastre y condiciones ambientales) para evaluar la dinámica de la estructura de las agrupaciones en el norte del mar Mediterráneo. Analizamos un conjunto de 18062 lances entre 10 y 800 m. de profundidad realizados anualmente, durante las últimas dos décadas, en 17 sub-áreas geográficas (GSAs) (Proyecto MEDITS). El análisis multi-tablas (STATICO-CoA) evidenció una fuerte estabilidad inter-GSAs en la organización de agrupaciones, con especificidades para algunas GSAs. Los factores estructurantes más estables se relacionaron con los gradientes combinados de clorofila a, fitoplancton, carbono y temperatura, correlacionados inversamente con los gradientes de profundidad, salinidad y nutrientes (eje 1 del STATICO-CoA, 93.74% de la variabilidad total). En la mayoría de las GSAs encontramos un patrón común que vincula la distribución de las especies a estos gradientes ambientales. La estimación de la presión pesquera de arrastre mostró un papel menor en la organización de las agrupaciones para la escala espacial y los años investigados (eje 2, 4.67%).

Palabras clave: composición específica; estabilidad; agrupaciones demersales; medio ambiente; presión pesquera; gran escala; análisis de co-inercia; STATICO-CoA.

Citation/Cómo citar este artículo: Mérigot B., Gaertner J.-C., Brind’Amour A., Carbonara P., Esteban A., Garcia-Ruiz C., Gristina M., Imzilen T., Jadaud A., Joksimovic A., Kavadas S., Kolitari J., Maina I., Maiorano P., Manfredi C., Micallef R., Peristeraki P., Relini G., Sbrana M., Spedicato M.T., Thasitis I., Vittori S., Vrgoc N. 2019. Stability of the relationships among demersal fish assemblages and environmental-trawling drivers at large spatio-temporal scales in the northern Mediterranean Sea. Sci. Mar. 83S1: 153-163. https://doi.org/10.3989/scimar.04954.30A

Editor: E. Massutí.

Received: February 27, 2018. Accepted: May 14, 2019. Published: September 2, 2019.

Copyright: © 2019 CSIC. This is an open-access article distributed under the terms of the Creative Commons Attribution 4.0 International (CC BY 4.0) License.

Contents

Summary
Resumen
Introduction
Materials and methods
Results
Discussion
Acknowledgements
References

INTRODUCTIONTop

Due to its geological evolution and environmental conditions, the Mediterranean Sea is characterized by relatively high biodiversity (Coll et al. 2010Coll M., Piroddi C., Steenbeek J., et al. 2010. The Biodiversity of the Mediterranean Sea: Estimates, Patterns, and Threats. PLoS ONE 5: e11842.). However, it is also one of the world’s most threatened marine regions due to environmental change and fishing activities (Ben Rais Lasram et al. 2010Ben Rais Lasram F., Guilhaumon F., Albouy C., et al. 2010. The Mediterranean Sea as a “cul-de-sac” for endemic fishes facing climate change. Glob. Chang. Biol. 16: 3233-3245., Coll et al. 2010Coll M., Piroddi C., Steenbeek J., et al. 2010. The Biodiversity of the Mediterranean Sea: Estimates, Patterns, and Threats. PLoS ONE 5: e11842., 2012Coll M., Piroddi C., Albouy C., et al. 2012. The Mediterranean Sea under siege: spatial overlap between marine biodiversity, cumulative threats and marine reserves. Glob. Ecol. Biogeogr. 21: 465-480.). These drivers might affect the structure and composition of species assemblages (Menge and Olson 1990Menge B.A., Olson A.M. 1990. Role of scale and environmental factors in the regulation of community structure. Trends Ecol. Evol. 5: 52-57., Greenstreet and Hall 1996Greenstreet S.P.R., Hall S.J. 1996. Fishing and the ground-fish assemblage structure in the north-western North Sea: an analysis of long-term and spatial trends. J. Anim. Ecol. 65: 577-598.), in particular those of demersal fish targeted by intensive commercial bottom trawling, such as in the northern Mediterranean Sea (Rochet et al. 2010Rochet M.J., Trenkel V.M., Carpentier A., et al. 2010. Do changes in environmental and fishing pressures impact marine communities? An empirical assessment. J. Appl. Ecol. 47: 741-750.). However, our knowledge on the impact of trawling and environmental forcing on demersal assemblages in the Mediterranean is still largely incomplete. On the one hand, previous studies analyzing fish composition at large scale in this sea have provided insights based on gridded species presence-absence maps drawn on the basis of expert-based knowledge (e.g. Mouillot et al. 2011Mouillot D., Albouy C., Guilhaumon F., et al. 2011. Protected and threatened components of fish biodiversity in the Mediterranean Sea. Curr. Biol. 21: 1044-1050., Coll et al. 2012Coll M., Piroddi C., Albouy C., et al. 2012. The Mediterranean Sea under siege: spatial overlap between marine biodiversity, cumulative threats and marine reserves. Glob. Ecol. Biogeogr. 21: 465-480.). However, these works were not necessarily focused on demersal fish assemblages caught by trawling. In addition, they were mainly restricted to indirect presence-absence data estimated at coarse spatial grain and they could not investigate temporal trends. On the other hand, some previous studies investigated large spatio-temporal trends in the Mediterranean Sea on the basis of fish abundance data collected by standardized scientific surveys. However, these studies were mainly focused on diversity components, such as species number (Gaertner et al. 2007Gaertner J.C., Bertrand J., Relini G., et al. 2007. Spatial pattern in species richness of demersal fish assemblages on the continental shelf of the northern Mediterranean Sea: a multiscale analysis. Mar. Ecol. Prog. Ser. 341: 191-203.) evenness, taxonomic, phylogenetic and/or functional diversity (Gaertner et al. 2010Gaertner J.C., Mérigot B., Rélini G., et al. 2010. Reproducibility of the multi-component aspect of species diversity across different areas and scales: towards the constitution of a shortlist of complementary indices for monitoring fish diversity? Ecography 33: 1123-1135., Granger et al. 2015Granger V., Fromentin J.M., Bez N., et al. 2015. Large-scale spatio-temporal monitoring highlights hotspots of demersal fish diversity in the Mediterranean Sea. Prog. Oceanogr. 130: 65-74., Brind’Amour et al. 2016Brind’Amour A., Rochet M.J., Ordines F., et al. 2016. Environmental drivers explain regional variation of changes in fish and invertebrate functional groups across the Mediterranean Sea from 1994 to 2012. Mar. Ecol. Prog. Ser. 562: 19-35.), but they did not deal with the structure and composition of demersal fish assemblages per se. Indeed they mainly focused on diversity indices which thus summarize in several single and synthetic values some features of assemblages, while some multi-variate analyses allows to assess the identity and contribution of species to patterns.

Thus, the aim of this study is to complete our knowledge in this domain by assessing large-scale spatial and temporal patterns of the structure and composition of demersal fish assemblages sampled by trawling over the northern Mediterranean Sea. More specifically, the main objectives are to assess both the stable and variable parts of the spatio-temporal relationships between the species composition of demersal fish assemblages and the environmental-trawling variables at large scale in the Mediterranean Sea. The analyses were based on species abundance data collected from 1999 to 2015 within a standardized scientific protocol consisting in 18062 sampled hauls performed along more than 20000 km of coastline over a wide bathymetric zone (10-800 m), seagrass excluded (Posidonia oceanica meadows, etc.).

MATERIALS AND METHODSTop

Survey methods and faunal data collection

Data was collected from annual bottom trawl surveys performed in May-July from 1999 to 2015 over the continental shelf (10 to 200 m depth) and the continental slope (200 to 800 m) of the Mediterranean Sea within the MEDITS scientific program (Bertrand et al. 2002aBertrand J.A., De Sola L., Papaconstantinou C., et al. 2002a. The general specifications of the MEDITS surveys. Sci. Mar. 66: 9-17., MEDITS Working Group 2017MEDITS Working Group. 2017. MEDITS-Handbook, Version n. 9. MEDITS Working Group, 106 pp.). The study area ranged from 34.33°N to 45.67°N and 5.22°W to 34.09°E and was divided into 17 operative Geographical Sub-Areas (GSAs), for which boundaries were established by the General Fisheries Commission for the Mediterranean Sea (management areas, according to the resolution GFCM/33/2009/2, http://www.gfcm.org, see Fig. 1 for the correspondence between GSA numbers and their names, and Table S1 for number of hauls per GSA).

figure1

Full size image

Fig. 1. – Fig. 1. – Study area and sampling sites based on the MEDITS program protocol with the position of 18062 hauls sampled between 1999 and 2015 in 17 Geographical Sub-Areas (GSAs); each color corresponds to one of the 17 GSA defined by the General Fisheries Commission for the Mediterranean (http://www.gfcm.org).

The sampling procedures were standardized according to a common protocol over GSAs and years. The standard device is a bottom trawl GOC-73 with 20-mm cod-end mesh size. The average vertical opening of the gear is 2 m and its wing-span 18 m (Bertrand et al. 2002bBertrand J.A., Leonori I., Dremiere P.Y., et al. 2002b. Depth trajectory and performance of a trawl used for an international bottom trawl survey in the Mediterranean. Sci. Mar. 66: 169-182.). All the tows were performed during daylight hours. Speed on the ground and duration of the tows were standardized to 3 knots and 30 min for shelf stations. For the slope this length of time was doubled (60 min) in order to overcome potential difficulties in the deep hauls (see Bertrand et al. 2002aBertrand J.A., De Sola L., Papaconstantinou C., et al. 2002a. The general specifications of the MEDITS surveys. Sci. Mar. 66: 9-17 for full details). For each haul, we collected raw fish abundance data, bathymetry, swept area (i.e. trawled surface) and haul coordinates.

Information was recorded by an underwater Scanmar system, to control the trawl geometry (horizontal and vertical openings, contact with the bottom), as well as to record swept area and trawling time, which enabled us to select the tows that were properly carried out. Among a total of 394 fish species sampled during MEDITS surveys (Relini 2015Relini G. 2015. Fish biodiversity in MEDITS surveys. Biol. Mar. Mediterr. 22: 176-177.), we considered 154 species efficiently sampled by the bottom trawl system used, and properly identified by teams involved in this large scale program (Gaertner et al. 2010Gaertner J.C., Mérigot B., Rélini G., et al. 2010. Reproducibility of the multi-component aspect of species diversity across different areas and scales: towards the constitution of a shortlist of complementary indices for monitoring fish diversity? Ecography 33: 1123-1135., 2013Gaertner J.C., Maiorano P., Mérigot B., et al. 2013. Large-Scale Diversity of Slope Fishes: Pattern Inconsistency between Multiple Diversity Indices. PLoS ONE 8: e66753., Granger et al. 2015Granger V., Fromentin J.M., Bez N., et al. 2015. Large-scale spatio-temporal monitoring highlights hotspots of demersal fish diversity in the Mediterranean Sea. Prog. Oceanogr. 130: 65-74.). Indeed, i) only demersal fishes (i.e. benthic and bentho-pelagic fishes) were considered, ii) individuals for which the identification at species level was complex were gathered at genus level (it was only the case for individuals of Lophius spp., see Table S2). This subset was thus designed in such a way as to minimize potential bias caused by gear-selectivity (e.g. Lefcheck et al. 2014Lefcheck J.S., Buchheister A., Laumann K.M., et al. 2014. Dimensions of biodiversity in Chesapeake Bay demersal fishes: patterns and drivers through space and time. Ecosphere 5: 14.), and to strictly limit the risk of a variability of accuracy in sampling identification between the different teams. The list of 154 species considered was set according to those used in previous MEDITS data based studies (Gaertner et al. 2010Gaertner J.C., Mérigot B., Rélini G., et al. 2010. Reproducibility of the multi-component aspect of species diversity across different areas and scales: towards the constitution of a shortlist of complementary indices for monitoring fish diversity? Ecography 33: 1123-1135., 2013Gaertner J.C., Maiorano P., Mérigot B., et al. 2013. Large-Scale Diversity of Slope Fishes: Pattern Inconsistency between Multiple Diversity Indices. PLoS ONE 8: e66753., Granger et al. 2015Granger V., Fromentin J.M., Bez N., et al. 2015. Large-scale spatio-temporal monitoring highlights hotspots of demersal fish diversity in the Mediterranean Sea. Prog. Oceanogr. 130: 65-74.; see Table S2 for the list of species considered). Raw abundance data was normalized by the swept area to compute species densities.

Environmental and trawling data

A total of eight explanatory variables were considered: one related to trawling pressure and seven to environment. Among the seven environmental variables, the mean depth of operated haul during surveys was computed. The six other variables (i.e. temperature, salinity, chlorophyll-a concentration, phosphate, nitrate and phytoplankton carbon biomass) were extracted from the Copernicus portal (http://marine.copernicus.eu/) when available for the years of the MEDITS survey (i.e. 1999 to 2015, thus leading to the exclusion of hauls sampled from 1994 to 1998). According to the variable, two models were considered to deliver their values. Both models have a horizontal grid resolution of 1/16×1/16° with 72 unevenly spaced depth levels (Oddo et al. 2009Oddo P., Adani M., Pinardi N., et al. 2009. A nested Atlantic-Mediterranean Sea general circulation model for operational forecasting. Ocean Sci. 5: 461-473.). The values of the variables were extracted according to the depth level that embeds the mean depth of the haul performed. Temperature (°C) and salinity (psu) were derived from the 3D Mediterranean Sea Physics Reanalysis model, whereas chlorophyll-a concentration (mg m–3), nutrients (phosphate and nitrate (mmol m–3) and phytoplankton carbon biomass (mmol m–3) were computed from the 3D Mediterranean Sea Biogeochemistry Reanalysis model (details available at http://marine.copernicus.eu/). Two values for each variable were prepared for each haul: the monthly-mean value over the year before the haul had been performed, and its standard deviation (sd). During a preliminary analysis, we investigated the pairwise redundancy/complementarity of the eight explanatory variables (mean and sd), by using a Draftsman’s plot and Spearman’s correlation coefficient. Because we found that mean and sd values were correlated, only means were used as explanatory variables in the statistical analyses presented below.

Given that data on fishing vessel locations (e.g. VMS/AIS data) were not available for all spatial and temporal scales of the MEDITS survey, a Multi-Criteria Decision Analysis (MCDA) has been employed to estimate a trawling pressure index for commercial bottom trawlers. As a rule, MCDA provides successive complementary techniques and procedures for structuring decision problems and deriving a way for a compromise in a transparent process (Malczewski 2006Malczewski J. 2006. GIS-based multi criteria decision analysis: a survey of the literature. Int. J. Geogr. Inf. Sci. 20: 703-726.). In our case, a GIS-MCDA process that combines influential geographical components affecting bottom trawling in terms of its distribution and intensity (such as fishing capacity, bathymetry etc. see more details in Fig. S1), as well as expert judgement opinions (Table S3), have been used based on a methodological approach proposed by Kavadas et al. (2015)Kavadas S., Maina I., Damalas D., et al. 2015. Multi-Criteria Decision Analysis as a tool to extract fishing footprints and estimate fishing pressure: application to small scale coastal fisheries and implications for management in the context of the Maritime Spatial Planning Directive. Mediterr. Mar. Sci. 16: 294-304.. The method was applied to bottom trawlers classified by total length (LOA) segmentation: <12 m, 12-24 m and 24-40 m. The generated trawling pressure index was calculated as the average of trawling pressure indices per LOA on a spatial resolution of 0.01×0.01 decimal degrees. The index resulted from the fuzzy product of an activity index and a suitability index (see details in Fig. S1).

The activity index was based on trawling vessel activities estimated by registration port for 1994, 2004 and 2014, considering the vessel length, gross tonnage and the number of vessels (source: European Fleet Register). An Inverse Distance Weighted (IDW) interpolation method was used to spatialize this index. Then a standardization with a linear Fuzzy Membership function was applied to scale the index between 0-1.

The suitability index reflected the suitability of trawling activity by combining information on several criteria: bathymetry (source: EMODNET), distance from the coast (estimated by ArcGIS ESRI, 2011), annual chlorophyll a (https://oceancolor.gsfc.nasa.gov) and fisheries restricted areas (legislation, source: MEDISEH, Colloca et al. 2013Colloca F., Spedicato M.T., Massutí E., et al. 2013. Mapping of nursery and spawning grounds of demersal fish. Mediterranean Sensitive Habitats (MEDISEH) Final Report, DG MARE Specific Contract SI2.600741, Heraklion (Greece).). Each criterion was assigned a grading value by expert judgment (i.e. a rank of order of importance from 0 to 5, Table S3).

The setting of this suitability index from the investigated criteria proceeded in the following steps: (i) creation of spatial information and calibration of each criterion according to a scale of evaluation and formation of the hierarchical structure of the multiple criteria problem; (ii) implementation of the Analytic Hierarchy Process (AHP, Saaty 1980Saaty T.L. 1980. The analytical hierarchy process: planning, priority setting, resource allocation. McGraw-Hill, New York, 309 pp.) to estimate the relative importance of the evaluation criteria (Table S4); (iii) application of the Weighted Linear Combination method using the weights (priority vectors) to estimate the suitability index (Table S4); (iv) standardization on a scale from 0 to 1 with linear Fuzzy Membership function (see Kavadas et al. 2015Kavadas S., Maina I., Damalas D., et al. 2015. Multi-Criteria Decision Analysis as a tool to extract fishing footprints and estimate fishing pressure: application to small scale coastal fisheries and implications for management in the context of the Maritime Spatial Planning Directive. Mediterr. Mar. Sci. 16: 294-304. for details of the overall methodology).

Finally, a value of the trawling pressure index was assigned to each MEDITS operated haul according to its spatial location (defined as mean position) and year. Trawling pressure index values of 1994, 2004 and 2014 were assigned to hauls performed in 1994 to 2003, 2004 to 2013 and 2014 to 2015, respectively. It is thus assumed that the index remained constant over each of these periods. Note that even if generally the index value remained relatively similar for a given location between 1994, 2004 and 2014, it varied in time within some GSAs between these periods (Figs S2A and S2B). In addition, the index showed contrasted spatial patterns among GSAs (Fig. S2A).

Statistical analysis

We investigated the stable patterns and spatio-temporal variability of the relationships between species composition of demersal fish assemblages and the environmental/trawling variables over the period under study at the scale of the northern Mediterranean Sea. To that end, we carried out a multi-table analysis, STATICO-CoA (Simier et al. 1999Simier M., Blanc L., Pellegrin F. et al. 1999. Approche simultanée de k couples de tableaux: Application à l’étude des relations pathologie végétale-environnement. Rev. Statist. Appl. 47: 31-46., Thioulouse et al. 2004Thioulouse J., Simier M., Chessel D. 2004. Simultaneous analysis of a sequence of paired ecological tables. Ecology 85: 272-283., Thioulouse 2011Thioulouse J. 2011. Simultaneous analysis of a sequence of paired ecological tables: a comparison of several methods. Ann. App. Stat. 5: 2300-2325.) that enabled us to investigate both the stable patterns and the spatio-temporal changes of these relationships. STATICO is a multivariate, multi-table method that allows to access the relationships between two kinds of information. Here, data were organized in a series of pairs of species and environmental-trawling tables associated with each of the 17 GSAs. STATICO, a co-inertia analysis, was used to assess not solely the relationship between species distribution and explanatory variables but also the degree of spatial-temporal reproducibility of this relationship. It is thus a very useful method to investigate key aims in ecology with data of large spatio-temporal scale compared to traditional analysis to study this relationship (Mendes et al. 2012Mendes S., Fernández-Gómez M.J., Pereira M.J., et al. 2012. An empirical comparison of canonical correspondence analysis and STATICO in the identification of spatio-temporal ecolog- ical relationships. J. App. Stat. 39: 979-994., Le Fur et al. 2019Le Fur I., De Wit R., Plus M., et al. 2019. Re-oligotrophication trajectories of macrophyte assemblages in Mediterranean coastal lagoons based on 17-year time-series. Mar. Ecol. Prog. Ser. 608: 13-32.). STATICO has been successfully applied on fauna and flora of various ecosystems, such as demersal fishes (Kidé et al. 2015Kidé O., Manté C., Dubroca L., et al. 2015. Spatio-temporal dynamics of exploited ground-fish assemblages faced to environmental and fishing forcings: Insights from the Mauritanian Exclusive Economic Zone. Plos ONE 10: e0141566.), diatoms (Mendes et al. 2012Mendes S., Fernández-Gómez M.J., Pereira M.J., et al. 2012. An empirical comparison of canonical correspondence analysis and STATICO in the identification of spatio-temporal ecolog- ical relationships. J. App. Stat. 39: 979-994.), estuarine fishes (Simier et al. 2006Simier M., Laurent C., Ecoutin J.M., et al. 2006. The Gambia River estuary: A reference point for estuarine fish assemblages studies in West Africa. Estuar. Coast. Shelf. S. 69: 615-628.), heteroptera (Slimani et al. 2017Slimani N., Guilbert E., El Ayni F., et al. 2017. The use of STATICO and COSTATIS, two exploratory three-ways analysis methods. An application to the ecology of aquatic Heteroptera in the Medjerda watershed (Tunisia). Env. Ecol. Stat. 24: 269-295.), macrophytes (Le Fur et al. 2019Le Fur I., De Wit R., Plus M., et al. 2019. Re-oligotrophication trajectories of macrophyte assemblages in Mediterranean coastal lagoons based on 17-year time-series. Mar. Ecol. Prog. Ser. 608: 13-32.), pelagic fishes and top predators (Certain et al. 2011Certain G., Masse J., Van Canneyt O., et al. 2011. Investigating the coupling between small pelagic fish and marine top predators using data collected from ecosystem-based surveys. Mar. Ecol. Prog. Ser. 422: 23-39.), pelagic larval and juvenile fishes (Carassou and Ponton 2007Carassou L., Ponton D. 2007. Spatio-temporal structure of pelagic larval and juvenile fish assemblages in coastal areas of New Caledonia, southwest Pacific. Mar. Biol. 150: 697-711.), zooplanctons (Mazzocchi et al. 2012Mazzocchi M.G., Dubroca L., García-Comas C., et al. 2012. Stability and resilience in coastal copepod assemblages: The case of the Mediterranean long-term ecological research at Station MC (LTER-MC). Prog. Oceanogr. 97-100: 135-151.), as well as indicators of sustainability (Gallego-Alvarez et al. 2015Gallego-Álvarez I., Galindo-Villardón M., Rodríguez-Rosa M. 2015. Evolution of sustainability indicator worldwide: A study from the economic perspective based on the X-STATICO method. Ecol. Indic. 58: 139-151.).

More in details, as recalled by Kidé et al. (2015)Kidé O., Manté C., Dubroca L., et al. 2015. Spatio-temporal dynamics of exploited ground-fish assemblages faced to environmental and fishing forcings: Insights from the Mauritanian Exclusive Economic Zone. Plos ONE 10: e0141566., STATICO (Simier et al. 1999Simier M., Blanc L., Pellegrin F. et al. 1999. Approche simultanée de k couples de tableaux: Application à l’étude des relations pathologie végétale-environnement. Rev. Statist. Appl. 47: 31-46.) is an application of the STATIS method (Lavit et al. 1994Lavit C., Escoufier Y., Sabatier R., et al. 1994 The ACT (Statis method). Comput. Stat. Data Anal. 18: 97-119.) to co-inertia operators (Dolédec and Chessel 1994Dolédec S., Chessel D. 1994. Co-inertia analysis: an alternative method for studying species-environment relationships. Freshw. Biol. 31: 277-294., Dray et al. 2003Dray S., Chessel D., Thioulouse J. 2003. Co-inertia analysis and the linking of ecological data tables. Ecology 84: 3078-3089.). In other words, STATICO “combines the objectives of STATIS (finding the stable part of the structure of a series of tables) and the objectives of co-inertia analysis (finding the common structure of two data tables)”, as stated by Thioulouse et al. 2004Thioulouse J., Simier M., Chessel D. 2004. Simultaneous analysis of a sequence of paired ecological tables. Ecology 85: 272-283. (see also Kidé et al. 2015Kidé O., Manté C., Dubroca L., et al. 2015. Spatio-temporal dynamics of exploited ground-fish assemblages faced to environmental and fishing forcings: Insights from the Mauritanian Exclusive Economic Zone. Plos ONE 10: e0141566. for a reminder of the description of the vectorial approach of STATICO). In comparison with the classical version of STATICO, STATICO-CoA (Thioulouse et al. 2004Thioulouse J., Simier M., Chessel D. 2004. Simultaneous analysis of a sequence of paired ecological tables. Ecology 85: 272-283.) used in the present study is based on the logic of correspondence analysis (CoA, based on the Chi2 distance Hill 1973Hill M.O. 1973. Reciprocal averaging: an eigenvector method of ordination, J. Econ. 61: 237-249.) to investigate the organization of species assemblages (Le Fur et al. 2019Le Fur I., De Wit R., Plus M., et al. 2019. Re-oligotrophication trajectories of macrophyte assemblages in Mediterranean coastal lagoons based on 17-year time-series. Mar. Ecol. Prog. Ser. 608: 13-32., or Gaertner et al. 2002Gaertner J.C., Bertrand J.A., Souplet A. 2002. STATIS-CoA: A methodological solution to assess the spatio-temporal organization of species assemblages. Application to the demersal assemblages of the French Mediterranean Sea. Sci. Mar. 66: 221-232. for a STATIS-CoA).

The data analyzed are a sequence of pairs of tables with the environmental/trawling variables and, separately, the species density, for each GSA. Each table is first analyzed using a CoA for the species tables, and a PCA for the explanatory tables. Species densities (n) were log (n+1) transformed to reduce the influence of dominant species. Environment-trawling data were centered and reduced (i.e. normed) in order to consider their different units. Each pair of tables is then linked by a co-inertia analysis (Dolédec and Chessel 1994Dolédec S., Chessel D. 1994. Co-inertia analysis: an alternative method for studying species-environment relationships. Freshw. Biol. 31: 277-294., Dray et al. 2003Dray S., Chessel D., Thioulouse J. 2003. Co-inertia analysis and the linking of ecological data tables. Ecology 84: 3078-3089.). Co-inertia analysis is a two-table coupling method, which allows a cross-table to be computed between the variables of the two tables (here between species and environmental/trawling variables) for each GSA. The resulting series of species and environmental/trawling variables cross-tables is then analyzed by the STATIS method (Lavit et al. 1994Lavit C., Escoufier Y., Sabatier R., et al. 1994 The ACT (Statis method). Comput. Stat. Data Anal. 18: 97-119.). This analysis is based on three main steps, each of them offering different kinds of results (see a STATICO flow chart in Kidé et al. 2015Kidé O., Manté C., Dubroca L., et al. 2015. Spatio-temporal dynamics of exploited ground-fish assemblages faced to environmental and fishing forcings: Insights from the Mauritanian Exclusive Economic Zone. Plos ONE 10: e0141566.): (1) the inter-structure step identifying the proximity between each of the K tables (here the 17 tables assessing the link between species composition and potential drivers for each of the GSAs); this step aims to compute the weighs of each of these K tables for the construction of a mean table of maximum inertia (compromise table). In other words, weights are used to build a linear combination of the series of tables, called «the compromise». These weights are the components of the first eigenvector of the RV coefficients matrix (see Lavit et al. 1994Lavit C., Escoufier Y., Sabatier R., et al. 1994 The ACT (Statis method). Comput. Stat. Data Anal. 18: 97-119.). The fit of each of the tables to the compromise table is given by the Cos2 value (absolute value ranges between 0 and 1, 1 being the best fit); (2) the compromise analysis, that defines factorial axes which express the stable part of the structures of the K tables, giving an ordination of the reproducible part of the relationships between species composition and environment/trawling that are common for the 17 GSAs across the years. If a/some GSAs contribute less or poorly to the compromise, it means that environmental conditions/trawling pressure were different and/or acted differently on assemblage composition than in the other GSAs; (3) the trajectories step, where species and environmental/trawling variables for each year can be projected as additional elements on the compromise axes in order to summarize their variability around the common structure. For the sake of clarity, a clustering analysis was added at step (2) for identifying on the compromise assemblages of species (see Mazzocchi et al. 2012Mazzocchi M.G., Dubroca L., García-Comas C., et al. 2012. Stability and resilience in coastal copepod assemblages: The case of the Mediterranean long-term ecological research at Station MC (LTER-MC). Prog. Oceanogr. 97-100: 135-151. and Kidé et al. 2015Kidé O., Manté C., Dubroca L., et al. 2015. Spatio-temporal dynamics of exploited ground-fish assemblages faced to environmental and fishing forcings: Insights from the Mauritanian Exclusive Economic Zone. Plos ONE 10: e0141566. for a similar approach). In this way, groups of species were identified with a hierarchical classification based on Euclidean matrix of pairwise distances between species on the first factorial plane (i.e. 1st and 2nd axes) of the compromise analysis. The UPGMA aggregation criterion, chosen on the basis of an objective approach among the main available aggregation criterion (see Mérigot et al. 2010Mérigot B., Durbec J.P., Gaertner J.C. 2010. On goodness-of-fit measure for dendrogram-based analyses. Ecology 91: 1850-1859.), was used to compute this classification. The optimum number of species groups was identified with the R package NbClust that computes and compares thirty indices to identify the optimum number of groups (Charrad et al. 2014Charrad M., Ghazzali N., Boiteau V., et al. 2014. An R Package for Determining the Relevant Number of Clusters in a Data Set. J. Stat. Softw. 61: 1-36.). All statistical analyses were performed using R 3.2.2 software (R Core Team 2015R Core Team. 2015. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. Available online at https://www.R-project.org/).

RESULTSTop

Reproducible part in the relationships of assemblages and environmental/trawling variables among GSAs

The common structure associated with the compromise (first axis of the inter-structure, Fig. S3) explains 52.51% of the total variability. In the construction of the compromise table, greater importance is given to the tables which have similar structures and limited importance to the other tables. Here, we found that the contributions (weights) of the GSAs to the construction of the compromise were relatively similar, ranging mainly from 0.19 to 0.29, except GSA 17 (Northern Adriatic Sea) weighing 0.14 (Fig. 2), with a standard deviation of 0.04. It means that the nature of the link between species assemblages and environmental/trawling drivers for GSA 17 will be less taken into account in the construction of the compromise table. Moreover, the fit of each of the tables to the compromise (see cos2 values, Fig. 2) showed both a common and a specific part of the assemblage organization in each of the seventeen GSA, with cos2 ranging mainly between 0.61 and 0.87 (except GSA 17 and 25 (Cyprus), with 0.38 and 0.56, respectively), with a standard deviation of 0.12. These weaker values observed for GSA 17 and 25 indicated that the structure of these GSAs contributed less to the construction of the compromise table. These results are confirmed by the RV correlation coefficients values among GSAs (Table S5).

figure2

Full size image

Fig. 2. – STATICO-CoA analysis. Cos2 and weights of each Geographical Sub-Area (GSA) to the compromise. Weight is the contribution of each GSA table in the construction of the compromise, and Cos2 is the fit of each GSA table to the compromise. GSA codes are given in Figure 1.

The use of a separate Correspondence Analysis (CoA) carried out independently in each GSA and the projection of the first axes of these separate CoAs on the first two axes of the STATICO-CoA compromise confirmed these results (Fig. S6). Except for GSA 17, the first axis of the compromise appeared strongly correlated to the first axis of each of the separate CoAs (Fig. S6). This result showed that the main assemblage structure was strongly similar in each of the sixteen GSA. Projection of the axis 2 of each separate CoA showed a stronger variability in the organization of the assemblages in each GSA, especially for GSA 17, 22 and 25. In short, these results (Figs 2 and S6) suggested a strong inter-GSA stability in the relationships between demersal species composition and environmental/trawling variables in the northern Mediterranean Sea. It also showed that a limited part of this relationship was specific to some GSAs. Along these lines, GSA 17 (Northern Adriatic Sea) and 25 (Cyprus) appear to be the most specific of the seventeen studied GSAs.

Common pattern in the environmental/trawling variables among GSAs

The first two axes of the compromise accounted for 98.41% of the total variability of the common structure on environmental/trawling variables and species with an overwhelming 93.74% for axis 1 (Figs 3 and 4). Thus, the first axis alone very well summarizes the stable part of the link between species composition and environmental/trawling drivers for all studied GSAs (although it is more limited for GSA 17, see Fig. S4). It shows that the dominant pattern of the environmental variables associated with the demersal assemblages variability in the different studied areas consisted in a chlorophyll a, phytoplancton carbon biomass and temperature gradient, inversely correlated with depth, salinity and nutrient gradients (see axis 1 in Fig. 3; pairwise plot and correlation between variables are provided in Fig. S5). The estimate of trawling pressure index only played a role on axis 2 (Fig. 3) which took into account a very small part of the compromise variability (4.67%).

figure3

Full size image

Fig. 3. – Projections of environmental and trawling variables on the compromise of the STATICO-CoA analysis. It shows the stable part of the species-environment-trawling relationships. Axis 1 and 2 explained 93.74% and 4.67% of the total variability, respectively. Species names are provided in Table S2.

figure4

Full size image

Fig. 4. – Projections of demersal species variables and assemblage groups on the compromise of the STATICO-CoA analysis. It shows the stable part of the species and environment- trawling relationships. Axis 1 and 2 explained 93.74% and 4.67% of the total variability, respectively. On the first factor plane are shown the two main species assemblages (in blue and red) obtained by hierarchical clustering (UPGMA criterion) and optimum number of groups approach (see Materials and Methods section). Species codes are given in Table S2.

Common pattern in species composition among GSAs

The projection of the species onto the first two axes of the compromise allowed identification of the species mainly contributing to the common structure of the demersal assemblages at the scale of the northern Mediterranean Sea (Fig. 4). It gave a typology of the species showing stable trends over the years for the 17 GSAs. There was a clear gradient in the distribution of species (Fig. 4). Two main fish assemblages were identified on the basis of the hierarchical clustering method and optimum number of groups approach computed on the species coordinates on the two first axes of the compromise analysis (Fig. 4, see Materials and Methods section).

As mentioned above, most of the variability is caught in the first axis, representing 93.74% of the total variability. A first group of species was associated with negative values of coordinates on the left part of this axis (i.e. higher temperature, phytoplancton carbon biomass and chlorophyll a), while a second group was located on the positive and right part of axis 1 (i.e. higher depth, nutrients concentration and salinity, Figs 3 and 4). The relevance of analysing axis 2 which represented 4.67% of the variability is limited. Only some species of the two main assemblages could be distributed according to the estimate of trawling pressure (upper left part) and depth (bottom right part).

More precisely, species associated with higher chlorophyll a, phytoplancton carbon biomass and temperature (i.e. left bottom part of the factorial plane, Fig. 3) were mainly Bothus podas, Synodus saurus, Arnoglossus imperialis, Dasyatis pastinaca and Trachinus radiatus (Fig. 4). Species collected in zones where depth, salinity and nutrient concentration were higher (i.e. right bottom part of the factorial plane, Fig. 3) included Alepocephalus rostratus, Centrophorus uyato, Heptranchias perlo, Bellottia apoda and Paralepis coregonoides (Fig. 4).

For each GSA, the projection on the compromise axes of the explanatory variables and species is provided in detail in Figure S4. It shows which GSAs contributed the most to the stable part of the species-environment dynamics (Figs 3 and 4) in the compromise analysis. Except for GSA 17, the stable part of the environmental/trawling variables on axis 1 of the compromise was clearly found for each separate GSA (i.e. similar pattern of explanatory variables in Figures 3 and S4). In contrast, projections of species showed more variability according to the GSAs considered (Figs 3 and S6).

Trajectories of environmental-trawling variables and species assemblages

The main trends of the temporal variability of environmental-trawling conditions and demersal assemblages are provided by their respective trajectories on the first two axes of the compromise for each GSA (Fig. 5). A year point corresponds to the mean coordinates on each axis of hauls performed during that year. The interpretation of the trajectories plots is done by analyzing in each GSA the contraction, the stretching of the cloud of points, and the proximity or remoteness along the axes. All these elements provide information on the relationships between the environment (Fig. 3) and fish assemblage variables (Fig. 4) along the studied years in a given GSA (Fig. 5). Temporal similarity in the trajectory patterns of fish assemblages and environment/trawling suggests that changes in assemblages are related to the environmental-trawling conditions. Overall, the temporal trajectories between environment, trawling features and assemblages did not match along the years analyzed, i.e. there was an absence of a clear coupling in time between environment and trawling features and species trajectories (Fig. 5, Table S6). A link is highlighted from the Kendall correlation coefficient τ computed among the projections of environmental/trawling variables and species with τ value superior to 0.6 only on x-axis (93.74% of the total variability of the common structure) for GSA 1 (τx=0.75) and GSA 20 (τx=0.78) (Table S6).

figure5

Full size image

Fig. 5. – Trajectories coordinates of environmental-trawling and assemblage variables for sampling years on the first factorial plane of the STATICO-CoA compromise for each Geographical Sub-Area (GSA). Axis 1 and 2 explained 93.74% and 4.67% of the total variability, respectively. Trajectory plots of the environmental-trawling variables are in red, and demersal assemblages in blue. Correspondence between GSA numbers and their names is given in Figure 1.

DISCUSSIONTop

The results of the STATICO-CoA analysis suggested a strong inter-GSAs stability in the spatio-temporal relationships between fish assemblage composition and environmental conditions/trawling pressure of the demersal assemblages in the northern Mediterranean Sea (Figs 2 and S4). The analysis highlighted that the most stable structuring variables (explaining over 94% of the total variability) were linked to a combined gradient of chlorophyll a, phytoplancton, carbon biomass and temperature, inversely correlated with the depth, salinity and nutrient gradients. Estimate of trawling pressure contributed relatively poorly to the common organization of the assemblages at the spatial scale and in the years investigated (less than 5% of the total variability). A clear gradient in the distribution of species (Fig. 4) linked to these environmental gradients was evidenced as the stable part of the species environment relationships (note that similar results were obtained for analyses performed without rare species, i.e. those present in less than 5% of the hauls, Gaertner et al. 2002Gaertner J.C., Bertrand J.A., Souplet A. 2002. STATIS-CoA: A methodological solution to assess the spatio-temporal organization of species assemblages. Application to the demersal assemblages of the French Mediterranean Sea. Sci. Mar. 66: 221-232., Figs S7 and S8). However, none of the large majority of GSAs considered showed temporal similarity between the temporal trajectories of the environment (including trawling) and the trajectories of fish assemblages, suggesting an absence of a clear coupling in time between environment, trawling pressure and species trajectories (Fig. 5, Table S5). However, our results also showed specificities in the relationships between species composition and environmental/trawling variables in GSA 1 and 20 (Fig. 5, Table S6), as well as GSA 17 Northern Adriatic Sea and GSA 25 Cyprus (Figs 2 and S3, Table S5).

The relative importance of abiotic and biotic factors in shaping assemblage patterns is considered to vary with spatial scale (Menge and Olson 1990Menge B.A., Olson A.M. 1990. Role of scale and environmental factors in the regulation of community structure. Trends Ecol. Evol. 5: 52-57.). For instance, at a local scale, physical and biotic factors may interact to influence local patterns of assemblage structure. Larger spatial scales are considered to be associated with an increase in the relative influence of variation in environmental or climatic conditions (Menge and Olson 1990Menge B.A., Olson A.M. 1990. Role of scale and environmental factors in the regulation of community structure. Trends Ecol. Evol. 5: 52-57.). Our results highlight the importance of environmental conditions with regard to species assemblage composition and structure at the scale of the northern Mediterranean Sea. Using environmental and trawling pressure information describing each hauls for the different GSAs, we provided insights on the determinism of the assemblage structure at this large scale. To our knowledge, there are few quantitative studies either in the Mediterranean Sea or worldwide that investigate simultaneously the effects of both environmental and anthropogenic factors, with a focus on the species composition, structure and abundance of exploited demersal assemblages (both target and non-target species) at large spatio-temporal scales. Firstly, in the Mediterranean Sea, studies performed on large-scale spatial and/or temporal distribution of fish species assessed mainly species diversity (e.g. species richness, evenness, taxonomic, phylogenetic and/or functional diversity) (Mouillot et al. 2011Mouillot D., Albouy C., Guilhaumon F., et al. 2011. Protected and threatened components of fish biodiversity in the Mediterranean Sea. Curr. Biol. 21: 1044-1050., Granger et al. 2015Granger V., Fromentin J.M., Bez N., et al. 2015. Large-scale spatio-temporal monitoring highlights hotspots of demersal fish diversity in the Mediterranean Sea. Prog. Oceanogr. 130: 65-74., Brind’Amour et al. 2016Brind’Amour A., Rochet M.J., Ordines F., et al. 2016. Environmental drivers explain regional variation of changes in fish and invertebrate functional groups across the Mediterranean Sea from 1994 to 2012. Mar. Ecol. Prog. Ser. 562: 19-35.). In contrast, there has been less focus on assemblage structure and mainly at a lower spatial scale, i.e. between the Straits of Gibraltar and the Gulf of Lions (Gaertner et al. 2005Gaertner J.C., Bertrand J.A., De Sola L.G., et al. 2005. Large spatial scale variation of dem- ersal fish assemblage structure on the continental shelf of the NW Mediterranean Sea. Mar. Ecol. Prog. Ser. 297: 245-257.).

Previous studies based on MEDITS data focusing on species diversity did not find an effect of a longitudinal gradient on diversity (Gaertner et al. 2007Gaertner J.C., Bertrand J., Relini G., et al. 2007. Spatial pattern in species richness of demersal fish assemblages on the continental shelf of the northern Mediterranean Sea: a multiscale analysis. Mar. Ecol. Prog. Ser. 341: 191-203., Granger et al. 2015Granger V., Fromentin J.M., Bez N., et al. 2015. Large-scale spatio-temporal monitoring highlights hotspots of demersal fish diversity in the Mediterranean Sea. Prog. Oceanogr. 130: 65-74., Brind’Amour et al. 2016Brind’Amour A., Rochet M.J., Ordines F., et al. 2016. Environmental drivers explain regional variation of changes in fish and invertebrate functional groups across the Mediterranean Sea from 1994 to 2012. Mar. Ecol. Prog. Ser. 562: 19-35.), in contrast to results provided by studies based on other types of data (presence-absence, coarser spatial grain, etc. see Coll et al. 2010Coll M., Piroddi C., Steenbeek J., et al. 2010. The Biodiversity of the Mediterranean Sea: Estimates, Patterns, and Threats. PLoS ONE 5: e11842., 2012Coll M., Piroddi C., Albouy C., et al. 2012. The Mediterranean Sea under siege: spatial overlap between marine biodiversity, cumulative threats and marine reserves. Glob. Ecol. Biogeogr. 21: 465-480., Mouillot et al. 2011Mouillot D., Albouy C., Guilhaumon F., et al. 2011. Protected and threatened components of fish biodiversity in the Mediterranean Sea. Curr. Biol. 21: 1044-1050.). Notably, using standardized data, Granger et al. (2015)Granger V., Fromentin J.M., Bez N., et al. 2015. Large-scale spatio-temporal monitoring highlights hotspots of demersal fish diversity in the Mediterranean Sea. Prog. Oceanogr. 130: 65-74. found that species richness showed relatively low values in the Adriatic Sea and Cyprus at both local and regional scales (i.e. alpha and gamma diversity, respectively), similarly to specificities we evidenced in assemblage organization in these two GSAs. They also highlighted that diversity indices remained stable over the last two decades. This result is in line with our findings that trawling pressure, estimated by the MCDA, had a low effect on assemblage composition at the scale of our study. This does not mean that intensive trawling, such as that practised in numerous GSAs of the northern Mediterranean Sea for more than 50 years, had no effect on the benthic ecosystems and the associated fauna (e.g. D’Onghia et al. 2003D’Onghia G., Mastrototaro F., Matarrese A., et al. 2003. Biodiversity of the Upper Slope Demersal Community in the Eastern Mediterranean: Preliminary Comparison Between Two Areas With and Without Trawl Fishing. J. Northw. Atl. Fish. Sci. 31: 263-273., 2005D’Onghia G., Capezzuto F., Mytilineou Ch., et al. 2005. Comparison of the population structure and dynamics of Aristeus antennatus (Risso, 1816) between exploited and unexploited areas in the Mediterranean Sea. Fish. Res. 76: 22-38., Farriols et al. 2017Farriols M.T., Farriols F., Somerfield P.J., et al. 2017. Bottom trawl impacts on Mediterranean demersal fish diversity: Not so obvious or are we too late? Cont. Shelf Res. 137: 84-102.). It may suggest that most of the impact of trawling pressure had probably altered the demersal fish diversity and composition before the beginning of the studied period (1994 for the MEDIT surveys) (Granger et al. 2015Granger V., Fromentin J.M., Bez N., et al. 2015. Large-scale spatio-temporal monitoring highlights hotspots of demersal fish diversity in the Mediterranean Sea. Prog. Oceanogr. 130: 65-74., Farriols et al. 2017Farriols M.T., Farriols F., Somerfield P.J., et al. 2017. Bottom trawl impacts on Mediterranean demersal fish diversity: Not so obvious or are we too late? Cont. Shelf Res. 137: 84-102.). Otherwise, possible uncertainty in estimating the trawling pressure spatialization in our study, due to the lack of spatial fishing effort data (such as from VMS data, not available at the spatio-temporal scale of this study), cannot be ruled out. Thus, to be able to establish more explicit links with fishery management perspective, the availability of finer spatio-temporal fishing effort data, such as VMS data, should be broaden and increased for the different GSAs. It may provide more direct assessment of the relationship between demersal assemblages and trawling pressure (e.g. Farriols et al. 2017Farriols M.T., Farriols F., Somerfield P.J., et al. 2017. Bottom trawl impacts on Mediterranean demersal fish diversity: Not so obvious or are we too late? Cont. Shelf Res. 137: 84-102.), and thus facilitate suggestions towards management.

Secondly, when considering other regions of the world, a recent study was based on two decades of scientific trawl surveys focused on exploited demersal fish assemblages in Mauritania (eastern Atlantic Ocean) among different depth strata and latitudinal areas (Kidé et al. 2015Kidé O., Manté C., Dubroca L., et al. 2015. Spatio-temporal dynamics of exploited ground-fish assemblages faced to environmental and fishing forcings: Insights from the Mauritanian Exclusive Economic Zone. Plos ONE 10: e0141566.). Similarly to our results, they highlighted a stronger effect of environmental conditions than trawling on the structure and composition of demersal fish assemblages. They found that chlorophyll a and sea surface temperature mainly influenced assemblages in areas where upwellings occurred. Effects of trawling pressure on the structure and composition of assemblages were relatively low. Similarly to our results, the temporal trajectories between environmental and trawling conditions and assemblages did not match over the entire time series they analyzed (only for some specific years and areas).

However, it should be noted that in other regions, targeted species were affected more directly by fishing and led to changes in the structure of demersal assemblages (e.g. Ansari et al. 1995Ansari Z.A., Chatterji A., Ingole B.S., et al. 1995. Community Structure and Seasonal Variation of an Inshore Demersal Fish Community at Goa, West Coast of India. Estuar. Coast. Shelf. S. 41: 593-610., Greenstreet and Hall 1996Greenstreet S.P.R., Hall S.J. 1996. Fishing and the ground-fish assemblage structure in the north-western North Sea: an analysis of long-term and spatial trends. J. Anim. Ecol. 65: 577-598., Levin et al. 2006Levin P.S., Holmes E.E., Piner K.R., et al. 2006. Shifts in a Pacific Ocean Fish Assemblage: the Potential Influence of Exploitation. Conser. Biol. 20: 1181-1190.). Regarding temporal trajectories among the environmental conditions and the species, the fact that for the large majority of GSAs a mismatch was observed across years (Fig. 5, Table S6) could be linked to (i) lags in response of assemblages in the face of changes in environmental conditions, (ii) approximate way that environmental values were computed (i.e. monthly mean across the year before the haul was performed) and/or (iii) the fact that assemblages could be sustained despite environmental changes. It may also suggest that factors other than those we have investigated may act on fish assemblages, such as biotic interactions.

In conclusion, our study provides a state of the spatio-temporal relationships between organizational patterns of demersal fish assemblages and several environmental and trawling descriptors, on the basis of a large standardized data set of scientific trawl samplings. Our analysis represents a step towards quantitative studies focused on the drivers of demersal fish diversity in the Mediterranean Sea. It suggests a strong inter-zone stability in the organization of fish assemblages, with specificities for some zones. The most stable structuring factors were linked to a combined gradient of chlorophyll a, phytoplancton carbon biomass and temperature, inversely correlated with a depth, nutrient and salinity gradient. A clear gradient in the distribution of species linked to these environmental gradients was evidenced as the stable part of their relationships. Investigating explicitly potential biotic factors (such as species aggregation, competition or dispersal limitations) would extend our findings within the framework of complementary analyses dedicated to this aspect.

ACKNOWLEDGEMENTSTop

This study has been carried out as part of the MEDITS survey program within the Data Collection Framework (DCF). The European Union and the Member States involved in the DCF are gratefully acknowledged. We are grateful to all colleagues who have participated to the MEDITS trawl surveys since 1994. We thank three anonymous reviewers for their constructive comments. The manuscript has been proof-read and corrected for the English by Michael Paul.

REFERENCESTop

Ansari Z.A., Chatterji A., Ingole B.S., et al. 1995. Community Structure and Seasonal Variation of an Inshore Demersal Fish Community at Goa, West Coast of India. Estuar. Coast. Shelf Sci. 41: 593-610.
https://doi.org/10.1016/0272-7714(95)90029-2

Ben Rais Lasram F., Guilhaumon F., Albouy C., et al. 2010. The Mediterranean Sea as a “cul-de-sac” for endemic fishes facing climate change. Glob. Chang. Biol. 16: 3233-3245.
https://doi.org/10.1111/j.1365-2486.2010.02224.x

Bertrand J.A., De Sola L., Papaconstantinou C., et al. 2002a. The general specifications of the MEDITS surveys. Sci. Mar. 66: 9-17.
https://doi.org/10.3989/scimar.2002.66s29

Bertrand J.A., Leonori I., Dremiere P.Y., et al. 2002b. Depth trajectory and performance of a trawl used for an international bottom trawl survey in the Mediterranean. Sci. Mar. 66: 169-182.
https://doi.org/10.3989/scimar.2002.66s2169

Brind’Amour A., Rochet M.J., Ordines F., et al. 2016. Environmental drivers explain regional variation of changes in fish and invertebrate functional groups across the Mediterranean Sea from 1994 to 2012. Mar. Ecol. Prog. Ser. 562: 19-35.
https://doi.org/10.3354/meps11912

Carassou L., Ponton D. 2007. Spatio-temporal structure of pelagic larval and juvenile fish assemblages in coastal areas of New Caledonia, southwest Pacific. Mar. Biol. 150: 697-711.
https://doi.org/10.1007/s00227-006-0389-y

Certain G., Masse J., Van Canneyt O., et al. 2011. Investigating the coupling between small pelagic fish and marine top predators using data collected from ecosystem-based surveys. Mar. Ecol. Prog. Ser. 422: 23-39.
https://doi.org/10.3354/meps08932

Charrad M., Ghazzali N., Boiteau V., et al. 2014. An R Package for Determining the Relevant Number of Clusters in a Data Set. J. Stat. Softw. 61: 1-36.
https://doi.org/10.18637/jss.v061.i06

Coll M., Piroddi C., Steenbeek J., et al. 2010. The Biodiversity of the Mediterranean Sea: Estimates, Patterns, and Threats. PLoS ONE 5: e11842.
https://doi.org/10.1371/journal.pone.0011842

Coll M., Piroddi C., Albouy C., et al. 2012. The Mediterranean Sea under siege: spatial overlap between marine biodiversity, cumulative threats and marine reserves. Glob. Ecol. Biogeogr. 21: 465-480.
https://doi.org/10.1111/j.1466-8238.2011.00697.x

Colloca F., Spedicato M.T., Massutí E., et al. 2013. Mapping of nursery and spawning grounds of demersal fish. Mediterranean Sensitive Habitats (MEDISEH) Final Report, DG MARE Specific Contract SI2.600741, Heraklion (Greece).

D’Onghia G., Mastrototaro F., Matarrese A., et al. 2003. Biodiversity of the Upper Slope Demersal Community in the Eastern Mediterranean: Preliminary Comparison Between Two Areas With and Without Trawl Fishing. J. Northw. Atl. Fish. Sci. 31: 263-273.
https://doi.org/10.2960/J.v31.a20

D’Onghia G., Capezzuto F., Mytilineou Ch., et al. 2005. Comparison of the population structure and dynamics of Aristeus antennatus (Risso, 1816) between exploited and unexploited areas in the Mediterranean Sea. Fish. Res. 76: 22-38.
https://doi.org/10.1016/j.fishres.2005.05.007

Dolédec S., Chessel D. 1994. Co-inertia analysis: an alternative method for studying species-environment relationships. Freshw. Biol. 31: 277-294.
https://doi.org/10.1111/j.1365-2427.1994.tb01741.x

Dray S., Chessel D., Thioulouse J. 2003. Co-inertia analysis and the linking of ecological data tables. Ecology 84: 3078-3089.
https://doi.org/10.1890/03-0178

Farriols M.T., Farriols F., Somerfield P.J., et al. 2017. Bottom trawl impacts on Mediterranean demersal fish diversity: Not so obvious or are we too late? Cont. Shelf Res. 137: 84-102.
https://doi.org/10.1016/j.csr.2016.11.011

Gaertner J.C., Bertrand J.A., Souplet A. 2002. STATIS-CoA: A methodological solution to assess the spatio-temporal organization of species assemblages. Application to the demersal assemblages of the French Mediterranean Sea. Sci. Mar. 66: 221-232.
https://doi.org/10.3989/scimar.2002.66s2221

Gaertner J.C., Bertrand J.A., De Sola L.G., et al. 2005. Large spatial scale variation of dem- ersal fish assemblage structure on the continental shelf of the NW Mediterranean Sea. Mar. Ecol. Prog. Ser. 297: 245-257.
https://doi.org/10.3354/meps297245

Gaertner J.C., Bertrand J., Relini G., et al. 2007. Spatial pattern in species richness of demersal fish assemblages on the continental shelf of the northern Mediterranean Sea: a multiscale analysis. Mar. Ecol. Prog. Ser. 341: 191-203.
https://doi.org/10.3354/meps341191

Gaertner J.C., Mérigot B., Rélini G., et al. 2010. Reproducibility of the multi-component aspect of species diversity across different areas and scales: towards the constitution of a shortlist of complementary indices for monitoring fish diversity? Ecography 33: 1123-1135.
https://doi.org/10.1111/j.1600-0587.2009.06259.x

Gaertner J.C., Maiorano P., Mérigot B., et al. 2013. Large-Scale Diversity of Slope Fishes: Pattern Inconsistency between Multiple Diversity Indices. PLoS ONE 8: e66753.
https://doi.org/10.1371/journal.pone.0066753

Gallego-Álvarez I., Galindo-Villardón M., Rodríguez-Rosa M. 2015. Evolution of sustainability indicator worldwide: A study from the economic perspective based on the X-STATICO method. Ecol. Indic. 58: 139-151.
https://doi.org/10.1016/j.ecolind.2015.05.025

Granger V., Fromentin J.M., Bez N., et al. 2015. Large-scale spatio-temporal monitoring highlights hotspots of demersal fish diversity in the Mediterranean Sea. Prog. Oceanogr. 130: 65-74.
https://doi.org/10.1016/j.pocean.2014.10.002

Greenstreet S.P.R., Hall S.J. 1996. Fishing and the ground-fish assemblage structure in the north-western North Sea: an analysis of long-term and spatial trends. J. Anim. Ecol. 65: 577-598.
https://doi.org/10.2307/5738

Hill M.O. 1973. Reciprocal averaging: an eigenvector method of ordination, J. Econ. 61: 237-249.
https://doi.org/10.2307/2258931

Kavadas S., Maina I., Damalas D., et al. 2015. Multi-Criteria Decision Analysis as a tool to extract fishing footprints and estimate fishing pressure: application to small scale coastal fisheries and implications for management in the context of the Maritime Spatial Planning Directive. Mediterr. Mar. Sci. 16: 294-304.
https://doi.org/10.12681/mms.1087

Kidé O., Manté C., Dubroca L., et al. 2015. Spatio-temporal dynamics of exploited ground-fish assemblages faced to environmental and fishing forcings: Insights from the Mauritanian Exclusive Economic Zone. PloS ONE 10: e0141566.
https://doi.org/10.1371/journal.pone.0141566

Lavit C., Escoufier Y., Sabatier R., et al. 1994 The ACT (Statis method). Comput. Stat. Data Anal. 18: 97-119.
https://doi.org/10.1016/0167-9473(94)90134-1

Le Fur I., De Wit R., Plus M., et al. 2019. Re-oligotrophication trajectories of macrophyte assemblages in Mediterranean coastal lagoons based on 17-year time-series. Mar. Ecol. Prog. Ser. 608: 13-32.
https://doi.org/10.3354/meps12814

Lefcheck J.S., Buchheister A., Laumann K.M., et al. 2014. Dimensions of biodiversity in Chesapeake Bay demersal fishes: patterns and drivers through space and time. Ecosphere 5: 14.
https://doi.org/10.1890/ES13-00284.1

Levin P.S., Holmes E.E., Piner K.R., et al. 2006. Shifts in a Pacific Ocean Fish Assemblage: the Potential Influence of Exploitation. Conser. Biol. 20: 1181-1190.
https://doi.org/10.1111/j.1523-1739.2006.00400.x

Malczewski J. 2006. GIS-based multi criteria decision analysis: a survey of the literature. Int. J. Geogr. Inf. Sci. 20: 703-726.
https://doi.org/10.1080/13658810600661508

Mazzocchi M.G., Dubroca L., García-Comas C., et al. 2012. Stability and resilience in coastal copepod assemblages: The case of the Mediterranean long-term ecological research at Station MC (LTER-MC). Prog. Oceanogr. 97-100: 135-151.
https://doi.org/10.1016/j.pocean.2011.11.003

MEDITS Working Group. 2017. MEDITS-Handbook, Version n. 9. MEDITS Working Group, 106 pp.
http://www.sibm.it/MEDITS%202011/principaledownload.htm

Mendes S., Fernández-Gómez M.J., Pereira M.J., et al. 2012. An empirical comparison of canonical correspondence analysis and STATICO in the identification of spatio-temporal ecolog- ical relationships. J. App. Stat. 39: 979-994.
https://doi.org/10.1080/02664763.2011.634393

Menge B.A., Olson A.M. 1990. Role of scale and environmental factors in the regulation of community structure. Trends Ecol. Evol. 5: 52-57.
https://doi.org/10.1016/0169-5347(90)90048-I

Mérigot B., Durbec J.P., Gaertner J.C. 2010. On goodness-of-fit measure for dendrogram-based analyses. Ecology 91: 1850-1859.
https://doi.org/10.1890/09-1387.1

Mouillot D., Albouy C., Guilhaumon F., et al. 2011. Protected and threatened components of fish biodiversity in the Mediterranean Sea. Curr. Biol. 21: 1044-1050.
https://doi.org/10.1016/j.cub.2011.05.005

Oddo P., Adani M., Pinardi N., et al. 2009. A nested Atlantic-Mediterranean Sea general circulation model for operational forecasting. Ocean Sci. 5: 461-473.
https://doi.org/10.5194/os-5-461-2009

R Core Team. 2015. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. Available online at https://www.R-project.org/

Relini G. 2015. Fish biodiversity in MEDITS surveys. Biol. Mar. Mediterr. 22: 176-177.

Rochet M.J., Trenkel V.M., Carpentier A., et al. 2010. Do changes in environmental and fishing pressures impact marine communities? An empirical assessment. J. Appl. Ecol. 47: 741-750.
https://doi.org/10.1111/j.1365-2664.2010.01841.x

Saaty T.L. 1980. The analytical hierarchy process: planning, priority setting, resource allocation. McGraw-Hill, New York, 309 pp.

Simier M., Blanc L., Pellegrin F. et al. 1999. Approche simultanée de k couples de tableaux: Application à l’étude des relations pathologie végétale-environnement. Rev. Statist. Appl. 47: 31-46.

Simier M., Laurent C., Ecoutin J.M., et al. 2006. The Gambia River estuary: A reference point for estuarine fish assemblages studies in West Africa. Estuar. Coast. Shelf Sci. 69: 615-628.
https://doi.org/10.1016/j.ecss.2006.05.028

Slimani N., Guilbert E., El Ayni F., et al. 2017. The use of STATICO and COSTATIS, two exploratory three-ways analysis methods. An application to the ecology of aquatic Heteroptera in the Medjerda watershed (Tunisia). Env. Ecol. Stat. 24: 269-295.
https://doi.org/10.1007/s10651-017-0370-6

Thioulouse J. 2011. Simultaneous analysis of a sequence of paired ecological tables: a comparison of several methods. Ann. App. Stat. 5: 2300-2325.
https://doi.org/10.1214/10-AOAS372

Thioulouse J., Simier M., Chessel D. 2004. Simultaneous analysis of a sequence of paired ecological tables. Ecology 85: 272-283.
https://doi.org/10.1890/02-0605

SUPPLEMENTARY MATERIAL

The following supplementary material is available through the online version of this article and at the following link:
http://scimar.icm.csic.es/scimar/supplm/sm04954esm.pdf

Table S1. – Number (n) of hauls analyzed per Geographical Sub-Areas (GSA). N=18062 hauls sampled between 1999 and 2015 within the 17 GSAs (see Fig.1 for map of locations and GSA names).

Table S2. – List of the 154 species included in the study. Species codes used in figures corresponds to the first four letters of the genus name and the first three letters of the species name.

Fig. S1. – Flowchart of the Multi-Criteria Decision Analysis (MCDA) used to compute the Fishing Pressure Index (FPI). AHP, Analytic Hierarchy Process; WLC, Weighted Linear Combination; FM, linear Fuzzy Membership function; GT, Gross Tonnage.

Fig. S2. – Fishing Pressure Index (FPI). A, maps for years 1994, 2004 and 2014. B, pairwise relationships of FPI among the three years 1994, 2004 and 2014. FPI values were generated for haul locations provided in Fig. 1 of the article. The black line represents the first bisector for which x=y (i.e same FPI values among the two years considered). See the Material and Methods section of the article for more details about the FPI computation.

Fig. S3. – Projection of each Geographical Sub-Area (GSA) on the two first factorial axes of the STATICO-CoA inter-structure (axis 1: horizontal, 51.56%; axis 2: vertical, 7.51%), with barplot showing the eigenvalues of each axis. Correspondence between GSA numbers and their names is given in Figure 1. Note that for this kind of multi-tables analysis, only axis 1 of the inter-structure has a meaning for the construction of the compromise, and is thus used to determine the contribution of tables to the compromise. The other axes are not used for that purpose, and should not be interpreted (axis 2 is used for display purpose).

Fig. S4. – Projection of the four factorial axes (arrows) of the separate Correspondence Analysis of each Geographical Sub-Area (GSA) on the two first factorial axes of the STATICO-CoA compromise (axis 1: horizontal; axis 2: vertical). Correspondence between GSA numbers and their names is given in Figure 1.

Fig. S5. – Draftsman plot of the eight environmental and trawling variables (Spearman correlation coefficient).

Fig. S6. – Plots for each Geographical Sub-Area (GSA) of the STATICO-CoA analysis: in the left column, projection of the eight environmental and trawling variables on the first factorial plane, and in the right column, projection of the average positions of species. Correspondence between GSA numbers and their names is given in Figure 1.

Fig. S7. – Projections of environmental and trawling variables on the compromise of the STATICO-CoA analysis based on a data set of 71 species (without 83 rare species, i.e. those present in less than 5% of the hauls). It shows the stable part of the species-environment- trawling relationships. Axis 1 and 2 explained 94.63% and 4.43% of the total variability, respectively. Species names are provided in Table S2.

Fig. S8. – Projections of demersal species variables and assemblage groups on the compromise of the STATICO-CoA analysis based on data set of 71 species (without 83 rare species, i.e. those present in less than 5% of the hauls). It shows the stable part of the species and environment-fishing relationships. Axis 1 and 2 explained 94.63% and 4.43% of the total variability, respectively. On the first factor plane are shown the two main species assemblages (in blue and red) obtained by hierarchical clustering (UPGMA criterion) and optimum number of groups approach (see Materials and Methods section). Species codes are given in Table S2.