Combined climate change and nutrient load impacts on future habitats and eutrophication indicators in a eutrophic coastal sea

Eutrophication and climate change will affect habitats of species and more generally, the structure and functioning of ecosystems. We used a three-dimensional, coupled hydrodynamic-biogeochemical model to investigate potential future changes in size and location of potential habitats of marine species during the 21 st century in a large, eutrophicated brackish sea (the Baltic Sea, northern Europe). We conducted scenario projections under the combined impact of nutrient load and climate change. Possible future changes of the eutrophication state of this sea were also assessed through two policy-relevant indicators. The results imply a physiologically more stressful environment for marine species by the end of the 21 st century: volumes of higher salinity water become more hypoxic/anoxic and the volumes of low salinity, oxic water increase. For example, these results impact and reduce cod reproductive habitats. The decrease is mainly climate change induced in the Baltic basins less directly in ﬂ uenced by in ﬂ ows of saline, oxic water to the Baltic Sea (E Gotland and Gdansk Basins). In basins more directly in ﬂ uenced by such in ﬂ ows (Arkona and Bornholm Basins), the combined effect from climate change and nutrient loads is of importance. The results

Climate change has already had an impact on coastal ocean ecosystems and hydrography (Richardson et al. 2012;Rhein et al. 2013;IPCC 2019).Likewise, eutrophication has been reported as a growing problem around the world, creating hypoxic (oxygen [O 2 ] concentration of 0-2 mL L −1 ) and anoxic (O 2 concentration < 0 mL L −1 , i.e., "negative oxygen," which correspond to the amount of oxygen needed to oxidize the hydrogen sulfide) areas in coastal oceans (Breitburg et al. 2018), for example, in Chesapeake Bay, the Gulf of Mexico, the South China Sea, the Black Sea, and the Baltic Sea (Rabalais et al. 2007;Murphy et al. 2011;Strokal and Kroeze 2013;Su et al. 2017;Meier et al. 2018c;Reusch et al. 2018;Murray et al. 2019).
The decline in O 2 concentration will most likely be further exacerbated by climate induced changes (Altieri and Gedan 2015;IPCC 2019).Climate models predict that areas with high precipitation will become even wetter in the future, especially in the Northern Hemisphere (Holopainen et al. 2016;IPCC 2019).Higher net precipitation increases river discharge, supplying more nutrients and freshwater to the sea.Consequently, eutrophication is accelerated due to the increased nutrient supply and primary production.This results in higher O 2 consumption due to the decomposition of the increased biomass.Enhanced supply of freshwater decreases salinity, which may act to increase the stratification and could potentially decrease the vertical mixing of O 2 saturated surface waters to depth.In addition, O 2 dynamics will be affected by temperature-dependent biological rates and water solubility (Irby et al. 2018).
The changes due to climate change and eutrophication are also seen in the projections for the Baltic Sea, for example, increased temperature and enhanced river discharge as the net precipitation over its catchment area increases (Eilola et al. 2012;Meier et al. 2012a,b;Neumann et al. 2012;Ryabchenko et al. 2016;Saraiva et al. 2019a,b).Because the spatial distribution of productivity and biota is strongly controlled by hydrographic conditions such as salinity, temperature, and O 2 concentration, the combined effects of climate change and eutrophication can be expected to have major impacts on habitats, the distribution and biodiversity of biota, and the functioning of ecosystems and food webs (Hinrichsen et al. 2016;Breitburg et al. 2018;Saraiva et al. 2019a).For example, the food web structures will become increasingly dominated by smaller organisms (Suikkanen et al. 2013) and changes in the environment will affect the fish community and may demand modifications to existing fisheries management policies (MacKenzie et al. 2007;Niiranen et al. 2013).The biota living in the Baltic Sea will potentially become more vulnerable to further stress from eutrophication and global climate change as they originate from both freshwater and marine environments, and even under a stable climate often live close to the limits of their physiological tolerances (Ojaveer et al. 2010).
The Baltic Sea is a semi-enclosed coastal ocean with an estuary-like circulation with negligible tidal influence.The water exchange with adjacent seas is limited because of its narrow and shallow entrances.The supply of saltier, O 2 -rich water occurs through episodic inflows of smaller and major water masses, which penetrate to deeper layers below the permanent halocline (Matthäus et al. 2008).Smaller inflows occur more frequently, but do not penetrate as deeply or flow as far east as the major inflows (so-called Major Baltic Inflows).These bathymetric and hydrographic properties result in a strong, permanent halocline at a depth of 60-80 m, limiting the vertical exchange between well-oxygenated surface water and weakly ventilated, oxygenpoor deep water.Consequently, the Baltic Sea is one of the most eutrophic seas in the world (Kemp et al. 2009;HELCOM 2017) and parts of it suffer from permanent hypoxia while other parts are periodically hypoxic or anoxic.The occurrence of low-oxygen waters is due to the excess inputs of a large supply of nutrients, in combination with long residence times and frequent occurrences of large cyanobacteria blooms (Carstensen et al. 2014).Therefore, this sea is under high pressure from human induced activity (The BACC II Author Team 2015; Reusch et al. 2018) and eutrophication is thought to be the most severe threat to biodiversity (HELCOM 2009).To improve the health and water quality of the Baltic Sea, the Helsinki Commission's (HELCOM) Baltic Sea Action Plan (BSAP) was signed by the EU and all countries adjacent to the Baltic Sea in 2007, and has been regularly updated at ministerial meetings (HELCOM 2007b(HELCOM , 2018c)).Eutrophication core indicators have been developed as a tool to follow progress toward achieving good environmental status (GES) by 2021.
In this study, we assess the future ecosystem functioning of the Baltic Sea.We investigate how volumes and locations of water masses with different combinations of salinity and O 2 will change in the future.This combination of variables is particularly important for estuarine biota, which can be a combination of marine, estuarine, and freshwater species (Ojaveer et al. 2010).As an example, we analyze the future projected change in volume of cod reproductive habitat (known as cod reproductive volume [CRV]), which is a climate-sensitive integrated parameter based on the minimum salinity and O 2 conditions that allow for successful cod egg development in the Baltic Sea (Plikshs et al. 1993;MacKenzie et al. 2000;Köster et al. 2005).These conditions are usually found near or below the permanent halocline in this area.We use our model framework to simulate how this habitat is expected to change, and the approach improves our understanding of how climate change and eutrophication interact and affect habitat characteristics in a eutrophic sea.Furthermore, we estimate the effect on two eutrophication core indicators, the summer chlorophyll a (Chl a) concentration in the surface water and the annual O 2 debt (a measure of the undersaturation of O 2 ) below the halocline (HELCOM 2013b(HELCOM , 2018a,b),b).

Method
Climate change impact studies are associated with several sources of uncertainty, including those due to model complexity and the greenhouse gas emission scenarios employed (Payne et al. 2016;Meier et al. 2019).A combination of multiple model ensemble members, global general circulation models (GCMs), and an increasing number of different greenhouse gas scenarios is necessary to robustly minimize uncertainty.While the computational demand is too great for the full suite of possible model simulations, we here combine a number of available GCMs and socioeconomic greenhouse gas emission scenarios and possible nutrient supply scenarios to estimate the potential range in future outcomes for the Baltic Sea.

Model overview
We utilize the three-dimensional regional coupled hydrodynamic-biogeochemical-ice model, RCO-SCOBI, which consists of the physical Rossby Centre Ocean (RCO) model and the Swedish Coastal and Ocean Biogeochemical (SCOBI) model.The RCO-SCOBI model has previously been evaluated and utilized for extensive long-term climatological studies up to 2100 (Meier et al. 2011b(Meier et al. , 2012a(Meier et al. ,b, 2016;;Eilola et al. 2013) as well as historical reconstructions from 1850 onward (Meier et al. 2018b).A thorough comparison of hindcast simulation by the RCO-SCOBI driven with regionalized ERA-40 reanalysis data (Samuelsson et al. 2011) for the period 1970-2005 was performed by Eilola et al. (2011).The authors concluded that the model captures much of the observed variability and that the response to changing physical conditions is simulated realistically.The hindcast simulations were used in this study to evaluate the model's ability to capture the temporal evolution of CRV.
RCO has a maximum depth of 249 m and is divided into 83 vertical levels with a thickness of 3 m and a horizontal resolution of 2 nautical miles.The model is a Bryan-Cox-Semtner primitive equation circulation model with a free surface and open boundary conditions in the northern Kattegat, which is located north of the belt between Denmark and Sweden.The ocean circulation model is coupled to a Hibler-type sea ice model and the ocean's sub grid-scale mixing is parameterized using a k-ε turbulence closure scheme with flux boundary conditions.For further details, see Meier (2001Meier ( , 2007) ) and Meier et al. (2003).
SCOBI represents a simplified set of biogeochemistry dynamics using nine pelagic and two benthic variables (Fig. 1).The pelagic variables are: nitrate, ammonium, phosphate, phytoplankton, one bulk zooplankton, and detritus.In addition, the model calculates the O 2 concentration and the hydrogen sulfide concentrations, which are represented by "negative oxygen" equivalents (1 mL H 2 S L −1 = −2 mL O 2 L −1 ) (Fonselius and Valderrama 2003).The benthic part consists of nitrogen (N) and phosphorus (P) pools.
The phytoplankton pool is divided into three functional algal groups: cyanobacteria, diatoms, and others (including flagellates).The phytoplankton assimilates N and P according to the Redfield molar ratio (C : N : P = 106:16:1; Redfield et al. 1963), which is also used to calculate carbon (C).The biomass is represented by Chl a according to a constant carbon to chlorophyll mass ratio (C:Chl a = 50:1).One of the algal groups has the potential to fix nitrogen, that is, to utilize molecular nitrogen (N 2 ) as a nitrogen source and consequently may continue to grow if phosphate is available.Other processes described in the model are temperature-dependent primary productivity, decomposition of detritus as well as nitrification and denitrification, grazing of zooplankton and excretion of detritus and dissolved inorganic nitrogen and phosphorus.These processes, as well as the temperature and salinity dependent O 2 solubility in the water, are also dependent on and/or affect the O 2 concentration in the water.For further details, see Eilola et al. (2009) and Almroth-Rosell et al. (2011).

Future scenarios Global climate projections
We use four climate change simulations to form an ensemble of transient runs for the time period 1961-2100.The atmospheric forcing for the RCO-SCOBI was taken from two established global GCMs, which were regionally adapted via a dynamical downscaling approach using the high-resolution Rossby Centre Atmosphere Ocean regional climate model (Doescher et al. 2002).The atmospheric model delivered output fields at a resolution of 0.24 .The global GCMs used are the coupled ocean atmosphere models HadCM3 from the Hadley Centre in the UK and ECHAM5/MPI-OM (ECHAM5) from the Max Planck Institute for Meteorology, Germany which gave good results for Europe and the North Atlantic ( Gordon et al. 2000;Jungclaus et al. 2006;Roeckner et al. 2006).Two different IPCC greenhouse gas emission scenarios were applied: the moderate SRES scenario A1B and the high emission SRES A2 scenario.HadCM3 was forced with the A1B greenhouse gas emission scenario and ECHAM5 with both the A1B and A2 scenarios (Nakicenovic et al. 2000).In addition, the ECHAM5 A1B scenario was simulated using two different sets of initial conditions (A1B1 and A1B3).The initial fields were taken from different years of the ECHAM5 of the multicentennial p1-control simulation for the preindustrial period.Thus, they rely on the same mean preindustrial climate but represent different years.Such different global realizations are usually produced to assess how natural variability impacts the climate in the long-term evolution after the preindustrial period.Therefore, they represent realistic potential trajectories for the climate evolution.In this article, however, we do not aim to study the effect of natural variability but consider the scenarios A1B1 and A1B3 as two different possible realizations for the future climate.The resulting four climate projections applied are HadCM3_A1B, ECHAM5_A1B1, ECHAM5_A1B3 and ECHAM5_A2.
Monthly mean river runoff forcing were calculated from regional climate projections, using a statistical model based on the difference between precipitation and evaporation over land (Meier et al. 2012a).At the northern Kattegat, open boundary vertical profiles of temperature, salinity, and nutrients were relaxed to climatological mean observations for the year 1978-2007.These boundary conditions do not change with time (Meier et al. 2012c).In case of inflow, temperature, salinity, nutrient (phosphate, nitrate, ammonium), and detritus values are nudged toward observed climatological seasonal mean profiles from the southern Skagerrak.Oxygen, phytoplankton, and zooplankton concentrations in the inflowing water are, however, not predescribed but obtain model values from the boundary in the northern Kattegat.Sea level variations at the boundary were calculated with a statistical model and the effect of the global mean sea level rise was not studied.For further details, see Meier et al. (2011a,b).

Nutrient load scenario
The nutrient loads from the river runoff to the sea were calculated as the product of the concentration of nutrients and the river discharge following Eilola et al. (2009) and Meier et al. (2011a,b).In this study, three scenarios are considered, ranging from a more pessimistic (Business As Usual [BAU]) to a more optimistic (BSAP) scenario: • REFerence (REF): current  loads from rivers and current atmospheric deposition continues in the future (see Eilola et al. 2009) • BAU: increased future nutrient loads from rivers assuming an exponential growth of agriculture in all Baltic Sea countries as projected in HELCOM (2007a) and a continuation of current atmospheric deposition.• Baltic Sea Action Plan (BSAP): reduced future river loads following HELCOM (2007a) and 50% reduced atmospheric deposition.
For further details, see Meier et al. (2012b, fig. 6).In total, 12 projections were performed (Table 1).For each of the three nutrient load scenarios, an ensemble average of the four climate projections for the analyzed parameters were calculated and averaged over two 30-yr periods, the hereafter called present day climate  and the future climate (2070-2099).The calculations of the ensemble average follow earlier approaches based on an equal weighting of the two GCMs, that is, both GCMs are equally represented in the ensemble average (one HadCM3-driven and three ECHAM5-driven simulations) (Meier et al. 2012b).As an illustration of the uncertainty associated with the spread of the ensemble members, we show the minimum and maximum values of the four climate projections.

Analysis parameters
We calculated water volumes of O 2 concentration and salinity during the summer (May-August) in the southern Baltic Sea.As an example of habitats with specific requirements of O 2 concentration and salinity, we show the results for CRV (May-August).In addition, two HELCOM eutrophication core indicators, the annual O 2 debt and the summer (June-September) concentration of Chl a, were analyzed.The definition of time period and the calculated areas for the two indicators differs from the water volumes of salinity and O 2 , in order to be comparable with the HELCOM results based on observations.

Oxygen and salinity
The climate projection ensemble averages were used to study changes in summer averages (May-August) of oxic (O 2 concentration > 2 mL L −1 ), hypoxic (O 2 concentration 0-2 mL L −1 ), and anoxic (O 2 < 0 mL L −1 ; e.g., H 2 S) water volumes, in addition Table 1.A summary of the 12 projections for the three nutrient loads (REF, BAU, and BSAP) and the four climate projections performed.

HadCM3_A1B
to changes in O 2 ranges at different salinities.This was performed separately for selected sub-basins of the Baltic Sea (Fig. 2, left panel) between the two 30 yr periods.The basins' total volume is shown in Table 2.We analyzed the entire water column as well as the bottom water, which is defined as the deepest 3 m.

Cod reproductive volume
Given earlier documented impacts of salinity and O 2 concentration on cod reproductive success as detected in field studies (Plikshs et al. 1993;Köster et al. 2005) and laboratory studies (Nissling and Westin 1991a,b;Westin and Nissling 1991;Wieland et al. 1994), we estimated how the CRV (Plikshs et al. 1993;MacKenzie et al. 2000) could change due to climate change and eutrophication using the climate ensemble average for the three nutrient load scenarios.This volume is defined as the volume of water with salinity > 11 and O 2 concentration > 2 mL L −1 ; the threshold environmental values are based on experimental studies of the effects of salinity and O 2 concentration on the successful survival, fertilization, and development of Baltic cod eggs.Water with this combination of salinity and O 2 concentration represents cod spawning habitat in the Baltic Sea and is found in the Arkona, Bornholm, E Gotland, and Gdansk Basins (Fig. 2, right panel).The geographic boundaries for the in situ observations are different compared to the boundaries for the model (Fig. 2).Therefore, the modeled CRVs are estimated using the same geographic coordinates as those used for estimating CRV from in situ observations (Plikshs 2014), and are therefore directly comparable.
For validation of the model skill of CRV, two dimensionless skill metrics-the Pearson correlation coefficient (r) and the mean of cost function (C)-were calculated (Eqs. 1, 2).
where P is model value, O observation, i data number, n the total number of data points, and std(O) is the standard deviation of the observations.Table 2.The total water volumes (km 3 ) for the basins and for the bottom layer defined as the deepest 3 m of the water column (Fig. 2, left) and the water volume of the basins calculated for CRV (Fig. 2, right).The basins are: the Arkona (AB), Bornholm (BH) and E Gotland (EGB) basins, the NW Baltic Proper (NWBP), the Bothnian Sea (BS), the Bothnian Bay (BB), the Gulf of Finland (GoF), and the Gdansk Basin (GB).Good and acceptable levels of r are achieved when r is higher than two-thirds (0.66) and one-third (0.33), respectively (Edman and Omstedt 2013).If the average model results for C fall within one standard deviation of observations, they are regarded as good, and as acceptable within two standard deviations (Eilola et al. 2009).The function, C, normalizes the difference between model results and observations with the standard deviation of the observations.We calculated the average and the standard deviation for each basin's time series of CRV due to the lack of replicates of observations.In addition, we calculated the same statistics of CRV for the first and second halves of the time series.This approach has been used in earlier studies ( Eilola et al. 2009;Almroth-Rosell et al. 2016;Edman et al. 2018) and is based on methods by Oschlies et al. (2010).

Eutrophication indicators
The two calculated HELCOM eutrophication core indicators are the present and future summer (June-September) concentrations of Chl a, which are part of the criterion "Direct effects of eutrophication," and the annual O 2 debt that is part of the criterion "Indirect effects of eutrophication" (HELCOM 2013b).The model calculations for Chl a and O 2 debt (HELCOM 2013b) follow the method of Carstensen et al. (2014) but with a simplified calculation of the halocline.The halocline depth is in the present study defined as the shallowest depth at each location where 68% of the salinity difference between the basin-averaged surface salinity and the salinity of the deepest part of the basin is exceeded (B.Gustafsson pers.comm.).The simulated O 2 debt below the halocline is calculated as the volume-and timeaveraged difference between the saturated O 2 concentration and the modeled O 2 concentration below the halocline for each basin.The time average includes the complete year during present day  and future (2070-2099) time slices.The O 2 debt is calculated for the Bornholm and E Gotland Basins as well as the Gulf of Finland, where the E Gotland Basin is considered as the sum of the E Gotland Basin and the NW Baltic Proper in accordance with the HELCOM methodology (HELCOM 2013a).The Chl a concentration is calculated as the average of the upper 9 m for the summer season for the Arkona, Bornholm and E Gotland Basins (sum of the E Gotland Basin and the NW Baltic Proper), the Bothnian Sea, the Bothnian Bay, and the Gulf of Finland.The concentration of Chl a and O 2 debt measures the effects of nutrient enrichment and the eutrophication status of the investigated area.The concentration of Chl a indicates the amount of phytoplankton and is therefore an important parameter for assessing the degree of eutrophication as it usually increases with increased nutrient concentration.The O 2 debt is a measure of the deficiency of O 2 (undersaturation) relative to the O 2 saturation state, that is, the water contains less O 2 than if it was in equilibrium concentration, for example, if the water had recently been in contact with the atmosphere.The results are presented as the relative changes between present day  and future (2070-2099) time slices for the combined effects of eutrophication and climate change.

Changes in oxygen concentrations and the distribution with salinity Vertically integrated water volumes through Baltic basins
The results for the ensemble average in the REF and BAU scenarios (considering all salinities) during the summer (May-August) indicate an expansion of hypoxic and/or anoxic water volumes by the end of the 21 st century (Fig. 3, upper panel).Overall, the greatest changes are in the E Gotland Basin and the NW Baltic Proper, while the changes in the Arkona Basin are minor.There is a shift in the Arkona, Bornholm, and E Gotland Basins from oxygenated water to hypoxic/anoxic waters, except for the BAU scenario in the E Gotland Basin where oxygenated/hypoxic water shifts to anoxic water.The latter pattern is also seen in the NW Baltic Proper in both REF and BAU scenarios, where the anoxic water volume increased by 14% and 22%, respectively.In the E Gotland Basin (BAU), the anoxic water volume increased by 13% and in the REF scenario the volume with less than 2 mL L −1 (hypoxic and anoxic) increased by 9%.This indicates a tendency toward lower future O 2 concentration during the summer in the E Gotland Basin and the NW Baltic Proper in the entire water column.The same pattern is seen in the Bornholm Basin, but with lower relative changes.
By contrast, the results for the BSAP scenario point to increased summer O 2 concentration by the end of the 21 st century (Fig. 3).The water volumes with O 2 concentration higher than 2 mL L −1 increased while volumes with lower O 2 concentrations decreased in all the basins, although the changes in the Arkona and Bornholm Basins are small.For the E Gotland Basin and the NW Baltic Proper, the increase in oxic water is 2% and 3%, respectively.Note, however, that the estimates in the BSAP scenario have large spread relative to the averaged values between the different climate projections.
These analyses describe the changes in volumes of water with different O 2 concentrations throughout the entire water column and regardless of salinity.However, the changes in O 2 Fig. 4. The figure is the same as Fig. 3 but for the bottom water, defined as the lowest 3 m of the water column.
concentration in water masses with specific salinity ranges are also important for specific biological processes or species.These results (Fig. 3, lower panel) are similar in the different basins for the different nutrient load scenarios and indicate that the largest relative change is in water volume for the oxic water.In this case, the relative decrease in volume is mainly in the salinity range of 7-9 and the increase is mainly in the range of 4-7, except for the NW Baltic Proper where the increase is in slightly lower salinity (range 3-6).Hence, there appears to be a shift toward more oxic conditions in water with lower salinity and toward lower O 2 concentrations in higher salinity water.The net relative changes in hypoxic/ anoxic water in the Arkona and Bornholm Basins are small.However, in the E Gotland Basin and the NW Baltic Proper, the pattern for the hypoxic water is the same as for the oxic water but at higher salinity.There is increased hypoxic water in salinity ranges 7-10 and 7-9 for the E Gotland Basin and the NW Baltic Proper, respectively, and a decrease at higher salinity.Anoxic and hypoxic waters both increase in the same salinity range, but there is no sign of decreased anoxic water.The distributions of water volumes with lower O 2 concentrations are therefore primarily impacted in more saline waters.

Bottom water volumes in Baltic basins
The results for the bottom layer (i.e., lowest 3 m of water column) of Baltic basins indicate an even larger shift toward more severe hypoxic and anoxic conditions in a future climate (Fig. 4, upper panel).In the Arkona and Bornholm Basins (REF and BAU), the bottom water turns from oxygenated to more hypoxic/anoxic waters.In the E Gotland Basin and the NW Baltic Proper, the oxic/hypoxic water shifts toward more anoxic bottom water.The most affected basin is still the NW Baltic Proper with relative increases in volumes of anoxic water of 22% and 27% in the REF and BAU scenarios, respectively.However, the increase of anoxic water in the E Gotland Basin is almost as severe with values of 18% and nearly 27% for the same scenarios.In the Bornholm Basin, there is a relatively strong increase in anoxic bottom water with values of 8% and 13% in the REF and BAU scenarios, respectively.In the BSAP scenario, the E Gotland Basin and the NW Baltic Proper increased bottom water O 2 concentrations by 3% while the Arkona Basin and the Bornholm Basin are almost unaffected.
The relative changes in bottom water O 2 distributions for specific salinity ranges are comparable with the changes for the entire water mass, except that the hypoxic/anoxic water volume increases more in all three nutrient load scenarios and in all basins except for the Arkona Basin.

Cod reproductive volume Skill assessment of modeled representations of CRV
The simulated CRVs for the average of May and August compared to the corresponding values calculated from observations of salinity and O 2 concentration for the Bornholm, E Gotland, and Gdansk Basins show similar patterns in the three basins (Fig. 5).In the observational data set used in this study (Plikshs 2014), the CRVs in the Arkona Basin are not calculated and are therefore excluded in Fig. 5.The hindcast model results are generally in good agreement with the estimates based on observations, especially in the Bornholm and Gdansk Basins.A larger difference is noticed in the E Gotland Basin in the early 1970s, where the model output has a higher CRV compared to the observation data.This difference is a result of a larger inflow of water through the sills in the model than seen in the observations during 1970 and 1972, and this inflow increases the salinity and the O 2 concentration (not shown).
The assessment of model skill (expressed as 1−r and C) in the Bornholm and Gdansk Basins is in the range of good performance (Fig. 5b).For the E Gotland Basin, the deviation in the time series in the early 1970s (Fig. 5a) results in a poor correlation (Fig. 5b, empty triangle), just over the acceptable limit for the entire time series (Fig. 5b, pink triangle).However, the cost function is within one standard deviation of measured data and therefore of good performance.

Century-scale changes in CRV
The spatial extent of the modeled present day summer (May-August) CRV covers the Arkona Basin to the Gulf of Finland (Fig. 6).However, the results for the modeled ensemble average of the future Baltic Sea indicate a decrease in its spatial extent (Fig. 6, colored area) and magnitude (Fig. 6, color scale) for all three nutrient load scenarios.The decline in spatial extent is most severe in the REF and BAU scenarios, but is also seen in the BSAP scenario in all basins; however, the decline

Wåhlström et al.
Future changes in marine habitats and indicators in the BSAP scenario is less severe compared to the REF and BAU scenarios (Fig. 6).
The E Gotland Basin formerly had the greatest extent of CRV (Fig. 7a) but the second largest regarding the percentage of the total water volume (together with the Arkona Basin, Fig. 7b) and the Bornholm Basin the largest.However, the E Gotland Basin loses 100% of the CRV in the REF and BAU scenarios by the end of the 21 st century and as much as 95% in the BSAP scenario (Table 3).The pattern is the same for the Gdansk Basin, where most of the water with the correct prerequisites for the CRV disappeared, that is, the most severe effects of eutrophication and climate change on CRV are in the E Gotland and Gdansk Basins.However, the present day results for these basins have the largest spread.
The Arkona and Bornholm Basins are less affected in a future climate compared to the E Gotland and Gdansk Basins (Fig. 7).In the Arkona Basin, the decrease in the relative change in CRV for the three nutrient load scenarios varied between 53% and 62% (Table 3), the least affected in the BSAP scenario.The same pattern is seen in the Bornholm Basin, which has the lowest decrease in relative change (43%) in the BSAP scenario, that is, slightly more than half of the CRV is preserved in a future climate.However, in the BAU scenario, the CRV decreased by 74% and thus the Bornholm Basin had the largest spread in relative changes for the CRV between the nutrient load scenarios.
The decline in CRV by the end of the 21 st century is mainly the result of the projected reduction in salinity (Fig. 7c), that is, the climate change effect.However, the nutrient loads influence the CRV indirectly through the effect on the concentration of O 2 (Fig. 7d), which is apparent in the Bornholm Basin where the largest differences between the nutrient load scenarios are seen (Fig. 7b).Thus, if the Bornholm Basin was only affected by climate change, the decrease in CRV would be identical in the three nutrient load scenarios.For the E Gotland and Gdansk Basins, the decrease in salinity is the most important factor even though these basins also suffer from decreasing O 2 concentration.In the BSAP scenario, the nutrient abatements will lead to an improved O 2 situation, and consequently all basins' CRVs are larger than in the REF and   .The basins are: Arkona (AB CRV ), Bornholm (BH CRV ), E Gotland (EGB CRV ), and Gdansk (GB CRV ).

HELCOM eutrophication core indicators Changes in oxygen debt
The present day modeled ensemble average O 2 debts below the halocline differ between the basins: the largest O 2 debt was found in the E Gotland Basin and the NW Baltic Proper (hereinafter referred to as EGB + NWBP for the indicators) and the lowest O 2 debts in the Bornholm Basin (Fig. 8).The present day (i.e., 1970-1999) modeled ensemble average values differ slightly from observations (Table 4).The observed value in the Bornholm Basin is somewhat smaller and for EGB + NWBP and the Gulf of Finland somewhat higher compared to the ensemble average.
The BAU scenario had the greatest influence on the O 2 debt in a future climate, increasing by 33% in both the Bornholm Basin and EGB + NWBP, while the increase in the Gulf of Finland was half that value (Fig. 8, Table 4).Although the O 2 debt also increased in the REF scenario, the increase was less severe compared to the BAU scenario.The O 2 debt response from implementing the Baltic Sea Action Plan (BSAP) resulted in decreased O 2 debt in all basins.The largest decrease was in the Gulf of Finland (11%), which was also the only basin, of the three examined, that reached the GES threshold by the end of the 21 st century (Table 4).

Changes in Chl a
The present day average Chl a concentration in the upper 9 m during the summer (June-September) was highest in the    for the Bornholm Basin (BH), the E Gotland Basin and the NW Baltic Proper (EGB + NWBP) and the Gulf of Finland (GoF).Reported status for 2007-2011(observed, HELCOM 2013a;;Fleming-Lehtinen et al. 2015) and threshold values to reach GES (HELCOM 2018b) are given.Relative changes (%) in the ensemble average of O 2 debt below the halocline for three nutrient load scenarios in the different basins of the Baltic Sea.Negative values are decreased O 2 debt in the future climate compared to the present period .Parentheses after the percentages indicate whether the targets are reached ("Yes" or "No") in the different nutrient load scenarios and basins by the end of the 21 st century.Gulf of Finland and lowest in the Bothnian Sea (Fig. 9).In general, these values are slightly underestimated compared to observations and all concentrations are below the GES except for the Gulf of Finland (Table 5).However, this is also the case for the observations other than the Bornholm Basin, where the observations for the Chl a concentration are higher than the GES.

BH EGB + NWBP GoF
For the REF and BAU scenarios, the concentrations of Chl a increased in all basins in the future climate compared to the present day (Fig. 9).The relative change in Chl a concentration in the BAU scenario is more or less twice as high as in the REF scenario in all the basins, except for the Bothnian Bay where the difference is smaller (Table 5).The Arkona Basin had the largest relative increase in the concentration of Chl a in the REF and the BAU scenario of 119% and 213%, respectively.In these two nutrient load scenarios, the Bornholm and EGB + NWBP had similar relative increases in Chl a concentration.These changes are two (three) times larger than in the BAU (REF) for the Bothnian Bay.The Bothnian Sea, the Bothnian Bay, and the Gulf of Finland are less affected than the three basins mentioned above, even though the Gulf of Finland had the highest absolute value for the concentration of Chl a by the end of the 21 st century.
In contrast to the REF and BAU scenarios, the BSAP scenario shows an improved eutrophication status in three out of six basins.The most improved basin is the Bothnian Sea, where the Chl a concentration decreased by 19% by the end of the 21 st century.However, in the Arkona Basin and the Bothnian Bay, the Chl a concentrations increased by 17% and 23%, respectively, and in the Bornholm Basin the concentration of Chl a is almost unaffected.
Three out of six basins (the Arkona Basin, the Bothnian Sea, and the Bothnian Bay) are estimated to be at GES for present day Chl a concentrations, whereas three others are not (the Bornholm Basin, the Gulf of Finland, and EGB + NWBP, assuming the lower value for GES).If the BSAP is implemented, almost all basins reach the GES by the end of the 21 st century.The exception is the Gulf of Finland, in spite of the decrease in Chl a concentration (13%, Table 5).In the BAU scenario, the Bothnian Bay is the only basin that attains the GES target and remains under the threshold in the REF scenario, together with the Bothnian Sea.

Discussion
We have analyzed how volumes of water with specific levels of O 2 and salinity are likely to change in a future climate in a large estuary.These water volumes, both under present and future climates, affect the species' habitats based on their differing physiological tolerances to local hydrographic conditions such as salinity and O 2 concentrations.Our results suggest that the present biota of the Baltic Sea will experience major changes in habitats, and the sizes and locations of their future habitats will be highly dependent on societal decisions regarding nutrient loading and climate change mitigation.Increased water volumes of lower salinity, particularly at the surface, and in northern and coastal areas, together with stabilized O 2 conditions, benefits freshwater species and will cause major losses of habitat for marine species.In contrast, in the deeper areas, if no action is taken, the present low O 2 conditions will become more severe and widespread, affecting especially associated benthic habitats.Consequently, marine benthic species will be particularly vulnerable to the combined effects of climate change and eutrophication.The results imply that if the nutrient mitigation plan (BSAP) is fully implemented, the net O 2 concentration could increase throughout in the entire water column and in the bottom water in the E Gotland Basin and the NW Baltic Proper, which agrees with the results of Saraiva et al. (2019a).

Cod reproductive volume
The evaluation of estimated CRVs from model results and observations was in good agreement overall and therefore, the method of calculating CRVs from model results was used for future projections.However, the model overestimated CRV in Table 5. Present day modeled ensemble average concentration of Chl a (μg L −1 , 1970-1999) for the Arkona Basin (AB), the Bornholm Basin (BH), the E Gotland Basin and the NW Baltic Proper (EGB + NWBP), the Bothnian Sea (BS), the Bothnian Bay (BB), and the Gulf of Finland (GoF).The reported status for summer during the period 1972-1980(observed, HELCOM 2018a) ) to reach GES (HELCOM 2018a).The last three rows show relative changes (%) in the concentration of Chl a for the nutrient load scenarios in the different basins.Positive/negative values are increasing/decreasing concentration of Chl a in the future climate compared to the present day.Parentheses after the percentages indicate whether the targets are reached ("Yes" or "No") in the different nutrient load scenarios and basins by the end of the 21 st century.the E Gotland Basin in the early 1970s, implying that the changes in CRV may be smaller than estimated.Nevertheless, the main results are complete loss of CRV in the future in this basin (REF and BAU scenario).In the other three cod spawning basins, the CRV decreases for all nutrient load scenarios.These results are in contrast to the study by Saraiva et al. (2019a) as their results indicate a modest increase in CRV in the BSAP scenario.However, their future scenarios had different boundary forcing's, for example, climate forcing and different assumptions of nutrient load scenarios.

AB
The decrease in CRV is the combined effect of climate change and nutrient loads in the most western basins (Arkona, Bornholm) where inflows of saline, oxic water from the northern Kattegat occur.In the basins further away from the source of the inflows (E Gotland, Gdansk), the decrease is mainly climate change-induced.This finding is similar to another recent analysis of the effects of nutrient load reduction and climate change on Baltic fish communities (Bauer et al. 2018).Thus, the basins closest to water inflows have higher salinities and are therefore less affected by the climate change induced salinity decrease.Instead, the CRVs are stressed by the nutrient-supply fueled decrease in oxygen levels.The opposite effect is true for regions further away from the inflows that are already more affected by nutrient loads.The ranges of the spreads, representing uncertainties, of the CRV projections from the four climate ensemble members are larger between GCM's than between emission scenarios.
The parameter thresholds used in this study, and in many previous studies of Baltic cod recruitment dynamics and ecology (Sparholt 1996;Köster et al. 2005;Heikinheimo 2008;Margonski et al. 2010), represent the minimum conditions that could allow some cod eggs to be fertilized and hatched.These thresholds do not represent the minimum volume necessary for the successful hatching of all cod eggs as such volumes would be considerably smaller than those estimated here.Temperature could potentially also affect the CRV.However, a direct incorporation of a temperature-dependence on the CRV has not been considered necessary in this study.This is because egg buoyancy is salinitydependent and nearly insensitive to temperature (Sundby and Kristiansen 2015).Furthermore, temperature impacts on CRV are indirectly accommodated via their inclusion in parameterizations of processes affecting oxygen production and consumption in the ocean-biogeochemical models.Temperature can have direct physiological impacts on cod eggs (e.g., hatch success).However, hatch success of Baltic cod eggs is temperature-independent in the range 3-9 C, and declines by about 10% at 11 C (Nissling 2004).Temperatures in cod spawning areas and times are generally within this range (Wieland et al. 2000;Nissling 2004;Köster et al. 2017): for example, during 1950-2013, temperatures exceeded 9 C only twice.Furthermore, in warmer years, we would anticipate that Baltic cod would spawn earlier, given earlier spawning phenologies in warmer years for this and other cod stocks (Wieland et al. 2000;Neuheimer and MacKenzie 2014;Mcqueen and Marshall 2017).
The projected changes in salinity and O 2 concentrations will have major impacts on the productivity and eventually the biomass of cod in the Baltic Sea, even if exploitation of marine resources is maintained at the maximum sustainable or lower levels (Lindegren et al. 2010;Niiranen et al. 2013).For example, cod feeding and condition during the past 20 yr have been linked to the spatial extent of bottom layer hypoxia, and cod in poor condition often have lower survival and less economic value to fisheries (Svedang and Hornborg 2014;Casini et al. 2016a,b;Eero et al. 2016;Neuenfeldt et al. 2019).As cod is a key predator in the food webs of the central Baltic Sea (Casini et al. 2009;Niiranen et al. 2013), a decline in its biomass will have important impacts on other trophic levels (e.g., zooplanktivorous fish; benthic consumers) and the overall dynamics of food webs and fisheries.

HELCOM eutrophication core indicators
Our analyses suggest that the two indicators are not likely to reach GES by 2100 in most basins under nutrient scenarios REF and BAU.If the planned nutrient mitigations are implemented, most of the basins reach the Chl a concentration of GES but this is not sufficient to address the O 2 debt.The two eutrophication core indicators investigated are interlinked because increased primary productivity affects the O 2 demand for the degradation of organic matter associated with higher phytoplankton production and decomposition.Higher O 2 demands lead to a deficiency in O 2 and, in the worst case, to hypoxic or anoxic areas thereby reducing the sizes of suitable habitats for aerobic species.A shift from oxygenated to anoxic conditions can also impact many biogeochemical processes (e.g., the phosphorus and nitrogen cycles) (Mortimer 1941;Sundby et al. 1992;Viktorsson et al. 2012;Bonaglia et al. 2014), which can have effects on eutrophication.
The Chl : C ratios and the Redfield ratio used in the biogeochemical model SCOBI may not perfectly describe the detailed dynamics of phytoplankton and inorganic nutrients, as discussed in Meier et al. (2011b).The constant Chl : C ratio used in the present model (from Fasham et al. 1990) reflect the variations of biomass more accurately while the deviations in chlorophyll content relative to biomass are potentially less well described.However, the main features of the biogeochemical processes are described in the model and evaluations of the models' performance in comparison with observations and other Baltic Sea biogeochemical models has shown good results.The comparison of results in the present and future climate scenarios has been discussed in previous publications (Eilola et al. 2009(Eilola et al. , 2011;;Meier et al. 2012aMeier et al. , 2018aMeier et al. , 2019;;Saraiva et al. 2019a,b).Additionally, parameterization of nutrient and Chl : C ratios according to, for example, Geider and La Roche (2002) and Macintyre et al. (2000) will be considered in future work.It has been shown that models with constant molecular plankton ratios show poorer performance in the northern Baltic Sea (Eilola et al. 2011), while Fransner et al. (2018) showed that including non-Redfieldian dynamics in biogeochemical models may improve the representation of the carbon dynamics in these boreal waters.The contrasting results indicate that more research is needed to fully understand the causes for deviations and to implement new process knowledge in future scenarios.The scenarios include many uncertain components beside the differences between detailed process descriptions in the biogeochemical models.A comprehensive discussion about the uncertainties of scenario simulations of biogeochemical cycles in the Baltic Sea are discussed more comprehensively by Meier et al. (2019).
Our results imply that the implementation of the BSAP has the capacity to improve the eutrophication state of the area, which will increase the ability of many species to survive and reproduce.The future eutrophication status until 2200 was investigated for different nutrient load scenarios by Murray et al. (2019).Their results in terms of whether or not GES is reached by 2100 partly diverge from this study, but this may partly be due to the exclusion of potential impact from climate change by Murray et al. (2019).Our results imply that both climate change and nutrient loads need to be considered in future eutrophication assessments of the Baltic Sea.This suggestion was also emphasized by the 2018 HELCOM Ministerial Meeting (HELCOM 2018c).

Results in a wider context
Changes in key parameters such as salinity and O 2 concentration are important for the Baltic biota because such changes affect the metabolic capacities, productivities, and distributions of both marine and freshwater species (Diaz and Rosenberg 2008;Ojaveer et al. 2010;Holopainen et al. 2016;Paiva et al. 2018).Fisheries and aquaculture businesses could also be influenced as some targeted (e.g., cod, plaice) or cultured (blue mussels) species may become less productive, abundant, or widespread (Tedengren and Kautsky 1986;MacKenzie et al. 2007;Maar et al. 2015;Hedberg et al. 2018).
Even though changes in temperature do not directly affect the CRV, the O 2 concentration is indirectly affected by temperature-dependent physical and biogeochemical processes such as sea-air exchange and seawater solubility of O 2 , decomposition of organic material, nitrification, and phytoplankton growth (Eilola et al. 2009).In the latter biogeochemical processes, the temperature effect is exponential.The growth of phytoplankton is a source of O 2 but occurs mainly in the upper part of the water column, while the O 2consuming decomposition of organic matter occurs in the entire water column as well as in the sediment.Thus, the main O 2 supply to the deeper waters occurs through irregular inflows of denser and O 2 -rich water from the Kattegat area.These inflows impact the marine ecosystem as they change the physical and biogeochemical variables and processes as well as nutrient fluxes from the sediments.However, the inflows increase the stratification, contributing to hypoxic conditions.During the stagnation period in the late 1980s and the early 1990s, Wieland et al. (1994) suggested that O 2 depletion was the most important factor limiting the reproductive success of Baltic cod.
It is unclear how climate change will affect smaller and Major Baltic inflows (Grawe et al. 2013;Schimanke et al. 2014).However, Meier et al. (2016) concluded that the dominating drivers for the future water exchange are increased river runoff and the global mean sea level height.These drivers promote hypoxic bottom area due to increased stratification counteracting the positive impact from inflows of O 2 rich water on the marine ecosystem.In Chesapeake Bay, a rise in sea level resulted in increased stratification (Hong and Shen 2012).However, Irby et al. (2018) found that the combination of decreased nutrient loads and sea level rise by 2050 would reduce the volume of water in Chesapeake Bay with O 2 concentrations below 3 mg L −1 .Furthermore, Laurent et al. (2018) investigated the impact of climate change on the Northern Gulf of Mexico by the end of the 21 st century and found that increased stratification is an important driver for decreased bottom water O 2 content.

Conclusions
The approach and findings from the present study form a basis for understanding how climate change and eutrophication may affect habitat characteristics and thereby species distributions and biodiversity in a eutrophic regional sea.The results imply that both eutrophication and climate change have the potential to impact future essential ecosystem components and to alter the present characteristics of both pelagic and benthic habitats regarding salinity and the concentration of O 2 .The results strongly support the implementation of the planned nutrient mitigations (BSAP) as this has the potential to improve the habitat environment for species with specific O 2 requirements that affect their ability to survive and reproduce in a future climate.Without this implementation, the volumes with sufficient O 2 levels to serve as suitable habitats for many species will be reduced.This might result in translocation to other areas and/or an ecosystem shift toward more freshwater species.More specific results include: • The future combined impact of nutrient load level and climate change increase the volume of low salinity, oxic (O 2 concentration > 2 mL L −1 ) water volumes in the Baltic Sea.Higher salinity water volumes turn more hypoxic/anoxic by the end of the 21 st century.• In a future climate, the CRV decreases with the combination of climate change and nutrient load level due to decreases in both salinity and O 2 concentration.The relative importance of these two factors differs throughout the Baltic Sea; nutrient loads impact the Arkona and Bornholm Basins more than the E Gotland and Gdansk Basins where the impact of climate change is expected to be more important.
• In future eutrophication assessments we suggest including the impact from climate change in addition to nutrient load scenarios.
• Multidisciplinary, multimodel approaches, which combine physical/biogeochemical parameters, can be effective strategies for assessing how multiple impacts affect indicators of potential habitat changes.These approaches can be applied to both field measurements and model outputs.

Fig. 2 .
Fig. 2. The model domain of the Baltic Sea area, divided into the calculated basins in this study: the Arkona (AB), Bornholm (BH) and E Gotland basins (EGB), the NW Baltic Proper (NWBP), the Bothnian Sea (BS), the Bothnian Bay (BB), and the Gulf of Finland (GoF) (left panel).The basins used for calculating CRV are the Arkona, Bornholm, E Gotland, and Gdansk basins (GB) with a subscript to distinguish the area of CRV basins from the model basin (right panel).

Fig. 3 .
Fig. 3. Upper panel: The changes in water volumes relative to the basins' size (%) for the three different ranges of O 2 concentration: > 2 (turquoise), 0-2 (purple) and < 0 mL L −1 (red) between the two periods 2070-2099 and 1970-1999 during the summer (May-August) for the basins' entire water column.This is shown for the three nutrient load scenarios REF (left panel), BAU (middle), and BSAP (right) in the Arkona, Bornholm and E Gotland basins and the NW Baltic Proper.The error bars are the maximum and minimum values for the four climate projections.Lower panel: The relative changes (color scale in %) in water volumes with different concentrations of O 2 (y-axis) with salinity (x-axis).The ranges of O 2 are the same as in the upper panel and the salinity interval is one unit between 0 and 20.A positive (negative) value means an increased (decreased) volume in the 2070-2099 period compared to the 1970-1999 period.

Fig. 5 .
Fig. 5. (a) CRV (km 3 ) in the Bornholm, E Gotland and Gdansk basins for the average of May and August from 1962 to 2007 for data calculated from observations (blue) and hindcast simulation (red).The modeled CRVs are estimated using the same geographic coordinates as those used for estimating CRVs from in situ observations.Note the different scale on the y-axis.(b) Statistical validation of the complete time series of modeled vs. observed CRV for the Bornholm (yellow square), E Gotland (pink triangle) and Gdansk (green circle) basins, and for the first (empty marker) and last (gray marker) 23 yr.The plot shows the cost function (C) relative to the Pearson correlation coefficient (expressed as 1−r).

Fig. 6 .
Fig. 6.The spatial extent of the ensemble average of the CRV during the summer (colored area for May-August) for the 1970-1999 (upper panel) and 2070-2099 periods for the three nutrient load scenarios (REF, BAU, and BSAP, lower panel).The color scale represents the CRV (km 3 ) in each grid cell.

Fig. 7 .
Fig. 7. (a) CRV (km 3 ) during the summer (May-August) for 1970-1999 (blue) and 2070-2099 periods in the Arkona (AB CRV ), Bornholm (BH CRV ), E Gotland (EGB CRV ), and Gdansk (GB CRV ) basins for the three nutrient load scenarios REF (yellow), BAU (red), and BSAP (green).(b) The percentage of the basins' total water volume of CRV, (c) salinity > 11 and (d) O 2 concentration > 2 mL L −1 for the same features as in (a).The bars show the ensemble means from four climate projections and the error bars show the maximum and minimum values for the four climate projections.
scenarios.Thus, the decrease in CRV in the E Gotland and Gdansk Basins is mainly climate change induced, while in the Arkona and Bornholm Basins the combined effect of climate change and nutrient loads are of importance.

Fig. 8 .
Fig. 8.The modeled ensemble average O 2 debt (mg L −1 ) below the halocline in the Bornholm Basin and in the E Gotland Basin and the NW Baltic Proper (EGB + NWBP) as well as in the Gulf of Finland for the 1970-1999 period (blue) and for the three nutrient load scenarios REF (yellow), BAU (red), and BSAP (green) for the 2070-2099 period.The error bars are the maximum and minimum values for the four climate projections.

Fig. 9 .
Fig. 9.The concentration of Chl a during the summer (June-September) for the 1970-1999 (blue) and 2070-2099 periods in the Arkona Basin, the Bornholm Basin, the E Gotland Basin and the NW Baltic Proper (EGB + NWBP), the Bothnian Sea, the Bothnian Bay, and the Gulf of Finland for the three nutrient load scenarios REF (yellow), BAU (red), and BSAP (green).The error bars are the maximum and minimum values for the four climate projections.

Table 3 .
Relative changes (%) in the CRV for the nutrient load scenarios in the different basins.Negative values are decreased CRV in the future climate compared to the present period

Table 4 .
Present day modeled ensemble average of O 2 debt (mg L −1 ,