Portal de Periódicos da CAPES

Overfishing of top predators eroded the resilience of the Black Sea system regardless of the climate and anthropogenic conditions

2010; Wiley; Volume: 17; Issue: 3 Linguagem: Inglês

10.1111/j.1365-2486.2010.02331.x

1365-2486

Marcos Llope, Georgi Daskalov, Tristan Rouyer, Vesselina Mihneva, Kung‐Sik Chan, А. Н. Гришин, Nils Chr. Stenseth,

Food Industry and Aquatic Biology

Global Change BiologyVolume 17, Issue 3 p. 1251-1265 Open Access Overfishing of top predators eroded the resilience of the Black Sea system regardless of the climate and anthropogenic conditions MARCOS LLOPE, MARCOS LLOPE Centre for Ecological and Evolutionary Synthesis (CEES), Department of Biology, University of Oslo, PO Box 1066 Blindern, NO-0316 Oslo, Norway 1Present address: M. Llope, Instituto Español de Oceanografía (IEO), Centro Oceanográfico de Cádiz, Puerto Pesquero, Muelle de Levante s/n, PO Box 2609, E-11006 Cádiz, Andalucía, Spain.Search for more papers by this authorGEORGI M. DASKALOV, GEORGI M. DASKALOV CEFAS Lowestoft Laboratory, Pakefield Road, Lowestoft, Suffolk NR33 0HT, UK 2Present address: Laboratory of Marine Ecology, Institute of Biodiversity and Ecosystem Research (Bulgarian Academy of Sciences), 18 Makedonia Str 9002, Varna, Bulgaria.Search for more papers by this authorTRISTAN A. ROUYER, TRISTAN A. ROUYER Centre for Ecological and Evolutionary Synthesis (CEES), Department of Biology, University of Oslo, PO Box 1066 Blindern, NO-0316 Oslo, Norway Centre de Recherche Halieutique Méditerranéenne et Tropicale, Institut Français de Recherche pour l'Exploitation de la Mer (IFREMER), Avenue Jean Monnet, BP 171, 34203 Sète cedex, FranceSearch for more papers by this authorVESSELINA MIHNEVA, VESSELINA MIHNEVA Institute of Fisheries and Aquaculture, Varna, PO Box 72, Varna 9000, BulgariaSearch for more papers by this authorKUNG-SIK CHAN, KUNG-SIK CHAN Department of Statistics and Actuarial Science, University of Iowa, 263 Schaeffer Hall, Iowa City, IA 52242, USASearch for more papers by this authorALEXANDER N. GRISHIN, ALEXANDER N. GRISHIN Southern Scientific Research Institute of Marine Fisheries and Oceanography (YugNIRO), 2, Sverdlov Street, 98300 Kerch, Crimea, UkraineSearch for more papers by this authorNILS CHR. STENSETH, NILS CHR. STENSETH Centre for Ecological and Evolutionary Synthesis (CEES), Department of Biology, University of Oslo, PO Box 1066 Blindern, NO-0316 Oslo, Norway Flødevigen Marine Research Station, Institute of Marine Research (IMR), NO-4817 His, NorwaySearch for more papers by this author MARCOS LLOPE, MARCOS LLOPE Centre for Ecological and Evolutionary Synthesis (CEES), Department of Biology, University of Oslo, PO Box 1066 Blindern, NO-0316 Oslo, Norway 1Present address: M. Llope, Instituto Español de Oceanografía (IEO), Centro Oceanográfico de Cádiz, Puerto Pesquero, Muelle de Levante s/n, PO Box 2609, E-11006 Cádiz, Andalucía, Spain.Search for more papers by this authorGEORGI M. DASKALOV, GEORGI M. DASKALOV CEFAS Lowestoft Laboratory, Pakefield Road, Lowestoft, Suffolk NR33 0HT, UK 2Present address: Laboratory of Marine Ecology, Institute of Biodiversity and Ecosystem Research (Bulgarian Academy of Sciences), 18 Makedonia Str 9002, Varna, Bulgaria.Search for more papers by this authorTRISTAN A. ROUYER, TRISTAN A. ROUYER Centre for Ecological and Evolutionary Synthesis (CEES), Department of Biology, University of Oslo, PO Box 1066 Blindern, NO-0316 Oslo, Norway Centre de Recherche Halieutique Méditerranéenne et Tropicale, Institut Français de Recherche pour l'Exploitation de la Mer (IFREMER), Avenue Jean Monnet, BP 171, 34203 Sète cedex, FranceSearch for more papers by this authorVESSELINA MIHNEVA, VESSELINA MIHNEVA Institute of Fisheries and Aquaculture, Varna, PO Box 72, Varna 9000, BulgariaSearch for more papers by this authorKUNG-SIK CHAN, KUNG-SIK CHAN Department of Statistics and Actuarial Science, University of Iowa, 263 Schaeffer Hall, Iowa City, IA 52242, USASearch for more papers by this authorALEXANDER N. GRISHIN, ALEXANDER N. GRISHIN Southern Scientific Research Institute of Marine Fisheries and Oceanography (YugNIRO), 2, Sverdlov Street, 98300 Kerch, Crimea, UkraineSearch for more papers by this authorNILS CHR. STENSETH, NILS CHR. STENSETH Centre for Ecological and Evolutionary Synthesis (CEES), Department of Biology, University of Oslo, PO Box 1066 Blindern, NO-0316 Oslo, Norway Flødevigen Marine Research Station, Institute of Marine Research (IMR), NO-4817 His, NorwaySearch for more papers by this author First published: 07 September 2010 https://doi.org/10.1111/j.1365-2486.2010.02331.xCitations: 75 Nils Chr. Stenseth, Department of Biology, Centre for Ecological and Evolutionary Synthesis (CEES), University of Oslo, PO Box 1066 Blindern, NO-0316 Oslo, Norway, tel. +47 22 854584, fax +47 22 854001, e-mail: [email protected] Re-use of this article is permitted in accordance with the Terms and Conditions set out at http://wileyonlinelibrary.com/onlineopen#OnlineOpen_Terms AboutSectionsPDF ToolsRequest permissionExport citationAdd to favoritesTrack citation ShareShare Give accessShare full text accessShare full-text accessPlease review our Terms and Conditions of Use and check box below to share full-text version of article.I have read and accept the Wiley Online Library Terms and Conditions of UseShareable LinkUse the link below to share a full-text version of this article with your friends and colleagues. Learn more.Copy URL Abstract It is well known that human activities, such as harvesting, have had major direct effects on marine ecosystems. However, it is far less acknowledged that human activities in the surroundings might have important effects on marine systems. There is growing evidence suggesting that major reorganization (i.e., a regime shift) is a common feature in the temporal evolution of a marine system. Here we show, and quantify, the interaction of human activities (nutrient upload) with a favourable climate (run-off) and its contribution to the eutrophication of the Black Sea in the 1980s. Based on virtual analysis of the bottom-up (eutrophication) vs. top-down (trophic cascades) effects, we found that an earlier onset of eutrophication could have counteracted the restructuring of the trophic regulation at the base of the food web that resulted from the depletion of top predators in the 1970s. These enhanced bottom-up effects would, however, not propagate upwards in the food web beyond the zooplankton level. Our simulations identified the removal of apex predators as a key element in terms of loss of resilience that inevitably leads to a reorganization. Once the food web has been truncated, the type and magnitude of interventions on the group replacing the apex predator as the new upper trophic level have no effect in preventing the trophic cascade. By characterizing the tipping point at which increased bottom-up forcing exactly counteracts the top-down cascading effects, our results emphasize the importance of a comprehensive analysis that take into account all structuring forces at play (including those beyond the marine system) at a given time. Introduction The Black Sea is a deep, mostly land-locked, basin in Eastern Europe. It is linked to the Mediterranean by the narrow straits of Bosporus and Dardanelles (Fig. 1a). The surrounding land area entertains intensive human activities and has experienced profound economical and societal changes in the formerly communist countries. That the Black Sea has undergone dramatic environmental changes in recent decades underlies its importance as a 'natural laboratory' for studying marine ecosystem dynamics (Mee et al., 2005; Daskalov et al., 2007; Oguz & Gilbert, 2007). Figure 1Open in figure viewerPowerPoint Black Sea and biological series. (a) Map showing the location of the Black Sea in Europe and the mouth of the Danube River. (b–e) Observations and predictions [as estimated from the individual generalized additive models, Eqns (1)–(4)] for phytoplankton, zooplankton, jellyfish, and planktivorous fish. The Black Sea is the world's largest meromictic basin consisting of a two-layer system separated by a permanent pycnocline (Sorokin, 2002). This density boundary effectively limits the vertical exchange between the oxygenated upper layer-influenced by the atmospheric and fluvial processes – and the almost completely isolated anoxic deep water. Despite its >2000 m depth, most of the biological activity (apart from bacteria) is hosted within the upper 100–150 m. The Black Sea is characterized by a positive water balance that results in a net outflow into the Mediterranean. With a drainage basin five times more extensive than the sea area (Ludwig et al., 2009) it works as a virtually isolated ecosystem, and is sensitive to distant anthropogenic activities. This terrestrial influence, together with a contrasting bathymetry and a cyclonic Rim Current (Stanev, 1990), contributes to the Black Sea horizontal zonation (Ragueneau et al., 2002). Two distinct regions can be recognized: the wide and shallow Northwest Shelf ( 1000 m). The latter is mostly isolated from the riverine inflow, which is known to be a key driver on the shelf. Although hydrographic processes, such as mesoscale eddies, filaments, and jets, effectively link these two subsystems together (Zatsepin et al., 2003), they have been seen to show biological differences (McQuatters-Gollop et al., 2008). Productivity of the shelf system appears to be primarily phosphorus limited whereas the open sea system would appear to be nitrogen limited and much more dependent on mixing processes for nutrient supply (Garnier et al., 2002). Climate affects the Black Sea via atmospheric transfer and riverine inflow. The latter has been demonstrated as a significant factor for the overall water balance and basin-scale circulation (Oguz et al., 1995), as well as nutrient loading from human activities in surrounding land. The Danube River provides about 70% of the freshwater inflow. Thirty-three and 56% of the phosphorus emissions are estimated to be derived from agriculture and urban settlements, respectively; only 8% is considered to be of natural origin (Kroiss et al., 2005). During the 1980s, the Black Sea underwent severe eutrophication caused by economical and lifestyle changes in the surrounding countries, including intensive animal farming and increasing use of agrochemicals and phosphate detergents. The physical environment of the Black Sea has a major influence across the food web at different time scales (Daskalov, 2003) and has been shown to be influenced by the Atlantic climate through cross-Europe atmospheric teleconnections (Polonsky et al., 1997; Oguz et al., 2006). The food web in the Black Sea is relatively simple and effects of both resource (bottom-up) and predation (top-down) have been identified. Major effects of predators at top and middle trophic levels have been found to drive system-wide trophic cascades (Daskalov et al., 2007). The overfishing of pelagic top predators in the 1970s, of planktivorous fish in the 1990s, and the unintentional introduction with ships' ballast water of the ctenophore Mnemiopsis leidyi (Konsulov & Kamburska, 1998) resulted in alternating changes in the abundance of the phytoplankton and zooplankton populations, which disturbed the structure and functioning of the entire pelagic food web (Kideys, 2002; Murray, 2005). The Black Sea have been populated, exploited, and explored by humans since the antiquity, but major anthropogenic changes such as fish stock collapses, cultural eutrophication, and invasions by alien species have occurred since the 1980s. Initially most changes were attributed solely to cultural eutrophication (Zaitsev, 1993; Bologa et al., 1995). More recently other factors, such as hydroclimate (Daskalov, 2003; Oguz et al., 2006), predation effects, and fishing (Bilio & Niermann, 2004; Daskalov et al., 2007) have been recognized as contributing to the changes. As put forward above, the recent history of the Black Sea is a combination of abrupt ecological events of great interest to the scientific community. Therefore, this system is an excellent location to study how the marine food web responds to various perturbations that, to varying degrees, occur in the world's oceans. Human activities affect ecological processes in a variety of ways. Harvesting and climate change (Stenseth et al., 2002), for instance, are known to have broad ecological consequences. It is less appreciated that activities in one ecological biome might affect the ecology of another biome. The sensitivity of the Black Sea to human-induced changes in the Danube watershed makes this system an ideal test basin to investigate the effect of socio-economical transformations on the marine biome. In this study we first address the dynamics of the Black Sea food ladder by estimating an individual model for each of the trophic levels: phytoplankton, zooplankton, gelatinous plankton, and fish. This set of models allows us to empirically study how the terrestrial, climatic, and marine (environmental and trophic regulation) effects influence the food web. The model formulation is tailored to detect and quantify the ecological thresholds at which a given covariate changes its effect on the response variable. In the second part we combine the previous empirically deduced relationships in one single statistical model. On this basis, the new model reproduces the observed biomasses based only on external drivers and the estimated relationships amongst trophic levels. With the focus on the trophic architecture of the food web, this global food-web model is run under hypothetical scenarios. Making use of a novel methodology, the present study aims to provide insight on how the marine food web restructures to accommodate changes in the intensity of different pressures (e.g., fishing or eutrophication) and by doing so assess the resilience of the Black Sea ecosystem as its capacity to buffer and withstand disturbance (Holling, 1973; Folke, 2006). Material and methods Trophic levels Previous work has established that cascading trophic interactions can explain the main patterns in the Black Sea time series (Daskalov, 2002, 2003; Daskalov et al., 2007). These interactions are detected across trophic levels and characterize the dominant flows of biomass in the food web. In the present study the system's food web complexity is compressed into five components, corresponding to four trophic levels: primary producers (phytoplankton), primary consumers (zooplankton), secondary consumers (planktivorous fish, jellyfish), and top predators (piscivorous fish). Although both planktivorous fish and gelatinous plankton feed on zooplankton, they are considered separately due to their different ecosystem functioning and management implications. Gelatinous plankton comprises Aurelia aurita and M. leidyi while planktivorous fish includes anchovy (Engraulis encrasicolus), sprat (Sprattus sprattus), and horse mackerel (Trachurus mediterraneus ponticus). Diet spectrum and trophic flow arguments (Shlyakhov & Daskalov, 2008) are at the base of such aggregation, see also Ecopath model in Daskalov (2002). Data series We used annual time series accounting for the various trophic levels (Fig. 1b–e) and environmental variables. The total database consisted of the biomass of phytoplankton (PHY), the biomass of zooplankton (ZOO), the gelatinous plankton biomass (GEL), the planktivorous fish biomass (FIS), fishing mortality (F), the predatory fish biomass (PRE), the sea surface temperature (SST), the North Atlantic Oscillation (NAO) index, and the total inorganic phosphorus loading in the Danube delta (P). The biological time series were compiled based on data from long-term monitoring collected by the Institute of Fisheries and Aquaculture, Varna (Bulgaria), and the Southern Scientific Research Institute of Marine Fisheries and Oceanography (YugNIRO), Kerch (Ukraine). Data were standardized to zero mean and unit variance, see details in supporting information material in Daskalov et al. (2007). This dataset is intended to be representative of the whole Black Sea. This is particularly valid for fish stocks which are estimated using population models applied to data from all Black Sea fisheries (Prodanov et al., 1997; Daskalov et al., 2008). All plankton components however, might be biased, giving the Northwest Shelf dynamics a proportionally larger weight than the open sea because of the higher productivity and intensity of processes as well as more accurate and frequent sampling along the shelf areas. Fishing mortality (F) was estimated as the ratio of total catch to biomass of the three dominant species of planktivorous fish in terms of biomass and catches (Prodanov et al., 1997). This index is meant to account for the cumulative 'trophic' effect of the fisheries on FIS and through them on other groups such as jellyfish and plankton. Predatory fish biomass (PRE) includes bonito, bluefish and mackerel, all pelagic fish predators mainly feeding on FIS. The SST time series consists of annual mean values over the whole Black Sea area extracted from the ICOADS dataset published in http://ingrid.ldeo.columbia.edu/SOURCES/.NOAA/.NCDC/.ERSST/.version2/.SST/ The NAO index corresponds to the difference in normalized sea level pressures between Lisbon (Portugal) and Reykjavik (Iceland) over the winter season and was extracted from http://www.cgd.ucar.edu/cas/jhurrell/indices.html Total inorganic phosphorus loading (P, tonnes) were measured at the Vilkovo station of Kilya branch of the Danube River. Data were compiled and analysed by Daskalov (2003) based on Juravleva & Grubrina (1993) and Weber (1993). Variations in phosphorus loading reflect well human activities in the catchment area (Kroiss et al., 2005) and as such can be considered a proxy for the anthropogenic forcing in the Black Sea system. Statistical analysis The annual averages of the trophic levels' biomasses were used as the response variable and regressed against the various biotic (i.e., the other trophic levels) and environmental variables in the year before. The regression analysis was performed using generalized additive models (GAM) (Hastie & Tibshirani, 1999). Model estimation. To avoid model over-fitting, the number of knots used in each of the GAM splines were kept to a maximum of four. As we were interested in characterizing nonadditive responses in relation only to the relative abundance of the various trophic levels and the environmental conditions, time ('years') was not used as a predictor. These precautionary measures and the common model selection procedures (see below) ensure the parsimony of the models and that the simulations are only based on the dynamic structure of the system. Model selection. Model selection was based on a step-wise approach, aimed at removing covariates with a P-value >0.05 and minimizing the generalized cross validation (GCV) criterion of the model (Wood, 2000). The GCV is a proxy for the model's out-of-sample predictive performance and it is analogous to Akaike's Information Criterion (Akaike, 1974). The residuals of the models appeared to be uncorrelated over time and followed a normal and homoscedastic distribution in all cases (Fig. S1), except for the gelatinous plankton model (Fig. S1c). To test the effect of these two outliers on the fitted model, we refitted the model by including two dummy variables accounting for the outliers (see text and Fig. S3 in SI). Threshold GAM. Several regime shifts have been reported in the system (Daskalov et al., 2007) indicating that the food-web interactions and their relationship with the environment might be nonadditive (i.e., different across regimes). To account for regime-dependent relationships, we used a modified GAM formulation, the threshold generalized additive model (TGAM). This GAM formulation allows for nonadditive effects of the explanatory variables below and above a certain value of a threshold variable (or a combination of variables), i.e., the regression structure is allowed to switch between two GAMs. The threshold is estimated from the data. Detection of regime-dependent dynamics. To compare threshold models (TGAM) with the fully additive model (GAM) formulations (i.e., without threshold) it is necessary to account for the additional parameter used for the threshold search (Ciannelli et al., 2004). The above-mentioned GCV is only a (good) approximation of the real CV and it does not take into account the fact that a grid search has been put in place to find the value of the threshold. Thus, we used the genuine CV to compare models (Table S1), which equals the average squared leave-one-out prediction errors; the leave-one-out prediction is obtained by removing one data case at a time from the model fitting and predicting its value from the resulting model. Sensitivity analysis. CV was also used to assess the predictive performance of the final set of models (see details in supporting information, Figs S4–S8). Simulations. The fitted models were used to simulate the observed dynamics after linking the trophic levels together. Specifically, we used the observations at time t to predict the various trophic levels at time t+1. Once we got the first prediction for the different trophic levels (at t=2), the latter were input as biological variables in the various models to predict the subsequent values at t+2, t+3, …, t+n. The covariates were fixed at their observed values. By doing this, we let the food web interact according to the estimated models. Also, as we did not used 'time' as a predictor (see above), the simulations are only based on the dynamic structure of the system. Noise was added to the biotic variables by sampling (with replacement) the model residuals. To preserve the contemporaneous correlation of errors, a whole vector of errors for the four trophic levels corresponding to a randomly sampled year was used at a time. One thousand Monte Carlo simulations were run for each trophic level from which the mean and the 95% prediction bands were calculated. Scenario construction. This skeletal food-web model was afterwards used to investigate the evolution of the system under different conditions (i.e., scenarios). The procedure consisted of three steps: (a) we defined scenarios where some variables were either increased or decreased by a percentage of the mean (e.g., −25%, −15%, +15%, +25%), (b) these modified variables (or scenarios) were input to the various models and, (c) the 'simulated system' (i.e., the biomass for the different trophic levels under a given scenario) was investigated with reference to the prevailing food-web control (bottom-up vs. top-down) using phase space plots. All the models were coded in r (v 2.5.1) (R Development Core Team, 2007) using the TGAM library (created by K.-S. Chan) that relies on the mgcv library (Wood, 2006). All the plots (except 1, 3) were made with r. Figure 3Open in figure viewerPowerPoint Food-web regulations. Schematic representation of the main trophic interactions under high (a) and low (b) biomass of planktivorous fish, which roughly coincides with the opposite state for gelatinous plankton. Arrows pointing upwards represent resource control (positive effect between consecutive trophic levels). Arrows pointing downwards represent predator control (negative effect). Cascading effects are represented by dashed lines crossing through a trophic compartment. The threshold effect of phosphorus on the phytoplankton dynamics is represented by an oblique dashed line. Results Black Sea ecological dynamic structure The most appropriate model structure found for each trophic level is shown below [Eqns (1)–(4)], where denote nonparametric smooth functions (natural cubic splines) with the first argument enclosed in the parentheses being the covariate and the second argument the estimated degrees of freedom of the splines. The threshold variables and the threshold values delineating the regimes are also given. In the case of bivariate threshold, the regimes are delineated by a line estimated from the data (see Fig. 2d). The residuals showed no serial auto-correlation (Fig. S1) indicating that the following set of models captured most of the system's variability (an average of 70% of explained variance, see details in Table 1 and observations vs. predictions from these models in Fig. 1): (1) (2) (3) (4) Figure 2Open in figure viewerPowerPoint Statistical models. Threshold estimation (first column) and partial plots showing the main biotic and abiotic effects for each of the trophic levels: phytoplankton (a), zooplankton (b), gelatinous plankton (c), and planktivorous fish (d). For the univariate thresholds, phosphorus (a) and fish (b–c), the threshold estimation (generalized cross validation minimization) and threshold value (θ) defining the low (black) and high (red) regime are indicated. For the fish model (d) the blue line (θ) corresponds to the bivariate threshold that assigns the space made by the two variables (zooplankton and jellyfish) to the low (black dots) and high (red dots) regimes. The individual effects are referred either to the low (black) or the high (red) regime of the threshold variables. Those effects acting throughout the whole range of the threshold variable are shown in blue. The y-axis indicates the partial additive effect that the term on the x-axis has on the response variable. The numbers in parentheses on the y-axis indicate the estimated degrees of freedom, which also appear in Table 1. Residuals check (independence, normality, and homoscedasticity) and regime assignation of the actual levels of the threshold variable are shown in Fig. S1. Table 1. Generalized additive models (GAM) models results PHY ZOO Estimate P-value Estimate P-value Intercept 0.114 0.235 Intercept −0.178 0.009 Threshold (θ) 47.52 Threshold (θ) 0.634 Regime Covariate edf P-value Regime Covariate edf P-value P≤θ ZOO 1.00 0.101 FIS≤θ GEL 1.54 θ NAO 2.23 <0.001 FIS≤θ SST 1.00 θ PRE 2.90 0.017 R 2 (adj)=0.710 R 2 (adj)=0.816 GEL FIS Estimate P-value Estimate P-value Intercept 0.033 0.697 Intercept≤θ 0.700 θ 1.378 <0.001 Regime Covariate edf P-value Regime Covariate edf P-value FIS≤θ ZOO 1.66 θ ZOO 1.00 0.023 Z/G>θ NAO 2.88 θ PHY 1.00 θ ZOO 1.16 θ, respectively. Note that for the fish model the threshold is defined by a line (intercept: 1.07, slope: 2.18) and not a single value (Fig. 2d). All regimes are defined in terms of the lag 1 of the threshold variables. These results support that the various trophic levels relate nonadditively to the environment and other trophic levels because the models including thresholds are preferred to their fully additive equivalents, based on CV (Table S1). The nonadditivity consists of the responses switching between two distinct regression functions upon crossing a level (threshold) given by a threshold variable(s) that could be either environmental [Eqn. (1)], biological [Eqns (2) and (3)], or a combination of two biological variables [Eqn. (4)]. GAMs are relatively complex regression techniques in terms of the mathematical formulas behind the smoothers, but are very intuitive when presented pictorially by plotting the graphs of its component functions (Fig. 2). GAMs also enjoy the advantage of being nonparametric (i.e., there is no need to a priori specify the functional forms between the response and the explanatory variables). This characteristic gives great flexibility as we let the data tell us what these functional forms look like. Phytoplankton. Phytoplankton showed a nonadditive response corresponding to different levels of the first lag of phosphorus load (Fig. 2a and Table 1). When this was low, the biomass of zooplankton had a slightly negative effect, suggesting that the latter were able to efficiently graze on phytoplankton. When the levels of phosphorus were high, negative NAO values had a strong effect indicating enhanced climate-driven primary productivity. A positive winter NAO index is associated with cold and dry air masses in southern Europe and the Black Sea region because the westerly winds take a more northwards direction. Conversely, a negative NAO index implies milder winters, with warmer air temperatures and less dry/more wet atmospheric conditions over the Black Sea due to the more direct effect of the Westerlies over the region (Oguz, 2005). Negative NAO years are therefore associated with greater run-off and higher temperatures (Polonsky et al., 1997; Konsulov & Kamburska, 1998; Oguz et al., 2006). The combination of favourable atmospheric conditions (i.e., negative NAO) and high phosphorus emissions results in increased phytoplankton biomass. For the whole range of phosphorus emissions, we found a positive effect of gelatinous plankton on phytoplankton suggesting a cascading effect through predation on zooplankton. Overall the model explained 71% of the variance (see R2 in Table 1) and the predictions matched very well the observations, not only for the low frequency oscillations but also for the high frequencies. See how the predictions are able to capture most of the observed peaks in Fig. 1b (see also the out-of-sample prediction performance in Fig. S4a). The emissions of phosphorus over the years and residuals check are shown in Fig. S1a. Zooplankton. The dynamics of zooplankton were found to shift between two regimes delineated by the level of the lag 1 of planktivorous fish abundance (Fig. 2b and Table