Spatio-temporal trends in diversity of demersal fish assemblages in the Mediterranean

1 Intituto Español de Oceanografía, Centre Oceanogràfic de les Balears, Moll de Ponent s/n, 07015 Palma de Mallorca, Illes Baleares, Spain. (MTF) (corresponding author) E-mail: mt.farriols@ieo.es. ORCID iD: https://orcid.org/0000-0002-7704-6504 (FO) E-mail: xisco.ordinas@ieo.es. ORCID iD: https://orcid.org/0000-0002-2456-2214 (EM) E-mail: enric.massuti@ieo.es. ORCID iD: https://orcid.org/0000-0002-9524-5873 2 COISPA Tecnologia and 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 (LC) E-mail: casciaro@coispa.eu. ORCID iD: https://orcid.org/0000-0003-4876-9874 3 CNR-IAMC, Mazara del Vallo, Via Vaccara 61, Mazara del Vallo, TP, Italy. (MDL) E-mail: manfredi.dilorenzo@libero.it. ORCID iD: https://orcid.org/0000-0003-3786-5772 4 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 5 Dipartimento di Biologia Animale ed Ecologia, Università di Cagliari, Cagliari, Italy. (CF) E-mail: follesac@unica.it. ORCID iD: https://orcid.org/0000-0001-8320-9974 6 Instituto Español de Oceanografía (IEO), Centro Oceanográfico de Málaga, Fuengirola, Málaga, Spain. (CG-R) E-mail: cristina.garcia@ieo.es. ORCID iD: https://orcid.org/0000-0003-2767-4200 7 Institute of Oceanography and Fisheries Split, Šetalište Ivana Meštrovića 63, 21000 Split, Croatia. (II) E-mail: igor@izor.hr. ORCID iD: https://orcid.org/0000-0001-7101-9575 8 Institut Français de Recherche pour l’Exploitation de la Mer (Ifremer), UMR 212 Ecosystèmes Marins Exploités (EME), Centre de Recherche halieutique Méditerranéenne et Tropicale, 34203 Sète, France. (AJ) E-mail: Angelique.Jadaud@ifremer.fr. ORCID iD: https://orcid.org/0000-0001-6858-3570 9 Centro Interuniversitario di Biologia Marina ed Ecologia Applicata, Viale Nazario Sauro 4, 57128 Leghorn, Italy. (AL) E-mail: ligas@cibm.it. ORCID iD: https://orcid.org/0000-0003-1036-3553 10 Laboratorio di Biologia Marina e Pesca di Fano, Dip.to B.E.S., Università di Bologna, Fano, Italy. (CM) E-mail: chiara.manfredi3@unibo.it. ORCID iD: https://orcid.org/0000-0002-2852-4856 11 Fishery Research Institute of Slovenia, Ljubljana-Smartno, Slovenia. (BM) E-mail: bojan.Marceta@zzrs.si. ORCID iD: https://orcid.org/0000-0001-8137-4474 12 Hellenic Center for Marine Research, Iraklion, Crete, Greece. 13 University of Crete, Biology Department, Stavrakia, Heraklion, Crete. (PP) E-mail: notap@hcmr.gr. ORCID iD: https://orcid.org/0000-0002-8608-078X 14 Institute of Oceanography and Fisheries, Split, Croatia. (NV) E-mail: vrgoc@izor.hr. ORCID iD: https://orcid.org/0000-0002-528-4512


INTRODUCTION
The Mediterranean is considered one of the Large Marine Ecosystems of the world, owing to its bathymetry, hydrography, productivity and trophic webs (Duda and Sherman 2002). It is a semi-enclosed sea connected to the Atlantic Ocean through the Gibraltar Strait, to the Black Sea through the Dardanelles Strait and to the Red Sea through the artificial Suez Channel (Fig. 1). It acts as a concentration basin, and evaporation is higher in its eastern basin, causing the water level to decrease and salinity to increase from west to east (Coll et al. 2010). While temperature also increases eastwards (Coll et al. 2010), surface productivity, organic matter availability at the seafloor and the biomass of megabenthic fauna of deep ecosystems decrease eastwards (Bosc et al. 2004, Danovaro et al. 1999, Tecchio et al. 2011. The Mediterranean has narrow continental shelves and a large area of open sea. In fact, the continental shelf covers about 20% of the Mediterranean bottoms, whereas the slope covers about 60% (Sardà et al. 2004). Therefore, a large part of this basin can be classified as deep sea (Coll et al. 2010).
The high species richness, coupled with a high proportion of endemism, makes the Mediterranean one of the world's 'biodiversity hotspots' , Coll et al. 2010, Lejeusne et al. 2010. Environmental variables such as temperature, productivity and distance from the Strait of Gibraltar have been shown to be causes of fish species richness distribution (Ben Rais Lasram et al. 2009, Meléndez et al. 2017. However, this high biodiversity is presently threatened by the combined action of anthropogenic impacts, introduction of alien species and climate change (Bianchi et al. 2012). Among human activities, fisheries are one of the most important factors affecting marine resources and ecosystems.
It is well known that fisheries have profoundly modified the structure of marine ecosystems (Dayton et al. 1995, Hall 1999, Kaiser and de Groot 2000. Effects of fishing on marine ecosystems include shifts in the food-web structure due to changes in predator-prey relationships (Kaiser et al. 2002); changes in size structure due to vulnerability and selection of fishing for large individuals (Gislason 2002, Jennings and Dulvy 2005, Daan et al. 2005; genetic selection of species with particular life-history traits, such as a higher growth rate and earlier age-at-maturity Fonteneau 2001, Jørgensen et al. 2007); changes in the spatial distribution of target species (e.g. Ciannelli et al. 2013); effects on the population of non-target species (Pranovi et al. 2001, Ordines et al. 2014; and decrease of habitat complexity and changes on the benthic community structure (e.g. Callaway et al. 2002).
The natural resources of the Mediterranean have been subject to human exploitation since ancient times, when coastal communities started to use different fishing gears, some of which are still in use (Farrugio et al. 1993). Dramatic long-term changes in marine communities took place before the industrialization of fisheries that occurred in the 1950s, and have already been documented in some areas, such as the Adriatic Sea (Fortibuoni et al. 2010). Until 1950, the exploitation of Mediterranean resources was limited to fishing areas shallower than 200 m depth. In the last few decades, with the decline of stocks on the continental shelf, increasing market demand and the introduction of new technologies, trawl fisheries have expanded offshore and towards the deeper waters of the continental slope (Roberts 2002, Morato et al. 2006 to target valuable resources such as red shrimps (e.g. Demestre and Martín 1993, Guijarro et al. 2008, Masnadi et al. 2018. In this area, fisheries are assessed within the framework of the General Fisheries Commission for the Mediterranean (GFCM), the regional fisheries management organization of the Mediterranean. Of the 27 Mediterranean stocks of fishing target species assessed by the GFCM in its last report, about 80% were considered overexploited (GFCM 2016). The presence of a high diversity of species and the absence of large monospecific stocks comparable to those inhabiting some wide areas of the open oceans are characteristic desembocado en una sobrexplotación de sus principales stocks comerciales. A través de datos independientes de las pesquerías recogidos en el marco de las campañas MEDITS desarrolladas durante las dos últimas décadas se ha estudiado el patrón de diversidad de peces demersales en el Mediterráneo a través de largas escalas temporales y espaciales para evaluar si este patrón se ve afectado por el estado general de sobrexplotación de sus recursos demersales. A continuación las tendencias detectadas en la diversidad han sido comparadas a la variación espacio-temporal del esfuerzo de la pesca de arrastre a través del Mediterráneo. Nuestros resultados muestran una estabilidad e incluso recuperación de la diversidad de peces demersales en el Mediterráneo junto a valores altos de diversidad en las plataformas continentales de las Islas Baleares, Cerdeña, Sicília y el mar Egeo. La alta diversidad de las asociaciones de peces demersales a escala tanto espacial como temporal está asociada a una reducción del esfuerzo pesquero. La inclusión de especies distintas a las objetivo a través de índices de diversidad es relevante en la implementación de la aproximación ecosistémica a la gestión de las pesquerías. of the Mediterranean demersal fisheries (Farrugio et al. 1993). Assessment at a community level is therefore crucial, particularly due to the multispecies nature of the bottom trawl fishery and also because a decline in the diversity of demersal assemblages has been reported due to fishing exploitation (e.g. Ungaro et al. 1998, Sabatini et al. 2013, Farriols et al. 2017). Assessment at a community level is also a requirement for the implementation of an ecosystem-based management of fisheries (Browman and Stergiou 2004).
The aim of this work is to study the demersal fish diversity pattern in the Mediterranean at a large spatial and temporal scale and to assess whether this pattern is being affected by the general fishing overexploitation of demersal resources in the area. To do so, we used fishery-independent data collected under the framework of the MEDITS trawl surveys carried out during the last 20 years. The detected trends were compared with the spatio-temporal variation in bottom trawl fishing effort in the Mediterranean Sea.

Data
Demersal fish were collected during MEDITS bottom trawl surveys conducted from 1994 to 2015 in 14 geographical sub-areas (GSAs) along the European coasts of the Mediterranean Sea. Some GSAs have gaps in their sampling years: i) GSA 5 started sampling in 2001; ii) there are no data in 2002 for GSA 8 (technical problem of the boat); and iii) there are no data for 2002 Sampling was performed during spring-summer in daylight hours using the GOC73 experimental gear, whose efficiency for catching demersal species has been tested by Fiorentini et al. (1999) and Dremière et al. (1999). For more details about the sampling strategy and protocol see Bertrand et al. (2002) and the MEDITS hand-book, instruction manual version 9 (MEDITS Working Group 2017).
A total of 17540 hauls performed between 46 and 866 m depth were analysed (Table 1, Fig. 1). Hauls shallower than 46 m depth were excluded from the analysis because they could not be found for all GSAs. The catch of each sample was sorted, identified to species level, counted, weighed and standardized to square km by using the horizontal opening of the net and the distance covered in each haul. Species with a pelagic or mesopelagic behaviour, like most species of the families Myctophydae (e.g. Ceratoscopelus maderensis) and Cupleidae (e.g. Engraulis encrasicolus), were excluded from the analysis. A species accumulation curve for each GSA was performed and we confirmed that differences in number of species were not due to differences in the number of hauls considered for each GSA (Table 1, Fig. S1).

Fish assemblages and diversity
Cluster analysis was used to analyse the structure of demersal fish assemblages and to identify different assemblages according to depth strata in each GSA. Rela-  tionships among samples were detected by hierarchical agglomerative clustering with group-average linkage after a forth root transformation of the data. The distance used to make groups was the Bray-Curtis similarity. These analyses were performed using PRIMER 7 (Clarke et al. 2014). The calculus of diversity indices explained below was made taking into account the groups of samples obtained from the cluster analysis. The N 90 diversity index was calculated following the procedure described by Farriols et al. (2015). It is the mean number of species contributing up to 90% of within-group similarity calculated from abundance data expressed as N km -2 and assigned a priori to groups. The calculation of N 90 starts with the calculation of the contribution of each species to the withingroup similarity using the Bray-Curtis similarity index (Bray and Curtis 1957), as proposed by Clarke (1993): where y ij is the abundance of the species i at the sampling site j, y ik is the abundance of the species i at the sampling site k, p is the total number of species in j and k, and min (y ij , y ik ) is the minimum value of the abundance of species i between the sampling sites j and k, also considering zeros.
The contribution of each species i to the total similarity of the group S i is the mean value of S jk (i) for the assigned group, and the total similarity in a group (Sim) is the addition of S i for all the species in the group: Then, the contribution of S i is calculated as a percentage of Sim. Species contributions are calculated for each re-sampling in a jack-knife routine, which removes a number of samples each time, producing lists of contribution to similarity by species in each re-sampling. Because the groups of samples for each GSA, strata and year were large, we removed 10% of the samples in each re-sampling with a 50% replacement. That is, 50% of samples removed in a re-sampling must be different from previous ones. In this way, we obtain values of deviation for N 90 other than 0 for groups with a large number of observations. The N 90 diversity index is the mean number of species which accumulate up to 90% of within-group similarity in all the re-samplings. N 90 was calculated using R scripts, version 3.1.1 (R Core Team 2014). Similarity percentage analysis (SIMPER) for each group of samples was also undertaken to see their species composition.
Diversity indices, such as species richness (S) and Pielou evenness (J'), which have shown some kind of response to fishing impact for demersal fish assemblages in the Mediterranean (Farriols et al. 2017), were also included in this work. These traditional indices are also helpful for comparison with previous works. S is the raw number of species in each haul and J' was calculated as follows: where p i is the proportion of all individuals belonging to species i and S is the total number of species in the sample.

Fishing effort
Information on annual fishing effort was collected from the working group reports of the GFCM (http:// www.fao.org/gfcm/data/safs/en/) and the Scientific, Technical and Economic Committee for Fisheries (STECF, https://stecf.jrc.ec.europa.eu/reports/medbs). Fishing effort data were compiled by trawl fleet targeting different species. The units vary between the different reports, being mainly provided in terms of number of vessels, kilowatts per days at sea and gross tonnage per days at sea (see Table S1).
To estimate fishing effort in each depth stratum obtained from cluster analysis, the strata were associated with the main target species of the fleets. Because target species varied between GSAs, we considered i) Mullus barbatus or Mullus surmuletus for the continental shelf; ii) Nephrops norvegicus or Parapenaeus longirostris for the shelf break/upper slope; and iii) Aristeus antennatus or Aristaeomorpha foliacea for the lower slope.
To compare temporal trends in fishing effort and demersal fish diversity, the longest series of fishing effort available for each GSA and depth stratum regardless of the kind of units were selected. When we found no values of fishing effort for a certain GSA, experts were contacted to obtain a trend in number of vessels in that area.

Temporal and spatial analysis
In order to analyse temporal trends in diversity, linear regressions were fitted to the mean values of S, J' and N 90 for each year, GSA and depth stratum. Linear regression analyses with the annual values of fishing effort in each GSA (see Fishing effort section in Materials and Methods) and depth stratum were also performed. The exploration of the scatter plots of the time series together with the comparison of Pearson (assuming linear pattern) and Spearman (suitable also for other monotonic patterns than the linear) correlation coefficients were done. The values of both correlation coefficients were similar, indicating that the detected trends could be fitted using a simple linear model. Thus, the linear regression and the Pearson coefficient of correlation were presented along with the coefficient of determination (i.e. variance explained). These analyses were carried out with R, version 3.1.1 (R Core Team 2014).
In order to observe spatial differences in diversity by GSA, time series of mean values and standard deviation of each diversity index (see Data section for years included in each GSA) were plotted. For those series with a significant temporal trend, the diversity values at the beginning and the end of the time series were plotted instead of mean values and standard deviation.
SIMPER analysis for each group of samples from N 90 was also performed to see differences in species composition in each GSA. The percentage of contribution of each species to within-group similarity was calculated as the mean value of species contributions to similarity, taking all groups of observations by year and stratum for each GSA into account.

Community structure
Results from cluster analysis detecting main fish assemblages for each GSA are shown in Figure 2. Three groups of samples were selected from most GSAs, corresponding to a level of similarity of between 30% and 40%. Maximum, minimum and mean depths of each cluster group per GSA were obtained. According to these depth values, samples were grouped in three different depth strata: shelf, shelf break/upper slope and lower slope (Table 2,Fig. 2). For GSAs 17 and 23 only two groups were selected. GSA 17 has no samples below 350 m, and GSA 23 has a negligible sample number over 496 m (Table 2). For GSAs 7 and 20 samples in the lower slope group were not enough to calculate the N 90 throughout the time series, so both lower slope groups were omitted from the temporal and spatial analysis.
In 9 out of the 12 GSAs showing lower slope samples, samples from the shelf break/upper slope clustered with samples from the lower slope. The exceptions were GSAs 6, 11 and 22, where samples from the shelf break/ upper slope clustered with those from the continental shelf. Minimum, maximum and mean depths for each group of samples from cluster analysis are shown in Table 2. Mean depth of continental shelf samples ranged from 76 m in GSA 1 to 125 m in GSA 10, while for the shelf break/upper slope they ranged from 180 m in GSA 17 to 421 m in GSA 7, and for the lower slope between 496 m in GSA 11 and 699 m in GSA 7.

Temporal trends
Although the analysis of temporal evolution for N 90 , S and J' showed no significance in most GSAs and depth strata (Table 3, (Table 3).
When quantitative analysis in temporal evolution of fishing effort could be made, the detected significance mainly showed a decreasing trend (Table 4,Fig. 3). This is the case of the continental shelf in GSAs 1, 5, 6 and 7, the shelf break/upper slope in GSAs 1, 11, 17 and 18 and the lower slope in GSAs 5 and 11. It increased only on the continental shelf and the lower slope of GSA 18. Expert knowledge suggested increasing trends in fishing effort for the lower slope in GSAs 20, 22 and 23 and decreasing trends on the continental shelf of GSAs 8,9,16,20,22 and 23 and on the lower slope of GSAs 1, 7 and 8 ( Table 4, Fig. 3).
In 6 of the 7 cases in which an increment of N 90 was detected, it coincided with a decrease in fishing effort (Tables 3-4, Fig. 3). Of the 11 cases showing increases in species richness, S, only in 3 cases was the increase in S coupled with a decrease in fishing effort. In 5 cases there was no trend in fishing effort, while in only one case the increase in S was coupled with an increase in fishing effort. In 2 cases, no information on the temporal evolution of fishing effort was available (Tables 3-4, Fig. 3). 10 GSAs showed significant trends in Pielou eveness, J'; in 2 of them J' increased and fishing effort decreased, while 3 GSAs showed a decrease in both J' and fishing  Decreasing Increasing effort. Two GSAs showed a decrease in J' coupled with no trend in fishing effort, while no information on fishing effort trend was available in 3 cases (Tables 3-4, Fig. 3).

Spatial patterns
Mean values of N 90 , S and J' showed differences between GSAs and depth strata (Fig. 4). Regarding   . Note that in some cases the trend of time series does not match the arrow's direction (see Table 2).  (Table 5). The species with the highest percentage contribution to within-group similarity were T. minutus on the continental shelf and G. argenteus on the shelf break/upper slope in GSA 7; Lepidotrigla cavillone on the continental shelf and C. agassizi on the shelf break/upper slope in GSA 20; and L. cavillone on the continental shelf and A. sphyraena on the shelf break/upper slope in GSA 23. In GSA 17, M. merluccius was the species with the highest contribution on both the continental shelf and the shelf break/ upper slope (Table 5).
Some species showed a high percentage contribution to within-group similarity in most years of the time series and for most of the GSAs (Table 6). On the continental shelf, these species were S. hepatus, L. cavillone, M. barbatus and M. merluccius. On the shelf break/upper slope, only G. argenteus was present in all GSAs, while on the lower slope those species were P. blennoides, G. melastomus and Etmopterus spinax.

DISCUSSION
The results have confirmed that demersal fish assemblages are highly structured in the Mediterranean. In fact, we were able to identify three common assemblages in most GSAs corresponding to the continental shelf, shelf break/upper slope and lower slope strata of each area. There were only two GSAs, the northern Adriatic Sea and Crete, that did not present lower slope assemblages, due to the shallower depth surveyed in these areas compared with the rest of GSAs. Although the number of samples was not enough to follow their temporal series, the Gulf of Lions and the eastern Ionian Sea also followed this depth structure. The results confirm the findings of previous works on the structure of demersal assemblages in the Mediterranean, showing that for fishes (Ungaro et al. 1999, Labropoulou and Papaconstantinou 2004, García-Ruiz et al. 2015 and other taxonomic groups (Tserpes et al. 1999, Colloca et al. 2003, Massutí and Reñones 2005 they are strongly organized along a depth gradient.
Despite the similar bathymetric gradient along the Mediterranean, the results showed differences in the bathymetric limitations and composition of demersal fish assemblages between GSAs. This is not surprising considering that oceanographic conditions vary between GSAs, and bathymetric distributions of communities respond according to these variations. In fact, we incorporated cluster analysis to escape from the assumption that communities are structured according to MEDITS strata, and we made an analysis based on real assemblages for each GSA. Therefore, the analysis of demersal fish diversity based on cluster analysis for each particular area of the Mediterranean is more accurate than the assignation of a depth stratum to the samples analysed for the whole Mediterranean, the method that has been generally used up to now.
Our results show a stability and even recovery of demersal fish diversity in the Mediterranean. Of the 114 temporal series analysed, only 27% showed a significant trend, with an increasing pattern in 71% of the cases showing significant trends. N 90 and species richness (S) showed increasing trends in most cases (87.5% and 84.6%, respectively), while Pielou evenness (J') was the indicator that showed the highest proportion of decreasing trends (60%). This stability was also shown in the only study analysing long temporal series from bottom trawl survey data (1994)(1995)(1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012) for the whole Mediterranean (Granger et al. 2015). These authors took into account three scales of analysis corresponding to 18 GSAs, 7 biogeographical zones and 2 basins at a depth ranging from 10 to 800 m. The assemblages by depth were not considered, which could explain why they did not detect any recovery.
The continuous increase in fisheries in the last few decades has led to the overexploitation of the main commercial stocks in most Mediterranean areas (Colloca et al. 2013, Sartor et al. 2014). However, bottom trawl fisheries in the Mediterranean have decreased recently, due to the economic losses of this activity (Quetglas et al. 2017;Table 4) and the implementation of additional management measures, such as the prohibition of bot- Table 5. -Similarity Percentage analysis (SIMPER) summary table of species appearing in the 90% cut-off of within-group similarity. A is the mean abundance (individuals km -2 ) of each species, and %C is the mean value of the percentage contribution of each species to within-group similarity, taking into account each SIMPER made by group of GSA, depth strata and year. Depth strata are shelf, shelf break/upper slope (SB/US) and lower slope.
GSA  tom trawling within 1.5 nautical miles of the coast (EC Regulation 1967). However, this recent measure has possibly displaced part of the trawl fishing effort from the shelf to deeper bathymetric zones (Tserpes et al. 2011). Our results show that increasing trends in N 90 and S and decreasing trends in J' coincide in some cases with decreasing trends in bottom trawl fishing effort. There could therefore be a cause and effect relation, because it is in accordance with the expected effect of fishing on biodiversity. The increasing trend in N 90 with decreasing fishing effort reinforces the previous results that confirm the usefulness of this index for detecting the effects of fishing on demersal fish diversity. The increase of evenness with increasing fishing effort has been suggested by some authors (Murawski 2000, Zhou et al. 2010 due to the reduction of dominant species by fishing (Cury et al. 2000, Rice 2000 and has been confirmed by the study of the effects of fishing on evenness indices (D'Onghia et al. 2003, Farriols et al. 2017). However, the expected increase in S and decrease in J' with decreasing fishing effort are not always observed in our results. There are some differences in the aspects of diversity that each of these indices capture. Increasing N 90 values with decreasing fishing effort indicate an increase in the frequency of occurrence and the evenness of the distribution of species abundances due to expan-sion to areas with the most favourable environmental conditions. On the other hand, an increase in S and a decrease in J' with decreasing fishing effort implies an increase in the number of species and an increase in the dominance of some species, respectively. Although both number of species and evenness are also affecting N 90 , the calculation of each of these indices is extremely different. N 90 takes into account the homogeneity or heterogeneity of all the samples of a stratum and year for each GSA in its calculation and involves the most frequent and abundant species in the group without losing species identity through the comparison among all the samples in the group. In contrast, S and J' in a group are calculated from their mean values and consequently species identity is lost. This may explain why extreme values of fishing effort were needed to detect the effects of fishing in S and J' in previous works (Farriols et al. 2017). In some cases, N 90 showed no trend when there was a trend in fishing effort and viceversa. This finding could be due to several causes. It is either too early to detect the effects of decreasing fishing effort on demersal fish diversity or the decrease is not sufficiently important to change the diversity trend. Similarly, increasing trends in fishing effort could not result in a decrease in fish diversity due to the adaptation of demersal fish communities to fishing exploitation.  ). -Similarity Percentage analysis (SIMPER) summary table of species appearing in the 90% cut-off of within-group similarity. A is the mean abundance (individuals km -2 ) of each species, and %C is the mean value of the percentage contribution of each species to within-group similarity, taking into account each SIMPER made by group of GSA, depth strata and year. Depth strata are shelf, shelf break/ upper slope (SB/US) and lower slope. Table 6. -Number of years that each species contributed to the 90% cut-off of within-group similarity, taking into account each similarity percentage analysis (SIMPER) by GSA, depth stratum and year during the time series. Depth strata are shelf, shelf break/upper slope (SB/ US) and lower slope. Species  GSA  1  5  6  7  8  9  10 11 16 17 18 20 22 23   Shelf  Argentina sphyraena  ----3  4  1  9  13  -2  9  6  -Arnoglossus imperialis  --------4  -----Arnoglossus laterna  8  -3  6  -12 10  -10  4 12 11

Stratum
It must also be considered that there is a high complexity in the evaluation of fishing effort in the whole area. Available temporal series used to analyse fishing effort do not cover the whole time series of demersal fish diversity in all cases (see Table S1), and the inclusion of more years of the temporal series to the analysis could lead to different trends of fishing effort. In addition, as the nominal spatio-temporal pattern of fishing effort on a Mediterranean-wide level is not available, the use of different effort estimates in the areas may increase the uncertainty of the model. Moreover, number of fishing vessels is a poor proxy for effort, because does not account for other capacity changes (e.g. length overall, or kilowatts), because it does not account for technological creep and or temporal and spatial changes of fishing operations (Anticamara et al. 2011). For instance, though the regulation decreasing fishing capacity has been in place since 1991, the gross tonnage of the fleets may be increasing because boats are decreasing in number (decommissioning) but increasing in size over time (e.g. Fortibuoni et al. 2017) or because vessels have increased trawling time. This issue is enormously relevant for the Mediterranean Sea, where fisheries are managed by effort control and technical measures in contrast to quotas (northern EU seas; see Cardinale and Scarcella 2017) and should be considered when the results are interpreted. However, in the Spanish and French Mediterranean, restrictions on hours of trawling would not permit an unlimited increase in fishing effort with a decreasing number of vessels (REAL DECRETO 1440/1999, de 10 de septiembre; Arrêté n° 99-162 du 10 juin 1999).      In any case, a more appropriate indicator than number of vessels should be used for fishing effort whenever possible.
Regarding spatial patterns, we did not find the expected longitudinal decreasing west-east pattern in species richness observed in previous works on fish communities (Quignard andTomasini 2000, Coll et al. 2010). Nor was this trend observed for N 90 and J' in any depth stratum. The absence of a western/eastern decreasing trend further suggests that primary production or temperature regime are possibly not the major factor explaining large-scale patterns of diversity in demersal fish assemblages, as suggested by Gaertner et al. (2007). However, it is difficult to compare our results with diversity values obtained with non-standardized data mainly collected from fish inventories from other works. Moreover, due to the limited sampling approach (i.e. data concerning only one guild of fishes or limited to specific depths, gear or habitat), some of regional inventories are useless for comparative studies (Psomadakis et al. 2012). Recent studies based mainly on standardized time series data also question the previously considered west-east decreasing diversity trend in the Mediterranean (Gaertner et al. 2013, Granger et al. 2015, Peristeraki et al. 2017. The highest diversity values were found on the continental shelf of insular areas, such as the Balearic Islands, Sardinia, Sicily and the Aegean Sea. This higher diversity can be explained by taking into account the peculiarities of the distinct biogeographic sectors within the Mediterranean (Lejeusne et al. 2010), which can be characterized by the shallow water biota (Bianchi et al. 2012). In the Strait of Sicily, for example, the meeting of western and eastern Mediterranean species produces a peak in fish species richness in the central Mediterranean (Ben Rais Lasram et al. 2009, Garofalo et al. 2007). The greater sampling effort of the present work compared with previous ones (Morri et al. 1999, Koukouras et al. 2001) could affect the unexpectedly high diversity values found in the Aegean Sea. The presence of algae facies deeper than 50 m around the Balearic Islands is likely to enhance demersal fish diversity in this area. Coralligenous and maerl communities are very characteristic of the Mallorca-Menorca continental shelf up to 85-90 m depth (Canals andBallesteros 1997, Ordines and, and this has been pointed out as a plausible reason for the differences observed between the coastal demersal resources of the Balearic Islands and the adjacent Iberian Peninsula (Massutí and Reñones 2005). In fact, habitat type has been shown to affect the distribution of demersal species, most of them being more abundant and showing a better condition in maerl and Peyssonnelia beds , which have also shown high diversity of fish.
The results of SIMPER analyses reinforce the idea of maerl and Peyssonnelia beds causing high diversity values also on the continental shelves of Sicily, Sardinia and the Aegean Sea. The species Serranus cabrilla, Scyliorhinus canicula and Mullus surmuletus, which in the Balearic Islands have shown to be more abundant in these habitats  contribute to N 90 mainly in this archipelago, Sardinia, Sicily and the Aegean Sea. Similar habitats to those found on the Balearic shelf have been reported in some of these areas, like the Aegean Sea (Georgiadis et al. 2009). The presence of a higher number of vulnerable species like demersal chondrichthyans in the Balearic Islands, Sardinia, Sicily and the Aegean Sea (Bertrand et al. 2000, Damalas and Vassilopoulou 2011, Ramírez-Amaro et al. 2015 compared with adjacent areas could also contribute to the higher fish diversity values found there.
The spatial distribution of the bottom trawl fishing effort by GSA shows that the number of vessels per km 2 is low on the continental shelf of the Balearic Islands, Sardinia and the Aegean Sea (Colloca et al. 2017). The coincidence of areas with a low fishing effort with areas with a high diversity is in accordance with previous works, in which higher values of N 90 and S and lower values of J' were associated with areas with a low fishing effort (Farriols et al. 2017). The lower fishing effort exerted by the relatively smaller bottom trawl fleets in these areas could have preserved, at least to some extent, their fish diversity along with a better conservation of their sensitive and essential habitats, such as maerl and Peyssonnelia beds. These habitats are precisely those most affected by the low selectivity and damaging collateral effects of bottom trawling on seabed communities, which decrease the presence of biogenic habitats, leading to a reduction in the biodiversity on exploited bottoms (e.g. Norse and Watling 1999, Smith et al. 2000, Hiddink et al. 2006. Spatial patterns of demersal fish diversity on the shelf-break/upper slope and lower slope of the Mediterranean are different to those detected on the continental shelf. Areas with the highest diversity values on the continental shelf do not coincide with areas with the highest diversity values on the shelf break/ upper slope and lower slope. Although the assignment of depth strata was different in previous works and the comparison is not straightforward, a different pattern on shelf and slope areas was also observed for species richness (Gaertner et al. 2007(Gaertner et al. , 2013. This is likely due to differences in the distribution of cumulative threats to marine biodiversity, which are mainly concentrated in coastal areas and on the continental shelf of the Mediterranean (Coll et al. 2012), and to the presence of particular habitats on the shelf break and slope bottoms, which may represent potential hot spots of biodiversity (Danovaro et al. 2010). Although the distribution of deep-sea diversity is different to that on the continental shelf, it is affected by similar factors: changes in spatial distribution of fishing effort together with habitat type. For example, higher N 90 values on the slope of northern Spain could be related to the presence of submarine canyons in the area where high values of biodiversity have been reported (see Fernández-Arcaya et al. 2017 for a review). However, the description of deep-sea habitats has just been implemented for some particular areas of the Mediterranean, and this information is not exhaustive at all (Danovaro et al. 2010). Moreover, an intensive habitat mapping based on MEDITS samples would be useful to relate demersal fish assemblages to their corresponding habitats, as has been done in some continental shelf areas (e.g. . The outcomes of the present study show that at large temporal and spatial scales bottom trawl fisheries have reduced the diversity of demersal assemblages in the Mediterranean. However, in recent decades a generally stable scenario or even a slight recovery trend have been highlighted. This result would have not been expected if the alarming overexploitation status of Mediterranean stocks were taken into account, which underlines the importance of using diversity indices to study the effects of fishing on demersal assemblages. Therefore, a change from the assessment of demersal resources based on exploited monospecfic stocks to one based on the study of whole demersal fish assemblages is needed due to the high multispecificity of the bottom trawl fishery in the Mediterranean (Caddy 1993, Lleonart andMaynou 2003). The inclusion of species other than target ones made in this work through diversity indices is therefore important for the implementation of an ecosystem-based fisheries management (Browman and Stergiou 2004).

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/sm04977esm.pdf