Elasmobranch spatial segregation in the western Mediterranean

1 Instituto Español de Oceanografía, Centre Oceanogràfic de les Balears, Moll de Ponent s/n, 07015 Palma, Spain. E-mail: joan.moranta@ba.ieo.es 2 University of Glasgow, Division of Ecology and Evolutionary Biology, Graham Kerr Building, Glasgow G12 8QQ, United Kingdom. 3 Centre for Ecological and Evolutionary Synthesis (CEES), Department of Biology, University of Oslo, PO Box 1066 Blindern, 0316 Oslo, Norway.


INTRODUCTION
The ecological importance of elasmobranchs is becoming increasingly clear (Stevens et al., 2000;Bascompte et al., 2005).Most elasmobranchs are predators at or near the top of marine food chains; thus, their removal affects the structure and function of marine ecosystems.These effects result in changes at species level in both abundance and size structure, which can lead to changes in the trophic interactions at the community level due to selective removal of predator or prey species, removal of competitors and species replacement (Stevens et al., 2000).As a result, understanding the ecological role, biology, distribution, and life histories of elasmobranchs is an important first step towards the preservation of marine biodiversity and ecosystem functioning.At a global scale, there is current concern over the effects of increased fishing intensity on the status of the elasmobranch populations (Stevens et al., 2000).Their low recovery rates are due to sensitive life-history stages (Cailliet et al., 2005), late maturity and low fecundity, which make elasmobranchs highly vulnerable to overexploitation (Stevens et al., 2000).This is of particular concern in the highly exploited Mediterranean, where sharks and rays make up a significant proportion of the by-catch (Carbonell et al., 2003) and are targeted at various bathymetric ranges, which has led to their drastic decline in this area over the last two centuries (Ferretti et al., 2008).Protecting these fish requires detailed information about their distribution and biology; however, this information is not available for many Mediterranean species and regions.
Knowledge on the spatiotemporal variation in species distribution and habitat requirements is essential in order to implement the ecosystem based management approach to fisheries, as well as to carry out long-term monitoring programs and ad hoc studies (e.g.Gaertner et al., 2005;Ciannelli et al., 2008;Moranta et al., 2008b;Ordines et al., 2009).However, further indepth studies are required in order to gain more insight into several harvested species in a community context.Previous studies suggested that elasmobranch species in the area of the Balearic Islands are segregated in relation to trophic resources and space (Moranta et al., 2008a), although this has never been thoroughly investigated.The present study focused on understanding the spatial distribution in relation to biological, environmental and anthropogenic factors in the areas along the continental shelf and slope off the coast of the Balearic Islands.Our study aimed to describe the depth related trends and the geographic patterns shaping the elasmobranch community and the distributions of the most abundant species.To do so, we analysed elasmobranch community descriptors and species-specific information.The specific objectives of this study were: (i) to establish whether there is a spatial pattern (referring to bathymetric and geographic differences in community descriptors) in the elasmobranch community, and (ii) to investigate the effect of fishing exploitation on the spatial variability of species-specific indicators (i.e.abundance and mean length).

Study area
There are large geographical barriers separating the Balearic Islands from the mainland (Fig. 1).The depth between the islands and the mainland ranges between 800 m and 2000 m and the minimum distance to the mainland is 90 km.The continental shelf around Mallorca and Menorca is narrow and is characterised by rocky coasts and the absence of submarine canyons (Massutí and Reñones, 2005).Its morphology is mainly conditioned by emergent geological structures rather than by the input of sediment from the continental shelf.The bottom of the continental shelf is covered with gravel and coarse sand, with the exceptions of north and south Mallorca, where muddy sediment is present in certain bays and on the slope.In these areas, the channels between Mallorca and Menorca, which have a maximum depth of 100 m, broaden the continental shelf.Sediments are mainly of biogenic origin as there is no river runoff, which reduces the presence of mobile sediments (Canals and Ballesteros, 1997).
The area of Menorca and the northeast coast of Mallorca are less exploited by fisheries than the south and northwest coast of Mallorca, which is characterised by a larger trawl fishery, both in terms of the number of boats and their size (Moranta et al., 2008b;Ordines and Massutí, 2009).
trawl fishing surveys at depths of 36 to 755 m around the continental shelf and the upper slope of Mallorca and Menorca (Balearic Islands) in the western Mediterranean (Fig. 1).The surveys were conducted on board the R/V Francisco de Paula Navarro and the R/V Cornide de Saavedra.The surveys were based on a stratified random sampling design and the sampling gear used was the GOC73, which has been commonly used in bottom trawl surveys throughout the Mediterranean since 1994 (Bertrand et al., 2002, Abelló et al., 2002).
The effective trawling time ranged from 20 to 60 min depending on the depth (20 min at <100 m, 30 min between 100 and 200m and 60 min at >200 m) with an average towing speed of 3 knots.The arrival and departure of the net to the bottom, as well as the horizontal and vertical openings of the net (16.4 and 2.8 m respectively), were measured using a SCANMAR system.This swept area information was used to standardise the catch (abundance and biomass).For each trawl, the individual fish species were classified by total number, size, sex (in the case of elasmobranchs) and weight.Further information on the sampling design and on the characteristics of the gear is available in Bertrand et al. (2002).

Spatial analysis
Analyses of abundance and biomass matrices of species composition were performed by means of the PRIMER 6.0 analytical software (Clarke and Gorley, 2001).Hauls containing a single species as well as those species recorded in less than 5% of the samples were not included.The Bray-Curtis index was selected as the similarity coefficient (Clifford and Stephenson, 1975).Different assemblages were identified using non-metric Multi-Dimensional Scaling (MDS).
Generalized additive modelling (GAM) (Hastie and Tibshirani, 1990) was used to fit spatial estimates of community descriptors (B, biomass; A, abundance; NS, number of species; H, diversitycalculated using the Shannon-Wiener Index; MW, mean weight) and abundance and mean length of the main elasmobranch species.GAMs are non-parametric regressions with the main advantage of not requiring a priori specification of underlying functional forms between dependent and independent variables.Our objective using this technique was to detect the spatial pattern once the effect of depth had been taken out.To model the community descriptors our additive formulation was: (l,j) ] + e t, (l,j) where GD is the general descriptor (B, A, NS, H and MW; log-scale if the response variable was not normally distributed) at longitude l and latitude j at time t (year), a t is a year-specific constant, s is a two-dimensional non-parametric smoothing function describing the geographic effect (thin plate regression spline, Wood, 2000), with a maximum of 27 degrees of freedom [28 knots]), g is a one-dimensional nonparametric smoothing function (cubic spline, Wood, 2000, with a maximum of 3 degrees of freedom [4 knots]), D is depth (m) and e is a normally distributed error term.
In order to model the abundance and the mean length of the main species (S. canicula, G. melastomus, R. clavata and R. miraletus) we used a similar model.Since we aimed to assess the main bathymetric and geographic pattern, we did not account for the interannual variation which could depend on internal demographic processes and the population dynamics.The species-specific formulation was: D t,(l,j) ] + e t, (l,j) where SP in this case is the species-specific descriptor (A, abundance and ML, mean length; both log-scale if the response variable was not normally distributed) and β the intercept.Finally, to model species-specific abundance we also included zero data; thus, the data showed a zero-inflated distribution.To tackle this data distribution, we used a model approach to analyse zero-inflated spatiotemporal data (Liu and Chan, 2010).This approach assumes that the response follows a distribution from the zero-inflated 1-parameter exponential family, with the further assumption that the probability of zero inflation is a monotone function of the mean response function (Liu and Chan, 2010).This technique bypasses the problems for analysing zero-inflated data that either focus only on the non-zero data or model the presence-absence data and the non-zero data separately (Liu and Chan, 2010;Liu et al., 2010).More details of this methodology are available in Liu and Chan (2010), which describes the R library used (COZIGAM, COnstrained Zero-Inflated Generalized Additive Model).The model selection criterion (MSC) was based on the minimisa-tion of logarithmic marginal likelihoods and generalized cross validation for COZIGAM (Liu and Chan, 2010) and regular GAM (Wood, 2000).MSC values indicate the best model compromise between model complexity and fit to the observed data.Regular GAMs for modelling the GD of the community and the ML of the four species analysed were performed using the R library mgcv (Wood, 2006).All the models, calculations and data pre-processing were coded in R (version 2.7.1;R development core team, 2007).

Species assemblages
During the surveys, a total of 51753 specimens (ca 7.8 tons of biomass) belonging to 27 chondrichthyan species and 8 families was collected (Table 1).
The MDS analysis for the samples taken below 50 m depth identified 4 main groups with 60% similarity.These included sample assemblages grouped at the following depth intervals: 55-185 m, 190-290 m, 330-460 m and >500 m (Fig. 2).
The results of the SIMPER analysis identified the relative importance of several elasmobranch species corresponding to different bathymetric assemblages (Table 2).
In shallow waters, at depths between 55 and 185 m, Scyliorhinus canicula and Raja miraletus were the most important species in terms of frequency of occur- rence.In the depth range between 190 and 290 m, S. canicula remained important, but was closely followed by Raja clavata in terms of abundance and frequency of occurrence.The depth range between 330 and 460 m was dominated by Dipturus oxyrinchus and Galeus melastomus, which made up 31% and 22% respectively of the total catch abundance.Finally, the samples taken deeper than 500 m were dominated by the sharks G. melastomus and Etmopterus spinax, which made up 46% and 38% respectively of the total catch abundance.

Spatiotemporal variability of the community descriptors
Data collected during 8 years, coming from 10 bottom trawl fishing surveys totalling 605 hauls, were included in the GAM models.Table 3 summarises the best model for each community descriptor.
In all cases, inter-annual variability (temporal factor) and bathymetric and geographic effect (smoothing terms) had a significant effect on the spatial distribution of the community descriptors.The percentages of deviance explained by the models ranged from 19.3% for diversity to 48.3% for mean weight.Figure 3 shows the inter-annual variability and the functional form of the temporal factor and smoothing terms.
The abundance and biomass (Fig. 3a, b) showed similar temporal patterns, and the average values of 2002 and 2007 were higher in comparison to other years.The bathymetric distribution showed different maxima at 200 m and 300 m for biomass and abundance respectively.The remaining variance explained   by the geographic effect showed both higher abundance and biomass around Menorca and decreased abundance around Mallorca with the lowest values on the continental shelf southwest of Mallorca.Mean fish weight (Fig. 3c) showed a general increase during the first five years, then a slight decrease in the following years.The effect of depth on mean fish weight decreased down to 400 m, and then increased.The variance of mean fish weight explained by the geographic effect was generally larger in the area around Mallorca.
Both the interannual variability and bathymetric effect showed consistent results for the number of species and diversity (Fig. 3d, e), with an increase in recent years and greater values in shallow waters with peaks around 200 m.However, the variance left by the other covariates and explained by the spatial effect showed more species in the area around Menorca, while the diversity was higher on the shelf southwest of Mallorca and lower around Menorca.

Geographical and bathymetric effect on individual species
GAMs were also used to analyse data on the four main species in terms of abundance and frequency of occurrence: S. canicula, R. clavata, R. miraletus and G. melastomus (Table 1).Table 4 shows a summary of the results of the best models for abundance (Fig. 4) and mean length (Fig. 5) for the selected species.
Our results show the effect of depth on the spatial pattern of the two species-specific indicators (Table 4 and Figs. 4 and 5).Once the bathymetric effect had been accounted for, a geographic pattern for abundance was only observed in two species: S. canicula and R. clavata.Only depth significantly affected the abundance of G. melastomus (Fig. 4a) and R. miraletus (Fig. 4b), and consequently their spatial distribution was exclusively depth-dependent.The highest abundance found for G. melastomus was between 300 and 400 m.In contrast, R. miraletus occupied shallow waters with a maximum around 75 m.Both the geographic variability and the depth significantly affected the abundance of S. canicula (Fig. 4c) and R. clavata (Fig. 4d).The maximum abundance of S. canicula was found between 200 and 300 m, with a lower abundance around Mallorca compared to Menorca.R. clavata showed the highest abundance around a depth of 200 m in two areas: the shelf west of Mallorca and the area northeast of Menorca.The lowest abundance was in the Menorca channel.
The mean length analysis showed consistent and complementary results to those obtained by the species abundance analysis.In contrast to the abundance analysis, the best models always contained both bathymetric and geographic effects (Table 4).Shelf distrib-Table 4. -Model selection for the generalized additive model of abundance and mean length of the four species studied.Estimated intercepts (b) and estimated degrees of freedom for the two nonparametric terms are shown.The base model denotes the complete model including both bathymetric, s(depth), and geographic, s(l,j), terms.Asterisks denote significance levels: ***: p<0.001;**: p<0.01; *: p<0.05 and ns: not significant or p≥0.1.Note that the model selection criterion (M.S.C.) applied for the abundance models (COZIGAM) was the logarithmic marginal likelihood, while GCV (Generalized Cross Validation) was applied for regular GAMs.Figures showing the most parsimonious models are indicated.DE refers to the percentage of variance explained by regular GAMs (i.e.mean length); note that DE is not shown for abundance analyses due to the methodological constraints of COZIGAM (see Liu and Chan, 2010 for details).Note that the G. melastomus mean length intercept was considerably lower compared to the other mean length analyses due to the log-transformation prior to the analyses.uted species, S. canicula, R. miraletus and R. clavata, showed a high abundance of large individuals in shallow waters (Fig. 5b, c, d; although for S. canicula and R. clavata large individuals were also found deeper within their distribution range), while the slope species G. melastomus showed an opposing pattern, with large individuals in deeper waters (Fig. 5a).Some interesting features can also be identified within the geographic effect.Of special interest were the low mean length values of G. melastomus, which were only observed on the north slope of Mallorca (Fig. 5a).It is also important to note that, although S. canicula had a lower abundance around Mallorca compared to Menorca (Fig. 5c), the individuals off Mallorca were larger  compared to those off Menorca (Fig. 5b).Only one (shelf off Menorca) of the two areas of high abundance of R. clavata showed large individuals (Fig. 5c).The highest numbers of large individuals of R. miraletus were found on the shelf west of Mallorca and the area northeast of Menorca (Fig. 5d).

DISCUSSION
The four elasmobranch assemblages detected in the present study of the Balearic Islands were found to correspond to the bathymetric regions of the surveyed area: the first assemblage had an approximate depth range of 55 to 185 m, which corresponds to the continental shelf; the second had a depth range of 190 to 290 m, which corresponds to the shelf break; the third had a range of 330 to 460 m, which corresponds to the upper slope; and finally the fourth assemblage started at 500 m and went down to the final sampling depth (755 m), which corresponds to the middle slope.The depth is the main factor for establishing these elasmobranch communities, which has also been reported in earlier studies in this area (Massutí and Moranta, 2003), as well as studies from other parts of the Mediterranean (Demestre et al., 2000a).Previous studies conducted off the Balearic Islands on the whole fish community revealed similar bathymetric assemblages (Moranta et al., 1998;Massutí and Reñones, 2005).
The results obtained for abundance, biomass and number of species of the elasmobranch community in relation to geographic position (i.e.bathymetry independent) clearly evidence the direct effect fishing has on these community descriptors in the studied area.For all three, the areas with low values coincided with the areas that have high fishing intensity; likewise, the areas with high values of these community descriptors in our study are closely associated with regions where fishing intensity is lower (Fig. 4a, b, d).
At species level, the bathymetric trend was evident in the distribution of the four elasmobranch species studied.In S. canicula the abundance peak at 250 m was most likely associated with juveniles, usually found at depths around 200 m, while adults had a wider bathymetric distribution, and commonly occurred in the shallows (Sánchez, 1993;Sánchez et al., 2002;Rodríguez-Cabello et al., 2004).This pattern of bathymetric distribution in S. canicula was further confirmed by the results of the mean length of this species with adults (high values of mean length) mainly distributed in the shallows, and juveniles (low values of mean length) distributed around 250 m.R. clavata showed a similar bathymetric distribution pattern to S. canicula, with adults below a depth of 100 m, and a peak in abundance of juveniles in the depth range between 200 and 300 m.According to this, juveniles of R. clavata showed a bathymetric overlap with S. canicula juveniles.However, once the effect of depth had been eliminated, R. clavata was found to be the most abundant on the shelf west off Mallorca.Thus, even though they overlap bathymetrically, the juveniles of R. clavata and S. canicula do not overlap in space.This indicates that there could be intraspecific competition for space and trophic resources, which is known to be among the most important factors influencing segregation of marine species (Haedrich, 1997;Moranta et al., 2008b).Nevertheless, other ecologically similar species in the area (e.g.other non-elasmobranch fish) could also be competing for the same resources.
The ray R. miraletus was found in shallower waters in comparison with the other species studied, with a bathymetric range down to 150 m and a peak in abundance between 50 and 100 m.These findings are similar to those of the previous study conducted in the central Mediterranean, where R. miraletus was found between 50 and 150 m (Relini et al., 1999), and off Tunisia, where the distribution of this species was reported down to a depth of 200 m (Capapé and Quignard, 1974).The abundance peak of G. melastomus between 300 and 400 m was also associated with a larger number of juveniles, in this case found in the shallower part of its bathymetric range; thus, showing bathymetric overlap with juveniles of S. canicula.According to our results, the adults of G. melastomus were generally distributed deeper within their bathymetric range, which was also confirmed by the mean length results for this species, which clearly showed increasing values with increasing depth.The study by Tursi et al. (1993) in the Ionian Sea reports G. melastomus juveniles almost exclusively between 200 and 300 m, whereas in the depth range between 400 and 650 m the population was mixed, and included juveniles, maturing individuals and adults.G. melastomus expressed a "bigger-deeper" distribution trend, with the relative size increasing with depth (Fig. 5a), unlike the other three species whose adults were mainly found in shallow waters (Fig. 5b, c, d).The possible ecological reason for the "bigger-deeper" phenomenon in G. melastomus can be assigned to the biomass peak of decapod crustaceans between 400 m and 700 m (Maynou and Cartes, 2000), which are an important part of the diet of this species (Valls et al., 2011).Crustaceans make up the majority of the diet of R. miraletus, S. canicula and R. clavata (Valls et al., 2011).This may originate a phenomenon of competitive exclusion due to the exploitation of similar resources, in this case in favour of G. melastomus at a depth interval between 400 m and 700 m.Based on the competitive exclusion principle, other species are forced to retreat to the bathymetric and/or geographic range to which they are most highly adapted in relation to the other potential inhabitants.According to this phenomenon it is possible that a strong selection pressure for diversification in the trophic strategies of these four species is the main explanation for their geographic and bathymetric segregation.Although the reproductive cycle of the species can have different seasonal patterns (i.e.recruitment events or adult spawning aggregations) which can affect the bathymetric and geographic distribution of the species, it is difficult to address this issue due to the lack of seasonal surveys.
The four elasmobranch species studied also showed evidence of the direct effect of fishing on their spatial distribution and demography (i.e.mean length), which is consistent with the elasmobranch community descriptors.It is well known that fishing effort in the Balearic Islands is displaced towards greater depth (Moranta et al., 2008b;Guijarro et al., in press), and thus has different effects on the populations of the two most common sharks in the area, S. canicula and G. melastomus.In our study, S. canicula expressed the opposite spatial pattern of abundance and mean lengthlow abundance around Mallorca but with higher mean length in the same area.This suggests that there are different mortality rates for juveniles (found deeper) and adults (shallow distribution) (see above).The increase in adult populations around Mallorca could be related to the scavenging and opportunistic behaviour of this species, which has been observed to be enhanced due to the discards of fishing activity both in the Atlantic (Olaso et al., 1998) and the Mediterranean (Demestre et al., 2000b).In contrast, G. melastomus was expected to display a homogeneous depth-dependent spatial distribution.However, the smallest mean length was observed in the area northwest Mallorca, which suggests there was a reduction in both the spatial distribution and their demography caused by the higher fishing effort in this region.Combined with the "bigger-deeper" distribution trend this species expresses, fishing in the area results in high selectivity for larger individuals.The outcome is an age/size truncation often referred to as an 'age truncation effect' (see review in Hsieh et al, 2010), which may have consequences for the population and community dynamics.
The large scale removals of certain species of predators reduces the environmental heterogeneity created by these organisms and can even impact organisms in the water column (Thrush and Dayton, 2010).For instance, the impact of declines in cod stocks and other predatory fish on the Scotian Shelf had a detrimental effect on the populations of shrimp, crabs and other benthic organisms as well as the nutrient dynamics (Frank et al., 2005).Therefore, overexploiting elasmobranchs in the area of the Balearic Islands can have a detrimental effect on the whole ecosystem by disturbing the complex equilibrium these predators provide.Our study suggests that these ecosystem alterations could be spatially heterogeneous due to clear differences in the species-specific consequences of the fishing exploitation.
The results from this study increase our knowledge of elasmobranch biology and distribution within the investigated areas.Combining the bathymetric and geographic distribution, our study shows a complete spatial segregation of the main elasmobranch species despite the apparent similarity of their diets and the geographical availability of resources.Our results also highlight that fishing is most likely heterogeneously affecting the spatial structure and dynamics of the elasmobranch community, and at the same time, the entire ecosystem due to both the depletion of the relative abundances (around Mallorca) and the demographic structure of certain species (e.g.juveniles of S. canicula around Mallorca, and adults of G. melastomus northwest of Mallorca).Thus, the elasmobranchs in the area of the Balearic Islands could potentially be a very valuable indicator of the fishing impact on the ecosystem, and thus monitoring future changes in the elasmobranch community should be given greater importance.With this in mind, it is important that the fishery around the Balearic Islands develops technical measures in order to protect the predatory fish populations in this area, including the elasmobranchs, which would subsequently ensure a healthy and sustainable ecosystem.

Fig. 1 .
Fig. 1. -Maps of the studied area off Mallorca and Menorca Islands (Balearic Islands, western Mediterranean).Black dots represent the locations of the samples collected during the MEDITS and MIGJORN cruises, the contours around the islands represent depths of 200 m and 800 m respectively.

Fig. 2
Fig. 2. -MDS ordination indicating the groupings with 60% similarity in different depth ranges, acquired through cluster analysis of elasmobranch samples obtained during MEDITS and MIGJORN trawl surveys.

Fig. 3 .
Fig. 3. -Partial effects of year (left column), bottom depth (middle column), and position (right column) on elasmobranch community abundance (a), biomass (b), mean fish weight (c), number of species (d) and diversity (e) estimated from the spatially explicit generalized additive model.For the position plots, orange indicates high predicted values, green medium values and blue indicates low values.Shaded areas on the bottom depth plots and lines going through the dots on the year plots are 95% confidence intervals, and tick marks on the x-axes of the depth plots indicate sampling intensity.

Fig. 4 .
Fig. 4. -The partial effects of bottom depth and geographic position on G. melastomus (a), R. miraletus (b), S. canicula (c) and R. clavata (d) abundance, estimated from the generalized additive model.For the position plot, orange indicates high predicted values, green medium values and blue indicates low values.Shaded areas on the depth plot are 95% confidence intervals, and tick marks on the x-axes of the depth plots indicate sampling intensity.Only those maps resulting from significant models were included.Note that the x-axes for the plots representing the bathymetric distribution are different in order to take into account the distribution ranges of the species.

Fig. 5 .
Fig. 5. -Partial effects of bottom depth and geographic position on G. melastomus (a), S. canicula (b), R. clavata (c) and R. miraletus (d) mean length, estimated from the spatially explicit variable coefficient generalized additive model.For the position plot, orange indicates high predicted values, green medium values and blue indicates low values.Shaded areas on the depth plot are 95% confidence intervals, and tick marks on the x-axes of the depth plots indicate sampling intensity.Note that the x-axes for the plots representing the bathymetric distribution are different in order to take into account the distribution ranges of the species.

Table 1 .
-Species caught off the Balearic Islands (western Mediterranean) between depths of 36 m and 755 m.Total abundance (A: number of specimens) and biomass (B: kg), frequency of occurrence (F), depth range (D: m) and size range (S: mm) are shown by species for each surveyed area.

Table 2
. -The top 5 ranking species for each depth interval by abundance identified by SIMPER analysis (the abundance is expressed as a percentage of the total catch (%) for each bathymetric range and frequency of occurrence (f) as the number of samples in which the species was caught in relation to the number of samples taken at each depth stratum).

Table 3
. -Results of the generalized additive model for the community descriptors analysed (A, abundance; B, biomass; NS, number of species; H, diversity; MW, mean weight).Note that the selection criteria are not shown because the most complete model including both bathymetric, s(depth), and geographic, s(l, j), terms was always the most parsimonious.For each model the intercept is year-specific and is shown in Figure 3. Asterisks denote significance levels: ***, p<0.001; **, p<0.01; *, p<0.05.DE refers to the percentage of variance explained.