Potential impacts of climate change on the primary production of regional seas: A comparative analysis of five European seas

Regional seas are potentially highly vulnerable to climate change, yet are the most directly societally important regions of the marine environment. The combination of widely varying conditions of mixing, forcing, geography (coastline and bathymetry) and exposure to the open-ocean makes these seas subject to a wide range of physical processes that mediates how large scale climate change impacts on these seas’ ecosystems. In this paper we explore the response of five regional sea areas to potential future climate change, acting via atmospheric, oceanic and terrestrial vectors. These include the Barents Sea, Black Sea, Baltic Sea, North Sea, Celtic Seas, and are contrasted with a region of the Northeast Atlantic. Our aim is to elucidate the controlling dynamical processes and how these vary between and within these seas. We focus on primary production and consider the potential climatic impacts on: long term changes in elemental budgets, seasonal and mesoscale processes that control phytoplankton’s exposure to light and nutrients, and briefly direct temperature response. We draw examples from the MEECE FP7 project and five regional model systems each using a common global Earth System Model as forcing. We consider a common analysis approach, and additional sensitivity experiments.

Comparing projections for the end of the 21st century with mean present day conditions, these simulations generally show an increase in seasonal and permanent stratification (where present). However, the first order (low- and mid-latitude) effect in the open ocean projections of increased permanent stratification leading to reduced nutrient levels, and so to reduced primary production, is largely absent, except in the NE Atlantic. Even in the two highly stratified, deep water seas we consider (Black and Baltic Seas) the increase in stratification is not seen as a first order control on primary production. Instead, results show a highly heterogeneous picture of positive and negative change arising from complex combinations of multiple physical drivers, including changes in mixing, circulation and temperature, which act both locally and non-locally through advection.


Background and motivation
Regional seas are the areas where society interacts most directly with the marine environment and as such have substantial socio-economic importance. For example, it is here that the large majority of the extraction of Living Marine Resources is concentrated (Stock et al., 2011;Watson and Pauly, 2001) and that the need to identify and ensure 'Good Environmental Status' is most pressing. How global scale climate change might impact regional, coastal and shelf seas is therefore one of the key issues currently facing environmental science. It is well established that the Ocean-Atmosphere General Circulation Models (OAGCMs) used in the CMIP 1 and IPCC 2 processes are primarily designed to provide reliable information at an ocean-basin to global and decadal to centennial scale. The participating climate models in CMIP5 that include a representation of the marine ecosystem generally have a resolution of ~1 o or coarser in the ocean (Bopp et al., 2013;Taylor et al., 2012). This is substantially inadequate for resolving regional sea processes, and so some form of downscaling is required. Given the non-linear and interconnected nature of the system, this generally requires dynamical rather than statistical approaches. The ultimate goal is to provide reliable projections into the future to, for example, aid policy decisions or inform the public debate on the need for mitigation action. It is not, however, just an issue of resolution: a suite of specific physical and ecosystem process act in regional seas, which along with their particular geographic setting act to shape the climatic impacts and lead to responses that may be very different from the wider global ocean. In this paper we explore the response of five regional sea areas to potential future climate change, acting via atmospheric, oceanic and terrestrial vectors and focusing on primary production as the engine that drives the marine ecosystem. Drawing on the experience of the MEECE project 3 , we contrast five very Dynamic downscaling is increasingly used to explore the impacts of climate change in regional seas in both physics only (e.g. Adlandsvik, 2008;Olbert et al., 2012) and coupled physics-ecosystems (e.g. Holt et al., 2012;Neumann, 2010;Omstedt et al., 2012) studies. An alternative approach is to develop fine scale and multi-scale approaches to global modelling.
For example, Gröger et al. (2013) rotate the pole of an otherwise coarse resolution ocean component of an OAGCM to give higher resolution in European seas. However, such an approach limits the potential to utilize regionally adapted models. Given the need to run multiple process experiments in this uncertain and evolving field, and the need for multiple regions of interest, the use of global models as the general tool for regional seas climate impact studies is still many years off (Holt et al., 2013 The impact of climate change in regional seas is largely a boundary value problem, with the large scale climate system impacting on the smaller scale regional sea. There are feedbacks, for example between the ocean and atmosphere physical system (e.g. Schrum et al., 2003) and between the regional seas and global physical and biogeochemical cycles.
However, on the decadal/regional time/space scales considered here we presume these are not of first order importance, and we focus on regional seas as 'driven' systems. There is a wellrecognised need to explore the uncertainty in climate change impact studies. This is largely rooted in the numerical experiment design. This has several facets, including the future scenario, the treatment of natural variability, the choice of driving OAGCMs, the approach to forcing the regional model and the structure and parameters of the regional model itself. The treatment of uncertainty arising from these facets would ideally take a probabilistic approach, with multiple simulations conducted to span the uncertainty space. (e.g. Howard et al., 2010;Lowe et al., 2009;Tinker et al., In review). However, here the dimensionality of this space is large and exploring this uncertainty just at a 'minimum-maximum' level is exceptionally challenging, let alone attempting to define a Probability Density Function. Hence, we must turn to process understanding to identify the nature of impacts in a semi-quantitative fashion.
The identification of the significant pressures on the system, the processes that mediate these and the resulting sign of change is an important first step.
The remainder of this introduction focuses on some principles of biophysical interactions that mediate climate change impacts in the regional seas in question, focusing on forcing and response. Sections 2 briefly describes the models and the experiment design, this is developed more fully in Appendix 1. Section 3 presents a comparison of results between these five regions and more detailed analysis of each region. The common themes in the results, specific short comings and ways forward are discussed in section 4.

Climatic forcing in regional seas
Future climate change scenarios are prescribed in terms of the anthropogenic Greenhouse Gas (GHG) emissions. These perturb the radiative forcing, modifying the climate system, and this impacts regional seas via three vectors: atmospheric, oceanic and terrestrial.
Two robust changes in atmospheric conditions under global warming are an increase in air temperature and in lower troposphere water vapour (Held and Soden, 2006). These would be expected to lead to a decrease in both sensible and latent heat flux from the ocean until the ocean reaches a dynamic thermal equilibrium, although the former is not necessarily clear in the IPCC Assessment Report (AR) 4 models (Held and Soden, 2006). For shallow shelf seas this equilibration is rapid (seasonal), while for the deep basins it is very slow, being determined largely by the over-turning circulation. Changes in wind speed/direction and shortwave radiation are also projected by the OAGCMs, and crucial for forcing regional seas.
However, accurately simulating the regional atmospheric circulation, such as the Northeast Atlantic storm track (Lowe et al., 2009), and the cloud amount/type often challenges these models.
Changes to the distribution of open ocean nutrients play a key role in changes to oceanshelf nutrient fluxes, and consequently for shelf sea primary production (Gröger et al., 2013;Holt et al., 2012). Permanent stratification in the open ocean is necessarily increasing due to global warming, since the open-ocean-atmosphere system is not in balance with changes in lower atmosphere heat and water vapour content (Capotondi et al., 2012). While changes in temperature stratification is a globally applicable effect, changes in haline stratification are expected to play an important role at a basin scale, generally following the amplification of the hydrological cycle (Held and Soden, 2006). For example, Capotondi et al (2012) suggest an increase in haline stratification is important in the North Atlantic and Arctic regions pertinent to this study. Increasing permanent stratification reduces winter mixing, decreases winter mixed layer depths, and so reduces the total amount of nutrients entrained into the winter mixed layer. Hence the amount of nutrients available for phytoplankton growth decreases, a well-documented impact of climate change in the open ocean (Bopp et al., 2013;Steinacher et al., 2010). This is shown schematically in Figure 2, and simplistically this amount decreases linearly with decreases in winter mixed layer depth arising from climate change. It is important to note that this change is driven by permanent stratification and winter mixing; the role of seasonal stratification is considered below. On a global scale winter deepening of the mixed layer is a universal phenomenon, except close to the equator or in ice covered regions (see de Boyer Montégut et al., 2004). Hence, it is not surprising this is seen as a leading order effect globally.
Terrestrial impacts are in the form of changes in riverine inputs of freshwater, nutrients and optically active constituents. Owing to the coarse resolution of the simulations in this study (see Appendix 1) compared with that needed to simulate the details of river plumes (<1 km), we limit our consideration to terrestrial impacts driving changes in large scale nutrient budgets, rather than immediate/local effects. Here these changes are limited to changes in river flows in two of the regions (see section 2), while changes to riverine nutrient concentrations are not considered.
Alongside GHG emission induced climate change are many modes of natural variability (e.g. Grossmann and Klotzbach, 2009) on many time scales (seasonal to decadal scales are our focus here), and one of the important challenges is to distinguish this variability from longer term change. This 'signal to noise' ratio (e.g. Deser et al., 2012;Hawkins and Sutton, 2009) dictates the earliest forecast horizon at which the climate change signal can be detected, which has important implications for adaptation measures. The Forcing-Response view needs also to be considered in light of the time scales needed for the system to adjust to changing external conditions. The required adjustment times can be long (multi-year to multicentennial), particularly when processes related to benthic recycling and deep-basin to surface exchange are considered, so the distinction between natural variability and trend can be blurred by this signal propagation as much as by long term modes of external variability.
For example, the deep water exchange times for the Black Sea are ~400 years (Murray et al., 1991) and for the Baltic Sea ~30 year .

Biophysical response to climatic forcing in regional seas
Most physical processes active in regional seas are potentially impacted by global change resulting from GHG emissions, and these impacts will in turn affect the biogeochemistry and lower trophic levels of the ecosystem. Detailed descriptions of the processes can be found in Robinson and Brink (1998) and Holt et al. (In Press). The response of the coupled physical-lower trophic level system to climatic forcing depends on three paradigms of biophysical interaction: i. Transport processes that set the overall elemental (of carbon, nitrogen etc) and chemical energy budget available for biological activity of a particular region, largely on longer time scales than the biological processes. ii. The seasonal and mesoscale processes that mediate the phyoplankton's exposure to light and nutrients, on similar time scales to the biological response. iii. Direct physiological response to the environment (e.g temperature). These need to be considered in the appropriate oceanographic context, notably the transport and mixing regime, the presence/absence of sea ice and the time scale of exposure to wider oceanic conditions. For transport process (i) we consider changes over the flushing periods of surface waters resulting from changes in circulation, ventilation and exterior (e.g. deep ocean) conditions, including long term/large scale terrestrial and oceanic influence. We only briefly consider the direct temperature effects on growth and reaction rates (iii), as this is largely a chemical and biological effect, with little relation to the hydrodynamics once the changes in the temperature field are established.
In relation to (i), the regional scale transport of water and its constituent properties from the open ocean and from deep waters, alongside benthic efflux, river input, atmospheric deposition and biogeochemical losses (e.g. denitrification), sets the overall elemental budget of N, P and Si in a region. In this context, the regional seas naturally divide depending on their exposure to open-ocean exchange and the natural contrast here is between seas connected with the open ocean on time scales much shorter than those of the climate signal under consideration (here the Barents, North and Celtic Seas) and those nearly enclosed basins (Black Sea, Baltic Sea) where the timescales of oceanic-exchange are similar to or much longer than the climatic time scales under consideration. Similarly there are distinct variations in terms of riverine input, with the northwest Black Sea shelf and southern North Sea showing high levels of riverine inputs, and the rest comparatively moderate or low levels (Artioli et al., 2008). A reduction in upper ocean nutrient levels as described above would naturally be expected to impact on adjoining shelf/regional seas that receive a significant fraction of their nutrients from surface oceanic waters, as is the case for downwelling shelves such as the Northwest European continental shelf (Hydes et al., 2004;Vermaat et al., 2008).
The impact of changes in oceanic nutrient transport has been seen in future scenario simulations (Gröger et al., 2013;Holt et al., 2012). However, this impact is highly uncertain, depending as it does on the fidelity of the global ocean biogeochemical model, the details of the vertical mixing (Sinha et al., 2010;Steinacher et al., 2010) and the ocean-shelf transport processes.
Once these background budgets are set, seasonality (see examples in Figure 3) is a primary determinant of the resulting primary production. The seasonality of phytoplankton growth is dependent on the interplay of light and nutrients, and how physical conditions control the phytoplankton's exposure to these. These physical processes particularly relate to vertical transport (mixing, upwelling), which in turn is controlled by stratification and the driving forces of wind, tide and buoyancy. This seasonality has been characterised by Longhurst (1995), and updated with more recent observational (D'Ortenzio et al., 2012) and statistical approaches (Reygondeau et al., 2013). In mid-to high latitudes the triggers of spring and autumn blooms and their relation to mixing are reasonably well established, (Chiswell, 2011;Huisman et al., 1999;Sverdrup, 1953;Taylor and Ferrari, 2011). Seasonally ice covered regions like the north-eastern Barents Sea and the northern Baltic Sea form additional provinces within these regimes, which are characterized by an extended light limited period and reduced productivity compared to ice free regions at similar latitudes.
Since a decrease in sea ice is a very robust feature in future climate projections (Overland and Wang, 2007), increasing productivity is expected to be a first order response to future climate change in seasonally ice covered regions. The potential strength of this change, however, is modulated by the local biogeochemical conditions. In highly stratified, more tropical regional seas, much is still uncertain about this seasonality (as discussed below in relation to the Black Sea), owing to a less predictable and weaker seasonal signal, and often a paucity of observations. This overview of forcing and response sets the background for the specific results from the five regional models, which are now considered.

The models
Three model systems are considered here (Table 1 and Appendix 1): POLCOMS-ERSEM (northwest European Continental shelf; NWS), ECOSMO (North Sea-Baltic Sea and Barents Sea), and BIMS (Black Sea), which use generally similar approaches. The hydrodynamic models use regular, quadrilateral grids in the horizontal and either geopotential levels (ECOSMO) or terrain following coordinates (POLCOMS, POM) in the vertical. The primitive equations of motion and of tracer (temperature, salinity and ecosystem variables) evolution are solved on these with a time-stepping approach and a dynamic (1 or 2-equation) turbulence closure scheme to describe vertical mixing. The physics models are coupled (through transport, mixing, temperature and light climate) to lower trophic level ecosystem models that divide the ecosystem into several nutrient, producer and consumer boxes, and cycle one or more elements among these. The models differ in details of the numerical solution of the equations of motion and how the ecosystem is partitioned. Each has been extensively validated against regional observations using both reanalysis and climate model forced simulations. References for the models are provided in Table 1 and more details are provided in Appendix 1. In the case of the northwest European Continental shelf and Black Sea some aspects of the climate change simulation considered here have been previously published (Cannaby et al., 2015;Chust et al., 2014;Holt et al., 2012).

Experiment design
There are two classes of dynamically downscaled experiments: transient and timeslice simulations. The former drives the downscaled model with lateral and surface boundary conditions taken from the OAGCM, starting from the present-day and running for typically many decades into the future (e.g. Olbert et al., 2012;Tinker et al., In review). After an initial adjustment period during which the model evolves from the initial conditions (the 'spinup'), the simulation can be analysed for the full range of variability and trends. In the timeslice approach, as used in this and many other studies (e.g. Adlandsvik, 2008;Holt et al., 2010) the model is driven by both future and present day conditions in two separate experiments. After a spinup period in both, the two can be compared for statistically significant differences. This approach is significantly more flexible than the transient approach, with several options for manipulating the driving data, for example: climate delta approaches, reconstructions and bias corrections (MEECE, 2011).
The timeslice approach comes with two specific issues. First is the adjustment to future conditions: the approach assumes that, after the spinup period, the model is not sensitive to the initial conditions, or else this sensitivity will manifest itself as a false climate change response. This can, to some extent, be ameliorated by deriving future initial conditions from the driving model. However, the OAGCMs rarely have a good representation of the benthic system, so the issue of benthic spinup and adjustment is an ongoing concern with the timeslice approach. The second issue is the relationship between the longer period modes of natural variability and the difference between the timeslices can be difficult to assess, since only a 'snap shot' of the variability is available from the simulations. Due to the length of these simulations some contamination by natural variability is always a possibility. While the full spectrum of variability is not available in these simulations, we can go some way to assessing this by testing the statistical significance of the difference between the future and present time slice compared with the inter-annual variability. Here we use a Kruskal-Wallis test (Kruskal and Wallis, 1952) to test whether the annual mean netPP in the two timeslices are sampled from distributions with the same median value.
For the experiments considered here, future climate forcing is provided by the IPSL-CM4 OAGCM (Marti et al., 2006), chosen because of ready access to high frequency forcing data that includes a biogeochemical component. The quality of the regional model output will always be limited by that of the driving model and here the direct use of a global ESM will have some bearing on the results, particularly on the wind field. There are approaches that may mitigate this, such as regionally downscaled atmospheric models, which have been shown to improve wind fields in coastal regions (Feser et al., 2011). These were not available for this work, but a statistical downscaling approach is used in the Black Sea case. Here each model uses a timeslice approach to compare potential conditions at the end of the century (2080-2099, identified as A1B) under the SRES A1B 'business as usual' scenario (Nakicenovic and Swart, 2000) with typical present day conditions . The POLCOMS-ERSEM and BIMS simulations include changes to the riverine freshwater forcing (simply scaled by the local precipitation changes), the ECOSMO models do not. The riverine nutrient concentrations are not changed between the timeslices. Two approaches are employed: a direct forcing approach in which future and present conditions are both taken from the OAGCM or a delta change approach whereby a reanalysis forced simulation is used for the present conditions and a mean monthly change is imposed on this for the future conditions. Further details are provided in Appendix 1. Figure 4 provides some examples from the timeslices of this driving data. It shows distinct spatial patterns relevant to the regional seas under consideration here. Notably, an increase in air temperature that is much stronger in Arctic regions and, to a lesser extent, over the European continent than over the Northeast Atlantic. The change in wind stress shows a very varied picture, with a notable increase in the Arctic and Norwegian coast, otherwise a decrease. The shortwave radiation shows a modest increase except in the Arctic and western Norwegian sea, where it strongly decreases. This figure also emphasises how few grid cells of the atmospheric model cover each region.

Analysis and sensitivity studies
A particular issue faced when considering multiple regional simulations is the wealth of possible analyses and model diagnostics available from the model data sets to meet the requirements of each region. To compare different models and regions, a common approach is required, which is necessarily more simplistic than would be considered for a single model or region. Here we use a simple heuristic framework to provide a comparison between the models' response to climate change, accepting that this combines regional differences in response with differences in internal model dynamics, and is not equally well suited to all the regions considered. We divide the annual phytoplankton growth into three periods by identifying: bloom start time (T1) : netPP is greater than 20% of the annual maximum (following Platt et al., 2009); bloom stop time (T2) when surface nitrate is less than 20% of winter maximum; and growth end time (T3) when netPP is less than 20% of the annual maximum. The three periods (see Figure 3) are then pre-bloom (PB; T3 to T1), spring bloom (SB; T1 to T2) and summer growth (SG; T2 to T3). This approach aims to separate out the production that occurs during the transition from winter to summer mixing and light conditions, from that occurring during the summer conditions (e.g. seasonally stratified or shallower mixed layers). Hence, T2 is chosen to mark the time when seasonal surface nutrient depletion becomes apparent. Because the seasonal cycle is much less marked in the Black Sea, these thresholds are set to 60% for this region.
By exploring properties averaged over these three stages the response to climate forcing can be seen in more detail, distilling a substantial volume of information into a tractable form. The factors we consider that may affect the net primary production are changes in the length and timing of the growing season, the wind stress and the potential energy anomaly (PEA, an integral measure of stratification; Simpson and Bowers, 1981). The PEA is defined by: where h is the water depth (here the integration is limited to 200m), g is the gravitational acceleration, ρ the density and z the (positive up-wards) vertical coordinate. An overbar indicates an average over the same depth as the integration. This represents the energy (per depth) needed to mix the water column. It is a commonly used metric for stratification as it is an integral quantity that does not relate to a particular vertical structure or threshold, and connects with simple theories of stratification evolution (Sharples and Simpson, 1996;Simpson and Hunter, 1974).
These properties are divided temporally between these three stages and spatially between regions showing positive and negative overall netPP change ( Figure 5). Tables 2 shows values for the CNTRL experiment and relative or absolute change as appropriate. To explore the relation between netPP changes and PEA, wind stress and length of growing season, we also consider scatter plots of these quantities for each region (Figures 6,7,8).
To explore the relative importance of various external drivers, we present a further series of POLCOMS-ERSEM driver sensitivity experiments. These experiments provide an illustration of the potential range of impacts from key drivers of shelf sea production.
However, the importance of an individual external driver will vary depending on the region and, to some extent, the model system in question. Unfortunately, similar simulations are not available for the other model systems considered here. The approach we adopt is to start with the self-consistent set of future forcing (from the IPSL-CM4 model) and systematically remove aspects of the climate change signal by running the POLCOMS-ERSEM model with individual future forcing variables replaced by present day values, taken from a random year in the control period (1981)(1982)(1983)(1984)(1985)(1986)(1987)(1988)(1989)(1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999)(2000). The difference between the standard run and this experiment then quantifies the climate change effect attributable to this forcing variable. This differs from a more conventional sensitivity experiment (e.g. Skogen et al., 2011) in the treatment of non-linearities. Here we identify the effects of an aspect of the forcing and all non-linearities associated with it under future conditions. We consider five forcing experiments, each an 18-year simulation following 5-years of spinup, for: wind (W), short wave radiation (L), air temperature and relative humidity (A), boundary nutrients (B) and precipitation (P). We also consider an experiment to investigate the temperature dependence of the ecosystem model (T). This experiment is similar to that used in a global context by Taucher and Oschlies (2011); the change due to physiological/chemical temperature effects are quantified using a pair of timeslice experiments with the ecosystem model rates fixed to their values at 10 o C.

Contrasting the climate change impacts on primary production in five regional models
Here we present results from the five regional models. First an overview of the change in primary production and a common analysis across the five regions is presented. This is followed by a regionally specific analysis. Figure 1 shows the average net primary production (netPP) for the five models and the driving OAGCM. The difference in detail between the regional models and the global model shows the potential for 'added value' by regional models, i.e. they provide information that is not available from the global model alone. Quantifying the extent to which this is achieved requires a more detailed, case by case assessment of the mean state and variability than can be included here. However, these models have each been assessed in reanalysis forced simulations and their control simulations validated, in a mean sense, against observations (Table 1), showing that the general spatial patterns such as these are well represented and the method of forcing used here does not substantially compromise the simulations. It is important to note that the POLCOMS-ERSEM NWS and ECOSMO North Sea models show good qualitative agreement where they overlap, with POLCOMS-ERSEM giving higher values than ECOSMO. Differences such as this are common in comparisons of ecosystem models of primary production (Lenhart et al., 2010;Skogen and Moll, 2000) since the netPP values are sensitive to details of the cycling of nutrients through the various model compartments and to the treatment of the light climate. Moreover the observational constraints are weak, as discussed in the context of the NW European shelf by Holt et al (2012). A noticeable difference between these two models is in the Norwegian Trench and Skagerrak regions where the POLCOMS-ERSEM model shows a minimum in production not seen in ECOSMO. This may well be due to an inadequate treatment of Baltic inflow in the POLCOMS-ERSEM model. In contrast ECOSMO is a coupled Baltic-North Sea model, so this is well modelled.

Change in net primary production: overview and common analysis
The overall climate change effect on netPP is shown in Figure 5 as the absolute difference in mean values between the two timeslices. The regional model results in Figure 5 are shaded to show whether the climate change signal is statistically significant, at the 90% level, compared with inter-annual variability. We see this is the case for most of the Baltic Sea and the Black Sea, except in the latter case in the central transition zone between regions with positive and negative change. In the North Sea (for both POLCOMS and ECOSMO), in the region with negative change, that change is generally statistically significant, whereas in the region with positive change it is not. This is also the case in the Celtic Seas (in POLCOMS), except for a region of the central English Channel and eastern Irish Sea. A lack of statistically significant change in a region implies that either the drivers are weak or that multiple drivers are compensating. This is explored for the NWS in more detail below.
The global model shows a general decrease in primary production at mid latitudes and an increase at high latitudes, reflecting the regimes suggested by Steinacher et al (2010).
However, it should be born in mind that this general picture overlies a high degree of interregion and inter model variability in the individual global models. The regional models each show a mixture of both positive and negative change at lower latitudes, where this particular global model is uniformly negative. This arises because the regional models are able to produce a more detailed process response than the global, with multiple processes acting depending on the detailed regional conditions. Particularly apparent is the absence of regions of positive change at lower latitudes in the global model. Whereas the regional models all show regions of positive change. Hence, the global model cannot simply be seen as providing a spatial average of the regions under consideration; they differ qualitatively from the global model on the scales being considered here. Table 2  Explaining these differences between regions and sub-regions in terms of the processes that mediate biophysical interaction is a central focus of this contribution.
The mean annual cycle of netPP averaged over sub-regions for each model (Figure 3) shows the seasonal cycle and its change between CNTRL and A1B. The details of the seasonality underlying this spatial response ( Figure 5) show a temporal picture that is similarly diverse. All regions except the Black Sea show a distinct growing season, modulated by solar forcing. This shows small overall variations (Table 2) except in the SNE Atlantic, which shows a substantial increase in growing season (according to this definition), as the maximum netPP is significantly reduced and growth starts earlier. The seasonality of the Black Sea is a special case that is considered in more detail below.
These seasonal cycles provide a guide to partitioning the growing season as described above and we now considered the change in netPP divided into the three stages and between Together these results suggest a general decoupling of the netPP response from the stratification response. In the spring bloom case stratification relaxing light limitation may be dominating over its effect on nutrient resupply in many cases. While in the summer case stratification clearly has an effect in regions where the netPP is decresing, This reflects a basic response to increased seasonal stratification (Figure 6) of reducing the thermocline production.
The production in the pre-bloom (PB) period, and its change, is small in all cases in shelf seas, partly reflecting the construction of the stages, but also showing that the 'tails' of the annual cycle are not greatly changed (Figure 3), except in the Black Sea (NW shelf) and Barents Sea (NE). Hence the general increase in PEA during the PB period does not affect the netPP during this period, accepting that the netPP can also change with the timing of the stages. In contrast the Black Sea, Baltic Sea and Barents Sea show a substantial increase in PB netPP, while the SNE North Atlantic shows a substantial decrease in PB netPP.
In terms of wind stress effects (Table 2; Figure 7) we might expect a decrease in wind in the SB to lead to an earlier bloom and increased production. In all regions except the Baltic, SNE Atlantic and Black Sea, the regions with positive Δ netPP experience a smaller increased (or a decrease) in wind stress than the corresponding negative region. During the summer growth period, an increase in wind speed might be expected to increase netPP. here the picture is less clear.
We now consider a more detailed exploration of these results with a brief analysis of each region. We group together the tidally active shelves, the isolated deep basins, and the Arctic sea we consider. They are contrasted with the POLCOMS-ERSEM behaviour in the NE Atlantic. This shows a substantial increase in PEA and a consequent reduction in netPP, particularly in the subtropical gyre (SNE), generally following the results of the global ESMs.

Tidally active shelf seas: Celtic Seas and North Sea
The northwest European Continental Shelf is characterised by being well mixed for, at least, several months during the year and so (except for the Norwegian Trench) lacks permanent stratification. Seasonal stratification exists across most of the shelf, except for the shallower regions of the southern North Sea, English Channel, and central and eastern Irish Sea, which remain generally well mixed throughout the year, although experience episodic, weak stratification (see Holt and Umlauf, 2008;Pingree and Griffiths, 1977). The seasonally The lack of permanent stratification in these regions implies the first order open-ocean response identified above is absent (as shown schematically in Figure 2). In this case changes in the seasonality of the heat flux and wind forcing are the primary effects changing shelf sea seasonal stratification. Considering the seasonally stratified case as a simple two layer system with only weak diapycnal mixing, the summer stratification is essentially set by the difference between spring and summer temperatures, and so relative changes in these under climate warming conditions will determine the changes in thermal stratification: summer warming faster (slower) than winters leads to an increase (decrease) in stratification. This is illustrated in Figure 4, which shows that in this case summer air temperatures are warming faster than winter, and so the tendency is to increase seasonal stratification, and this is seen in the generally increasing PEA during the summer growth period. This change in stratification is a detail of the forcing that is not necessarily robust and a different OAGCM may equally give a change of opposite sign. However, the non-linear equation of state (EoS; here we use UNESCO (1981)) provides for an increase in density stratification that is robust to changes in OAGCM: at warmer temperatures the same temperature stratification leads to greater density stratification. For example, in the weakly stratified cases during spring bloom initiation or intermittent production in well-mixed waters, a 1 o C surface to bottom temperature difference corresponds to a 11.4% larger density difference at an SST of 12 o than at 10 o C. This means the same variations in tidal and wind mixing conditions are more effective at initiating stratification (and hence blooms) under increased temperatures. This is sufficient to account for the small but consistent hastening of spring bloom times seen in Table 2, Figure 8. The effect enhances the feedback in the stratification initiation process and so would increase seasonal stratification changes above that expected from the seasonality of the forcing alone. This effect has not been fully explored but model simulations of future scenarios have tended to show an increase in stratification (as here and Adlandsvik, 2008;Holt et al., 2010).
The efficiency of a bloom, in terms of total primary production achieved, in seasonally stratified conditions is dependent on how much of the available nutrients can be used before stratification starts limiting the vertical diffusive flux needed to replenish the nutrients in the euphotic zone (shown schematically by the shaded area in Figure 2). A bloom that starts earlier than the on-set of stratification sufficiently strong to inhibit this vertical flux would be expected to be more efficient. Hence, this efficiency would be sensitive to the details of the stratification and mixing conditions during the bloom, and how they change. In this case stratification can act both to alleviate light limitation at the start of the bloom and introduce nutrient limitation if it is sustained. In the well-mixed cases more frequent and stronger intermittent stratification will initiate more blooms throughout the growing season, these will each last until mixing increases again and light limitation is reintroduced.
To help explain the regional difference in netPP response in terms of changes in forcing we turn to the forcing-response sensitivity experiments. The change in net primary production associated with each of these experiments is shown in Figure 9. These results are shaded where the change in netPP attributed to a particular driver is less than 90% significant, compared to the inter-annual variability in the total netPP, noting that small changes of multiple drivers may still add to give a significant change. They show boundary nutrients have a predominantly negative effect in this case, as identified by Holt et al (2012), and that these impact the whole region with large areas of statistically significant change across the Similarly precipitation shows a generally weak, positive effect, suggesting a reduction in salinity stratification, but is generally not significant. This is at odds with the ensemble mean picture given by Capotondi et al (2012), but simply reflects the single choice of forcing model.
The direct temperature effect on netPP through growth rates in this model is a modest increase; it is only statistically significant in some small isolated regions, such as close to the coast in the German Bight. This effect results from a more rapid growth during the spring bloom, and enhanced recycling, countered by enhanced grazing. This is accompanied by a more substantial reduction in annual mean phytoplankton biomass (not shown). Here, however, we find the air temperature response is mainly due to the various stratification effects noted above.
Hence, to describe the regional differences in the context of the analysis in section 3.1, we first note all regions (both positive and negative Δ netPP) in the POLCOMS-ERSEM model show an earlier and increased spring bloom, in contrast to the spring blooms in ECOSMO in the North Sea, which follows the sign of the overall netPP change (but this model has a less rapid spring increase in netPP; Figure 3). In the North Sea the negative Δ netPP regions are predominantly seasonally stratified. In POLCOMS-ERSEM, these experienced earlier and more efficient spring blooms (Table 2; Figure 8), but this is counted by reduced oceanic nutrient inflow and increased summer stratification. The positive Δ netPP regions of the North Sea are predominantly well mixed. Here the only air temperature effects available are the equation-of-state effect on intermittent stratification and recycling effects (generally small); changes to seasonality of stratification are not relevant here, but the PEA in this nominally well-mixed regions increases by 4-18% (but remains small). So we attribute the net increase in netPP here to the combination of intermittent stratification being more frequent and effective, alongside wind, light and recycling effects, all overcoming the negative effect of reduced ocean nutrient input.
In the Celtic Seas, the region of negative Δ netPP change is around the shelf break and the change arises from the same factors as in the stratified North Sea. The region of positive Δ netPP includes both seasonally stratified and well mixed waters, making this analysis a little less clear. The spring blooms are earlier and significant more efficient, and there is a notable increase in the pre-bloom PEA (38%). Summer PEA increases by ~13% and we can assume the seasonally stratified regions are subject to the same netPP reduction as in these areas of the North Sea. However, this is smaller than the increase during the spring bloom, leading to these areas being classified at net positive Δ netPP. Similarly this netPP reduction in the summer growth period is counted by continued increase in the unstratified region (due to the same processes identified in the well-mixed North Sea case) to give a negligible area mean Δ netPP for this region (Table 2).
Hence the picture in the NWS is of a large scale negative effect arising from the reduction in the transport of oceanic nutrients onto the shelf augmented by the effects of increasing stratification reducing the summer growth in seasonally stratified areas. This is mitigated by several processes inducing a positive Δ netPP. These are predominantly earlier spring blooms and enhanced bloom efficiency (in the stratified Celtic Sea this is sufficient to counter the negative effects, while in the stratified North Sea it is not) and increased weak intermittent stratification in nominally well mixed areas, relaxing light limitation, making blooms more frequent and productive. These effects are combined with more minor, positive, wind, light and recycling effects. Hence, the regions of positive change that are not statistically significant in Figure 5 are not necessarily areas of weak impact, but more likely to be areas where multiple drivers compensate leaving a weak change according to this metric. We might then expect changes to be more apparent in other aspects of the ecosystem, e.g. the functional group composition, where the compensation is less effective.
To give a first assessment of the robustness of these results, Figure 10 shows results from a small (4-member) ensemble in the North Sea. The two simulations described here are augmented with two further POLCOMS-ERSEM simulations, using an alternative OAGCM for forcing (HADCM3; Gordon et al., 2000;Pope et al., 2000) and an alternative forcing approach (the delta-change method, see appendix 1). The mean netPP change is shown along with an indication of how well the model's agree in sign. This is quantified at each grid cell by counting the number permutations of pairs of models that agree in sign (a value of 6 indicates all models agree in sign).This identifies three regions where all four models agree: the regions of negative change following the inflow of water from the Northeast Atlantic into the North Sea and around Dogger bank, and the region of positive change along the coast of continental Europe (although this is patchy near coast). The ensemble mean response is substantially weaker in most parts of the North Sea due to contrasting signals in the single ensemble members. We note that the model using an alternative OAGCM is the only one that shows a different sign of change (positive) in the central North Sea.

Isolated regional seas: Baltic Sea and Black Sea
The Baltic and Black Seas have restricted exchange with the open ocean and are both highly stratified basins, in part due to high riverine freshwater input. The consequent twolayer circulation brings denser, more saline water into both seas, from the Mediterranean (Black Sea) and North Sea-North Atlantic (Baltic Sea). In the deeper Black Sea (maximum depth > 2000 m) a permanent halocline at 150-200 m depth prevents deep winter convection and leads to anoxic conditions over much of the water mass. A similar halocline between 60-100 m is present in the Baltic Sea that also limits winter convection. The low salinities in the Baltic furthermore limit mixing throughout the winter and favour the development of a winter thermocline and ultimately sea ice (e.g. Backhaus, 1996;Schrum and Backhaus, 1999).
The seasonal characteristics of phytoplankton growth in the Black Sea are characterized by high inter-annual and multiannual variability, and comparatively weak seasonality (Vantrepotte and Melin, 2009). This is reflected in Figure 3 with the two deep basin regions showing weak seasonal variation. The NW shelf shows an apparent seasonal cycle, but also much intra-annual variability. This should not necessarily be interpreted as a genuine seasonal cycle, but rather a vestige of the strong inter-annual variability. Specifically, a change in bloom timing in the late 1980s following the introduction of Mnemiopsis leidyi (introduced in the model at the start of 1988), as well as inter-annual changes in river input that drive production on the shelf. Hence, a 20-year mean is insufficient to average-out this variability. The detailed picture of the climate change signal has, therefore to be treated with some caution here. The low seasonality arises because the shallow mixed layer depths ( ~5 m in summer to ~70-140 m in winter; Oguz, 2008) are not so great as to provide a strong constraint on phytoplankton growth during the winter months. Nezlin et al. (2002), in an analysis of the remotely sensed data record, suggests that, away from the northwest shelf, the Black Sea has a seasonal cycle of phytoplankton biomass akin to that typical of subtropical regions (e.g. Longhurst et al., 1995). They suggests that surface phytoplankton biomass peaks during September and October in the open Black Sea and production continues throughout the winter months, followed by a smaller bloom during spring with the onset of seasonal stratification Nezlin, 2008;Nezlin et al., 2002).
Spring blooms (e.g. as reported by Demidov, 2008;Oguz et al., 2001a) may be triggered by reduced mixing, but the sustained winter production limits the build up of surface nutrients and so also the bloom strength. These results generally concur with this picture. Surface nutrients are resupplied from depleted summer levels by convective and wind driven mixing in autumn and winter (Oguz et al., 1996), triggering the strong autumn bloom.
The picture provided by the common analysis is a complex one in the Black Sea, partly due to the approach not being so well suited to a region with such limited seasonality. The change in growing season length (Table 2)  nutrients within the surface mixed-layer and the regional scale response of the Black Sea to changes in the wind driven circulation, which in turn influences the distribution of Danube plume waters. This arises because of changes in the basin wide cyclonic circulation, particularly in the intense 'rim current', a cyclonic geostrophic current that approximately follows the 200m contour (Oguz et al., 1992). This mean circulation is driven by the predominantly positive wind stress curl (Korotaev et al., 2011), which is substantially modified in these experiments. This gives a response entirely different to open ocean regions at a similar latitude.
The seasonal characteristics of phytoplankton growth in the Baltic Sea are highly subbasin dependent. Different regions exhibit a combination of spring-bloom-, polar-bloom-, and eutrophic coastal system characteristics, but with a strong seasonal cycle owing to the high latitude (in contrast to the Black Sea) (e.g. Daewel and Schrum, 2013). Figure 3 shows substantially different seasonal behaviour between the Baltic Proper and the Gulf of Bothnia, which is dominated by seasonal ice cover. Both show a positive change, and the statistics calculated for Table 2 are limited to the Baltic Proper. In the Baltic Sea the thermocline is generally shallow (compared to the euphotic depth), so nutrients below the summer thermocline are generally accessible and production tends not to be limited by stratification.
As with the Black Sea, the substantial increase in stratification does not translate to a reduction in netPP. In fact changes in the Baltic Sea are uniformly positive and dominated by an increase during the SG period. Sediment water exchange, oxygen deficiency, and denitrification and nitrogen fixation are all particularly important for the Baltic Sea ecosystem and modulate phytoplankton production by modulating the N/P ratio in the water column, and thereby the limiting processes for phytoplankton production (Rodhe et al., 2006). Moreover, the retreat of sea ice, the maximum ice extent and the different length of the ice free period, which varies greatly locally, play key roles in structuring the seasonal phytoplankton dynamics in the Baltic Sea (discussed in more detail in 3.2.3).
The central Baltic Sea is thought to be forced by an interaction between bottom topography and buoyancy (Sarkisyan et al., 1975). It is characterised by a general cyclonic gyre (as in the Black Sea), while anticyclonic gyres can be found in the Bornholm Basin, the Gdansk Basin and north of the Gotland Basin (Lehmann et al. 2002). The circulation is highly variable and strongly dominated by changes in the large scale atmospheric forcing and associated changes in the regional wind field (Lehmann et al., 2002). During positive North Atlantic Oscillation (NAO + ) conditions, strong Ekman currents drive increased up and downwelling along the coasts, while under NAOconditions ventilation is strongly reduced (Lehmann and Myrberg, 2008). Upwelling in the Baltic Sea is highly relevant for ecosystem dynamics through replenishing depleted nutrients in the surface layer, and while it acts locally it affects the entire basin. Since Baltic Sea deep waters are often anoxic and nitrate depleted the mechanism changes the N/P ratio in the surface waters and hence favours the production of nitrogen fixing cyanobacteria (e.g. Daewel and Schrum, 2013;Janssen et al., 2004). The wind stress in this region overall decreases slightly in this future climate scenario but shows a significant increase in westerly winds during the autumn and winter, which leads to enhanced upwelling and increased winter nutrient levels for the subsequent growing season. These effects are explored in more detail by Pushpadas et al (2015).
Hence, it is apparent that alongside changes to the permanent stratification, these regions are highly susceptible to changes to the wind strength, gradients and direction, and consequent impacts on circulation and upwelling, in a way that is not applicable to the open ocean.

Seasonally ice covered Arctic Seas: Barents Sea and Baltic Sea
The Barents Sea is a large broad part of the Arctic shelves. It is strongly controlled by advection: inflowing warm, salty and nutrient rich Atlantic water is subject to water mass transformation and leaves the Barents Sea towards the Arctic Ocean as dense and cold water (e.g. Årthun et al., 2011). The north-eastern part of the Barents Sea is dominated by the presence of Arctic water and sea ice. A pronounced frontal system separates the seasonally stratified ice free and highly productive south-western region from the year round stratified, seasonally ice covered and less productive north-eastern region. The position of the polar front, and so the area occupied by these two regimes, is advectively controlled by the overall amount of incoming Atlantic Water. Hence this sea would be expected to be sensitive to changes in inflowing water in marked contrast to the 'well-mixed box' view of ocean transport influencing properties (Gordon et al., 1996;Holt et al., 2012).
The Atlantic water region of the Barent Sea is more nutrient rich and characterized by ice free conditions throughout the entire winter (e.g. Årthun and Schrum, 2010). This supports generally higher production compared with the coastal, northern and north-eastern parts of the Barents Sea. The latter regions are subject to a substantially reduced growing season due to the light limiting effect of sea ice (Wassmann et al., 2010). The minimum sea ice cover in the Barents Sea is not reached until September and large parts are ice covered when solar radiation is at its seasonal maximum. This mismatch between seasonal radiation and ice-free conditions reduces the growth period substantially in the eastern parts of the Barents Sea, (Figure 3). A similar effect is apparent in the Gulf of Bothnia (the northern most part of the Baltic Sea). Seasonal ice cover, in the current climate, occurs here from November to the beginning of May and causes light limitation during the spring bloom in wide parts of the Gulf of Bothnia and a considerable delay of the bloom compared to the Baltic Proper ( Figure 3). The major signal in the climate change projection in both seas is a significant increase in the earliest production in the previously ice-covered parts of the seas. This can be explained by the reduced sea ice development and earlier retreat of sea ice in large parts of the seas.
The sea ice extent in both seas exhibits strong inter-annual variability. The Baltic shows larger variations in the present day climate, with sea ice extending into the Baltic Proper and even the entire Baltic Sea being ice covered in severe winter conditions. The variations in Baltic Sea ice extent are largely controlled by variations in air temperature (e.g. Omstedt and Nyberg, 1996;Tinz, 1996) and show little sensitivity to wind changes (Schrum and Backhaus, 1999). In the future climate simulation the air temperature increase (Figure 4) can, therefore, be considered the major driver for the simulated decrease in Baltic Sea ice conditions (not shown here), which is in accordance to earlier climate change downscaling simulations for the Baltic (Meier, 2006). In contrast, the dynamics of sea ice variations in the Barents Sea are intimately connected with the wind-driven inflow of Atlantic Water and linked to changes in air temperature (Årthun and Schrum, 2010;Goosse and Holland, 2005).
Hence, the direct temperature effect is complemented by a pronounced impact of wind driven circulation changes. A warmer Barents Sea is generally characterised by more inflow of Atlantic Water, less sea ice and increased heat loss from the ocean (Årthun and Schrum, 2010). Earlier studies identified a connection with reduced air pressure and increased cyclonic circulation (A dlandsvik and Loeng, 1991;Bengtsson et al., 2004). Similar atmospheric changes also arise from the climate change scenario considered here ( Figure 4) and are responsible for decreasing sea ice and increasing inflow of Atlantic water. The latter is partly replacing the stratified polar water regime and is reflected in an overall substantial decrease in PEA in the Barents Sea ( Figure 6). Wind stress changes are particularly strong during the pre-bloom phase (Figure 7). The changes in netPP in the Barents Sea indicate two different regimes. The northern and north-eastern Barents Sea show an increase in production ( Figure 5), caused by the larger extent of nutrient -rich Atlantic Water mass, reduced sea ice cover (and hence reduced light limitation). The latter is also reflected by an earlier bloom start ( Table 2). The coastal southern Barents Sea shows a decrease in production ( Figure 3 and Figure 5), most like due to reduced nutrients in inflowing Atlantic water dominating over the effects of a reduction in sea ice. The distinction between the coastal and northern Barents Sea is also crudely present in the global model ( Figure 5), but detailed regional patterns, such as local extrema along Novaya Zemlya, Bear Island and Svalbard are missing.

Discussion and Conclusions
The impacts of climate change in regional seas are far from straightforward. A myriad of physical processes can potentially act as vectors transferring the larger scale oceanic, atmospheric and terrestrial variability and change to regional sea physics, biogeochemistry, lower trophic level ecosystems, and so up the foodweb. These processes act on a wide range of time scales and are strongly dependent on the prevailing conditions of an individual regional sea basin. Here we have explored some of the physical mechanisms driving this interaction, drawing on a set of regional model simulations to provide illustrations. Based on the discussion of forcing and responses provided above, we can hypothesize the causes of change in each region, whilst noting that conclusively attributing a particular causal relationship is very challenging in these non-linearly connected systems. This is summarised in Table 3, separating each model domain into the regions that shows positive and negative netPP change. While some regions show a relatively straightforward response, such as the change in nutrient inputs and increased stratification in the northern North Sea, other regions show multiple drivers, such as the complex combination of factors that lead to a projected increase in netPP in the Celtic Seas. The response of the Black Sea, similar to the central Baltic Sea, is an important special case, as unlike the other regions a change in the circulation and upwelling is seen as a primary driver of change here (Cannaby et al., 2015). For the Barents Sea and the northern Baltic Sea, sea ice changes are considered to be the primary driver of the changes in primary production.
It should be emphasised that each of these simulations only represent a single possible view of future conditions, and no quantitative assessment of the likelihood of this occurring has been made. Moreover, differences in modelling approaches are combined with the climate change signal. So in comparing regions, while the qualitative differences in forcing-response are robust, the quantitative differences are less so. The comparison of the four model realisations in the North Sea ( Figure 10) provides some assessment of this. Moreover, the timeslice approach always risks 'contamination' by long-term variability. We provide a first assessment of the significance of the change seen here, but are limited by the variability provided by these comparatively short simulations. Full transient simulations (Tinker et al., In review;Tinker et al., 2015) are required to properly assess the relation between climate change and variability in this context, and without these we are unable to accurately assess the biases resulting from the timeslice approach. Spin-up is also an issue with the timeslice approach, whereby the timeslices may still be adjusting to the initial conditions, this is ameliorated to some extent by using initial conditions from the driving OAGCM. However, these conditions will still have to adjust, over some time period, to the dynamics of the regional model. Again, this issue would largely be resolved by using a transient simulations.
A general conclusion of this work, in terms of experiment design, is that there are many benefits from transient simulations that outweigh their computational cost.
Considerable focus has been given in the global context to the impacts of changes in permanent stratification on vertical nutrient resupply (Capotondi et al., 2012) leading to a decrease in primary production (Bopp et al., 2013;Steinacher et al., 2010) as a general pattern that emerges from the regional variations in the global models. While, there has been some acknowledgements of direct temperature effects and changes to growing season, these have received less attention (Sarmiento et al., 2004;Taucher and Oschlies, 2011). These relate to less robust aspects of ecosystem models (in the case of direct temperature effects) and climate models, in the cases of growing season effects. In regional seas, however, these permanent stratification effects are less clear and we must consider all potential processes, with priority being dictated by regionally specific conditions. In seas that are shallow compared with turbulent boundary layer thicknesses and so are well mixed for part of the year this leading order effect is absent. Hence this property may shelter these seas from some direct impacts of climate change. Indeed they show a mixture of positive and negative response, which is generally weaker than in the open ocean. In the permanently stratified seas we consider here (the Black Sea and Baltic Sea) the impact of climate change on stratification is highly modulated by changes in circulation and overturning/mixing in its effect on nutrient resupply and primary production. These more dynamic processes are seen to be the leading order effects in these cases. Hence, the global view of climate change causing a decrease in netPP (except in high latitudes) must be reconsidered at a regional scale, particularly in the light of all of these regional models showing regions of increasing netPP (for various reasons). The interplay of the seasonal cycle and sea ice cover forms another important driver for the Arctic and sea ice covered regions, the Barents Sea and the northern Baltic Sea.
The view we have provided here only considers a subset of the possible range of processes mediating climate impacts in regional seas. This is due to limitations of model resolution excluding some process, only considering some regional seas (e.g. Eastern Boundary upwelling shelves have not been considered), and the depth of the analysis presented not being able to fully elucidated the role of some processes. We can see where the processes treated here fit in the broader picture, by speculating on a wider range of possible effects and identifying those we have considered here and those excluded. This is summarised in Figure 11, along with a hypothesised sign of change in netPP. This particularly identifies the climatic impacts on mesoscale (e.g. shelf break internal tides and fronts) and near coastal (e.g. river plumes) processes have been neglected by this study (and other similar ones). Other neglected processes include changes to optical properties (Lohmann and Wiltshire, 2012), effects of sea level rise on tides (Pickering et al., 2012) and changes to the nutrient concentrations in upwelling water (Rykaczewski and Dunne, 2010).
The resolution considered in these models (and those used for similar studies) is marginal for the consideration of many processes. This is particularly the case in near coastal, frontal and shelf-slope regions of the ocean margins; in fact where we would expect the primary production to be highest. Several physical processes (e.g. Rossby adjustment scales) show a ~h 0.5 relationship between horizontal scale and water depth, i.e. a ~10 fold decrease in scale from 4000m to 40m. Hence, a1/10 o regional model in shallow water is in some sense comparable to a 1 o global model. The models domains used here are well-established and so have not fully taken advantage of the continued growth in computer power over recent years, in terms of their grid resolution and there is now an opportunity to address these issues, accepting that increases in resolution must be tensioned against the need for multiple process experiments, longer simulations and ensembles. As far as the ecosystem models are concerned, the parameterisation of temperature effects is obvious area for potential development, e.g. in ERSEM all groups share the same parameter value. Hence, all aspects of the system are changed by an equal amount and tend to cancel in terms of the annual production, but not the biomass. A more sophisticated approach, for example that differentiates between autotrophic and heterotrophic processes (e.g. Wohlers et al., 2009), might be expected to give a larger direct temperature response.
In terms of forcing, changes in the large scale ocean circulation and details of the gyre structure potentially have important implications for the characteristics of the water transported on-shelf that is not well treated by the coarse resolution ocean model used for forcing here (e.g. Holt et al., 2014;Molinari et al., 2008). Moreover, improved atmospheric resolution either in the global model or through downscaling and fully coupled systems has the opportunity to represent the changes in drivers with more fidelity, particularly in wind driven effects such as those in the Black Sea, Baltic Sea and the North Atlantic storm track.
The enclosed basins' connections with the wider ocean are at very narrow straights, which require high resolution coupled basin models (Black Sea -Mediterranean and Baltic Sea -North Sea) and high resolution wind forcing, to accurate simulating the processes driving these exchanges, and so capture events such as the 'Major Baltic Inflows' (Gustafsson, 1997;Omstedt et al., 2004), and how these might change.
While the view of several competing processes acting with both positive and negative sign ( Figure 11) really needs to be considered on a case-by-case basis, some general principle can be identified. When there are multiple effects of different sign these will tend to mitigate the climate change impact, suggesting some regional seas will generally be less vulnerable to climate change effects than the open ocean. This can act both locally and spatially, i.e. advective and diffusive transport will tend to reduce effects across gradients of negative and positive impact. This will not be the case in enclosed regional seas, where a single dominant effect can have a large impact that is not mitigated by exchange with neighbouring regions of a different sign. Hence, we might expect the enclosed regional seas to be more highly impacted. This is indeed the case in these experiments for the Black sea and Baltic Sea in contrast to the response of the North Sea and Celtic Seas. In both the Black and Baltic Seas this impact is driven by changes in wind effects leading to changes in circulation patterns in the former and upwelling rates in the later, i.e. both highly dependent on the detailed conditions in the basin.
Another consequence of multiple, competing processes is that uncertainties are enhanced. Simplistically, uncertainties for uncorrelated processes add in quadrature even if the effects are of different sign and cancel, i.e. fractional uncertainty can substantially increase. This is compounded by the fact that many of the processes considered here relate to less well modelled, regionally specific, aspects of the OAGCM forcing such as details of the wind fields. Hence, the climate change signal in these is substantially less certain than, for example, changes in air temperature. Assessing and addressing this uncertainty remains an ongoing challenge for future work.

Appendix I Model descriptions and experiment design POLCOMS ERSEM; Northwest European Continental shelf
The Model ERSEM ( Figure A1) is a well established, generic lower-trophic level/biogeochemical cycling model. Eight plankton functional types are represented, including phyto-, zooplankton and bacteria, along with the cycling of C, N, P, Si through pelagic (Blackford et al., 2004) and benthic (Blackford, 1997) ecosystems; the latter being critical for nutrient cycling in shelf seas. The model equations can be found in these two papers. The implementation of ERSEM considered here essentially matches that described in Blackford et al. (2004) with the treatment of abiotic (SPM, CDOM) absorption described by Wakelin et al. (2012). The parameter set matches that used by Blackford et al. (2004), except here we limit the carbon to chlorophyll ratio to better match observations Geider et al., 1997). In the current model we also included a resuspension flux of particulate organic material driven by tidal and residual bottom currents, following Wakelin et al. (2012).
POLCOMS (Holt and James 2001) is a three-dimensional hydrodynamic model using a quasi finite-volume approach, discretised on a B-grid in spherical-polar-terrain following coordinates. The Atlantic Margin Model configuration considered here has a resolution of 1/9 o latitude by 1/6 o longitude grid (~12 km) with 42 s-coordinate levels (Song and Haidvogel, 1994) in the vertical. This configuration is further described by Wakelin et al (2009).

Experiment Design
Here we adopt a direct forced time-slice approach whereby mean conditions in a future period are compared with mean conditions in a present day reference period to give a measure of the climate change signal, on the assumption that conditions in both time-slices are approximately stationary. The CNTRL simulation is a 23-year present day simulation forced by the IPSL-CM4 20C model for the nominal present day period 1980-1999 (1980 is repeated three times before this period is simulated). The OAGCM provides 6-hourly atmospheric air temperature, winds, pressure and relative humidity, and daily precipitation and short-wave radiation (the latter is modulated by the diurnal cycle). Surface fluxes are calculated by COARE v3 bulk formulae (Fairall et al., 2003). Monthly ocean currents, sea level, temperature and salinity are provided at the open boundaries. Tidal forcing is taken from a North Atlantic Tidal model (providing 15 constituents;Flather, 1981). Elevation and depth mean current boundary conditions (tidal and 5 day mean residuals) are applied using a flux/radiation scheme. Lateral boundary conditions for ERSEM use monthly values from the World Ocean Atlas (WOA; Garcia et al., 2006) for nitrate, silicate and phosphate, imposed with an up-wind advection boundary condition. Other variables use a 'zero-gradient' boundary condition. An exception is the detrital organic material fluxes, which are set to zero inflow concentration to avoid numerical instability. For freshwater fluxes, daily discharge data for 250 rivers are used from the Global River Discharge Data Base (Vörösmarty et al., 2000) and from data prepared by the Centre for Ecology and Hydrology as used by Young and Holt (2007). River nutrient loading matches that used by Lenhart et al. (2010), with raw data for the UK, Northern Ireland, Ireland, France, Norway, Denmark and the Baltic processed by van Leeuwen (CEFAS, UK) and raw data for Germany and the Netherlands was processed by Pätsch and Lenhart (2004). The Baltic exchange at the Belts is treated crudely as an inflow source using a mean annual cycle of depth averaged transport, salinity and nutrients. A constant spatial field of atmospheric nitrogen deposition (oxidized and reduced) is provided by EMEP (Cooperative Programme for Monitoring and Evolution of the Longrange transmission of Air Pollutants in Europe).
A1B is a future climate scenario representative of possible conditions in 2080-2100 under a business as usual emissions scenario: SRES A1B. The forcing matches CNTRL using the same OAGCM simulation run forward to this period. This OAGCM simulation includes the PISCES ecosystem model (Aumont et al., 2003) and we perturb the open-boundary nutrient values (nitrate, silicate and phosphate) by the fractional change between this timeslice and CNTRL. River flows are perturbed by changes in regional precipitation from the OAGCM, but are not related to changes over specific catchments. Riverine nutrient loads and atmospheric inputs are unchanged.
In all experiments we treat the first 5 years as 'spin-up' to allow the model to adjust to its lateral boundary and surface forcing conditions, so the results presented here are means for 18 years for CNTRL and A1B.

ECOSMO: North Sea -Baltic Sea, Barents Sea
The Model ECOSMO is a coupled physical-biogeochemical model, with the hydrodynamics based on the HAMSOM (HAMburg Shelf Ocean Model, Schrum and Backhaus, 1999), a 3D baroclinic coupled sea-ice model. The ecosystem component ( Figure A1) solves for 16 state variables that are divided into 3 phytoplankton functional groups (including cyanobacteria), 2 zooplankton functional groups, 3 detritus (including DOM), 3 sediment and 4 nutrient groups plus oxygen. Primary production in ECOSMO is limited by either nutrients or light, and the 3 major nutrient cycles (nitrogen, phosphorus and silicate), important for simulating nutrient limitation processes in the Nordic marine ecosystems, are resolved.
The online-coupled model system is applied to the combined North Sea -Baltic Sea system, using a spherical grid with a horizontal resolution of 6´ x 10´ and 20 vertical z-levels.
It has a free surface and allows for variable thickness in the last model layer, thereby resolving a realistic bathymetry. A second order Total Variation Diminishing (TVD) advection scheme is used. The non-diffusive properties of this scheme are especially important for preserving the permanent halocline in the Baltic Sea, providing better representation of water mass characteristics and nutrient distribution in the Baltic Sea.
The model system is also applied to the Barents Sea, here with a horizontal resolution of 7km and 16 vertical levels (Schrum et al., 2005;Årthun et al., 2011). The Barents Sea model setup matches the North Sea-Baltic Sea setup, except in this case the ECOSMO model is used offline in the Barents Sea: the ecosystem model uses daily mean fields from the hydrodynamic model.

Experiment Design
The reference simulations are hindcast simulations using atmospheric boundary conditions as provided by the NCEP/NCAR re-analysis (Kalnay et al. 1996), for the simulation period 1948-2008 (for the North Sea-Baltic Sea: Daewel & Schrum, 2013, for the Barents Sea physics: Årthun et al., 2011). Additionally, monthly means of land-based freshwater runoff and nutrient loads are included. For the Baltic Sea runoff, two available datasets (1950-1970: Bergström and Carlsson 19941971 (Lammers and Shiklomanov, 2000). At the open boundaries to the North Atlantic Ocean sea surface variations are prescribed from a coarser diagnostic model (North Sea) (Backhaus 1987) and from coarser North Atlantic model (Barents Sea) (Sandø et al., 2010 [1980][1981][1982][1983][1984][1985][1986][1987][1988][1989][1990][1991][1992][1993][1994][1995][1996][1997][1998][1999], and future scenario period: 2080-2099 are considered. The Delta change method is applied to both the atmospheric and oceanic boundary conditions (including nutrients). For the Barents Sea, a delta change for ocean temperature, salinity and sea level is not applied since these are extremely biased in the IPSL-CM4 model in the present day by unrealistic sea ice cover in the Barents Sea. The delta change method is applied to the atmospheric and oceanic nutrient boundary conditions. River runoff and nutrient loads were not available from the OAGCM and are kept unchanged for the climate change experiments.

BIMS: Black Sea
The Model BIMS (Black Sea Integrated Modelling System) is a physical-biogeochemical modelling system with the hydrodynamics based on the POM2K (Princeton Ocean Model) using the Mellor-Yamada level 2.5 turbulence parameterization (Blumberg and Mellor, 1987). The The ecosystem component, BIMS_ECO ( Figure A1) is based on the 1D model developed by Oguz et al. (2001b) and has 12 state variables consisting of two phytoplankton groups, micro-and mesozooplankton, bacterioplankton, labile pelagic detritus, nitrogen, nitrate, ammonium, the opportunistic heterotrophic dinoflagellate Noctiluca scintillans, and the gelatinous carnivores; Aurelia aurita and Mnemiopsis leidyi. M. leidyi is introduced into the model at the start of 1988.

Experiments design
A direct forcing approach, similar to the POLCOMS-ERSEM simulations, is adopted here. Surface atmospheric forcing are prescribed using 6-hourly fields of wind, fresh water fluxes (evaporation, convective precipitation and large-scale precipitation), and radiation fields (surface shortwave radiation, surface long-wave radiation, evaporative heat flux, and convective heat flux). A statistical downscaling approach is applied to the wind field to correct the bias in this (see Cannaby et al., 2015). Exchanges through the Bosphorous and Kerch strait and inflows of the eight largest rivers flowing into the Black Sea are defined by monthly mean volume flux, temperature, and nitrate flux climatologies (Ludwig et al., 2010;Ludwig et al., 2009). Nitrate flux from sediment suspension was parameterized on the North-Western Shelf at a constant rate throughout all simulations. For the future simulation (A1B) the volume flux of river water entering the basin was scaled by the difference in precipitation between the CNTRL and A1B simulation periods. Initial conditions are defined using monthly mean climatology for January from the WOA, since the IPSL-CM4 model has an unrealistic representation of conditions in the Black Sea. Both time slice simulations are 20 years long (CNTRL: 1980(CNTRL: -2000(CNTRL: and A1B: 2080(CNTRL: -2100 and augmented by a 5 year spin-up period to allow the biogeochemistry and ecosystem to adjust to the prevalent conditions. Table 1: List of models, references, and basic experiment properties Table 2: Summary of response from the common analysis. At each grid cell the model values are integrated temporally over the periods described in section 3. These are then averaged spatially over the areas showing an overall increase and an overall decreased in netPP. Variables shown are: percentage area of positive and negative change; netPP and fraction change between timeslices; growing season and change in days; netPP, PEA and wind stress for the three periods, positive and negative Δ netPP areas and whole region. Also shown is the change in timing of the three stages.  Table 3 Sign of change of key forcing variables from IPSL-CM4 (A1B compared with CNTRL), and hypothesis for cause of positive (+ve) and negative (-ve) change in netPP. Figure 1 Average annual net primary production (netPP) for the global OAGCM and the five regional models over a nominal 20-year present-day period (1981)(1982)(1983)(1984)(1985)(1986)(1987)(1988)(1989)(1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999)(2000). A: global; B: Barents Sea (ECOSMO); C: Northwest European Shelf (POLCOM-ERSEM); D: North Sea (ECOSMO); E: Baltic Sea (ECOSMO); F: Black Sea (BIMS_ECO).

Figure 2
Conceptual view of the difference in seasonal cycles between A: open-ocean/deep regional sea and B: seasonally stratified sea. 'Leakage' primarily reflects the long-term overturning circulation. As the mixed layer deepens from summer values (h s ) to winter values (h w ) a quantity of nutrients ~(h w -h s ) is entrained into the surface layer to provide the basis for the following season's production

Figure 3
Mean annual cycle of depth mean netPP from the five regional models. Solid lines are from the present day timeslice (CNTRL), dashed from future timeslice (A1B). Timeseries are for area mean over subregions of each model domain. The first panel also illustrates the three stages used in the common analysis.

Figure 4
Mean values from the atmospheric component of IPSL-CM4 for 1980 and change between mean for 2080-2100 (A1B) and this, absolute difference for air temperature and otherwise fraction change (A1B/CNTRL-1). Variables shown are 2m air temperature (Ta (2m) , °C) , surface wind stress (τ, m 2 s -2 ) and surface short wave radiation (SWR, Wm -2 ). Also shown is the difference between Summer (July-August) and Spring (March-April) temperatures and the change in wind stress for these periods.

Figure 5
Change in mean netPP between A1B and CNTRL for the five regional models, and the global OAGCM. All the regional models are forced by the IPSL-CM4 model and use a common timeslice, but the forcing methodology differs between models. The regional model results are shaded where the difference between the two timeslices is less than 90% significant given the inter-annual variability.

Figure 9
Forcing experiments with POLCOMS-ERSEM. Fractional change in netPP associated with five external drivers and their non-linear interactions: boundary nutrients (B); wind (W); short wave radiation (L); air temperature and relative humidity (A); precipitation (P), and the direct effects of temperature on growth rates (T). The model results are shaded where the change in netPP associated with each driver is less than 90% significant given the inter-annual variability of the total netPP.

Figure 11
Summary of physical processes mediating climate change impacts. These are organised according to the three paradigms of biophysical interaction identified in section 1, general process and whether the sign of change in netPP is expected to be positive or negative or either depending on the sign of forcing. This is identified as unknown if the sign of the effect is not straightforward given the sign of forcing. Processes in red are not specifically or satisfactorily considered here.

Figure A1
Schematics of the three ecosystem models used in this study. Highlights: • Barents, Black, Baltic, North and Celtic Seas, and the Northeast Atlantic are considered • The link between stratification and primary production climate impacts seen in much of the global ocean, is not so apparent here. • Many physic processes are identified as mediating climate impacts on primary production • Primary production is seen to increase in many coastal regions and deep-isolated basins