Probabilistic Ecosystem Model for Predicting the Nutrient

Policy Analysis pubs.acs.org/est

Probabilistic Ecosystem Model for Predicting the Nutrient Concentrations in the Gulf of Finland under Diverse Management Actions Jarno P. Vanhatalo,*,† Laura M. Tuomi,‡,§ Arto T. Inkala,∥ S. Inari Helle,† and J. Heikki Pitkan̈ en§ †

Fisheries and Environmental Management group (FEM), Department of Environmental Sciences, University of Helsinki, P.O. Box 65, FI-00014 University of Helsinki, Finland ‡ Finnish Meteorological Institute/Marine Research, P.O. Box 503, FI-00101 Helsinki, Finland § Finnish Environment Institute/Marine Research Centre, P.O. Box 140, FI-00251 Helsinki, Finland ∥ Environmental Impact Assessment Centre of Finland Ltd., Sinimäentie 10 A, FI-02630 Espoo, Finland ABSTRACT: Many countries define legislative targets for the ecological status of aquatic ecosystems. Fulfilling these legally binding targets requires often large scale and expensive management actions. The expected benefits from alternative actions are commonly compared with deterministic ecosystem models. However, from a practical management point of view the uncertainty in model predictions and the probability to achieve the targets are as essential as the point estimates provided by the deterministic models. For this reason, we extend a deterministic ecosystem model into a probabilistic form. We use the model for predicting the probability to achieve the targets set by EU’s Water Framework Directive (WFD) in Finnish coastal waters in the Gulf of Finland, one of the most eutrophicated areas of the Baltic Sea, under alternative management scenarios. Our results show that the probability to reach the WFD objectives for total phosphorus is generally less than or equal to 0.51 in all areas. However, for total nitrogen the probability varies substantially as it is practically zero in the western areas but almost 0.80 or higher in the eastern areas. It seems that especially with phosphorus, international co-operation is needed in order for Finland to fulfill the objectives of the WFD.

1. INTRODUCTION

It has been estimated that reaching the targets of the BSAP would require a reduction of 15 250 tonnes of phosphorus (P) and 135 000 tonnes of nitrogen (N),5 and that the costs of reaching the BSAP objectives in the whole BS would be 1620− 3 000 million euros per year.10,11 Managers are thus faced with a problem characterized by a strong need for high-cost abatement measures that would be carried out within a complex natural system with high inherent uncertainty. In order to provide quantitative tools to support the decision making, the effects of alternative management actions are commonly studied with ecosystem models that simulate the nutrient fluxes and algae concentrations as functions of nutrient loads to the sea.12,13 These models are commonly deterministic. However, when making management decisions, the uncertainty in the model prediction should be taken explicitly into account since it, for example, affects the expected utility of a decision and allows managers to explicitly weight between the possible benefits and risk of failure.14,15 Often it is not feasible to build a stochastic ecosystem model which would assess the uncertainty in prediction by construction. Instead, we need to rely on available deterministic

During the past decades, excessive nutrient loads have resulted in major changes in marine and freshwater ecosystems around the world. 1−4 The effects of nutrient enrichment are pronounced in regional seas, the catchments of which are highly altered and densely populated, and whose physical characteristics, like limited water exchange, extend the residence time of the nutrients.5 This is the case also in the Baltic Sea (BS), one of the largest brackish water areas in the world, where eutrophication is “perhaps the single greatest risk to the BS environment”.5 As eutrophication has many adverse effects, the society endeavors to find effective ways to combat it. Today, the protection of the BS is prescribed in the legislation of the European Union (EU), that is, the Water Framework Directive (WFD)6 and the Marine Strategy Framework Directive (MSFD),7 international conventions (e.g., Helsinki Convention),8 programs based on these conventions (e.g., the Baltic Sea Action Plan, BSAP),9 and national legislations of the coastal states. The BSAP aims for restoring the good ecological status of the Baltic marine environment by 2021, and the WFD binds the EU member states legally to strive for good ecological status in inland and coastal waters already by 2015. Both the WFD and BSAP define the ecological status using indicators such as the nutrient concentration in the water. © 2012 American Chemical Society

Received: Revised: Accepted: Published: 334

June 19, 2012 October 11, 2012 November 28, 2012 November 28, 2012 dx.doi.org/10.1021/es302475v | Environ. Sci. Technol. 2013, 47, 334−341

Environmental Science & Technology

Policy Analysis

Figure 1. The BS with a close up of our focus area the GoF. The top right map shows the measurement locations of the training data of the PEM. The lower right map shows the WFD areas in the Finnish coastal zone (WI arch. with yellow, EI arch. with green, WO arch. with red and EO arch. with blue) and the Längden measurement station in the WO arch. with black star.

models and assess the uncertainty in their predictions. In this work, we demonstrate how to build a probabilistic ecosystem model (PEM) upon a deterministic EIA-SYKE 3D ecosystem model.16−18 We use the PEM to calculate the probability to achieve Finland’s WFD target on the nutrient concentrations in the Finnish coastal waters in the Gulf of Finland (GoF) (see Figure 1) under alternative management measures. The PEM is used also to verify the accuracy of the EIA-SYKE model both spatially and temporally.

To implement the WFD in Finland, the Finnish coastal waters in the GoF have been divided into four areas:25 the southwestern inner archipelago (WI arch.), the GoF inner archipelago (EI arch.), the southwestern outer archipelago (WO arch.) and the GoF outer archipelago (EO arch.) (Figure 1). There are several indicators for the ecological status and each area has its own targets for them. In this work, we concentrate on the total N (totN) and total P (totP) concentration which are monitored from January to March. As a case example, we analyze how probably three alternative policies fulfill the objectives set by the WFD in these areas. The policies compared are Baseline (no additional actions compared to the situation in early 2000s), BSAP9 and the Finnish national high reduction scenario (FINH) devised according to the outlines of the Finnish water protection policy.26 The main assumptions concerning the nutrient loads under these scenarios are summarized in Table 1. 2.3. The 3-Dimensional Ecosystem Model. In order to assess how likely the alternative management actions fulfill the objectives set by the WFD we need a model description from the external loads to the water quality in the GoF. Such predictions are provided by the EIA-SYKE model, which is a deterministic 3-dimensional hydrodynamic-biogeochemical model with a simple nutrient cycle (inorganic nutrients) and two phytoplankton groups. The model has been developed mainly to describe the conditions in the GoF, but the model includes also the Gulf of Bothnia, the Gulf of Riga and the Baltic Proper. The flow model is 3D hydrostatic, z-coordinate model in Arakawa E-grid.27−29 The transports are calculated with TVD-superbee scheme and for vertical turbulence k-model is used. Data-assimilation was used in the Baltic proper for salinity and temperature layers located under halocline. The ice conditions have been taken into account based on the ice charts

2. MATERIALS AND METHODS 2.1. Study Area and Alternative Management Actions. The GoF is an elongated estuary in the northeastern BS surrounded by Finland, Estonia, and Russia (Figure 1). It has suffered from severe eutrophication during the past decades as nutrient loads relative to surface area have been 2−3 times the average in the BS.19 Although P loads have been decreasing since the 1990s especially in Russia,20,21 the GoF is still one of the most eutrophicated areas in the BS.5 The process is intensified due to the properties of the catchment as well as the gulf itself. The catchment maintains a population of over 12 million people22 and harbors St. Petersburg, the largest city within the BS. The sensitivity of the GoF for nutrient loads is high due to salinity stratification, which hampers vertical mixing, and a fairly long renewal time of water.23 The hindered vertical mixing of water leads to oxygen depletion, and extensive bottom areas are classified hypoxic or anoxic yearly.24 In anoxic conditions, the release of P from the sediments to the overlaying water is escalated, resulting in so-called “internal load” due to which the P concentrations in the GoF have remained high despite the decreasing trend in external nutrient loads.19 335

dx.doi.org/10.1021/es302475v | Environ. Sci. Technol. 2013, 47, 334−341


Policy Analysis

Table 1. External Load Scenarios for Finland, Estonia and Russiaa baseline Finland Change Estonia Change Russia Change Total Change

FINH

N

P

N

P

N

245

10 018

278

8313

2575

42 943

3098

61274

177 −28 217 −22 1194 −54 1588 −49

9043 −10 8236 −1 38 471 −10 55750 −9

129 −47 278 0 2575 0 2982 −4

6481 −35 8313 0 42 943 0 57737 −6

(%) (%) (%) (%)

BSAP

P

a

The numbers represent the biologically available P and N in tons per year. Figure 3. A time series from the Längden measurement station (shown in Figure 1) comparing the DIP and DIN data (blue dots) with the EIA-SYKE model (red curve) and the PEM (the predictive mean with a black curve and central 95% probability region with a shaded area). The scenario prediction after 2005 (1y-5y) corresponds to BSAP loads (see Section 2.1).

produced by the Finnish Ice Service. The biogeochemical model has a simple cycling for the dissolved inorganic P (DIP) and N (DIN) from algal cells via detritus back to DIN and DIP, and a separate module for describing nutrient flows from sediment to water under oxygenated and reduced conditions.17,30 Figure 2 shows a schematic description of the model

The input from the EIA-SYKE model to the PEM (Section 2.3) are the weekly average DIP and DIN concentrations over 0−20 m depth, which represents the water column where primary production takes place. In scenario predictions, which are five years ahead from 2004, only the external loads are altered and the other inputs to the EIA-SYKE model, such as the meteorological forcing and river runoffs (Figure 2), are replicated from the calibration time. 2.3. Probabilistic Ecosystem Model (PEM). In general, several sources of uncertainty can be recognized when using computer models to predict a real system. These include, for example, the parameter uncertainty, model inadequacy, residual variation, and observation error in the data.31 In our application, the parameters can be divided into context specific model parameters, such as the meteorological forcing, river runoffs and the growth rates (see Figure 2), and to the nutrient loads which change between the scenarios. In principle, we should set a prior distribution for the context specific parameters and by using the Bayes theorem32 infer their posterior distribution which would then be used in scenario predictions.31,33 In practice, the prior elicitation and solving the posterior distribution would not be feasible in our case, for which reasons we will treat the parameters of the EIA-SYKE model as known and fix them to the point values determined by the manual calibration (Section 2.2). Thus, we account for the uncertainty in the model prediction only through the model inadequacy, residual variation and observation error which is a similar approach to, for example, the spatiotemporal downscaler of Berrocal et al.18 The observation error represents the residual between the true nutrient concentration and its measurement. Since nutrient concentrations are non-negative, exhibit large annual changes and the variability in the observations increase with their magnitude we use a Gaussian distribution with square root transformation to model the observation error. That is, let f(s,t) denote the square root of the true average concentration, in the grid cell s on day t and relate it to the square root of the measured concentration through y(s,t) = f(s,t) + ε, where ε ∼ N(0,σ2) represents the measurement error. We relate the square root of the EIA-SYKE model output, x(s,Tt), to the expectation of f(s,t) through a linear model E[f(s,t)] = α0(s,t) +

Figure 2. Structure and the main processes of the biogeochemical module of the EIA-SYKE model. Calculated variables are DIN, DIP, two algae groups (cyanobacteria and other algae), detritus and P−Fe = Fe(III) oxide bound phosphorus.

and references provide a detailed description of the hydrodynamic27−29 and biogeochemical16,17 components. The model uses spatial resolution of 5 km and the vertical dimension is divided into 17 layers that are 2 m thick until the depth of 20 m and irregularly spaced after that. Due to the complexity of the biogeochemical processes and the limited amount of measured data, some of the process parameters, such as the growth rate parameters of the algae and coefficients and limiting factors in the CO2 fluxes, need to be calibrated in the model. This was done for years 2000−2004 using DIN, DIP and chlorophyll-a measurements from two intensive monitoring stations from the Finnish coast of the GoF. Different values for the parameters were tested based on the measurements and values presented in the literature. The values that resulted to the best match with calibration data set were further evaluated in validation against independent data set, including one coastal and three offshore stations.16 Figure 3 shows the time series from one of the calibration stations. 336



Policy Analysis

β(s,t)x(s,Tt), where Tt denotes the week containing day t. The multiplicative bias β(s,t) corresponds to a time and location specific regression coefficient and the additive bias α0(s,t) to the model inadequacy that is independent of the EIA-SYKE model output.31 This model is supposed to predict the expected value of the real process under conditions specified by the input parameters for the EIA-SYKE model. However, the process may not always take the same value if the inputs were repeated and this variation in the real process is called residual variability.31 We model the residual variability with additive term α1(s,t) so that our final model is

km spatial grid (the same as in the EIA-SYKE model) and over 0−20m depth. 2.5. Inference and Prediction. We give weakly informative Student’s-t prior distributions38 for the parameters of CFs and approximate their posterior distribution by a maximum a posterior (MAP)39 estimate found by gradient based optimization.40 After this, the posterior distribution of f(s,t) conditional on the parameters can be solved analytically35 and used to predict for new location-time pairs. The posterior distribution of an integrated variable, such as the average concentration ∫ f(s,t)ds/|S| on the sea area S is approximated via Monte Carlo (MC).41 Since the EIA-SYKE model is a highly nonlinear function of its inputs, its output with fixed parameter values is, in general, different from what the expected value of its output would be with respect to uncertain parameters. However, this difference is presumably at most moderate since the earlier experience from the regular use of the EIA-SYKE model (e.g., the yearly predictions for the summertime cyanobacteria blooms) suggests that its predictive accuracy is not sensitive to different loads and initial conditions. In general, by fixing the parameters we restrict the form of the EIA-SYKE model output and presumably decrease its predictive accuracy.31 The discrepancy between the EIA-SYKE model output and the training data is taken up by the bias terms. Due to finite temporal length-scales, they contribute the less to the mean and the more to the uncertainty the further into future we predict. After the correlations between the historical and the future bias terms have decreased close to zero, the predictive uncertainty of PEM saturates at level that is defined by their variation over the training time period. See Figure 3 for illustration. Thus, by fixing the context specific parameters of EIA-SYKE model we presumably decrease the predictive accuracy of PEM but at the same time increase the uncertainty about its future predictions as would be desired. The conversion from the dissolved inorganic nutrients predicted by PEM into total nutrients, for which the WFD targets are defined, is done using linear regression so that, for example, totP = a+b × DIP+ε, where ε ∼ N(0,λ2). The weights, a and b, and the variance of the error term, λ2, are given the conjugate normal-Inverse-Gamma prior39 after which their posterior is solved using measurements made on R/V Aranda in the GoF in 2000−2004. 2.6. Model Assessment. We assess the predictive performance of the EIA-SYKE model and the PEM using the leave-one-out cross-validation (LOO−CV).42 The LOO−CV is based on calculating for each training data point the posterior predictive distribution given all the other data, after which the predictive distribution and the observation are compared with some test statistics. Due to this procedure the model predictions are compared to data that have not been used in the training set. The LOO−CV predictive distributions conditional to the CF parameters can be solved analytically43 but the CF parameters would need to be inferred numerically with each training set separately. Since this is prohibitive from the time consumption point of view, the CF parameters of PEM are fixed to their MAP estimate conditional on all data. This should be a reasonable approximation as a single observation has negligible effect on the posterior of the CF parameters with a large data set.38 The test statistics used in the comparison are the root-meansquare error (RMSE), (1/nΣni = 1(E[ỹi|y\i]−yi)2)1/2, mean log predictive density (MLPD), 1/nΣni = 1log p(yi|y\i), and the

y(s , t ) = α0(s , t ) + α1(s , t ) + β(s , t )x(s , Tt ) + ε

which will be called the probabilistic ecosystem model (PEM). We model the bias terms with a Gaussian process (GP) which is a stochastic process that defines a probability distribution over functions.34,35 This implies a GP prior distribution for f(s,t) and by conditioning on the spatiotemporal monitoring data (Section 2.4) we can calculate its posterior distribution which, unlike the deterministic EIA-SYKE model, provides uncertainty estimates for the prediction. The GP is defined by a mean function, for example, mβ(s,t) = E[β(s,t)], and a covariance function (CF), for example, kβ([s,t],[s′,t′]) = E[[β(s,t)−mβ(s,t)][β(s′,t′)−mβ(s′,t′)]], which determine the properties, such as the smoothness and variability, of the process. Here, we assume a zero prior mean for α0 and α1, and a unit prior mean for β, which imply that the prior mean of f(s,t) equals the EIA-SYKE model output. The model inadequacy α0 and the residual variability α1 are not strictly separable but the assumption made here is that the former should, together with β, reflect the “global” discrepancy between the EIA-SYKE model prediction and observations, and the latter should reflect correlated noise. By global we mean discrepancies which are similar over large areas and long periods of time. These assumptions are encoded into the model by the choice of the CFs which commonly follow the separable form k([s,t],[s′,t′]) = ks(s,s′)kt(t,t′) for spatiotemporal processes.36 The smoothness assumption concerning α0 can be reflected by a Gaussian CF35 kG(v,v′) = σ2Gexp(−r2/2) in both spatial and temporal domain. Here r2 = Σmi = 1(vi′ − vi)2/l2i is the scaled euclidean distance between length m vectors v and v′, σ2G is the scale parameter and l2i is the length-scale in direction i. The CF of α1 is a product of a Matérn CF with 3/2 degrees of freedom km(s,s′) = σ2m(1 + √3r)exp(−√3r) for the spatial and an exponential CF ke(t,t′) = σ2e exp(−r) for the time domain that give rise to rather rough spatiotemporal processes.35 The CF of β is a product of a Matérn CF with 3/2 degrees of freedom for the spatial and a delta function with respect to a week for the time domain, km(s,s′)δ(Tt,Tt′), which corresponds to linear weights that are spatially correlated but independent across weeks.18 These covariance structures were the best among several alternatives that were compared as described in Section 2.6. 2.4. Training Data for the PEM. The training data for the PEM comprise of 3377 DIN and 3647 DIP observations from years 2000−2004 whose spatial coverage is illustrated in Figure 1. The data are obtained from the database of the International Commission for the Exploration of the Seas (ICES),37 which collects and distributes the Helsinki Commission (HELCOM) monitoring data collected by different institutes around the Baltic Sea. The data were averaged to one day resolution, 5 × 5 337



Policy Analysis

probability integral transform (PIT) values,44,45 ∫ yi0 p(ỹi|y\i)dỹi, where y\i denotes all the other data points expect yi. The RMSE assesses the predictive mean whereas MLPD takes into account also the uncertainty in the predictive distribution. The distribution of PIT values can be used to analyze how well the predictive distributions are calibrated as discussed by Gneiting et al.45 In order to include the EIA-SYKE model into the MLPD and PIT value comparison, we assign it a Gaussian noise model whose variance is inferred in a similar manner as the CF parameters in PEM.

is too intensive vertical mixing in the hydrodynamic model reported e.g. by Myrberg et al.29 During summers, the EIASYKE model overestimates DIP concentrations everywhere in the GoF and during springs it predicts too high DIN concentrations in the west. Otherwise there are no systematic differences between the EIA-SYKE model and the PEM predictions for DIN. 3.2. Expected Effects of Nutrient Load Reductions. Figure 5 shows the probability distributions of the average

3. RESULTS AND DISCUSSION 3.1. Model Assessment. The PEM shows improved behavior compared to the EIA-SYKE model in the RMSE and MLPD statistics (Table 2). The distribution of the Table 2. RMSE and MLPD Statistics of the EIA-SYKE Model and PEM Predictions for DIP and DIN Concentrations MLPD

RMSE

model/variable

DIP

DIN

DIP

DIN

EIA-SYKE PEM

−1.89 −1.08

−2.90 −2.03

11.24 4.88

74.28 34.75

residuals and PIT values (not shown here) show that the EIA-SYKE model predicts in average too high concentrations and occasional extreme residuals increase its RMSE statistics. According to the distribution of the PIT values, the predictive distributions of the PEM are well calibrated. The PEM can also aid in analyzing the performance of the EIA-SYKE model by indicating the simulation errors in the spatial dimension. The verification of deterministic ecosystem models is spatially limited, since the comparisons are typically done against measurement stations that have high temporal coverage, which are few in the Baltic Sea. Methods to verify the accuracy of deterministic models spatially are not commonly used, even though they give useful information of the model behavior. Our analysis shows that the EIA-SYKE model can simulate with reasonable accuracy the winter-time DIN and DIP concentration. However, it overestimates in general the DIP concentrations in the westernmost part of the study area which is illustrated in Figure 4. The evident reason behind this

Figure 5. The posterior predictive distributions of the average concentration of totP and totN in the Finnish WFD areas (see Figure 1) in January-March over the last three prediction years. The highlighted area shows the probability (Pr) to achieve the good water quality and black cross shows the EIA-SYKE model point estimate. The totP (totN) concentration is classified to be in good level if it is below 29 μg/l (488 μg/l) in the WI arch., (585 μg/l) in the EI arch., 27 μg/l (345 μg/l) in the WO arch., and 29 μg/l (510 μg/l) in the EO arch.46 Note that WFD objectives are not defined for totP in the EI arch.

January−March totP and totN concentrations in the Finnish coastal areas over the three last prediction years. The probability to achieve the WFD objectives for the totP concentration is the highest with the BSAP. The BSAP reduces the totP concentration more in the east than in the west. Even though Finland’s P load reduction under the FINH is higher than under the BSAP (see Table 1), the totP concentrations remain higher under the former than the latter. These results

Figure 4. DIP concentrations on the second week of January 2002 predicted by the EIA-SYKE model and the PEM (the posterior mean). The black crosses in the lower figure show the observations from that week. 338



Policy Analysis

are reasonable since in early 2000s Russia has been responsible for 83% of the estimated bioavailable P loads to the GoF and the BSAP would reduce the P load mostly in the east since Russia’s share of the reduction is also the greatest. Both policy options reduce the totN concentration compared to the baseline (Figure 5). The FINH reduces the concentration the most which is reasonable since the BSAP would reduce the N load only moderately more than the FINH (Table 1) and Finland’s reductions affect more directly the Finnish coastal waters than reductions by Russia and Estonia. The probabilities to reach the WFD objectives are higher in the east than in the west where it is practically impossible to achieve the objectives. In general, the probability to reach the WFD objectives for totP is below or equal to 0.51. With totN there is a clear difference between east and west as the probability to reach the target is practically zero in the west, but nearly 0.80 or higher in the east. Moreover, it seems unlikely that Finland could alone fulfill its WFD objectives for totP but, in the east, there is over 95% probability that Finland could fulfill the WFD objectives for totN with FINH. One reason for the large difference in the probability to reach the totN targets between the east and west is that the objective is much tighter in the west than in the east. The ecological status boundaries for nutrient concentrations have not been intercalibrated between EU member states and in Finland they have been defined by analyzing the frequency distributions of monitoring data from years 1962−2005.46 It has been suggested that also summertime nutrients should be included in the assessment work in the Finnish coastal waters, as the link between winter nutrient concentrations and summer phytoplankton growth is not straightforward.47 As our scenario predictions are thought to cover the subsequent 5 years from 2004, it is interesting to compare the results with today’s situation. Since 2005, Russia has practically fulfilled its BSAP obligations by implementing a chemical P removal in the wastewater treatment plants of St. Petersburg.20 Our results suggest that this should be seen in the P concentrations in the eastern GoF. Monitoring results from resent years seem to support this conclusion. However, due to release of P from the sediments to the overlaying water in anoxic conditions,19 the interannual variations in the P concentration are strong and final conclusions only few years after the load reductions cannot be made. Further, as our scenarios do not match with the actual development seen in the GoF in 2005−2009 (e.g., in Finland and Estonia the fulfillment of the BSAP obligations is still underway), a more elaborate comparison between our results and the reality is not meaningful. Figure 5 illustrates also why it is important to quantify the uncertainty in the forecast. E.g., when EIA-SYKE model predicts that the target is not met, the probability of success may vary from 0.00 to 0.51. It is important that managers are aware of this kind of vagueness of the results. Further, Burgman15 has noted that point estimates oblige decisionmakers to be risk-neutral, although they seldom are like it. Instead, risk-averse people tend to favor actions that lead to more certain outcomes, whose utilities may however be low, whereas risk-loving people aim at higher but more uncertain utilities.48 However, point estimates do not provide any measure of uncertainty related to different decisions, and thus deterministic models may results in dissatisfaction with decision-making process.15

Our treatment for the context specific parameters of the EIASYKE model follows the common approach to fix them according to the calibration of the computer model.18,41,49 As discussed in Section 2 this is not optimal and the correct approach would be to infer their posterior distribution. However, this is generally computationally demanding and we are aware of only few examples for such treatment with large scale computer models.33,50 Another source of parameter uncertainty left out from our treatment is related to the external loads which may play a major role in the goodness of the model prediction. For example, just recently it was found that the bioavailable P loads into the GoF may be up to 1000 tons higher than the previous estimate given in Table 1.51 This about 30% increase to the previous estimates comes from a Russian Fertilizer Plant area alongside river Luga, whose real discharges have not been available for the annual data gatherings of HELCOM until fall 2011.51 In decision making the uncertainty about the future loads should also be accounted for in the form of the implementation uncertainty. We have assumed the loading reductions to become realized with certainty, but in reality people will not necessarily act as anticipated. The relevance of implementation uncertainty has been demonstrated, for example, in the case of the Baltic salmon fisheries management.52 If the uncertainty related to the implementation of management actions is not acknowledged, the society’s capability to manage environmental risks can be overestimated. Furthermore, as our model covers only a 5-year period, the effects of climate change are not taken into account. However, within the Baltic Sea, climate change is predicted to increase precipitation and runoff that may cause intense eutrophication in coastal areas.53 From the management point of view, this means that in the future it may be even more unlikely to reach the water quality targets of the WFD. To conclude, assessing the effects of nutrient load reductions in a probabilistic context is essential for decision-makers, who have to make decisions that can cost substantially but whose outputs may be highly uncertain. Our analysis is focused on historical values in the sense that the reference years are in the early 2000s and we look at scenarios that are partly fulfilled already. However, here we have assumed that the reductions take place immediately whereas in reality they are implemented gradually. Thus, our analysis gives still valuable insight for the coming years. In addition, it seems that the WFD targets will not be met in the Finnish coastal areas by 2015, and thus more effort should be made to reach the targets by 2021, that is, by the end of the next six-year management cycle defined by the WFD. The methodology presented here could be used when planning these management actions.

■

AUTHOR INFORMATION

Corresponding Author

*Phone: +358 44 3765458; e-mail: jarno.vanhatalo@helsinki.fi. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work was prepared in the project “Integrated Bayesian risk analysis of ecosystem management in the Gulf of Finland” (IBAM), funded by the European Community’s Seventh Framework Programme under grant agreement no 217246 made with the joint Baltic Sea research and development 339



Policy Analysis

(20) Pitkänen, H.; Knuuttila, S.; Lehtoranta, J.; Peltonen, H. Eutrophication - Present status and future perspectives. In The Gulf of Finland, Finnish-Russian-Estonian Cooperation to Protect the Marine Environment, History and Prospects for the Future; Rintala, J-M. and Myrberg, K., Eds.; Ministry of the Environment of Finland: Helsinki, 2009. (21) HELCOM. The Fifth Baltic Sea Pollution Load Compilation (PLC-5). In Baltic Sea Environment Proceedings, No. 128; Helsinki Commission: Helsinki, 2011. (22) Julie, S.; Sindre, L.; Carl, F. Land cover and population density in the Baltic Sea drainage basin: A GIS database. Ambio 1996, 25, 191−198. (23) Andrejev, O.; Myrberg, K.; Alenius, P.; Lundberg, P. A. Mean circulation and water exchange in the Gulf of FinlandA study based on three-dimensional modelling. Boreal Environ. Res. 2004, 9, 1−16. (24) Hansson, M.; Andersson, L.; Axe, P. Areal extent and volume of anoxia and hypoxia in the Baltic Sea, 1960−2011. Report Oceanogr. 2011, 42, 1−76. (25) Vuori, K-M.; Bäck, S.; Hellsten, S.; Karjalainen, S. M.; Kauppila, P.; Lax, H-G.; Lepistö, L.; Londesborough, S.; Mitikka, S.; Niemelä, P.; Niemi, J.; Perus, J.; Pietiläinen, O-P.; Pilke, A.; Riihimäki, J.; Rissanen, J.; Tammi, J.; Tolonen, K.; Vehanen, T.; Vuoristo, H.; Westberg, V. Suomen pintavesien tyypittelyn ja ekologisen luokittelujärjestelmän perusteet (The basis for typology and ecological classification of water bodies in Finland). Finnish Environ. 2006, 807, 1−15. (26) Rekolainen, S.; Vuoristo, H.; Kauppi, L.; Bäck, S.; Eerola, M.; Jouttijärvi, T.; Kaukoranta, E.; Kenttämies, K.; Mitikka, S.; Pitkänen, H.; Polso, A.; Puustinen, M.; Rautio, L. M.; Räike, A.; Räsänen, J.; Santala, E.; Silvo, K.; Tattari, S. Rehevö ittävän kuormituksen vähentäminen Taustaselvitys osa I, Vesiensuojelun suuntaviivat vuoteen 2015 (Reducing nutrient loads to the surface waters Background study part I, Guidelines for Water Protection to 2015). Rep. Finn. Environ. Inst. 2006, 22. (27) Koponen, J.; Alasaarela, E.; Lehtinen, K.; Sarkkula, J.; Simbierowicz, P.; Vepsä, H.; Virtanen, M. Modelling the dynamics of a large sea area. Public Water Environ. Res. Inst. 1992, 7, 1−91. (28) Inkala, A.; Myrberg, K. Comparison of hydrodynamical models of the Gulf of Finland in 1995A case study. Environ. Modell. Software 2002, 17, 237−250. (29) Myrberg, K.; Ryabchenko, V.; Isaev, A.; Vankevich, R.; Andrejev, O.; Bendtsen, J.; Erichsen, A.; Funkquist, L.; Inkala, A.; Neelov, I.; Rasmus, K.; Rodriguez Medina, M.; Raudsepp, U.; Passenko, J.; Söderkvist, J.; Sokolov, A.; Kuosa, H.; Anderson, T. R.; Lehmann, A.; Skogen, M. D. Validation of three-dimensional hydrodynamic models of the Gulf of Finland. Boreal Environ. Res. 2010, 15, 453−479. (30) Korpinen, P.; Kiirikki, M.; Rantanen, P.; Inkala, A.; Sarkkula, J. High resolution 3D-ecosystem model for the Neva Bay and estuary Model validation and future scenarios. Oceanologia 2003, 45, 67−80. (31) Kennedy, M. C.; O’Hagan, A. Bayesian calibration of computer models. J. R. Stat. Soc., Ser. B 2001, 63, 425−464. (32) Bernardo, J. M. and Smith, A. F. M. Bayesian Theory; John Wiley & Sons, Ltd: Chichester, 2000. (33) Kennedy, M.; Anderson, C.; O’Hagan, T.; Lomas, M.; Woodward, I.; Gosling, J. P. Quantifying uncertainty in the biospheric carbon flux England Wales. J. R. Stat. Soc., Ser. A 2008, 171, 109−135. (34) O’Hagan, A. Curve Fitting and Optimal Design for Prediction. J. R. Stat. Soc., Ser. B 1978, 40, 1−42. (35) Rasmussen, C. E. and Williams, C. K. I. Gaussian Processes for Machine Learning; The MIT Press: Cambridge MA, 2006; p 248. (36) Finkenstädt, B.; Held, L.; Isham, V. Statistical Methods for SpatioTemporal Systems; Chapman & Hall/CRC: Boca Raton FL, 2007. (37) ICES. Dataset on Ocean Hydrography; The International Council for the Exploration of the Sea: Copenhagen, 2009. (38) Vanhatalo, J.; Pietiläinen, V.; Vehtari, A. Approximate inference for disease mapping with sparse Gaussian processes. Stat. Med. 2010, 29, 1580−1607. (39) Gelman, A.; Carlin, J. B.; Stern, H. S.; Rubin, D. B. Bayesian Data Analysis; Chapman & Hall/CRC: Boca Raton FL, 2004.

programme BONUS, and the Academy of Finland. We also thank the anonymous reviewers for comments that helped us to improve this work.

■

REFERENCES

(1) Smith, V. H.; Tilman, G. D.; Nekola, J. C. Eutrophication: Impacts of excess nutrient inputs on freshwater, marine, and terrestrial ecosystems. Environ. Pollut. 1999, 100, 179−196. (2) Smith, V. H. Eutrophication of freshwater and coastal marine systems: A global problem. Environ. Sci. Pollut. Res. 2003, 10, 126−139. (3) Smith, V. H.; Schindler, D. W. Eutrophication science: Where do we go from here? Trends Ecol. Evol. 2009, 24, 201−207. (4) Millennium Ecosystem Assessment. Ecosystems and Human WellBeing: Biodiversity Synthesis; World Resources Institute: Washington DC, 2005; p 100. (5) HELCOM. Eutrophication in the Baltic SeaAn integrated thematic assessment of the effects of nutrient enrichement and eutrophication in the Baltic Sea region. In Baltic Sea Environment Proceedings, No. 115B; Helsinki Commission: Helsinki, 2009. (6) European Parliament and Council. Directive 2000/60/EC of the European Parliament and of the Council establishing a framework for Community action in the field of water policy. Off. J. Eur. Communities: Legis. 2000, 327, 1−72. (7) European Parliament and Council. Directive 2008/56/EC of the European Parliament and of the Council establishing a framework for community action in the field of marine environmental policy. Off. J. Eur. Communities: Legis. 2008, 164, 1−22. (8) HELCOM. Convention on the Protection of the Marine Environment of the Baltic Sea Area; Helsinki Commission: Helsinki, 1992; p 45. (9) HELCOM. HELCOM Baltic Sea Action Plan. HELCOM Extraordinary Ministerial Meeting,Krakow, Poland, 15 November 2007; HELCOM: Helsinki, 2007; p 101. (10) HELCOM and NEFCO. Economic Analysis of the BSAP with Focus on Eutrophication. Final Report; COWI, 2007; p 112. (11) Gren, I. Costs and Benefits from Nutrient Reductions to the Baltic Sea. Swedish Environmental Protection Agency: Stockholm, 2008; p 69. (12) Pitkänen, H.; Kiirikki, M.; Savchuk, O. P.; Raike, A.; Korpinen, P.; Wulff, F. Searching efficient protection strategies for the eutrophied Gulf of Finland: The combined use of 1D and 3D modeling in assessing long-term state scenarios with high spatial resolution. Ambio 2007, 36, 272−279. (13) Wulff, F.; Savchuk, O. P.; Sokolov, A.; Humborg, C.; Morth, C. M. Management options and effects on a marine ecosystem: Assessing the future of the Baltic. Ambio 2007, 36, 243−249. (14) Kuikka, S.; Hildén, M.; Gislason, H.; Hansson, S.; Sparholt, H.; Varis, O. Modelling environmentally driven uncertainties in Baltic Cod management by bayesian influence diagrams. Can. J. Fish. Aquat. Sci. 1999, 56, 629−641. (15) Burgman, M. Risks and Decisions for Conservation and Environmental Management; Cambridge University Press: Cambridge, 2005; p 504. (16) Kiirikki, M.; Inkala, A.; Kuosa, H.; Pitkänen, H.; Kuusisto, M.; Sarkkula, J. Evaluating the effects of nutrient load reductions on the biomass of toxic nitrogenfixing cyanobacteria in the Gulf of Finland Baltic Sea. Boreal Environ. Res. 2001, 6, 131−146. (17) Kiirikki, M.; Lehtoranta, J.; Inkala, A.; Pitkänen, H.; Hietanen, S.; Hall, P. O. J.; Tengberg, A.; Koponen, J.; Sarkkula, J. A simple sediment process description suitable for 3D-ecosystem modelling Development and testing in the Gulf of Finland. J. Mar. Syst. 2006, 61, 55−66. (18) Berrocal, V.; Gelfand, A. E.; Holland, D. M. A spatio-temporal downscaler for output from numerical models. J. Agric. Biol. Environ. Stat. 2009, 15, 176−197. (19) Pitkänen, H.; Lehtoranta, J.; Räike, A. Internal nutrient fluxes counteract decreasing in external load: The case of estuarial Eastern Gulf of Finland, Baltic Sea. Ambio 2001, 30, 195−201. 340



Policy Analysis

(40) Vanhatalo, J.; Riihimäki, J.; Hartikainen, J.; Jylänki, P.; Tolvanen, V.; Vehtari, A. Bayesian Modeling with Gaussian Processes using the MATLAB Toolbox GPstuff (v3.4). arXiv 2012, 1206.5754v2. (41) Fuentes, M.; Raftery, A. E. Model evaluation and spatial interpolation by Bayesian combination of observations with outputs from numerical models. Biometrics 2005, 66, 36−45. (42) Vehtari, A.; Lampinen, J. Bayesian model assessment and comparison using cross-validation predictive densities. Neural Comput. 2002, 14, 2439−2468. (43) Sundararajan, S.; Keerthi, S. S. Predictive approaches for choosing hyperparameters in Gaussian processes. Neural Comput. 2001, 13, 1103−1118. (44) Dawid, A. P. Statistical theory - The Prequential Approach. J. R. Stat. Soc., Ser. B 1984, 147, 278−292. (45) Gneiting, T.; Balabdaoui, F.; Raftery, A. E. Probabilistic forecasts, calibration and sharpness. J. Roy. Stat., Soc. B 2007, 69, 243−268. (46) Vuori, K-M.; Mitikka, S.; Vuoristo, H., Eds.; Pintavesien ekologisen tilan luokittelu (Guidance on ecological classification of surface waters in Finland − Part 1: Reference conditions and classification criteria, Part 2: Environmental impact assessment). Environmental Administration Guidelines 3/2009, Finnish Environmental Institute: Helsinki, 2009. (47) Fernandes, J.; Kauppila, A.; Uusitalo, P.; Fleming-Lehtinen, L.; Kuikka, V.; Pitkänen, S. H. Evaluation of reaching the targets of the Water Framework Directive in the Gulf of Finland. Environ. Sci. Technol. 2012, 46, 8220−8228. (48) Clemen, R. T. Making Hard Decisions: An Introduction to Decision Analysis; Duxbury Press: Belmont CA, 1996. (49) Berrocal, V. J.; Raftery, A. E.; Gneiting, T.; Steed, R. C. Probabilistic Weather Forecasting for Winter Road Maintenance. J. Am. Stat. Assoc. 2010, 105, 522−537. (50) Dietzel, A.; Reichert, P. Calibration of computationally demanding and structurally uncertain models with an application to a lake water quality model. Environ. Modell. Software 2012, 38, 129− 146. (51) HELCOM. Significant Source of Phosphorus Found in Kingisepp, Russia, http://www.helcom.fi/projects/on_going/ balthazar/en_GB/river_Luga/ (accessed January 6, 2012). (52) Haapasaari, P.; Karjalainen, T. P. Formalizing expert knowledge to compare alternative management plans: Sociological perspective to the future management of Baltic salmon stocks. Mar. Policy 2010, 34, 477−486. (53) BACC Author Team Assessment of Climate Change for the Baltic Sea Basin. Regional Climate Studies; Springer-Verlag: Berlin Heidelberg, 2008.

341


Probabilistic Ecosystem Model for Predicting the Nutrient

Recommend Documents