Pelagic primary production in the coastal Mediterranean Sea: variability, trends, and contribution to basin-scale budgets

We estimated pelagic primary production (PP) in the coastal (<200 m depth) Mediterranean Sea from satellite-borne data, its contribution to basin-scale carbon fixation, its variability, and long-term trends during the period 2002–2016. Annual coastal PP was estimated at 0.041 Gt C, which approximately represents 12 % of total carbon fixation in the Mediterranean Sea. About 51 % of this production occurs in the eastern basin, whereas the western and Adriatic shelves contribute with ∼25 % each of total coastal production. Strong regional variability is revealed in coastal PP, from high-production areas (>300 g C m⁻²) associated with major river discharges to less productive provinces (<50 g C m⁻²) located in the southeastern Mediterranean. PP variability in the Mediterranean Sea is dominated by interannual variations, but a notable basin-scale decline (17 %) has been observed since 2012 concurring with a period of increasing sea surface temperatures in the Mediterranean Sea and positive North Atlantic Oscillation and Mediterranean Oscillation climate indices. Long-term trends in PP reveal slight declines in most coastal areas (−0.05 to −0.1 g C m⁻² per decade) except in the Adriatic where PP increases at +0.1 g C m⁻² per decade. Regionalization of coastal waters based on PP seasonal patterns reveals the importance of river effluents in determining PP in coastal waters that can regionally increase up to 5-fold. Our study provides insight into the contribution of coastal waters to basin-scale carbon balances in the Mediterranean Sea while highlighting the importance of the different temporal and spatial scales of variability.

1 Introduction

Coastal ocean waters (i.e., <200 m depth) are an important link between the land and the open ocean. They act as a buffer between terrestrial and human influences and the open ocean (Liu et al., 2000). Despite their relatively reduced extension (∼7 % of ocean surface area; Gattuso et al., 1998), they behold some of the most productive habitats on the planet. Therefore, they have a disproportionate importance in many basin-scale biogeochemical and ecological processes, including carbon and nitrogen cycling, and in the maintenance of marine diversity (Cebrian, 2002; Coll et al., 2010; Dunne et al., 2007). Besides, biological production of continental shelves supports over 90 % of global fish catches (Pauly et al., 2002; Pauly and Christensen, 1995).

Coastal seawaters support high primary production (PP), contributing to some 10 % of global ocean PP and up to 30 % if estuarine and benthic production is considered (Ducklow et al., 2001; Muller-Karger et al., 2005). These high rates of organic productivity occur in the coastal oceans due to the rapid turnover of the large inputs of nutrients and organic carbon from land. PP drives a significant carbon sink in the ocean (Field et al., 1998; Laws et al., 2000) and is a key regulator of ecological processes such as elemental cycling, trophic structure variabilities, and climate change (Bauer et al., 2013; Chavez et al., 2011). In coastal waters, physical and biological processes enhance the carbon transport out of the continental margins into the deep layers of the oceans, thus connecting terrestrial with deep oceanic systems (Cai, 2011; Carlson et al., 2001; Cole et al., 2007). The productivity of coastal sea areas is also of strategic socio-economic importance for many countries considering that PP constrains the amount of fish and invertebrates available to expanding fisheries, a primary resource for many coastal human communities (Chassot et al., 2010). The estimation and understanding of PP evolution and trends in the coastal seas is therefore essential to improve our knowledge of the oceanic carbon cycle.

Scaling up local measurements to estimate the contribution of coastal regions to global carbon fluxes has been hindered by the high spatial and temporal heterogeneity of these waters. Global models of oceanic systems produce carbon fixation estimates with a high degree of uncertainty in coastal regions (Muller-Karger et al., 2005). Coastal waters are complex because of the tight connection between terrestrial and oceanic systems. Terrestrial uploads of nutrients and organic matter originating from groundwater discharges, flash floods, or river runoff as well as exchanges with the seafloor strongly control the productivity of these waters (Basterretxea et al., 2010; Woodson and Litvin, 2014). The amplitude of seasonal variation in surface chlorophyll (chl) and surface temperature is often higher in coastal waters compared to the open ocean (Cloern and Jassby, 2008). Furthermore, coastal topography and its interaction with winds, waves, and currents generates a high variety of physicochemical niches for phytoplankton growth. Likewise, benthic–pelagic coupling allows the remineralization of nutrients present in shelf sediments during most intense storms. These episodic variations may constitute an important contribution to the overall productivity of shelf waters. Because of the high spatiotemporal heterogeneity in the main coastal subsystems and the concomitant lack of data, most estimated carbon fluxes in these subsystems have relatively high uncertainties (Bauer et al., 2013). In addition, direct human activities and climate change lead to a long-term variation in terrestrial fluxes and coastal biogeochemistry that can potentially have important consequences for the global carbon cycle (Gregg et al., 2003).

In the Mediterranean Sea, coastal and shelf areas represent about 21 % of the global basin (259 000 km²), which is a higher contribution than for the global ocean (Pinardi et al., 2006). Although the Mediterranean Sea is included amongst the most oligotrophic areas of the world oceans, it can display marked spatial productivity variations related to the variety of regional climate and oceanographic conditions as well as related to the multiple land-derived fluxes that locally fertilize the coastal waters (Goffart et al., 2002). Nutrient-rich effluents from human activities on the coast (domestic wastewater, fertilizers, industry, etc.) and natural river discharges affect continental shelf productivity in this sea, sustaining locally enhanced pelagic and benthic biomass. Nevertheless, the influence of some river flows has been notably reduced by damming affecting water chemistry and sediment loads and, thereby, the productivity of coastal waters at local and regional scales (Ludwig et al., 2009; Tovar-Sánchez et al., 2016). Moreover, intensive agricultural practices and urbanization have brought unprecedented use and contamination of coastal groundwater (Basterretxea et al., 2010; Tovar-Sánchez et al., 2014). For example, the use of fertilizers has resulted in higher nutrient flow into the Adriatic and in the lagoons of the Nile River, which has led to eutrophication (Turley, 1999). However, the impact of this anthropogenic nutrient enrichment may vary between regions, and modeling projections suggest spatial variations in PP as a result of climate change.

Accurate quantification of the coastal PP is fundamental for assessment of global carbon cycling in the Mediterranean Sea. Changes in PP have important effects on fish stocks that are socially relevant because of the economical dependency of many Mediterranean coastal communities on marine food products. Several studies have assessed PP at the scale of the entire Mediterranean Sea from satellite remote sensing data (Bosc et al., 2004; Bricaud et al., 2002; Lazzari et al., 2012). However, coastal areas were generally ignored in such studies, so that their contribution to basin-scale budgets is still largely unknown. Most coastal studies have a focus on specific regions and/or times (Estrada, 1996; Marty et al., 2002; Moutin and Raimbault, 2002; Rahav et al., 2013). Observed rates of climate change in the Mediterranean basin exceed global trends (Cramer et al., 2018) and future warming in the Mediterranean region is expected to be above global rates by 25 % (Lionello and Scarascia, 2018). Long-term responses of PP in coastal areas to climate forcing remain uncertain however because of the scarcity of adequate field datasets (Gasol et al., 2016).

In this study, we present major characteristics of pelagic PP in Mediterranean coastal waters based on satellite-borne observations for the period 2002–2016. First, we provide global estimations of PP in coastal waters, and we assess their contribution to basin-scale PP, their interannual variability, and long-term trends. Then, we regionalize the coastal waters based on their temporal patterns of pelagic PP using self-organizing maps (SOMs), and we analyze the contribution of each region to total coastal PP.

2 Materials and methods

2.1 Remote sensing data

We used the Mediterranean Sea Level-3 reprocessed surface chlorophyll concentration product (Chl L3) obtained from the EU Copernicus Marine Environment Monitoring Service (CMEMS). This product merges multi-satellite observations, and it is available at 1 d and 1 km resolution (https://resources.marine.copernicus.eu/product-detail/OCEANCOLOUR_MED_CHL_L3_REP_OBSERVATIONS_009_073/, last access: 2 August 2019). Specifically, the dataset used is “dataset-oc-med-chl-multi-l3-chl_1km_daily-rep-v02”, and the variable name used is “mass_concentration_of_chlorophyll_a_in_sea_water (Chl)” obtainable in a NetCDF-4 file format. This Chl L3 dataset is derived with an updated version of the regional algorithm MedOC4 (Mediterranean Ocean-Colour 4 bands MedOC4; Volpe et al., 2019) for deep pelagic Case-1 waters and the AD4 algorithm (ADriatic 4 band; Berthon and Zibordi, 2004; D'Alimonte and Zibordi, 2003) for Case-2 coastal waters (generally shallow and turbid waters).

Level-2 sea surface temperature (SST, ^∘C) at 1 d and 1 km was obtained from every available orbit from the Moderate Resolution Imaging Spectroradiometer (MODIS) aboard the Terra and Aqua satellites. Data were downloaded from the National Aeronautics and Space Administration (NASA) archive website (http://oceancolor.gsfc.nasa.gov/, last access: 18 July 2017). Only nighttime orbits were selected to avoid problems with skin temperature during daylight. Orbits with quality flags 0 (best), 1 (good), and 2 (questionable) in SST were included after checking their validity and accuracy in order to have a more complete dataset. Daily (24 h averaged) photosynthetically active radiation (PAR, in E m⁻²) was obtained as a Level-3 product at 9 km and 1 d resolution. This is the best available resolution at the NASA archive of MODIS and Medium Resolution Imaging Spectrometer (MERIS) data (https://oceancolor.gsfc.nasa.gov/l3/, last access: 30 December 2017).

All satellite-derived variables were remapped onto a regular 1 km spatial grid over the study area, by averaging all available pixels within each grid cell. For each parameter, outliers were removed whenever they exceeded 3 times the mean ± SD (standard deviation) of the time series. For the purpose of this study, coastal areas were defined as the waters lying between 5 and 200 m depth. Only values at depths exceeding 5 m depth were considered in order to reduce the possible influence of seafloor vegetation reflectance in chlorophyll concentration values (chl, mg m⁻³) at shallow waters. The analyzed time series covers the period from January 2002 to December 2016 for the Mediterranean Sea (30 to 46^∘ N and 6^∘ W to 37^∘ E, Fig. 1).

Figure 1 Map of the Mediterranean Sea showing the main basins, sea regions, surrounding countries and major rivers. Bathymetric data were obtained from ETOPO1 (Amante and Eakins, 2009). The black contour indicates the 200 m isobath, the limit of coastal waters as defined in the present study.

2.2 Primary production estimates

PP was estimated from satellite-derived chl, SST, and PAR values using the time, depth, and wavelength-resolved light-photosynthesis model of Morel (1991). This model was previously used for estimating PP in the Mediterranean Sea (Antoine and André, 1995) and at global scale (Antoine et al., 1996; Antoine and Morel, 1996) and performs well when compared to in situ measurements (Campbell et al., 2002; Friedrichs et al., 2009) or when compared to other similar algorithms designed for use with satellite observations (Carr et al., 2006; Saba et al., 2011). Instantaneous production at depth z (m) of the water column, time t of the day, and for absorption of irradiance at wavelength λ, $P(\mathit{\lambda },z,t)$ is calculated as

$$\begin{array}{}\text{(1)}& \begin{array}{rl}& P(\mathit{\lambda },z,t)=E(\mathit{\lambda },z,t)\text{chl}\left(z\right)a^{\ast }(\mathit{\lambda },z)\mathit{\Phi }\\ & \left(\mathrm{g}\mathrm{C}{\mathrm{m}}^{-\mathrm{3}}{\mathrm{s}}^{-\mathrm{1}}\right),\end{array}\end{array}$$

where $E(\mathit{\lambda },z,t)$ is the spectral scalar irradiance for wavelength λ, depth z, and time t of the day ($\mathrm{mol}\mathrm{photons}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{s}}^{-\mathrm{1}}$); $a^{\ast }(\mathit{\lambda },z)$ is the spectral chlorophyll-specific absorption coefficient of phytoplankton (m² mg chl⁻¹); and Φ is the quantum yield of photosynthesis for carbon fixation (mol C mol per photon; its possible spectral changes are ignored). Note that chl, a^∗, and Φ are not made variable with time.

The triple integration of (1) with respect to wavelength, depth, and time gives the daily column-integrated primary production, PP:

$$\begin{array}{}\text{(2)}& \begin{array}{rl}& \text{daily\hspace{0.17em}PP}=\mathrm{12}\underset{\mathrm{0}}{\overset{D}{\int }}\underset{\mathrm{0}}{\overset{min(Z_{\mathrm{p}}/Z_{\mathrm{b}})}{\int }}\underset{\mathrm{400}}{\overset{\mathrm{700}}{\int }}\mathrm{P}(\mathit{\lambda },z,t)\mathrm{d}\mathit{\lambda }\mathrm{d}z\mathrm{d}t\\ & \left(\mathrm{g}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}\right),\end{array}\end{array}$$

where the factor 12 is the conversion from moles to grams of carbon, D is the day length or hours of daylight (h), Z_p (m) is the depth where the photosynthetically available radiation (PAR) falls to 0.1 % of its value just below the sea surface (so approximately 1.5 times the euphotic depth), and Z_b is the bottom depth taken from the ETOPO1 database (Amante and Eakins, 2009). The time integration used intervals equal to $\mathrm{1}/\mathrm{30}$ of the day length (about 20 to 30 min depending on season). The depth integration used intervals equal to $\mathrm{1}/\mathrm{50}$ of Z_p and goes down to whichever is shallower between Z_p and Z_b. The spectral integration is performed over the visible range (400 to 700 nm) with a 5 nm resolution.

The spectral irradiance at a given depth z, $E(\mathit{\lambda },z,t)$, is calculated as (starting from just below the sea surface)

$$\begin{array}{}\text{(3)}& E(\mathit{\lambda },z,t)=E(\mathit{\lambda },z-\mathrm{d}z,t)e^{[-K_{\mathrm{d}}(\mathit{\lambda },z\left)\mathrm{d}z\right]},\end{array}$$

where the diffuse attenuation coefficient for downward irradiance, K_d(λ,z) (m⁻¹), is computed as a function of chlorophyll following Morel and Maritorena (2001):

$$\begin{array}{}\text{(4)}& K_{\mathrm{d}}(\mathit{\lambda },z)=K_{\mathrm{w}}(\mathit{\lambda },z)+\mathit{\chi }\left(\mathit{\lambda }\right)\text{chl}(z{)}^{e\left(\mathit{\lambda }\right)}.\end{array}$$

Details about how values are assigned to the parameters a^∗ and Φ, their dependence on temperature, and other features of this net primary production (NPP) model are to be found in Morel (1991) and Morel et al. (1996).

Daily PP calculation was performed every 7 d (starting at day 4), using the averaged chl for a 7 d window from day −3 to day +3 and for one grid point out of three. Therefore, the PP model was run with 8 d resolution and for one pixel out of three pixels. Later, daily PP data were interpolated to a 1 d by 1 km grid. Monthly PP data and their anomalies were derived from the abovementioned dataset.

The model was operated both for clear-sky conditions and for the actual MODIS PAR values, in which case a reduction of the clear-sky irradiance is uniformly applied across the entire day, as being the ratio of the daily satellite to clear-sky PAR values. Chl is assumed to be uniformly distributed with depth and equal to the satellite-derived value. This simplification was considered more appropriate for the generally shallow and well-mixed waters of coastal areas than the use of global parameterization of the shape of the vertical profile as a function of the surface chl value (e.g., Morel and Berthon, 1989; Uitz et al., 2006), whose validity outside of open ocean waters is not established.

From PP estimates, new (PP_new) and regenerated (PP_reg) production were calculated using the ratio of export production (PP_exp) to total production (PP) (i.e., ef-ratio) (Laws et al., 2000, 2011). Indeed, assuming a steady state, the export production must equal the new production fueled by new nutrients brought to the surface layers. The ef-ratio as a function of satellite-derived temperature and production can be obtained from the empirical relationship obtained by Laws et al. (2011):

$$\begin{array}{}\text{(5)}& {\displaystyle }\mathit{ef}=& {\displaystyle }{\displaystyle \frac{(\mathrm{0.5857}-\mathrm{0.0165}T)\text{PP}}{(\mathrm{51.7}+\text{PP})}},\text{(6)}& {\displaystyle }{\text{PP}}_{\text{exp}}=& {\displaystyle }{\text{PP}}_{\text{new}}=\text{PP}\times \mathit{ef},\text{(7)}& {\displaystyle }{\text{PP}}_{\text{reg}}=& {\displaystyle }\text{PP}-{\text{PP}}_{\text{new}},\end{array}$$

where T is temperature in degrees Celsius (^∘C) and PP is the daily production (mg C m⁻²).

As shown in Table 1, we report PP as vertically integrated values (PP, in g C m⁻²), spatially integrated estimates for certain basins or regions (ΣPP, in Gt C), or mean volumetric values (PP_VOL, in g C m⁻³). The coefficient of variation (CV) has been estimated for PP as the ratio of the standard deviation to the mean. While some authors include the Adriatic in the eastern basin (e.g., Bosc et al., 2004), we treated this region separately because its peculiarities (i.e., bathymetry, influence of rivers, eutrophic character) differentiate it from the rest of the Mediterranean Sea (Cushman-Roisin et al., 2001). Most of the Adriatic Sea has a shallow (<200 m) bathymetry, and it collects some 30 % of the freshwater flowing into the Mediterranean, acting as a dilution basin for the nutrients discharged by the Po and other Adriatic rivers and becoming one of the most human-impacted regions of the Mediterranean Sea (Ludwig et al., 2009; Micheli et al., 2013; Raicich et al., 2013a).

Table 1 Primary production acronyms used in this study, their units, and their definitions.

2.3 Coastal regionalization

We used a two-step classification procedure to define coastal regions along the Mediterranean based on their temporal PP patterns. First, nine regions (R1 to R9) were identified using a classification technique based on an unsupervised learning neural network (self-organizing maps or SOMs; Kohonen, 1982, 2001). Then, 18 alongshore marine ecoregions were obtained considering the most relevant cross-shore limits of the SOM-derived regions (Z1 to Z18).

SOM is an unsupervised neural network method that reduces the high-dimensional feature space of the input data to a lower-dimensional network of units called neurons. SOM is especially suited to extract patterns in large datasets of satellite data (Ben Mustapha et al., 2014; Charantonis et al., 2015; Farikou et al., 2015). Unlike other classification methods, like k-means, SOM tends to preserve data topology (i.e., preserves neighboring regions), and, therefore, it is particularly suited for pattern recognition (Liu and Weisberg, 2005). It allows adequate classification of areas with high spatial complexity and strong gradients. Similar neurons are mapped close together on the network, facilitating the visualization of patterns and a topological ordination of the classified areas, and the relative distance among neurons is obtained as a result of the analysis.

For typical satellite imagery, SOM can be applied to both space and time domains. Here, we have addressed the analysis in the time domain of the datasets, which allows regionalization of the studied area on the basis of similitudes in the time variation in PP. We chose a map size of (3×3), with nine neurons (for further details, see Basterretxea et al., 2018). We used a hexagonal map lattice in order to have equidistant neighbors and to avoid introducing anisotropy artifacts. For the algorithm initialization, we opted for linear mode, batch training algorithm, and ep, or Epanechnikov function, type neighborhood function since this parameter configuration produces the lower quantitative and topological error and computational cost (Liu et al., 2006). This methodology was applied in Basterretxea et al. (2018) for satellite-derived chlorophyll time series. For further details of the methodology used, see Hernández-Carrasco and Orfila (2018). These SOM computations were performed using the MATLAB toolbox of SOM v.2.0 (Vesanto et al., 2000a, b) provided by the Helsinki University of Technology (http://www.cis.hut.fi/somtoolbox/, last access: 12  April 2019).

2.4 Climate data

To identify possible drivers of long-term PP variability, we searched for correlations with two climate indices, the North Atlantic Oscillation (NAO) index and the Mediterranean Oscillation Index (MOI). The corresponding data were downloaded from the Climate Research Unit at the University of East Anglia (https://crudata.uea.ac.uk/cru/data/, last access: 21 August 2018) in monthly resolution. Climate indices are defined either as anomalies of a climate variable, using the difference between two geographical points, or as principal components (Hurrell, 1995).

NAO is the central mode of climate variability of the Northern Hemisphere atmosphere. It is based on the pressure difference between the middle of the North Atlantic Ocean and Iceland, which affects winter conditions in the North Hemisphere (Hurrell and Van Loon, 1997; Marshall et al., 2001). Positive NAO results in a relatively dry winter in the Mediterranean but a warmer and wetter winter in northern Europe, and vice versa. Because of its influence on precipitation, Mediterranean river inflows are generally anti-correlated with the NAO (Trigo et al., 2006).

MOI is the most widely used teleconnection index for the Mediterranean basin. It reflects differences in temperature, precipitation, circulation, evaporation, and other parameters between the eastern and western basin. There are different versions depending on the points of reference (Criado-Aldeanueva and Soto-Navarro, 2013). We used the version obtained as the normalized pressure difference between Gibraltar and Israel (Palutikof, 2003). Positive MOI phases are associated with increased atmospheric pressure over the Mediterranean Sea that promotes a shift of the wind trajectories toward lower latitudes leading to milder winters (Criado-Aldeanueva and Soto-Navarro, 2013). Under these conditions, reduced precipitation is observed in the southeastern Mediterranean region (Törnros, 2013). With some regional differences, NAO and MOI express relatively similar climate patterns over the Mediterranean Sea. They are highly positively correlated in winter and weakly but still significantly correlated in summer (Efthymiadis et al., 2011; Martínez-Asensio et al., 2014).

2.5 Statistical analyses

Linear temporal trends in the PP series were calculated using Theil–Sen slope adjustment (Sen, 1968) of the residuals of the de-seasonalized series. Time series have been de-seasonalized by removing the 8 d climatological mean for the original time series. Time series with >12 % of missing values were excluded from the analysis, and the analysis was calculated pixel by pixel. A Mann–Kendall statistics test is applied to each pixel, and only pixels with a trend statistically significant at the 95 % level were considered (Salmi et al., 2002). The use of this non-parametrical test is suitable for non-normally distributed data and has been previously used in the trend examination of remote-sensing chl time series (Colella et al., 2016; Kahru et al., 2011; Salgado-Hernanz et al., 2019).

3 Results

3.1 Coastal primary production

Annual primary production in coastal waters of the Mediterranean Sea (ΣPP_Coast) is estimated to be 0.041±0.004 Gt C, which represents some 12 % of total carbon fixation in the Mediterranean Sea (see Tables 2 and 3). Approximately, 80 % of this ΣPP_Coast is sustained by recycling processes, and, the rest, PP_new, is exported to the seafloor or to nearby areas. Although average surface chl concentration is 3-fold higher in coastal areas (0.3 mg m⁻³) than in open areas (0.11 mg m⁻³), the annual carbon fixation per surface area (PP_annual) over the shelf is, on average, 26 % lower than in the open ocean (100±91 and 136±40 g C m⁻², respectively; see Tables 2 and 3). The lower PP_annual values on the coast are due to depth integration quite often being stopped at shallow depths, i.e., Z_b<Z_p.

Figure 2 reveals the differences in PP_annual between the more productive shelf waters in the western basin and those in the eastern basin (90±39 and 73±86 g C m⁻², p<0.001). The Adriatic shelf is by far the most productive, particularly if mean values are considered (124±76 g C m⁻², Table 2). Annual carbon fixation is 97 % higher in the eastern (ΣPP${}_{\text{east}}=\mathrm{0.021}\pm \mathrm{0.002}$ Gt C) than in the western shelf (ΣPP${}_{\text{west}}=\mathrm{0.011}\pm \mathrm{0.001}$ Gt C), which is due to greater eastern surface area (about twice the western surface area; Table 2). PP_annual varies spatially from 90 to 250 g C m⁻² in the western shelf and from 50 to 400 g C m⁻² in the eastern basin. Contrastingly, PP_annual exceeds 100 g C m⁻² in the Adriatic basin, reaching values above 400 g C m⁻² on the northwestern Adriatic coast (Fig. 2b). The most productive coastal regions (>150 g C m⁻²) are located in regions enriched by major river outflows. Indeed, the highest values of the coefficient of variation in primary production (CV_PP) are observed at the mouth of the Nile, Ebro, Rhône, Tiber, Po, Neretva, or Nestos–Evros rivers. Nevertheless, in some coastal regions of the eastern basin like the Gulf of Gabès primary production is also outstandingly high (>300 g C m⁻²). Along the western North African coast, PP_annual also displays values >150 g C m⁻²; however, since the shelf is narrow, its contribution to ΣPP_Coast is marginal (Figs. 1 and 2). The annual volumetric productivity follows a similar pattern than the annual integrated production, with most values varying between 1–3 g C m⁻³ but reaching up to >10 g C m⁻³ in the most productive coastal regions of the Adriatic Sea and the Gulf of Gabès (Fig. 3a).

Table 2 Surface area and the corresponding basin percentage to surface area, (chl) ± standard deviation (SD), ΣPP, the corresponding basin percentage to ΣPP_Coast (% ΣPP_Coast), PP_annual median ± SD (PP_annual mean), and PP_VOL median ± SD (PP_VOL mean) for the Mediterranean Sea, open ocean waters, and coastal waters during the period 2002–2016. For ΣPP and PP_annual, Mediterranean Sea values were obtained summing open ocean water values and coastal water values.

^a Mean surface chl values obtained by averaging the 8 d and 4 km resolution of surface satellite chl values obtained from CMEMS (Salgado-Hernanz et al., 2019). ^b PP estimated by averaging published satellite data shown in Table 3. ^c ΣPP estimated adding coastal water data from this study to open ocean water data obtained from Table 3.

Table 3 Compilation of published values of PP_annual and ΣPP for the different Mediterranean basins.

^a The estimates from Antoine et al. (1995) and Morel and André (1991) have been corrected by a factor of 1.25 as recommended by Morel et al. (1996). ^b From Colella et al. (2003), who estimated it using f-ratios (the ratio between new and total PP) obtained from Boldrin et al. (2002). ^c Daily PP ($\mathrm{mg}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{d}}^{-\mathrm{1}}$) converted to annual PP ($\mathrm{mg}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$) multiplied by 365. ^d Conversion to carbon unit using photosynthetic quotient PQ=1. ^∗ Data included in the Data Published for Earth and Environmental Science PANGAEA data collection.

Figure 2 Mean distribution of (a) chlorophyll (chl_mean, in mg m⁻³), (b) annual PP (PP_annual, in g C m⁻²), and (c) coefficient of variation in PP values (CV_PP).

Figure 3 Mean distribution of (a) volumetric PP (PP_VOL, in g C m⁻³) and (b) coefficient of variation in volumetric PP values (CV_PP).

3.2 Long-term variability and trends

As shown in Fig. 4, variability in annual PP is dominated by short-scale variations (i.e., sub-decadal). The interannual variability is indicated by the low-frequency signal of the monthly mean anomalies. The filtered low-frequency signal to the anomalies of ΣPP and chl has been calculated using the MATLAB smooth function and applying the sgolay filter, which uses the Savitzky–Golay method, with a polynomial span degree of 17. This degree was therefore filtering about 8 months before and after every time step (about 1.5 years), then showing interannual variability. ΣPP_Coastal exhibits moderate interannual variability (up to 25 %) whereas basin-scale interannual variations range from 26 % in the Adriatic basin up to 28 % in the western basin and 29 % in the eastern basin. This value of interannual variability was calculated subtracting the year with the minimum annual PP from the year with the maximum annual PP and then dividing this value by the mean annual PP. Positive anomalies in ΣPP_Coastal extended between 2004 and 2011 (mean ΣPP${}_{\text{Coastal}}=\mathrm{0.044}\pm \mathrm{0.001}$ Gt C yr⁻¹ for the whole basin between those years, Fig. 4a and Appendix Fig. A1a). Conversely, the year 2012 was particularly unproductive in all three basins (0.037 for the whole basin, 0.010 for the western basin, 0.019 for the eastern basin, and 0.009 Gt C yr⁻¹ for the Adriatic basin). This negative anomaly marked the beginning of a less productive period, particularly noticeable in the eastern basin (Figs. 4, A1 and A2). PP_annual (g C m⁻²) and ΣPP (Gt C) in every year are shown in Figs. A1 and A2.

Figure 4 PP variability and trends for coastal waters in (a) the whole Mediterranean Sea, (b) western basin, (c) eastern basin, and (d) Adriatic coast. Black solid lines indicate the original monthly ΣPP_Coastal anomalies, and the filtered low-frequency signal is overlaid in blue (16-month width sgolay filter). Green solid lines indicate the filtered low-frequency signal for chl anomalies (mg m⁻³). The red line indicates the PP trend during the analyzed period (2002–2016), and the gray band indicates the year 2012.

As revealed by Fig. 5a, some regionally coherent patches of significant trend in PP are observed along the coast. Most of these regions presented declining PP trends, particularly along the North African coast where SST is increasing at a higher rate (Fig. 5b). Typical PP trend magnitudes observed along the Spanish Mediterranean and the North African coast from the Gulf of Gabès range from −0.05 to +0.05 g C m⁻² per decade. Some positive PP trends, exceeding +0.1 g C m⁻² per decade, can be determined in some coastal regions north of the Adriatic Sea.

Figure 5 Trends in primary production and sea surface temperature. Values correspond to the change per decade. (a) Theil–Sen trend in pelagic primary production (PP trend) estimated from daily values for the 2002–2016 period. (b) Trend in sea surface temperature (SST trend). Only significant trends (p<0.05) are shown.

A significant negative correlation was observed between monthly coastal anomalies of ΣPP and SST at the whole Mediterranean basin scale ($r=-\mathrm{0.66}$, p<0.001; Fig. 6a, Table A1), revealing a decrease in phytoplankton biomass as the sea warms up. In addition, we observed evidence of an inverse relationship between PP variability and the phase of the climate indices NAO and MOI ($r=-\mathrm{0.46}$, p<0.001, and $r=-\mathrm{0.24}$, p<0.001, respectively; Fig. 6b and c, Table A1). A similar relationship is observed between monthly coastal chl anomalies with SST ($r=-\mathrm{0.64}$, p<0.001) and climate indices NAO and MOI ($r=-\mathrm{0.60}$, p<0.001, and $r=-\mathrm{0.35}$, p<0.001, respectively; see Appendix Fig. A3) at the whole Mediterranean basin scale.

Additionally, the response of PP to climate variations varied seasonally. Indeed, NAO influenced coastal PP in summer, both in the western and in the eastern basin (r=0.25, p=0.06; r=0.22, p=0.08, respectively) and MOI variations were better correlated with PP global variations in spring (r=0.28, p=0.04), showing a higher impact in the Adriatic basin (r=0.37, p=0.02). No significant correlation was found during the winter or fall season for any index (see Table A1 in Appendix). Correlation analyses were performed using the Pearson product moment correlation. Differences between means were tested using the Kolmogorov–Smirnov test (Massey et al., 1951).

Figure 6 Relationship between coastal pelagic primary production (ΣPP anomalies, blue lines) and (a) SST anomalies, (b) NAO index, and (c) MOI index (red lines).

3.3 Coastal regionalization

The nine characteristic temporal PP patterns, their corresponding spatial distribution obtained from the SOM analysis, and the 18 zones in which the coastal region was classified are shown in Fig. 7. Generally, wider shelves present higher spatial complexity manifested as a larger number of SOM patterns. About 78 % of the shelf waters include R1, R2, R3, and R4 patterns. In particular, R1 and R2 characterize areas of low production, with scarce seasonality (PP${}_{\text{annual}}=\mathrm{44}\pm \mathrm{17}$ and 69±22 g C m⁻², respectively; Table 4) typically occurring in the southern and eastern Mediterranean (12 % and 29 % of the total surface area). They are representative of the productivity patterns in vast shelf regions in the Gulf of Gabès and Sirte and in the central Aegean (Z17, Z16, Z14 and Z11). R3 and R4 correspond to higher production and a wider range of variation (PP${}_{\text{annual}}=\mathrm{90}\pm \mathrm{32}$ and 98±39 g C m⁻²; Table 4 and Fig. 7). While R3 (18.6 % of total coastal surface) is frequent in shelf regions of the western basin (Z2, Z3, Z5, and Z7), R4 extends over the deepest areas of the Adriatic Sea (Z9 and Z10). Both patterns are representative of variations observed in 36.9 % of the coastal waters.

Figure 7 Regionalization of the coastal waters in the Mediterranean Sea based on their temporal patterns of pelagic primary production. (a) Characteristic temporal patterns of PP obtained from SOM classification (R1 to R9) and (b) coastal regions defined from alongshore variations in the SOM – regions (Z1 to Z18).

Table 4 Extension and primary production for each of the SOM-defined regions (R1 and R9). Mean annual PP is estimated by averaging mean daily PP and then multiplying it by the number of days of the year, i.e., 365.

As shown in Table 4, areas of low production with seasonal patterns R1 to R4 contribute to more than 62 % of total pelagic carbon fixation in Mediterranean coastal areas. In contrast, systems of higher production (PP_annual>280 $\mathrm{g}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{d}}^{-\mathrm{1}}$) barely contribute to 17 % of total production. These regions of enhanced production are generally constrained to regions of freshwater influence (ROFIs; Simpson, 1997) where terrestrial nutrients fuel coastal production. Indeed, R8 is almost exclusively restricted to the river mouths, and it presents elevated PP_annual values (1.29±0.50 $\mathrm{g}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{d}}^{-\mathrm{1}}$) and a wide range of variation of 0.67 to 2.14 $\mathrm{g}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{d}}^{-\mathrm{1}}$. An exception is the R7 pattern, which is exclusively located in the shallowest inner shelf of the Gulf of Gabès, and it is bounded by R5, a transition region between the inner and outer shelf. Unlike the other regions, where PP peaks in late winter–spring, maximum PP in R7 occurs in fall. Finally, R5 and R6 patterns correspond to transition regions accounting for 20.4 % of the total production and 14.8 % of the Mediterranean coast. While R5 mainly occurs in deltas, R6 is characteristic of the western Mediterranean shelf, including the North African coast (0.36–0.08 $\mathrm{g}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{d}}^{-\mathrm{1}}$; Fig. 7).

The SOM-based regionalization reveals two groups of coastal waters: those with low cross-shore variability and including only one or two SOM regions (i.e., Z1, Z13, Z16, and Z18) and those with strong cross-shore gradients including several SOM regions (i.e., Z4, Z9, Z12, Z15). The first pattern is typically observed in narrow continental shelf areas with low influence of river inputs whereas the second group is found in regions with a wider continental shelf such as ROFIs (the Rhône delta, the north and western coastline of the Adriatic Sea, and the Nile Delta) and in the Gulf of Gabès. The western Adriatic (Z9) and the Gulf of Gabès (Z17) are the largest contributors to ΣPP_Coast, together contributing 35.9 % of shelf production in the Mediterranean Sea, but, in the case of Z17, it is mainly due to its large extension (Table 5). PP is also high in the northern Alboran Sea (Z1), Nile delta (Z15), the western Adriatic (Z9), and Gulf of Lion (Z4; Table 5). With the exception of Z1, influenced by the entrance of waters from the Atlantic Ocean and by local coastal upwelling, these zones receive important riverine fluxes (Q).

4 Discussion

4.1 Coastal primary production

In this study we focused on the contribution of coastal waters to the overall pelagic PP in the Mediterranean Sea. While the mean coastal values for the Mediterranean (100±91 g C m⁻²) are somewhat lower than the mean values over the continental shelves of the world ocean (160±40 g C m⁻²; Smith and Hollibaugh, 1993), the impact of coastal pelagic PP on total basin production (12 %) is in the high range of the estimations for other seas (Muller-Karger et al., 2005). This estimation is subject to the uncertainties inherent to using satellite ocean color, which is limited to the upper ocean (down to 20 m at best in clear waters) and has poor performance in some areas (i.e., Case-2 waters). It nevertheless provides an assessment of net rates of carbon fixation in coastal areas that is consistent with global estimations of the contribution of coastal areas to oceanic production (Gattuso et al., 1998; Ducklow and McCallister, 2004; Muller-Karger et al., 2005). Bias in coastal chl estimations is mainly due to the presence of non-phytoplankton components such as colored dissolved organic matter (CDOM) or other terrestrial substances (Morel et al., 2006). These compounds originate from coastal erosion, resuspension in shallow areas, river inputs, or anthropogenic effluents. Likewise, they affect the propagation of photosynthetic radiation through the water column (Morel, 1991). However, the possible uncertainties and biases caused by chl estimation through satellite data might very weakly alter our estimation of coastal PP. Indeed, previous studies agreed that Case-1 waters are largely predominant in the coastal Mediterranean regions whereas Case-2 waters are reduced to less than 5 % of the whole basin (Antoine and André, 1995; Bosc et al., 2004; Bricaud et al., 2002). This constitutes some 23 % of the coastal waters with prevalence in the north Adriatic Sea, Gulf of Gabès, and around the Nile delta. In particular, they are confined to the north Adriatic Sea, Gulf of Gabès, and around the Nile delta where our PP estimations may present larger uncertainties (Antoine and André, 1995). However, PP_annual values off the Nile River delta, >100 g C m⁻² estimated here, are only slightly higher than those reported by Antoine et al. (1995) (80–100 g C m⁻²). The highest values have been reported for this region (>300 g C m⁻²), but, as shown in Fig. 2, they are restricted to a narrow coastal band. In the case of the Adriatic Sea, Umani (1996) reported values of PP from 50 to 200 $\mathrm{g}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$, while Zoppini et al. (1995) estimated PP rates from 210 to 260 $\mathrm{g}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$ in the northern coastal areas. Our estimations range between 100 and >350 g C m⁻² (with mean values of 123±106 g C m⁻²).

Table 5 Surface, river discharge flow (Q), annual mean PP (PP_annual), annual integrated PP (ΣPP), and its contribution with respect to the total coastal Mediterranean Sea PP for each of the 18 alongshore zones characterized in the Mediterranean Sea. Mean and standard deviation (SD) are calculated from 15-year averages (2002–2016). River discharge flow was extracted from Tockner et al. (2008).

Contrary to what we would have expected, we observed that the eastern basin contributes more than the western basin to overall coastal production (51 % and 25 %, respectively; Table 2). Its great extension (twice higher than the western basin) and the increased productivity in regions like Gabès, the Nile, and the northern Aegean Sea may explain greater coastal PP in the eastern basin than in the western basin. Additionally, due to the lack of large shallow and productive areas in the western basin, we observed few volumetric PP values above 30 g C m⁻³ in the western Mediterranean whereas high volumetric PP is more frequent in the Adriatic Sea and in the eastern Mediterranean, in shallow waters of the Gulf of Gabès, and in the Nile delta. Furthermore, carbon fixation from the Adriatic Sea represents 24 % of total coastal production, which is significant considering the area of this sea (19 % of Mediterranean coastal waters). The relevance of the contribution of the Adriatic Sea in overall coastal PP lies in two main characteristics: (1) coastal waters (<200 m) constitute a large part of the Adriatic Sea, and (2) about one-third of the river discharge in the Mediterranean is concentrated in the Adriatic Sea (see Table 5). Indeed, patterns in the northern Adriatic Sea reflect a variation in the drivers of PP with respect to other regions. For example, while internal processes (i.e., vertical diffusion and mixing) and atmospheric deposition drive PP in most coastal waters, production in the north Adriatic is mainly driven by fluvial sources of carbon and regeneration through bacterial pathways (Powley et al., 2016; Rodellas et al., 2015; Umani et al., 2007). Moreover, distinctive dynamics in this sea is driven by the influence of river outflows on stratification and general circulation patterns (Djakovac et al., 2012; Giani et al., 2012).

Regionally, the Alboran Sea, the Gulf of Lion, the Adriatic Sea, the Aegean Sea, the Nile delta, and particularly the Gulf of Gabès (up to 23 % of total coastal production) are the most productive zones. While some of these productive (i.e., PP_VOL>3 g C m⁻³) coastal regions are located in areas affected by river outflows, this type of production rapidly decreases in the offshore direction (see Fig. 3). Exceptions are the Rhône and the Po and less so the Nile, the influence of which extends far along the shelf. Other processes like mesoscale circulation are also important in some of these regions (i.e., Gulf of Lion), as demonstrated by Macias et al. (2017). Additionally, the exchanges with the more productive Atlantic and Black Sea waters in Alboran and in the northern Aegean and local enrichment processes in the Gulf of Gabès are major contributors to coastal productivity. Indeed, the Gulf of Gabès constitutes an anomaly in the eastern Mediterranean basin. Its shallowness (<50 m depth at 110 km off the coast), unique tidal range (maximum >2 m) promoting vertical mixing, and the lack of summer nutrient exhaustion undoubtedly contribute to its high productivity (Béjaoui et al., 2019).

Notwithstanding the importance of land inputs in the production of coastal Mediterranean waters, from ef-ratios, we estimated that on average only 22 %±20 % of the production is new, and the rest is sustained via regenerated sources (Fig. 8a–c). This PP_new value (Fig. 8b) is comparable to the mean organic carbon that sinks to the sea floor (28 %) estimated from Muller-Karger et al. (2005) and Pace et al. (1987) but higher than PP_new estimations provided by Vidussi et al. (2001) for oceanic waters in the eastern basin (15 % of total production). Contrarily to Vidussi et al. (2001), who estimated PP_new in the eastern basin, the coastal PP_new average here includes both eastern and western basins of the Mediterranean, but also the highly productive areas in the northern Adriatic. This could explain why the PP_new observed here is higher than the one observed for oceanic waters in the eastern basins. Additionally, high ef-ratios (>0.3) are observed in our case in the areas where nutrient inputs from the Atlantic and river effluents significantly enhance PP_new (Fig. 8a). Furthermore, ef-ratios present significant seasonality, varying between 0.26±0.04 in the most productive winter–spring season and 0.15±0.02 in summer, when the water column is strongly stratified and the food web shifts to a more recycling-dominated system.

Figure 8 (a) ef-ratio in coastal waters (<200 m) of the Mediterranean Sea and estimated values of (b) new (PP_new) and (c) regenerated production (PP_reg). Mean values for the period 2002–2016.

4.2 Long-term variability and trends

Available satellite ocean color data span about 20 years, so that temporal trends derived from their analysis highly depend on decadal variability (Henson et al., 2010). Despite these limitations, satellite observations of ocean color over the past 2 decades suggest a relationship between warming and reduced productivity in permanently stratified areas (Behrenfeld et al., 2006). Since clear tendencies of warming are observed in the Mediterranean Sea (Nykjaer, 2009; Pastor et al., 2018), intensification of stratification would decrease nutrient supply to phytoplankton and, thus, decrease PP (Behrenfeld et al., 2006; Stambler, 2014).

Barale et al. (2008), using chl anomalies derived from SeaWiFS data, observed a general decrease in chl biomass in the Mediterranean Sea over the period 1998–2003. However, some coastal areas in their study displayed the opposite tendency. Macias et al. (2015) anticipated no future global changes of integrated PP in the Mediterranean Sea from modeling results. They predicted a tendency to oligotrophication in the western basin and increase in the productivity of the eastern basin. Our study reveals that ΣPP_Coastal in the Mediterranean Sea varies nonlinearly, and a reduction of carbon fixation rates has been observed since 2012 (Fig. 4a). Overall negative coastal PP trends are reported here in both the western and in the eastern basin (−10.70 and −25.39 TC per decade; Fig. 4b and c). A spatial analysis of the long-term decadal variability reveals weak but spatially coherent and significant tendencies (p<0.05; Fig. 5). In particular, PP declines along the coasts of Spain and Africa. Conversely, trends in some areas of the Adriatic Sea are markedly positive (>0.1 g C m⁻² per decade), mainly in the proximity of the Po River. While negative tendencies seem to fit with the assumed model of PP limitation associated with increasing temperatures, the origin of the positive trend in the Adriatic basin is more uncertain. A plausible explanation is the variation in the flux and loads of the northern Adriatic rivers. For example, Giani et al. (2012) observed an increase in the Po River flow with increasing phosphate and dissolved nitrogen concentrations in the Po delta and its surrounding shelf waters. Alternatively, the Bimodal Oscillating System (BiOS), i.e., the feedback mechanism between the Adriatic and Ionian (Civitarese et al., 2010) peaking between 2004 and 2006, could have affected mass and nutrient exchanges between the Adriatic and the north Ionian Sea (Font et al., 2007; Schroeder et al., 2008; Šolić et al., 2008; Viličić et al., 2012).

Long-term decadal variations in the eastern and western basins are mostly coupled, suggesting that they share the same PP drivers at this basin scale (Fig. 4b and c). A major feature in the interannual pattern is a global decrease in production in 2012 that is extended to the following years in the eastern basin. Durrieu de Madron et al. (2013) reported peculiar atmospheric conditions in the Mediterranean Sea during 2012 that triggered a massive formation of dense water on the continental shelf and in the deep basin of the Gulf of Lion. A similar anomaly was described in the Adriatic shelf where unprecedented dense water generation was preconditioned by a dry and warm year, resulting in a significant reduction of coastal freshwaters and basin-wide salinity increase (Mihanović et al., 2013; Raicich et al., 2013b). Additionally, Pastor et al. (2018) observed an anomalous temperature increase in the Mediterranean Sea during summer 2012. From our analysis, we infer that this climate-related event had a strong influence on the global coastal PP of the Mediterranean Sea.

Several studies have reported influence of climate variations on the coast (Belgrano et al., 2008; Cloern et al., 2007; Tiselius et al., 2016). In agreement, we observed an influence of climate-scale variability on coastal productivity as suggested by the inverse correlations between ΣPP and SST and, more loosely, by NAO and MOI (Fig. 6). While these correlations emphasize the pre-eminent role of climate variability in the regulation of interannual to decadal-scale coastal productivity, the pathways through which this control of the atmosphere over coastal productivity is exerted are complex and may regionally differ (Grbec et al., 2009). Climate can influence phytoplankton growth by the direct effect of temperature on algal metabolism, by changes in basin-scale circulation (including exchanges with adjacent seas), by regulating nutrient supply through variations in the thermocline intensity, by changes in wind patterns affecting mixing and dust deposition pathways, or through changes in precipitation that have direct influence on wet deposition and on river runoff. These effects are modulated by changes in the biota and in the interaction between organisms (e.g., Molinero et al., 2005). The relative importance of climate-driven processes relative to other productivity-enhancing processes depends on regional characteristics and may be seasonally varying. For example, variations in dust deposition, which may sustain up to 50 % of new production in the Levantine basin (Kress and Herut, 2001; Herut et al., 2002), are expected to be more important on the eastern and southern Mediterranean coasts because of their proximity to the Saharan dust sources. Likewise, variations in cooling and vertical mixing are expected to be more effective during late winter when PP peaks and when diatoms dominate in the Mediterranean Sea (Lacroix and Nival, 1998; Marty, 2002; Marty and Chiavérini, 2010).

Our results reveal that, in contrast to other regions like the North Sea (Capuzzo et al., 2018) or the Arctic Ocean (Gregg et al., 2003), the coastal Mediterranean Sea did not globally display a marked decline in PP during the last decades. We suggest that in some coastal areas, a decrease in vertical nutrient supply though the thermocline may be compensated for by other nutrient sources. Variations in atmospheric deposition, groundwater and river outflows, and the influence of human activities through changes in landscape use and domestic wastewater management are important sources of nutrients in the ecosystem and thus act as major drivers of PP in these waters (e.g., Paerl et al., 1999; Powley et al., 2016). As a consequence of human activities, both terrestrial and coastal ecosystems have experienced progressive nutrient enrichment (Conley et al., 2009; Deegan et al., 2012). However, while this effect is evidenced in shallow nearshore waters, its influence in the ocean is estimated to be minimal (Wang et al., 2018). In the Mediterranean Sea, high coast population growth rates and concomitant food demand have resulted in a dramatic increase in water demand for irrigation farming and fertilizer use (Ryan, 2008). Indeed, while the freshwater discharge of Mediterranean rivers has significantly reduced during recent decades (∼20 %), the corresponding total nitrogen inputs to coastal seas are estimated to have increased by a factor of up to 5, fueling PP in river-influenced areas (Ludwig et al., 2009). While the importance of groundwater in the Mediterranean Sea could be comparable to that of rivers (Rodellas et al., 2015) and generalized nitrification of Mediterranean coastal aquifers is acknowledged (EEA, 1995; Zalidis et al., 2002), general trends in groundwater discharges remain largely unknown.

4.3 Coastal regionalization

Coastal regionalization reveals marked differences in coastal water PP in the Mediterranean Sea. Annual values range from 215±124 g C m⁻² in the north Alboran Sea (Z1) to 48±17 g C m⁻² along the coasts of Egypt and Libya (Z16). These values are globally lower than published data; yet, literature values in coastal waters are highly variable depending on methodology, depth, and/or sampling date. For example, García-Gorriz and Carr (2001) estimated annual PP of 300–900 g C m⁻² for Z1, but Morán and Estrada (2001) narrowed this range to mean values between 121 and 366 g C m⁻² depending on distance from the coast. Pugnetti et al. (2008) reported mean values of 150 g C m⁻² that are almost 2 times higher than our values at Z8. In the lower range, Sournia (1973) estimated 30–60 g C m⁻² in Z16, which is in accordance with our values for this zone. In this sense, despite the limitations inherent to satellite data, the present work provides estimations based on a long data record (14 years) and a homogeneous methodology.

The nine characteristic temporal patterns obtained from the SOM analysis (Fig. 7) reveal small differences in PP among the different regions. Most variations are due to changes in the magnitude of annual carbon fixation, although seasonality varies little. Exceptions are R7, R8, and R9 which represent the dynamics in coastal regions regulated by terrestrial inputs. Likewise, interannual variations are highly coherent among regions, following the basin-scale pattern shown in Fig. 4, including the remarkable decline in productivity during 2012. Exceptions are R1, representing the dynamics in the Gulf of Sirte and R7 in the Gulf of Gabès where a different interannual variability suggests alternative sources of PP variability in this region. Indeed, the Gulf of Gabès is a region displaying consistently high chl and PP in most studies (e.g., Bosc et al., 2004; Barale et al., 2008). Drira et al. (2008) reported high biomass and toxic dinoflagellate blooms in the inner shelf of the Gulf of Gabès where surface nitrate concentration often exceeded 1 µM. This enrichment is associated with degradation of the water quality attributed to industrial and urban activities (Hamza-Chaffai et al., 1997; Zairi and Rouis, 1999). However, even though these waters may suffer from eutrophication, satellite-borne data overestimate chl within these waters, as revealed by Katlane et al. (2011), who observed constant high turbidity and suspended matter of industrial origin affecting these waters but also reflection from the bottom affecting MODIS data. This suggests that general chl algorithms may be particularly inaccurate in this region.

The magnitude of coastal PP has often been related to both shelf width and magnitude of river discharge (Liu et al., 2010). Our data do not display a general relationship between shelf width and PP_annual (Table 5 and Fig. 9). Indeed, wide shelves with important river discharge flux from the Po, the Rhône, and the Nile rivers display high productivity (Z4, Z9, and Z15 >170 g C m⁻²) whereas production is low in the narrowest shelves like Z2, Z5, and Z16. However, PP_annual in some regions with important river inflows, like Z10, are significantly lower (89±37 g C m⁻²). In other regions like Z1 and Z8 PP is high despite the lack of important freshwater sources.

The role of river discharges depends both on Q and on nutrient loads. Most Mediterranean rivers have lost their natural flows, and their discharges to the sea are strongly regulated by dams and water abstractions. Consequently, their outflow to the Mediterranean Sea is highly uncoupled from weather and climate variability. For instance, some rivers flowing into the Adriatic and Ionian seas, like the Acheloos, Nestos, or Aliákmon, today present high to maximum discharge in July due to peak hydropower production (Skoulikidis et al., 2009).

Figure 9 Relation between primary production, shelf width, and river discharge flow (Q). Bubble colors indicate the PP_annual (g C m⁻²) for each of the 18 defined zones (see Fig. 7 and Table 5).

Nutrient and organic matter loads have increased globally during the last century (Beusen et al., 2016). Concentrations exported to nearby seas depend on the combined effects of lithology, urban effluents, industry, and agriculture in catchment basins that are often difficult to quantify. The land use and land cover (LULC) data collection provides indices of the threat of potential development for setting land and water quality policies. Rivers like the Rhône and Po with important influence on coastal productivity flow through extensive areas and therefore accumulate the impact from anthropogenic activities. Agricultural practices and urban effluents can strongly determine the concentration and molar ratios of the nutrients flowing into coastal waters. For example, despite the flux of the Nile River having been drastically reduced after the operation of the Aswan Dam (from 47 to 17 km³ yr⁻¹, Ludwig et al., 2009), the coastal region is still highly productive. A remarkable increase in the concentrations of nitrate derived from fertilizers and sewage is responsible for this sustained productivity (Turley, 1999; Nixon, 2003, 2004). Conversely, pollution pressures in the western Balkan basins are relatively low, and the Neretva (Z10), running through a karstic region in Croatia, displays low nutrient levels (Ludwig et al., 2009; Skoulikidis et al., 2009).

Finally, other oceanographic processes determine the productivity of coastal regions. In particular, Z1 and Z18 in the Alboran Sea are comparatively more productive than other areas. The influence of winds and circulation patterns favoring subsurface water upwelling and higher productivity in the northern Alboran Sea were described by García-Gorriz and Carr (2001). Also, localized patterns of relatively high primary production were found in persistent deep water density fronts resulting from the interaction of Modified Atlantic Water (MAW) and Mediterranean water by Lohrenz et al. (1988).

5 Conclusions

In summary, pelagic PP in coastal shelves of the Mediterranean Sea during the period 2002–2016 was estimated in this study for the first time using available satellite ocean color products. We estimated that 12 % of PP of the Mediterranean Sea is attributable to coastal pelagic production, and from that about 80 % of this carbon fixation is sustained by regenerated pathways. High PP spatial variations were observed among the different regions, as mainly driven by major river effluents, exchanges with nearby seas (i.e., Black Sea and the Atlantic Ocean), and local processes. Our study shows that some coastal areas are indeed highly productive (>400 g C m⁻²) and sustain a large percentage of overall coastal production. Indeed, their temporal variability could be of paramount importance to understand variations in higher trophic levels (e.g., Piroddi et al., 2017). Despite that temporal variability being dominated by interannual and sub-decadal variations, our analysis reveals a weak negative global PP trend in the Mediterranean Sea related to climate-driven patterns (i.e., temperature increase). Nevertheless, long-term effects can be regionally variable (i.e., PP trends in the Adriatic Sea are positive), and variations in fluvial nutrient inputs, together with other processes such as ocean warming in coastal regions, including heat waves, deserve a closer look as longer ocean color databases become available. Finally, we identify 18 along-shelf zones based on their temporal PP patterns. Two main PP groups were observed: zones with strong cross-shore gradients, typically found in wider estuarine regions, and homogeneous zones within narrow continental shelf areas. These two types of coastal waters clearly characterize the coastal area of a sea where coastal waters are otherwise strongly influenced by ocean conditions.

Appendix A

Figure A1 Annual coastal variations in ΣPP (Gt C yr⁻¹) for (a) the whole Mediterranean Sea, (b) the western basin, (c) the eastern basin, and (d) the Adriatic basin. Dashed black lines indicate the mean ΣPP_Coastal for the period 2002–2016 at each study region. Red and blue bars indicate years above and below average values, respectively.

Figure A2 Annual mean coastal variations, PP_annual ($\mathrm{g}\mathrm{C}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{yr}}^{-\mathrm{1}}$), for (a) the whole Mediterranean Sea, (b) the western basin, (c) the eastern basin, and (d) the Adriatic basin. Dashed black lines indicate the coastal mean PP_annual at each study region. Red and blue bars indicate years above and below average values, respectively.

Figure A3 Relationship between coastal pelagic chlorophyll (chl anomalies, green lines) and (a) SST anomalies, (b) NAO index, and (c) MOI index (red lines).

Table A1 Correlations between yearly and seasonal total carbon fixation anomalies (ΣPP) for the coastal Mediterranean waters and its sub-basins with the climatic indices NAO and MOI.