Modelling of European hake nurseries in the Mediterranean Sea: An ecological niche approach

An ecological niche modelling (ENM) approach was developed to model the suitable habitat for the 0-group European hake, Merluccius merluccius L., 1758, in the Mediterranean Sea. The ENM was built combining knowledge on biological traits of hake recruits (e.g. growth, settlement, mobility and feeding strategy) with patterns of selected ecological variables (chlorophyll-a fronts and concentration, bottom depth, sea bottom current and temperature) to highlight favourable nursery habitats. The results show that hake nurseries require stable bottom temperature (11.8–15.0 (cid:2) C), low bottom currents (<0.034 m s (cid:2) 1 ) and a frequent occurrence of productive fronts in low chlorophyll-a areas (0.1– 0.9 mg m (cid:2) 3 ) to support a successful recruitment. These conditions mostly occur recurrently in outer shelf and shelf break areas. The prediction explains the relative balance between biotic and abiotic drivers of hake recruitment in the Mediterranean Sea and the primary role of unfavourable environmental conditions on low recruitment in speciﬁc years (i.e. 2011). The ENM outputs particularly agree spatially with biomass data of recruits, although processes such as ﬁshing and natural mortality are not accounted for. The seasonal mapping of suitable habitats provides information on potential nurseries and recruitment carrying capacity which are relevant for spatial ﬁsheries management of hake in the Mediterranean Sea. an open access the CC BY license (http://creativecommons.org/licenses/by/3.0/).


Introduction
Understanding spatial patterns in population dynamics is a necessary prerequisite to protecting critical habitats, and thus ultimately in ensuring sustainable management of fishery resources (Berkeley et al., 2004;Caddy, 2000). Reducing fishing effort on juveniles in particular is vital if populations are harvested at maximum sustainable yield, especially in areas where juveniles are vulnerable to unselective fishing gears (Caddy, 2009), as is often the case in the Mediterranean Sea (Colloca et al., 2013). Furthermore, information on critical reproductive habitats provides an insight into the likely spatial structure of population units or stocks, whilst an insight into environmental conditions required for successful recruitment allows scientists and managers to better depict the interaction between environmental parameters and stock recruitment relationships. Information on the spatial aspects of population ecology and interactions with relevant ecosystem components are needed to implement an Ecosystem Approach to http://dx.doi.org/10.1016/j.pocean.2014.11.005 0079-6611/Ó 2014 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/3.0/). Fisheries Management (Link, 2013;Pauly et al., 2011) that is required as part of the implementation of the Marine Strategy Framework Directive (European Commission, 2008) and is recognized as a fundamental principle underpinning the revised Common Fisheries Policy. As a matter of fact the Council Regulation (EC) 1967/2006, concerning management measures for the sustainable exploitation of fishery resources in the Mediterranean Sea, specifically requires the inclusion of spatial aspects such as the establishment of fishing protected areas in order to protect nurseries and/or spawning areas. Among the commercial species, the European hake (Merluccius merluccius, L. 1758) is one of the most important in the Mediterranean Sea with total landings of 22547 tons in 2011 (GFCM-FAO 1 ). All the available assessments in the Mediterranean Sea have underlined that the status of hake stocks is characterised by high fishing mortalities on juveniles (Colloca et al., 2013). Considering the large size that hake can reach (more than 100 cm Total Length), coupled with a low size at first capture of the Mediterranean fine-meshed trawling (Bethke, 2004), the protection of hake nurseries has been proposed as an effective measure to improve size composition of catches (Caddy, 1999).
In order to identify appropriate areas to be closed to fishing, many authors have regionally studied the spatial distribution of the European hake juveniles and identified the main nurseries as areas where the highest concentrations of young-of-the-year remain remarkably stable over the years (Carlucci et al., 2009;Colloca et al., 2009;Lembo et al., 2000;Fiorentino et al., 2003;Lleonart, 2001;Murenu et al., 2010;Tserpes et al., 2008). Although the stability of the nurseries over time implies the existence of favourable habitats, only a few regional research initiatives have focused on the identification of the ecological factors which make some areas more suitable compared to others for hosting high concentrations of 0-group hakes in the Mediterranean Sea. These factors such as wind mixing, temperature, currents, fronts and primary production were identified independently in different Mediterranean regions (Abella et al., 2008;Bartolino et al., 2008a;Hidalgo et al., 2008;Lleonart, 2001) and no clear explanation of their role in the nurseries' functioning could be provided. In this paper, the underlying assumption for the feeding habitats of hake recruits in the Mediterranean Sea relies on the importance of productive fronts (chlorophyll-a fronts), by means of their long lifetime, to efficiently transfer the flow of energy along the food chain up to top predators. It is well known that productive oceanic features (chlorophyll fronts) are key vectors of the oceans' productivity along the food chain (Belkin et al., 2009;Druon et al., 2012Druon et al., , 2011Kirby et al., 2000;Le Fèvre, 1986;Olson et al., 1994;Polovina et al., 2004Polovina et al., , 2001. Bakun (2006) highlights the importance of frontal systems as sub-seasonal meso-scale environmental processes that may often be critical to regulating population-scale reproductive success, as in the Strait of Sicily with the semi-permanent eddies and fronts produced by the Atlantic Ionian Stream (see e.g. Fortibuoni et al., 2010;Garofalo et al., 2011). Very recently Alemany et al. (2014) showed that marine fronts represent important fishing areas even for demersal resources, as the distribution of fishing fleets and fishing effort are positively associated with frontal zones. In the case of hake juveniles, feeding on vertically migrating preys is an ecological characteristic that is presumably linked to the occurrence of chlorophyll-a fronts. Feeding intensity of hake was significantly correlated with major phytoplankton bloom events with a delay from one (Cartes et al., 2004) to two months (Hidalgo et al., 2008) presuming that most hake prey were pelagic (euphausiids, clupeids) and they may reach high densities after exploiting local phytoplankton blooms (Cartes et al., 2004).
Ecological niche models (ENMs), also termed Species distribution modelling (SDM) (Elith and Leathwick, 2009;Guisan and Thuiller, 2005;Peterson and Soberón, 2012) have become increasingly popular tools in the study of marine species distribution (Bentlage et al., 2013;Friedlaender et al., 2011;Tyberghein et al., 2012;Wiley et al., 2003). ENMs are spatially-explicit methods for modelling the ecological requirements of a given species and predicting its potential distribution in geographical space. They encompass numerous conceptual approaches and analytical tools, all underpinned by the niche concept formalized by Hutchinson (1957), i.e. the n-dimensional hypervolume formed by the environmental conditions a species can tolerate and within which populations can survive (Hirzel et al., 2002). ENMs basically work by distinguishing between ecological and geographical space. Given a set of species occurrence (or abundance) data across geographical space coupled with variations of a set of environmental factors in the same geographical space, ENMs are used to (i) explore the relationship between observed species occurrence and environmental variables, (ii) define, in the ecological space, the environmental variables that govern or limit the species distributional potential in the geographical space, (iii) predict, by projection back onto geographical space, the species occurrence also in areas where the distribution is unknown. Expected output of ENMs are maps of suitable or unsuitable habitats for studied species (Hirzel et al., 2002), which have many challenging applications in ecological studies (see Guisan and Thuiller, 2005 for a review) and are currently recognized as powerful tools for supporting appropriate management and conservation plans of marine resources. Despite these potentialities, to date few ENMs application have specifically addressed at a regional spatial level the Mediterranean Sea (Azzellino et al., 2012;Azzurro et al., 2013;Druon et al., 2012Druon et al., , 2011Langer et al., 2012;Sarà et al., 2013) and those based on distributional information of demersal species are rare (Hattab et al., 2013).
Data collected of decadal scientific surveys in the Mediterranean Sea coupled with environmental data from remote sensing and circulation models were used to apply the ENM approach, with the aim of identifying the most suitable environmental conditions which could promote the aggregations of 0-group hake in nursery areas. It is important to note that the present modelling approach refers to potential -rather than effective -habitats since the identified environmental conditions of nurseries are projected back onto space and time. The distribution of the realized nurseries and the dynamics of recruitment should, at model level, include other factors such as spawning stock biomass, connectivity from spawning to nursery grounds, predation and/or fishery pressure. This work will nevertheless provide relevant information to explain and monitor the environmentally-driven variability of hake recruitment, as well as to identify priority protection areas of hake recruitment.

Materials and methods
The methodological approach used in our ENM is essentially composed of four main steps (Fig. 1), namely: (1) identify the main life-history and ecological traits of 0-group hake based on literature; (2) process biomass indices for hake recruits and environmental covariates; (3) identify a suite and relevant thresholds of environmental variables related to the recruits' lifetime to describe the nursery habitat characteristics and finally (4) develop a habitat model to classify on a daily basis the degree to which each portion of the study area (model grid cell) is either suitable or unsuitable for recruitment. All variables were projected on the finest horizontal grid of the satellite ocean colour data which was used (NASA MODIS-Aqua sensor), i.e. at the resolution of 1/24°(about 4.6 km).

Biological traits of 0-group hake
This first step of ENM consists in gathering the relevant ecological traits of hake recruits as regards to life stages and relation to their environment, starting with the identification of 0-group hake in our dataset. We used the main source of standardised information about distributions, abundances and size compositions of the demersal resources in the region: the MEDITS bottom trawl survey program (Bertrand et al., 2002). Specifically we use data collected from 1994 to 2011 in the different FAO-GFCM Geographic Sub-Areas (GSA) of the EU and bordering countries (Spain, France, Italy, Malta, Slovenia, Croatia, Montenegro, Albania, Greece and Cyprus) (Fig. 2).
We considered as recruits those specimens that have settled on the bottom, becoming available to the fishing gear in well-defined habitats at the end of their larval -pelagic stage and which remain in these habitats before dispersing or migrating (Bartolino et al., 2008b). Due to the lack of large-scale studies on the dispersal behaviour of hake juveniles we assumed that the 0-group of hake was composed of specimens recently settled on the bottom and sharing similar habitat preferences. The spatial and inter-annual variability of the first cohort total length (TL) (Fig. 3) suggests that spawning period and growth rate are variable, depending on regional trophic conditions in addition to differences in sampling time (mostly from June to July with minor sampling in May and August). On the basis of previous studies on hake recruits in the Mediterranean Sea (Abella et al., 2008;Bartolino et al., 2008b;Colloca et al., 2009;Fiorentino et al., 2003;Murenu et al., 2010), we selected a threshold of 15 cm TL to identify the portion of the 0-group hake to be included in the model. This threshold corresponds to the 90th percentile size value of the 0-group (i.e. hake below 20 cm) in all GSAs (Fig. 3).
An important aspect of the ENM is to define the most significant period prior to sampling during which recruits were bound to the seabed at the end of their planktonic life phase. The literature reports a wide variability of the modal growth rate of 0-group hake in the Mediterranean Sea from 0.8 to 2.53 cm month À1 (Orsi Relini et al., 1992;Morales-Nin and Aldebert, 1997;Lleonart, 2001; Species' ecological analysis    Belcari et al., 2006) also in relation to environmental variability (Mellon-Duval et al., 2010;Morales-Nin and Moranta, 2004;De Pontual et al., 2013). Using a mean growth estimate of 1.25 cm month À1 , which is the most commonly reported value for important nurseries, hakes from 8.5 to 13 cm TL (first and third quartiles of all GSAs from survey) collected mostly in June-July (74% of hauls) were born 6.8 to 10.4 months earlier, i.e. from July-August to November-December of the previous year (see Table 1). The duration of the pelagic stage duration and fish size when settling to the bottom also appears to be variable. While Bozzano et al., 2005) found a minimum size of 5 cm TL for hake in settlement areas, Palomera et al. (2005) stated that M. merluccius begin to settle on the bottom at a size between 1.1 and 1.6 cm TL which corresponds to an age of over one month. Overall, the duration of the pelagic stage is therefore likely to be of about 1.5-2 months. Based on this, we estimated that bottom settlement for hake sampled during MEDITS surveys started in September-October of the previous year for the bigger sampled recruits to last until March-April for the smaller ones (see Table 1). This is in agreement with a seasonal minimum size of TL > 3 cm observed during winter and spring in the north-west Mediterranean Sea (Lleonart, 2001). Finally, previous studies have shown that European hake juveniles undertake daily feeding migrations towards the sea surface at night (Bozzano et al., 2005;Carpentieri et al., 2008;De Pontual et al., 2012;Orsi Relini et al., 1997, 1989Papaconstantinou and Stergiou, 1995). The diurnal migration of fish above 5-7 cm TL (Bozzano et al., 2005;Orsi Relini et al., 1997) started two to three months later, i.e. in December-January for the largest sampled fish and May-June for the smaller ones. Considering the above elements on variable growth rates and duration of life stages, we integrated the preferential habitat from February to June in order to take into account the environmental conditions that were effectively experienced by most of the sampled recruits.
Regarding the horizontal mobility of hake juveniles in their first months of life it is well known that the first stage of hake's life is  Diurnal migration pelagic. However, it is still not clear if larvae and post-larvae are passively transported by currents or they are able to use them to increase their mobility (Staaterman and Paris, 2014). Whether they reach passively or not the nursery area, it is hypothesised that 0group hake is relatively static horizontally at the scale of the model resolution (ca. 5 km) after they settle to the sea ground. When the vertical migration starts at about the size of 5-7 cm (Bozzano et al., 2005;Orsi Relini et al., 1997), the increased swimming capability presumably allows hake recruits to cope with the horizontal current to remain in the preferred habitats. Such spatial stability is also supported by the recurrence of the main hake nurseries found in different Mediterranean regions (Carlucci et al., 2009;Colloca et al., 2009;Fiorentino et al., 2003;Lleonart, 2001;Murenu et al., 2010;Tserpes et al., 2008).
Stating that 0-group hake have limited-mobility, we explored the range of favourable environment conditions for most of the demersal life of recruits extracting on a monthly basis the environmental variables of the grid cell corresponding to the selected hauls for the period 0-5 months before sampling (see also Table 1).

Data processing
The second step of our framework focuses on the collection and suitable preparation of input data for the model.

Biomass indices of hake recruits
The computation of abundance of hake recruits was performed using MEDITS data. The biomass index (hereafter BI, kg km À2 ) of each haul was calculated by only considering specimens of the size classes below 15 cm TL, by converting size-class numbers per haul into size-class weight per haul using the GSA-specific lengthweight (LW) curve provided by the literature and MEDITS regional coordinators (see Table 2) and by normalising the total weight by the estimated swept surface. The biomass index was preferred to the density index since the modelling refers to the feeding habitat and to the trophic relationship between primary productivity and the growth of 0-group hake. The third quartile biomass values (i.e. above 75th percentile considering all GSAs from 1994 to 2011) were finally selected to encompass the optimal environmental conditions for the growth of hake recruits in the most important regional nurseries. The selection of high biomass levels was necessary since recruits are present at low levels in the majority of hauls.

Chlorophyll
Surface chlorophyll-a concentrations (CHL) and fronts are used as a proxy for food availability for hake nurseries (see introduction above). CHL data were used at a daily time scale from the MODIS-Aqua (http://modis.gsfc.nasa.gov) (July 2002-2011) and SeaWiFS (http://oceancolor.gsfc.nasa.gov/SeaWiFS/) (1998-2010) ocean colour sensors. The MODIS spatial resolution of 1/24°(about 4.6 km) is used to identify meso-scale CHL fronts. The SeaWiFS resolution is double and CHL data were interpolated to the MODIS-Aqua's grid (which was chosen as the model grid). We used one set of CHL data at a time to derive a daily habitat and the sensor-specific habitat maps were then merged. This process allowed a substantial gain of habitat coverage, similar to the gain obtained by the merging of CHL data ($20%, Maritorena and Siegel, 2005), in relation to differences in observing time and thus in cloud cover. Daily CHL data were also pre-processed using iterations of a median filter in order to recover missing data on the edge of valid data. The median filter and Gaussian smoothing procedure (see Druon et al., 2012 for details) additionally allowed for the recovery of ca. 8% of the CHL data. The relative gain in coverage is much higher after the gradient calculation (CHL fronts) with +38% for the CHL gradients and +57% for the habitat coverage. The front enhancement of daily CHL data was calculated with an edge-detection algorithm which was showed to perform better than the histogram methods in detecting horizontal gradients given clear viewing conditions (Ullman and Cornillon, 2000). Note that the daily time scale is required here to allow the identification of CHL fronts which would be blurred or would disappear if using time-integrated data. If the daily data was used for front computation, we extracted a 3day mean CHL value in case daily data was unavailable in order to substantially improve both the comparison with hake biomass and model coverage. The 24 h variability of CHL level was thus stated to be low. The minimum value of chlorophyll-a horizontal gradient (gradCHL) in a 10 km-radius was also computed to estimate the threshold above which chlorophyll-a fronts are defined as relevant for hake nurseries.

Current velocity and temperature
As mentioned in the introduction, currents and temperature are likely to play a key role for hake nurseries. While temperature may simply be a tolerance criteria for the habitat, the relationship between hake recruits and Sea Surface Currents (hereafter SSC) is likely to be complex and indirect since, on the one hand, SSC are likely to influence productivity (chlorophyll-a levels and fronts) and, on the other hand, 0-group hake are likely to avoid surface waters with high currents to remain in their preferred habitat. In contrast, young hakes are always observed close to the ground during the day (Arneri and Morales-Nin, 2000;Orsi Relini et al., 1997) so that Sea Bottom Current (hereafter SBC) is likely to directly influence their distribution. We therefore focused on SBC rather than SSC and investigated specifically if SBC may impact hake nurseries through processes such as settlement at seabed and food availability. Moreover, high bottom current velocities prevent from deposition of the particulate organic matter and may limit 0-group hake feeding at seabed. Hake recruits are indeed recognised to be distributed over seabed substrate with high organic content (Maynou et al., 2003) and, reversely, to avoid sediment with very low organic matter content (Lleonart, 2001).
The physical data were produced by MyOcean Consortium (http://www.Myocean.eu), a marine core service within the European Global Monitoring for Environment and Security (GMES) Program whose objective is to develop an integrated capacity for ocean monitoring and forecasting. Monthly mean data of temperature and current velocities for the period 1998-2012 were extracted from a hydrodynamic model (Mediterranean Forecasting System) which has 72 unevenly spaced vertical levels and includes Table 2 Length-weight (LW) relationships specific of geographic sub-areas used for estimating the biomass index in kg km À2 . Values of a and b coefficients for the LW equation cm were provided through regional estimations mostly using MEDITS data (W in g and L in cm). The weight (in g) is provided for L = 15 cm. a variational data assimilation scheme for temperature and salinity vertical profiles and satellite sea level anomaly (Oddo et al., 2009). Original data at the resolution of 1/16°(ca. 6-7 km) were interpolated on the MODIS-Aqua grid. Monthly data were linearly interpolated to daily values. Such monthly to daily interpolation is believed to be relevant for detecting the seasonal changes that most impact hake nurseries, especially for sea bottom variables for which the variability is low. The meso-scale features are also assumed to be well represented by the similar original resolution than CHL. Sea Surface Temperature (hereafter SST) is taken from the upper model layer (ca. 3 m) while SSC is taken as the mean of the four upper layers of the MyOcean model (ca. 13.5 m) in order to capture the transport of the mixed layer. Sea Bottom Temperature (hereafter SBT) and current velocity are taken from the deepest model layer. The current intensity was investigated in the habitat model regardless the direction.

Bottom depth and seabed substrate
The bottom depth from GEBCO at 1.8 km resolution (http:// www.gebco.net/data_and_products/gridded_bathymetry_data/) was interpolated to the model's grid ($4.6 km) and the depth of hauls at the cell centre was extracted for comparison with MEDITS values. Relatively little difference of depth was observed (quartiles are 5, 14 and 30 m, i.e. 5%, 12% and 24% in relative values respectively). The low arithmetic mean of depth difference (1.7%) indicates that the difference is mostly due to the absence of interpolation from the cell centre (GEBCO) to the haul position (MEDITS). GEBCO depth data was used in the model to project back on the map the preferred depth range of nurseries using MEDITS data.
The main seabed substrate types were analysed in the western Mediterranean Sea using data from EUSeaMap habitat classification (EMODNET Project, http://jncc.defra.gov.uk/page-5040) at the resolution of 0.0027°(ca. 0.3 km). The central position of the haul was used to extract the seabed substrate at the original resolution using the intermediate classification (seven categories). This information is not yet available for the eastern Mediterranean Sea.

Nursery environmental characteristics
As a third step of our ENM we (a) select a set of environmental variables linked with the hake recruitment dynamics and, in relation to the recruits biomass, we (b) explore their seasonal variability during the recruits' lifetime in order to identify relevant environmental threshold values. The biomass data used in the model are haul values above the 3rd quartile, i.e. >8.4 kg km À2 (n hauls = 3355). The analysis of the effect of environmental variables on recruit biomass was made for the period 1994-2011 for the static variables (bottom depth and seabed substrate) while the period 1998-2012 was used for the physical variables of Myocean and 1998Myocean and -2010Myocean and and 2002Myocean and -2011 was used for the CHL data of SeaWiFS and MODIS-Aqua sensors respectively.
In order to analyse further the link of each selected environmental variable on the nurseries' productivity for the period 0 to 5 months prior to sampling, we performed a cluster analysis following the procedures reported in Berthold et al. (2010) and Hartigan (1975). The variables used were hake BI levels, bottom depth and the mean and standard deviation of the other variables, i.e. CHL (log transformed, MODIS-Aqua sensor), horizontal gradient of CHL (gradCHL, log transformed, MODIS-Aqua sensor), SBC and SBT. To estimate the similarity of data points between clusters, we used the Euclidean distance due to the normally set goal of minimizing the within-cluster sum of square errors, and we computed the k-means as clustering method (MacQueen, 1967).
In k-means clustering, the number of clusters k is first chosen and the cluster centres are initialised randomly. Each data point is then assigned to the closest cluster based on a selected distance measure (similarity) and updated cluster centre. At each iteration step, the new cluster centres are computed as the mean vectors of the assigned data points. These two steps, data point assignment and cluster centre update, are repeated until the cluster centres do not change any more or until a sufficient number of iterations are performed. The Matlab's k-means function was used with 500 iterations and the Euclidian distance setting. We performed before clustering the z-score-transformation (Berthold et al., 2010) where each data variable is normalised to zero mean and unit variance to guarantee that each selected variable has equal influence for the minimization of the within-cluster sum of squares objective function.
We tested sequentially from two to five clusters to conclude that four clusters described best the dependency of hake nurseries to the selected environmental variables at the scale of the Mediterranean Sea. The 5th and 95th percentile values were mostly chosen for the habitat model because they represent extreme environmental boundaries and these values were consistent with the clustering analysis.

Formulation of the ecological niche model
Once the environmental variables are selected and the threshold values are set using the dependency with hake biomass (mostly the 5 and 95 percentile values of high 0-group BI) and the cluster analysis, the last step consists in defining the specific ecological niche of hake nurseries through the areas of favourable feeding conditions (represented by CHL concentrations and CHL gradient) and tolerances of seabed temperature and current velocity. In order to classify the degree to which each cell is either suitable or unsuitable for the recruits of hake, a function was formalised considering a daily favourable habitat in the range from 0 and 1. The areas meeting daily the biotic and abiotic requirements of the habitat model are then integrated during the demersal stage to represent at best the most favourable environmental conditions of the 0-group hake sampled by the MEDITS surveys. The preferred habitat of hake nurseries is therefore expressed as relative frequency of occurrence and relates to the environmentally-driven potential development of 0-group hake independently of mortality by predation and fishing. These differences between potential and realised habitat impeded a classical validation exercise which compares the latter with the species' biomass. We however compared the habitat occurrence for each quartile of BI as a way to evaluate the capacity of ENM to depict the observed pattern of hake nurseries with particular focus on the absence of false negative. Interannual trends of BI and potential habitat integrated over the basin are shown to illustrate the year-to-year robustness of the model. A qualitative validation is also proposed with maps of BI and potential habitat for specific years. Fig. 4 maps the mean biomass of hake below 15 cm TL collected by the annual MEDITS surveys from 1994 to 2011 (all BI levels were used in this map). Lower BI levels (9-25 kg km À2 ) are usually present on most shelves except in the northern Adriatic Sea, while high BI levels (>25 kg km À2 ) are mostly located in the vicinity of the shelf break (e.g. northern Catalan Sea, eastern Gulf of Lions, south-west of Sardinia, northern Tyrrhenian Sea, south-east of Sicily and Saronic Gulf in the southern Aegean Sea). Most of areas with mean BI levels above 35 kg km À2 show occasional values above 100 kg km À2 .

Biomass index and depth distribution
Hake nurseries were mostly found on the shelf and shelf break area in the range 38-312 m (5-95th percentile values), with a median value at 119 m (Fig. 5). However, highest BI values were globally found in the shelf break area (about 175-375 m, Fig. 5).

Environmental variables during the sampling period
During the sampling period, mostly from May to August, high 0group abundance occurs in areas where the daily surface chlorophyll-a (CHL) is relatively low, i.e. in the range 5-95th percentile of 0.09-0.56 mg m À3 for MODIS-Aqua sensor (n = 751) and of 0.10-0.78 mg m À3 for SeaWiFS sensor (n = 665) with values for the 50th percentile values being 0.16 and 0.19 mg m À3 , respectively. The 5th percentile value of chlorophyll-a gradient is 6.6 Â 10 À4 and 9.1 Â 10 À4 mg m À3 km À1 for MODIS-Aqua and SeaWiFS respectively. Both CHL and gradCHL at the day of sampling showed no correlation with 0-group hake BI of the upper quartile (r = 0.07, p = 0.01). These biotic limits, considered to be a first estimation of the model parameters during summer, primarily represent mesotrophic areas of the shelf and shelf break between the eutrophic waters under the influence of river plumes and the generally oligotrophic oceanic waters.
Since part of the 0-group hake are migrating vertically for feeding, we investigated the potential tolerance of both the surface and near bottom physical variables on the biomass. The comparison of SST/SBT during summer sampling for all the model grid cells from 38 to 324 m and the SST/SBT at the location of high recruit BI shows no particular difference for lower limits ($0.15°C, Table 3). On the opposite, the SST and SBT 95th percentile values of hake nurseries is lower by 1.2 and 2.8°C respectively for the same depths and months than the surrounding environment. SST and SBT at the day of sampling showed no correlation with high 0group hake BI. Hake nurseries thus appear to be strongly limited by high temperature during summer, and especially near the seabed, so as to focus in this study on SBT rather than on SST (also in agreement with the reasoning on SBC/SSC above).
SBC shows no correlation with hake BI during the sampling period (|r| < 0.02, p > 0.05), but 95% of hake nurseries during summer are in areas characterised by low SBC (under 0.026 m s À1 , Table 3).
Regarding seabed substrate, the observed upper BI quartile of 0group hake in the western basin showed the following distribution among the main categories: mud (42%), muddy sand (33%), sandy mud (13%), sand (10%), coarse sediment (2%), seagrass and rock (<1%). If the preferred substrate of 0-group hake contained a component of mud, the 10% of nurseries in the western Mediterranean Sea which are still on sand cannot be neglected. It is suggested that  inaccuracies in the substrate classification might have occurred in the original database due to a highly heterogeneous substrate distribution on the continental shelf compared to the limited records. At any rate, because the substrate type might not be accurate or discriminant enough and because the coverage is lacking for the eastern basin, this environmental characteristic was not selected for the model. Additionally, note that the low SBC and the chlorophyll-a fronts highlighted in the following sections are likely to emphasize areas of relatively organic-rich sediments, i.e. with a muddy component.

Environmental conditions of nursery habitats during months prior to sampling
Although some significant correlation was identified for CHL and SBC with hake BI for the period 0-5 months before sampling, the correlation still remains weak (|r| < 0.04). This result suggests that multiple and/or temporary unfavourable conditions are likely to affect the biomass.
Compared to summer, levels of preferred CHL for the 0-5 months period before sampling increase with the 5-95th percentile values of 0.10-0.91 mg m À3 and a median of 0.26 mg m À3 (MODIS-Aqua sensor). The range of SBC is also wider with 5-95th percentile values of 0.001-0.034 m s À1 while SBT values are globally lower with 5-95th percentile values of 11.8-15.0°C (Fig. 6). Note that extreme levels of BI (>100 kg km À2 , n hauls = 228) only occur if the physical conditions are stable from winter to summer. Biomass levels of hake recruits above 100 kg km À2 are indeed characterised by a narrow and stable range of SBT (13.8 o C ± 1°C) and by low SBC (below 0.032 m s À1 ). The higher CHL levels than of BI above 8.4 kg km À2 appear to be marginal compared to the physical limitations except for the minimum values which appear to be more sustained.
The cluster analysis resulted in a descriptive classification of hake nurseries at the scale of the Mediterranean Sea (see Figs. 7 and 8). Cluster 1 (blue star symbols) is the largest with 42% of hauls within the last BI quartile. This cluster is characterised by intermediate BI levels, relatively high depths (shelf break area), low CHL, Table 3 Comparison of lower and upper limits (5th and 95th percentile values) of sea surface temperature (SST), sea bottom temperature (SBT) and current (SBC) in high BI (Biomass Index) of 0-group hake (1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)  ⁄ The area is also delimited by the minimum and maximum of longitude and latitude of MEDITS data. Values in bold highlight substantial difference between high BI of 0-group hake and their surrounding environment, i.e. revealing a particular preference or tolerance. Fig. 6. 5th, 50th and 95th percentile values of (a) surface chlorophyll content (3-day mean, MODIS-Aqua sensor), (b) sea bottom temperature (SBT) and (c) sea bottom current (SBC) from 5 to 0 months prior to sampling for all hauls as regards to 0-group hake biomass (above the third quartile, solid line, and above 100 kg km À2 , plus dash line).
gradCHL and SBC levels as well as low variability of SBT around the optimum temperature (see Fig. 8 for median values). Cluster 2 (yellow circles) represents 26% of hauls and corresponds to the highest BI levels, i.e. about 70% higher than the other three clusters. Environmental conditions of cluster 2 are similar to cluster 1 (shelf break area) except for productivity levels (CHL and gradCHL) which are intermediate and more stable. Cluster 3 (red plus symbols) represents 21% of high-biomass hauls. While the hake BI of cluster 1 and 3 are similar, cluster 3 is characterised by lower depths (mid-shelf), higher productivity (higher CHL, gradCHL) and by the highest variability of SBT. Cluster 4 (green crosses) represents only 11% of high-biomass hauls. Cluster 4 is characterised by intermediate depths (outer shelf), productivity levels and SBT variability and by particularly high SBC levels.
In other words, highest concentrations of hake recruits (median of 33 kg km À2 , cluster 2) occur in relatively deep waters (shelf break area) where conditions of bottom temperature, bottom currents are seasonally stable, and where the plankton productivity is enhanced and stable (lower standard deviation values). In case of normal (though relatively low) plankton productivity over the shelf break area where the physical conditions are stable, BI values of recruits are intermediate-high (median of 20 kg km À2 , cluster 1). The same intermediate-high BI values (median of 19 kg km À2 , cluster 3 and 4) are encountered on the mid-and outer shelf where plankton productivity is higher together with the bottom temperature variability (cluster 3-4) and bottom current levels (cluster 4).

Habitat modelling and parameterization
The cluster analysis highlights the generally low trophic conditions of the deep -and often productive -hake nurseries. At depths greater than ca. 130 m, nurseries are equally or more productive and more persistent than at shallower depths due to the higher stability of physical conditions and despite equal or lower CHL levels. With the purpose of realistically representing highly different levels of trophic conditions in hake nurseries, we introduced in the habitat model a second level of trophic habitat (sub-optimal) over which the abiotic limitations apply. The optimal trophic habitat represents the larger frontal systems which, by their size and persistence, identify productive water masses with potentially well-developed food webs (most of clusters 2, 3 and 4). The sub-optimal trophic habitat refers to smaller -less productive -frontal systems which may be sufficient to sustain intermediate levels of 0-group hake biomass levels within an optimal physical environment (most of cluster 1). We defined three threshold values for CHL and two for gradCHL that delimit the optimal, sub-optimal and unfavourable trophic habitat with daily values of 1, 0.3 and 0 respectively (Fig. 9). The value of 0.3 was chosen as an ad-hoc value for the sub-optimal trophic habitat as it represents a substantially less favourable feeding habitat (about 3-fold) than the optimal trophic conditions (of value 1) and is markedly above 0.
As a result of the above considerations, the preferred habitat of 0-group hake has two levels of trophic proxies (small and large CHL gradient and content), a preferred range of bottom depth and bottom temperature and a maximum value of sea bottom current. The flowchart of the habitat model translates into the following equation with a resulting daily favourable habitat of value 0, 0.3 or 1 for each grid cell (refer to Fig. 9 for the trophic habitat): Hake Nursery Habitat Day;Cell ¼ Trophic Habitat 0=0:3=1 Ã Depth range 0=1 Ã SBT range 0=1 Ã SBC max 0=1 The model was parameterised using the preferred range of the selected environmental variables consistent with high productivity of hake recruits. The large 5th and 95th percentile values of the detected environmental conditions were chosen as boundary limits for suitable nursery habitat (Table 4) because (a) they are sensible values bearing in mind they correspond to the upper quartile BI values of hake recruits during their development period (0-5 months prior to sampling), (b) these values agree with the interpretation made out of the cluster classification Figs. 7 and 8. A slightly wider range (2.5-97.5th percentiles) was selected for bottom depth to account for differences between MEDITS and GEBCO data in shallow and deep environments. The intermediate thresholds for CHL and gradCHL were chosen using the cluster analysis and the differences between the intermediate-high and high hake biomass values in deep grounds (clusters 1 and 2).

Outputs of the habitat model
We present in this section the spatio-temporal distribution of model results and we evaluate them with observed biomass levels. Potential habitat for hake nurseries from February to June for the contrasted years 2008 and 2011 are shown with the corresponding distribution of recruit's biomass (Fig. 10). The biomasses in the main nurseries as described in Fig. 4 show medium to high values for 2008, while levels are substantially lower in 2011. The main nurseries and the overall substantial decrease in 2011 are generally predicted by the model. However restricted areas representing a false-positive (low hake biomass and high preferred habitat) are identified in the considered period (2008), e.g. in the north Aegean Sea (potential overestimation). Note that shallow nurseries appear to be more vulnerable to environmental change since merely all nurseries on the shelf show low BI levels and habitat occurrence in 2011 compared to 2008. Fig. 11 shows the annual occurrence of favourable habitat in the cases of biomass absence, below and above the third quartile biomass. A minimum of about 25% of favourable habitat is required for a biomass in the upper quartile (>8.4 kg km À2 ) and corresponding median values of favourable habitat are in the range 40-65%. False negative prediction is fairly restricted while false positive is more  9. Definition of the three trophic habitats (unfavourable, sub-optimal and optimal, of value 0, 0.3 and 1 respectively in the model) based on levels of surface chlorophyll content (CHL) and horizontal chlorophyll gradient (gradCHL). The suboptimal and optimal trophic habitats refer to small and large productive frontal systems respectively. common (as expected for a potential habitat) but median values of biomass absence show no favourable habitat. Median levels for each quartile are generally consistent with the occurrence of habitat for a given year (the higher biomass, the more frequent the favourable habitat) except in 2007 where the median values for the low and high BI levels are at the same level. Fig. 12 spatially details the mean occurrence of favourable 0group hake habitat from 1998 to 2011, and therefore provides an Table 4 Model parameters defining hake nurseries in the Mediterranean Sea. The preferred range corresponds to the 5th and 95th percentile values of high-biomass hake recruits (above the third quartile) in the period 0 to 5 months before sampling (except for depth where the 2.5-97.5th percentiles values were used).  estimate of the nursery habitat persistence. The comparison of the decadal biomass distribution (Fig. 4) or the year-to-year comparison ( Fig. 10) between observation and prediction shows an overall spatial agreement of nurseries. However, there is a potential overestimation of the prediction in the southern Adriatic and northern Aegean Seas, and an underestimation in south-west Sardinia and the Strait of Sicily. The monthly variability of favourable habitat expressed in mean relative surface (Fig. 13) indicates an increase of potential nurseries of about 60% in winter-spring compared to summer-autumn at the scale of the Mediterranean Sea with $1.7% and 1.1% of the basin respectively. The mean favourable habitat from February to April is more than twice more favourable than from September to November. Inter-annual variability is however high with differences of winter-nurseries occurrences up to ±50% (e.g. 2005-2007-2008 against 2010-2011). The last years of the study (2010)(2011)(2012), and particularly 2011, showed a substantial decrease of the annual favourable habitat overall in the Mediterranean Sea of about 28%. The annual trends of habitat and recruit BI medians (solid line and square solid line in Fig. 14) are in agreement for the period 2003 to 2011 except for 2009 for which an increase of favourable habitat is associated with a decrease of recruit biomass. The median biomass of juveniles and adults (square dash line, MEDITS data, Fig. 14) generally follows the trend of recruits with one or two years delay. Total hake landings of the Mediterranean (Europe area, FAO) peaked in the mid-1990s at 52,000 tons before to sharply decrease to an average of $25,000 tons in the period 2000. From 2009, the biomass of recruits and juveniles/adults as well as the nursery habitat and total landings show all a substantial negative trend.

Discussion
The availability of standardised indices of abundance by size of hake together with environmental data covering wide areas of the Mediterranean allows the investigation of relationships among 0group hake abundance and the ecological factors that affect (c) 2011 hake 0-15cm biomass (kg.km -2 , observed, n hauls = 564) (d) 2011 potential hake nurseries (predicted occurrence of favourable habitat) Fig. 10 (continued) recruitment dynamics in space and time. The use of the ecological niche model (ENM) approach for modelling potential hake nurseries in the Mediterranean Sea provides a synoptic view of the recruitment carrying capacity of this widely-spread species. The biomass distribution of hake recruits can be explained by the balance between food availability and hydrological tolerance over most of their lifetime near the sea bottom of the shelf and shelf break areas. The large agreement of potential hake nurseries with the observed distribution of recruit biomass in space and, to a lower degree, in time demonstrates the robustness of using a modelling approach that is largely driven by a species' ecology. Although no strict quantitative validation of the model is given, the quasi-absence of high 0-group biomass related to low potential habitat (model false-negative, Fig. 11) provides good confidence in the model results. The difficulty to strictly validate the results is a limitation of the approach, together with the horizontal resolution (ca 4.6 km) in archipelago areas (e.g. in the Aegean Sea and along the Croatian coasts) where the ocean model is unlikely to capture the level of required environmental variability.
The high recruitment success of hake and the ability of this species to sustain important fisheries may result from its biological and ecological flexibility. Hake spawning takes place over extended time periods, recruits are found over large areas of the continental shelf and the shelf break area, and hake is able to feed on a wide range of trophic resources (Lleonart, 2001). The habitat model identified biotic and abiotic limitations of the recruitment success of hake nurseries in agreement with the observation, and particularly regarding tolerance levels for temperature, water currents and proxy for food availability. Hake prefers cold waters and has a limited tolerance to warm temperature (De Pontual et al., 2013;Jolivet et al., 2012). Indeed, young hake were not observed in areas of the north-west Mediterranean Sea where bottom temperature is above 15°C (Lleonart, 2001). It is worth noting that hake eggs showed a temperature preference in the range of 10.5-12°C in the west of British Isles (Coombs and Mitchell, 1982). Similar ranges of temperature than highlighted by the model are found in the shelf break area of the Bay of Biscay where juvenile hake are abundant (11.3-13.9°C at 210 m, Casey and Pereiro, 1995). Our results suggest that about five months of stable biotic and abiotic conditions favour the success of hake recruitment. Besides the narrow range of preferred temperature that can be found in the outer shelf and shelf break area, favourable conditions include a high frequency of productive fronts that likely enhance prey availability in line with the findings of Alemany et al. (2014) who identified a strong correlation between semi-permanent fronts and demersal resources and even a stronger link for fish preys of low trophic level. The high resolution and integrated characteristics of biotic variables were in particular recognised to Fig. 11. Distribution of hake below 15 cm TL biomass (kg km À2 , MEDITS data) with the corresponding mean occurrence of potential habitat from February to June (in % of available habitat detection including SeaWiFS and MODIS-Aqua sensors) at haul location from 2003 to 2011. The same boundary values were used for each year (no biomass, below and above the third quartile biomass level of 8.4 kg km À2 ). favour sound predictive science in the field of demersal fish habitat determination (Johnson et al., 2013). Our results also support the hypothesis that low current velocity at sea bottom would favours the settlement of hake post-larvae after their pelagic stage and/ or the feeding of hake recruits in habitat favouring the deposition of particulate organic matter (Maynou et al., 2003). More generally, this work suggests that hake nurseries are preferably located where productive fronts frequently occur near the shelf break area. While in movement, these productive fronts may stand long enough, i.e. for several weeks or months, to sustain a full trophic chain and attract predators, including the 0-group hake. The frequency of occurrence of these meso-scale features in the shelf break area may well determine the carrying capacity of demersal nurseries. A particular case in the present model is represented by the Central Adriatic Sea, where the main nursery of hake is situated in the Pomo pit (Arneri and Morales-Nin, 2000). In such area the nutrient-rich waters generated in winter in the northern sector accumulate in this depression (Artegiani et al., 1997), making the Pomo pit a peculiar site in the Mediterranean of strong nutrient re-cycling processes.
The link of the 0-group mean biomass with total landings one year later since 2003 suggests a higher relative catch of age-1 specimens and increased fishing pressure after the sharp decrease of total landings after the mid-1990s (Fig. 14). In addition, hake land-  ings for Spain, Italy and France decreased by 28% in 2012 compared to 2010 (STECF, https://fishreg.jrc.ec.europa.eu/web/datadissemination/home) which correlates with the substantial decrease of 0-group biomass and preferred habitat in 2010 and 2011. This decrease in the 0-group recruitment could be due to an environmental pressure -a lower frequency of productive fronts likely in relation with the enhancement of the seasonal thermoclinewhich may have strongly limited the availability of food for recruits. Note that fishing mortality of recruits should have decreased since 2010 with the legally-binding use of a larger mesh size in EU countries. The MEDITS surveys represent an appropriate sampling strategy since they cover most of hake nurseries of the northern Mediterranean Sea with a standardised protocol and sampling generally occurs after the most favourable period for recruitment (see Fig. 13). The habitat model reveals in addition important potential hake nurseries with good spatial agreement in unsampled areas such as in the shelf and shelf break grounds off Morocco and Tunisia (Fig. 12) in agreement with other sampling data (see CopeMed II, 2012 for Morocco and Garofalo et al., 2008 for Tunisia). The prediction at a large spatial scale also allows foreseeing the potential connectivity of nurseries with spawning and juvenile grounds, which is an essential element for further understanding the ecology of that species and contributing to identification of stock boundaries of hake in the Mediterranean. Past studies based on trawl survey data have identified important hake nursery areas in the Mediterranean but the findings reflect the situation during the sampling season; thus only temporal ''snapshots'' were provided. The approach followed in the current study allowed to describe most of the variability of nurseries even if, as mentioned above, potential and effective habitats show differences notably through the mortality by predation and fishing. A potential second peak of recruitment is highlighted by the model as a result of a favourable habitat during summer and autumn (result not shown) mostly in the eastern Gulf of Lions, the northern Aegean Sea and specific areas of the eastern Ligurian and northern Tyrrhenian Seas. This is in line with past studies that have reported relatively high number of recruits in those areas during the summer-autumn months (Abella et al., 2008;Belcari et al., 2006Belcari et al., , 2001Orsi Relini et al., 2002;Tserpes et al., 2008). These seasonally-persistent nurseries are particularly important for hake recruitment since they appear to be favourable areas even under severe environmental pressure such as in 2011 (Fig. 10).
The convergence, front formation and near bottom stability at the shelf break can all be factors helping small sized hake to avoid dispersion towards unfavourable habitats after the bottom settlement. The limited horizontal mobility of 0-group hake of about few model cells of $4.6 km over several months thus infers that the seasonal integration of potential habitat represents what the sampled fish experienced in terms of environmental conditions. An analysis of the effective connectivity with the spawning population would also be likely required to locate the most productive areas within the optimal habitat (Nagelkerken et al., 2013). An attempt to depict the most temporally persistent nursery areas of hake in the Mediterranean Sea based on 0-group biomass has been carried out within the EU project MEDISEH (http://mareaproject.net/contracts/5/overview/). Besides the differences between potential and effective habitats, a preliminary comparison with the results obtained by this project showed that our modelling approach correctly located the potentially suitable and persistent hake recruitment habitats.
The estimation of the nursery carrying capacity due to environmental factors is essential information for fisheries management. The monitoring of strength of recruitment and nurseries of hake provides the basis for a preventive management strategy, where high inter-annual variability should be taken into account for precautionary management purposes (Lleonart, 2001). This approach also allows identifying optimal spatial protection measures to reduce fishing mortality of undersized recruits as specifically required by the Council Regulation (EC) 1967/2006 for the Mediterranean Sea. Such habitat mapping could indeed enhance fisher avoidance of these areas in the context of the landing obligation of the revised EU-Common Fisheries Policy. The spatial based approach to fishery management becomes a necessity since it (i) introduces a critical environmental dimension in stock assessments, (ii) will help in appraising the impact of spatial measures in management strategy evaluations (provided the spatial dimension is introduced in stock assessment models), (iii) will contribute to the inclusion of fishery considerations in the EU Maritime Spatial Planning (European Commission COM (2013) 133 final, n.d.) and (iv) will support the implementation of an Ecosystem Approach to Fisheries Management as recently advocated by the revised Common Fisheries Policy (Regulation EU 1380/2013, Part III, Title 1, Article 9).
In conclusion our results suggest that productive fronts in the shelf break area represent important nursery areas for European hake, a widely-distributed demersal species in the Mediterranean Sea. A stable abiotic environment for bottom temperature (SBT of about 13.8°C ± 1°C) and bottom currents (SBC below 0.034 m s À1 ) are also necessary conditions for a high biomass of hake recruits. Information on the essential habitat for hake reproduction in the Mediterranean Sea shall contribute to properly implement the spatial management of fisheries required by EU policies, and notably to limit recruits' mortality by fishing.