Looking beyond stratification : a model-based analysis of the biological drivers of oxygen deficiency in the North Sea

1University of Hamburg, Department of Informatics, Scientific Computing, Bundesstraße 45a, 20146 Hamburg, Germany 2Centre for Environment, Fisheries and Aquaculture Science (Cefas), Lowestoft, Suffolk, NR33 0HT, UK 3University of East Anglia, School of Environmental Sciences, Norwich, NR4 7TJ, UK 4University of Hamburg, Institute for Hydrobiology and Fisheries Science, Olbersweg 24, 22767 Hamburg, Germany 5University of Hamburg, CEN, Institute of Oceanography, Bundesstraße 53, 20146 Hamburg, Germany 6Federal Maritime and Hydrographic Agency, Bernhard-Nocht-Straße 78, 20359 Hamburg, Germany 7CNRS, UMR 7144, Equipe Chimie Marine, Station Biologique de Roscoff, Place Georges Teissier, 29680, Roscoff, France 8Dalhousie University, Department of Oceanography, 1355 Oxford Street, Halifax, Canada


Introduction
Low oxygen (O 2 ) conditions (concentrations < 6 mg O 2 L −1 ; OSPAR-Commission, 2003), often referred to as O 2 deficiency, occur regularly in the North Sea.A major process regulating the seasonal dynamics of bottom O 2 is the occurrence and duration of thermal stratification (e.g., Greenwood et al., 2010;O'Boyle and Nolan, 2010), which limits the vertical exchange of O 2 between the oxygenated surface layer and the deeper layers.In combination with events of enhanced primary production, and the subsequent degradation of organic matter, this favours the evolution of O 2 deficiency (e.g., Diaz and Rosenberg, 2008).Although the northern North Sea reveals strongest stratification, lowest O 2 concentrations occur in the central North Sea, which is shallower and where the duration of stratification is shorter and shows highest year-toyear variability.Thus, one can argue that stratification is an important prerequisite for O 2 deficiency, but its severity and duration is controlled by the complex interaction between the hydrodynamical condition and the biogeochemical processes involved.
The North Sea is a temperate, semi-enclosed shelf sea adjacent to the northeastern Atlantic ocean.It has an average depth of about 90 m (Ducrotoy et al., 2000) with northward increasing bottom depth.The North Sea circulation is characterised by a cyclonic pattern mainly driven by the southward Atlantic inflow across the shelf edge defining its northern boundary.Lenhart and Pohlmann (1997) showed that about 85 % of the incoming Atlantic water is recirculated north of the Dogger Bank, a shallow area with water depth less than 40 m and 300 km of zonal extent at about 55 • N, 2 • E ( Kröncke and Knust, 1995).The circulation south of the Dogger Bank is governed by the inflow through the English Channel and follows the continental coast.At the southern tip of Norway it joins the Norwegian coastal current leaving the North Sea at its northern boundary.
Stratification in the North Sea reveals some substantial regional differences.While the shallower southern parts are permanently well-mixed due to the strong influence of the M 2 tidal component (Otto et al., 1990), the deeper parts north of 54 • N reveal seasonal, mostly thermal stratification (e.g., Burt et al., 2014;Pingree et al., 1978;van Leeuwen et al., 2015).Seasonal haline stratification occurs to a lesser extent along the Norwegian coast.The transition between these permanently mixed and seasonally stratified regions occurs gradually (Pingree et al., 1978).In consequence, even areas relatively near to the coast, which are affected by high riverine nutrient run-off, often reveal stratified conditions at subseasonal timescales (e.g., Burt et al., 2014).
In the 1980s, events of O 2 deficiency reaching values below 3 mg O 2 L −1 occurred regularly in the stratified southeastern central North Sea and in the German Bight (Brockmann and Eberlein, 1986;Brockmann et al., 1990;Rachor and Albrecht, 1983).During that time, the problem of low O 2 conditions in the North Sea reached public awareness in rela-tion to eutrophication as demersal animals died across a large area due to these low O 2 concentrations (von Westernhagen et al., 1986).Eutrophication, or in other words, high anthropogenic nutrient loads mainly supplied by rivers (Brockmann et al., 1988;Jickells, 1998;Rabalais et al., 2010), may raise the ambient nutrient concentrations followed by an increase in biomass production.Under given physical conditions, eutrophication thus causes an enhanced supply of organic matter sinking into the subsurface layer and reinforces O 2 consumption near the sea floor due to bacterial remineralisation.
Even though the second International Conference on the Protection of the North Sea (INSC-2) prescribed a 50 % reduction of river nutrient loads (inorganic nitrogen and phosphorus) in order to mitigate the effects of eutrophication (de Jong, 2006), Fig. 1 shows that O 2 deficiency remains a persistent problem in the North Sea up to the present day.According to Kemp et al. (2009) these events can be classified as "persistent seasonal".Low bottom O 2 concentrations may cause death of benthic organisms or fish eggs as well as avoidance of the affected areas by benthic species.Therefore, low O 2 concentrations constitute a major indicator of eutrophication (category 3 indicator, i.e., "evidence of undesirable disturbance"; OSPAR-Commission, 2003) and concentrations lower than 6 mg O 2 L −1 result in the classification as "problem area" in terms of eutrophication within the OSPAR assessment (Claussen et al., 2009).In the present study the term "oxygen deficiency" is used in this OSPAR context rather than "hypoxia".While O 2 deficiency is clearly defined within OSPAR by the 6 mg O 2 L −1 threshold, hypoxia refers to the negative impact of low O 2 concentrations on organisms.An overview of the impact of hypoxia on marine biodiversity can be found in Vaquer-Sunyer and Duarte (2008).Further descriptions on the ecological disturbance of different levels of low O 2 concentrations are summarised by Friedrich et al. (2014) and Topcu et al. (2009).
Despite the relevance of the bottom O 2 concentrations for the assessment of the ecological status of an ecosystem, O 2 measurements are sparse and either temporally or spatially limited.In addition, it is difficult to place the measurement at the right time and location to obtain a comprehensive picture of the duration and spatial extent of summer O 2 deficiency (Friedrich et al., 2014).One way to address this problem is to analyse the representativeness of available data with respect to eutrophication assessment (Brockmann and Topcu, 2014).
Only in recent years continuous measurements for the North Sea have become available by, e.g., the SmartBuoy programme of Cefas (Centre for Environment, Fisheries and Aquaculture Science, UK; Greenwood et al., 2010) or the MARNET programme (MARine Monitoring NETwork in the North Sea and Baltic Sea) of the BSH (Federal Maritime and Hydrographic Agency, Germany).These monitoring programmes provide daily time series of O 2 and related parameters (e.g., temperature, salinity, chlorophyll) in different depths and allow for the analysis of the temporal evolu-  Topcu and Brockmann, 2015).tion of stratification and O 2 concentrations at the location of observation.Greenwood et al. (2010) published the first data from continuous measurements of bottom O 2 concentrations for two sites ("North Dogger" and "Oyster Grounds") in a European shelf sea.Using these measurements, the dynamic interaction between stratification and the evolution towards low bottom O 2 concentrations can be observed, as well as the rapid recovery to saturated O 2 conditions after the breakdown of stratification due to mixing in autumn.However, even these continuous measurements did not provide sufficient information to fully understand the processes which caused the observed O 2 evolution.Greenwood et al. (2010) and Queste et al. (2013), who extended the locally confined findings by Greenwood et al. (2010) to the spatial scale using survey data from August 2010 and ICES historical data, refer to "plausible mechanisms" like vertical mixing or advection when the measurements could not be explained in detail.In consequence, Greenwood et al. (2010) stated that the data provided insight into the processes affecting the O 2 dynamics but models are required to further elucidate the significance of the seasonal drivers.
Ecosystem models produce a temporally and spatially consistent picture on O 2 and can therefore provide insight into the balance between the physical and biological factors and processes governing the evolution of the bottom O 2 concentrations.Thus, they can help understand and interpret measurements of O 2 and related parameters and can further describe the history of events of low O 2 conditions.
In this study we use the three-dimensional physicalbiogeochemical model system HAMSOM-ECOHAM to pro-vide a detailed description of the current state of the North Sea in terms of its O 2 conditions, and the processes leading to low bottom O 2 concentrations.The interpretation of the model results will enable the following questions to be answered: what are the main drivers for the O 2 dynamics in the various subregions of the North Sea?Why are certain North Sea regions more susceptible to low O 2 conditions than others despite similar stratification patterns?
For this purpose, we first validate the simulated bottom O 2 concentrations with respect to their temporal evolution and spatial distribution in order to show that the model captures the main features.Subsequently, we present a regional characterisation of the parameters controlling the bottom O 2 dynamics and propose a simple O 2 deficiency index which extends this characterisation to the entire North Sea.Finally, we attribute the individual contributions of the governing processes to the temporal and spatial variability of the overall O 2 evolution, under particular consideration of the continuous O 2 measurements at North Dogger (Greenwood et al., 2010).

The ECOHAM model
Our study is based on the coupled physical-biogeochemical model system HAMSOM-ECOHAM.The physical model HAMSOM (HAMburg Shelf Ocean Model; Backhaus, 1985) is a baroclinic primitive equation model using the hydrostatic and Boussinesq approximation (Pohlmann, 1991).HAM-SOM provides the temperature (T ) and salinity (S) distribution, in addition to the advective flow fields and the vertical turbulent mixing coefficient, which are used as forcing for the biogeochemical model ECOHAM (ECOsystem model HAMburg).For a detailed description of HAMSOM the reader is referred to Pohlmann (1991).Further information on the application of HAMSOM can be found in Backhaus and Hainbucher (1987) and Pohlmann (1996).
The air-sea exchange of O 2 at the sea surface constitutes an important physical process besides the effects of advective transport and vertical diffusion in the interior water column.The air-sea flux of O 2 in the present application is parameterised according to Wanninkhof (1992).In relation to the biology, the O 2 cycle is linked to the C cycle by photosynthesis, zooplankton respiration and bacterial remineralisation.While photosynthesis is a source of O 2 , zooplankton respiration and bacterial remineralisation act as O 2 sinks.
A further sink of O 2 is nitrification, the bacterial transformation of ammonium to nitrate.Within ECOHAM, this process is light-dependent and links the O 2 cycle to the N cycle.Nitrification only occurs under aerobic conditions (i.e., concentrations > 0 mg O 2 L −1 ), which is a realistic constraint for the pelagic North Sea environment.It is light-dependent, being stronger under low light conditions.Pelagic denitrification is implemented, but is negligible as it only occurs under anaerobic conditions.Pelagic anaerobic ammonium oxidation (anammox) is not implemented, however, it can be neglected for the same reason.Except for primary production, the biological processes involved in the O 2 cycle are not temperature-dependent in the present model setup.
For the representation of the benthic remineralisation processes a simple sediment module is used.A layer of zero extent is defined below the deepest pelagic layer of each water column.There the deposited organic matter is collected and remineralised (Pätsch and Kühn, 2008).The benthic remineralisation of the organic matter is defined as a first-order process with relatively high remineralisation (C, N, P) and dissolution rates (Si; opal) preventing year-to-year accumulation of deposited matter.The released dissolved inorganic matter is returned directly into the pelagic bottom layer.Different rates are applied to organic C, N, P and Si resulting in different timescales for the release into the pelagic.In ECOHAM, the O 2 cycle is affected by the benthic remineralisation in a direct and indirect way.First, the remineralisation in the sediment is accompanied by the direct reduction of the O 2 concentrations in the pelagic bottom layer above.Second, inorganic nitrogen is released from the sediment in the form of ammonium, which can be nitrified within the water column under O 2 consumption.According to Seitzinger and Giblin (1996), who suggested a tight coupling between benthic nitrification and denitrification, benthic denitrification depends on the benthic O 2 consumption in our model.Direct benthic nitrification and benthic anammox are neglected as the sediment has zero vertical extent (Pätsch and Kühn, 2008).
For a more detailed description of the ECOHAM model, including the full set of the differential equations and parameter settings of ECOHAM, the reader is referred to Lorkowski et al. (2012).A detailed description and analysis of the O 2 module can be found in Müller (2008).

Model setup and forcing data
The model domain extends from 15.250 • W to 14.083 • E and from 47.583 to 63.983 • N and comprises the entire North Sea, large parts of the northwestern European continental shelf and parts of the adjacent northeastern Atlantic.The horizontal resolution is 1/5 • with 82 grid points in latitudinal direction and 1/3 • with 88 grid points in longitudinal direction.The horizontal grid of the model domain is shown in Fig. 2. The vertical dimension with a maximum depth of 4000 m is resolved by 31 z-layers with a surface layer of 10 m.The vertical has a resolution of 5 m between 10 and 50 m depth, which is relevant for the calculation of the MLD (Sect.2.2).Below 50 m, the layer thicknesses successively increase with depth.
The model system was run over the period 1977 to 2012.HAMSOM was initialised with a monthly-averaged climatology based on the World Ocean Atlas (WOA; Conkright et al., 2002).The meteorological forcing was derived from NCEP/NCAR reanalysis data (Kalnay et al., 1996;Kistler et al., 2001) and provides 6 hourly information for air temperature, cloud coverage, relative humidity, wind speed and direction.Short wave radiation was calculated from astronomic insulation and cloud coverage applying a correction factor of 0.9 (Lorkowski et al., 2012).The data were interpolated to the model grid and time step according to O'Driscoll et al. (2013) and Chen et al. (2014).Daily freshwater run-off data for 249 rivers were provided by Cefas and represent an updated data set of that used by Lenhart et al. (2010) covering the entire simulation period.The same data set encompassed nutrient loads used for the ECOHAM.
At open boundaries, surface elevation was prescribed as a fixed (Dirichlet) open boundary condition (OBC) according to the M2 tide, while for horizontal transport velocities radiation OBCs were applied.For tracers (T and S) radiation and radiative-nudging OBCs were used in the case of inflow and outflow, respectively.A detailed description of the OBCs is provided by Chen et al. (2013).The HAMSOM simulation was carried out with a 10 min time step.
ECOHAM was run off-line with a time step of 30 min using the 24 h averages of the hydrographic and hydrodynamic fields generated by HAMSOM.In the model setup used, short wave radiation is attributed to the first layer (surface) only and the specific effect of light attenuation due to SPM and planktonic self-shading on the thermal structure is not taken into account.A sensitivity study allowing for deeper light penetration and feedback on the thermal structure confirmed this effect to be only of minor importance (not shown).
For the biogeochemical state variables a climatology of depth-dependent monthly averages was prescribed at the boundaries and solely for DIC yearly changing data were provided (Lorkowski et al., 2012).To include the effect of SPM on the light climate, a daily climatology from Heath et al. (2002) was used.Data for atmospheric N deposition were compiled using a hybrid approach.This was required since the overall simulation period  exceeds the period of data available from the EMEP (Cooperative program for monitoring and evaluation of the long-range transmissions of air pollutants in Europe) model (1995)(1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012).First, the EMEP results for total deposition of oxidised (NO x ) and reduced nitrogen (NH 3 ) were interpolated to the model grid.Second, we calculated the average annual deposition rates for the NO x and NH 3 for each grid cell, based on the 1995-2012 EMEP data.The resulting spatially resolved arrays of average deposition rates were subsequently normalised by the spatial average of the entire domain to yield the spatially resolved anomaly fields.Finally, gridded deposition rates for individual years were obtained using (1) the gridded anomaly fields, (2) EMEP's spatially averaged (over our model domain) deposition rates for year 2005, and (3) long-term trends (normalised towards year 2005) for the temporal evolution of European emissions of NO x and NH 3 (Fig. 2 in Schöpp et al., 2003).The output of the biogeochemical simulation was stored as daily values (cumulative fluxes, state variable snapshots) for the entire domain and simulation period.

Extracting stratification parameters from model results
Stratification constitutes the prerequisite for the potential evolution of low O 2 conditions in the North Sea.In this study (1) its duration and (2) the mixed layer depth (MLD) are used to describe stratification.Seasonal stratification in the North Sea is mainly T -driven (Burt et al., 2014), except for the regions of haline stratification along the Norwegian coast.As observations do not cover the entire model domain and simulation period we determined the duration of stratification and the MLD from the simulation results.For this purpose we developed a simple 2-step algorithm based on T .First, the stratified period is determined using a temperature difference criterion: S strat is a switch defining if a water column at location (x, y) and time t is stratified (S strat = 1) or not (S strat = 0) depending on the temperature difference T between the surface and bottom depth H .The critical temperature difference T crit = 0.05 K was determined by evaluating different T crit against the temporal evolution of simulated bottom O 2 at different locations within the model domain.In addition, periods of stratified conditions are only considered as such, if they last for at least 5 days without any interruption.Otherwise bottom waters are considered to be ventilated again.
In the second step, in the case of stratification the MLD of a model water column is determined using the vertical T gradient T / z: T / z is calculated for each grid cell interface within the considered water column.T represents the T difference between two vertically adjacent model layers and z represents the distance between the centre points of these two grid cells.D is then defined as the depth level of the interface where T / z has its maximum.As the described stratification and MLD criterion differ significantly from common MLD criteria (e.g., Table 1 in Kara et al., 2000), an evaluation is provided in Appendix A.

Validation data
For the validation of the model results we used observation data from different sources.The data sets can be subdivided into two types: (1) temporally resolved, localised data and (2) spatially resolved "snapshots".The first type was used for the validation of the seasonal evolution of O 2 , whereas the second type was used to validate the general spatial patterns and year-to-year variability of bottom O 2 during late summer.

Localised, temporally resolved data -Cefas-SmartBuoy and MARNET
Cefas operates a network of SmartBuoys to provide autonomous in situ measurements of physical, chemical and biological parameters (Mills et al., 2005).A SmartBuoy was located directly north of the Dogger Bank ("North Dogger") at 55  (Greenwood et al., 2010).For validation purposes the O 2 data derived from the sensor at 85 m depth were used.North Dogger data are published and can be accessed according to Greenwood et al. (2016).
The BSH operates a continuous monitoring station at 54 • 10 N, 6 • 21 E (see Fig. 2, region 1; hereafter referred to as station "Ems").The O 2 saturation is measured hourly using opto-chemical sensors (optodes).Sensors are located in 6 and 30 m depth, respectively, and the bottom depth is 33 m.The applied sensors have a resolution of 0.03 mg O 2 L −1 and an accuracy better than 0.26 mg O 2 L −1 .Before deploying the sensors a 0-100 % calibration is conducted, and they are re-calibrated after operation to quantify any drift.In addition, a regular on-site validation takes place using a calibrated fast optode (accuracy of ±2 %) or by applying the Winkler titration (accuracy better than ±1 %).

Spatially resolved "snapshot" data -the North Sea programme
During the North Sea programme, carried out by the Royal Netherlands Institute for Sea Research (NIOZ) with support from the Dutch Science Foundation (NWO) and the European Union, the North Sea was sampled from 18 August to 13 September 2001, and from 17 August to 5 September 2005 and 2008.The North Sea was covered by an approximate 1 • × 1 • grid, sampling approximately 90 stations in each of the years (Bozec et al., 2005(Bozec et al., , 2006;;Salt et al., 2013).During each cruise, a total of 750 water samples were collected for dissolved O 2 .In 2001, the O 2 concentrations were determined by the Winkler titration using a potentiometric end-point determination with an accuracy of ±2 µmol O 2 kg −1 (less than ±0.07 mg O 2 L −1 depending on T and S).In 2005 and 2008, the O 2 concentrations were obtained applying the spectrophotometric Winkler approach with a precision of less than 0.03 mg O 2 L −1 .A detailed description of the measurement system used is given in Reinthaler et al. (2006).The data for the years 2001 and 2005 have been published and can be accessed according to Thomas et al. (2012) and Thomas and Borges (2012), respectively.
The data available were gridded to the model grid (Fig. 2).In the case of multiple measurements for the same model grid cell and date, the average of these measurements was used for validation.To compare our model results to these data, we calculated the averages and standard deviations of our simulation over the observation period of the corresponding year.

Deriving a regional O 2 characterisation of the
North Sea

Identification of the key parameters
For the development of a regional O 2 characteristic, potential controlling factors were analysed in relation to bottom O 2 .
Besides stratification, eutrophication is considered as a major driver for developing low O 2 conditions (e.g., Diaz and Rosenberg, 2008;Kemp et al., 2009).Thus, primary production within the mixed layer and the resulting organic matter export into the layers below the MLD must be considered to be the main source for degradable organic matter.In addition, organic matter can be advected from surrounding waters in the form of phyto-or zooplankton and detritus, subsequently sinking out of the mixed layer.
Another important criterion is the water volume below the thermocline (Druon et al., 2004).A smaller volume separated from the surface due to stratification holds a lower initial inventory of O 2 than a larger volume even though concentrations can be similar or even higher in the smaller volume.Thus, our set of O 2 -related characteristics consists of mixed layer primary production (PP mld ), horizontal advection of organic matter into and out of the mixed layer (ADH org,in and ADH org,out ; including phyto-/zooplankton and detritus), vertical organic matter export below the MLD (EXP org ; only detritus) and mixing of O 2 below the MLD (MIX O 2 ), and the sub-MLD volume V sub .
To detect regional characteristics within the North Sea area, we defined four different sub-domains encompassing 4 × 4 model water columns each (see Fig. 2, red boxes): (A) southern North Sea (SNS) under strong tidal influence, (B) southern central North Sea (SCNS) with high year-toyear variability in stratification, (C) northern central North Sea (NCNS) with a dominant summer stratification each year, and (D) northern North Sea (NNS) with a dominant summer stratification each year and a strong influence of the Atlantic.For all these regions, the parameters described above were calculated for the years 2000-2012 relative to a reference depth D ref , which is defined as the bottom depth of the model layer directly below the annual maximum MLD among all four regions.We decided to use a D ref > MLD to ensure that for the different regions all parameters were determined on a comparable level.This implies that the values for PP mld , ADH org,in and ADH org,out are integrated from the surface to D ref , whereas EXP org and MIX O 2 are the vertical fluxes through D ref .
The same D ref was applied to all regions, but year-to-year variations were allowed.
To determine the annual maximum MLD, we first calculated the stratification period for the 4 × 4 regions B-D using Eq.(1).Region A was excluded from this calculation as no persistent MLD developed due to tidal mixing.In this context, S strat (t) of a region is only 1 if S strat (x, y, t) = 1 for all 16 water columns within a 4 × 4 region.The daily MLD for each water column within a region was calculated by apply-ing S strat (t) to Eq. ( 2), and subsequently the daily MLD of the region is defined as the median of these 16 daily values.The annual MLD for each region was then determined as the median of this daily time series.Finally, the annual maximum MLD among all four regions is used to determine the reference depth D ref , which is defined as the bottom depth of the layer directly below this maximum MLD.
The values for these O 2 -related quantities were calculated for individual years relative to D ref and temporally integrated over the period from 1 April to 30 September (hereafter "summer").Consequently, the average values over the entire period 2000-2012 are calculated and presented in Table 1, additionally including the average O 2 concentrations at the beginning and end of the summer period as well as the average duration of stratification.

Development of a spatially resolved index for North Sea O 2 deficiency
In order to obtain a North Sea wide indicator for O 2 deficiency under stratified conditions, it is necessary to extend the regionally confined characteristic described in the previous section.For this purpose, we extract the key factors affecting O 2 from this regional information and combine them into a single index -the oxygen deficiency index (ODI).The ODI aims to represent the main spatial and temporal patterns of O 2 deficiency in the North Sea under stratified conditions, while being as simple as possible and incorporating only a very limited number of parameters.Stratification period, organic matter export and subthermocline volume are considered as the key parameters controlling the bottom O 2 dynamics.Surface primary production can be used as a proxy for organic matter export assuming that most of the exported organic matter is produced locally.Bottom depth can be used as an indicator for the sub-MLD volume assuming only minor fluctuations of the MLD during the summer stratified period.In addition, the bottom depth directly influences the amount of organic matter reaching the bottom layer relative to the amount being produced near the surface, due to the exposure of sinking matter to pelagic remineralisation.Thus, the following key factors are used for the calculation of this index: (longest continuous) stratification period (t strat ; in days), summer surface primary production (PP mld ; in g C m −2 ; 1 April to 30 September), and bottom depth (D bot ; in m).
First, individual dimensionless indices are calculated for each of these quantities.The individual indices range between 0 and 1, indicating conditions counteracting and supporting O 2 deficiency, respectively.The calculation of the stratification and production indices, I strat and I pp , is based on the work by Druon et al. (2004) and reads as I Q i (x, y) represents the index corresponding to the actual value of the quantity Q i (x, y) with its defined upper and lower thresholds, Q i,max and Q i,min .For t strat , Q i,max and Q i,min are set to 50 and 150 days, respectively.Stratification periods of less than 50 days are considered to be too short to facilitate the evolution of O 2 deficiency, while periods longer than 150 days are considered seasonally well-stratified.The lower threshold for PP mld was set to 120 g C m −2 as PP mld does not reach lower values in most parts of the North Sea.The upper threshold was set to 200 g C m −2 as such high values and even higher are simulated in the southeastern North Sea.
For the depth index, I D , a different definition was chosen as lowest O 2 concentrations occur in areas of intermediate depth, where seasonal stratification can develop and the O 2 inventory is limited due to a small volume below the thermocline.Therefore, we defined I D as follows: (4) D bot represents the actual bottom depth at location (x, y).D peak = 40 m is the bottom depth we found to be most favourable for O 2 deficiency in the North Sea.The lower threshold D min = 25 m corresponds to the maximum MLD we found for the shallower southern North Sea.The upper threshold D max = 90 m was chosen to exclude the areas where the initial O 2 inventory is sufficient to prevent O 2 deficiency due to the large volume below the thermocline.Finally, the ODI combines the three individual indices according to the following equation: (5) Here, I Q i and w Q i represent the index for a quantity and the related weight, respectively.The values for t strat are referred to by Q 1 and those for PP mld by Q 2 .The equation for ODI implies that it is zero in areas where I D = 0.The stronger weighting of PP mld implies that variations in the ODI between different years are more strongly affected by variations in summer surface productivity than by the duration of stratification.
The ODI ranges between 0 (low risk of O 2 deficiency) and 1 (high risk) and is calculated for each water column (x, y) within the model domain.By this we obtain a spatially resolved indicator for O 2 deficiency in the North Sea, which helps regionalise the North Sea in terms of O 2 conditions.

Quantification of driving processes: spatial and temporal variability, and data interpretation
In order to quantify the processes driving the O 2 dynamics in different regions, we calculated O 2 mass balances for three different regions encompassing 2 × 2 grid cells (see Fig. 2, regions 3-5).First, mass balances for the entire volume below the thermocline (hereafter "sub-MLD") in region 3 are compared with the corresponding bottom layer mass balances to identify differences between the bottom layer dynamics and the dynamics within the entire sub-MLD volume.This is done for 2 years, 2002 and 2010, to analyse variations between these years.Region 3 was chosen as it shows the lowest bottom O 2 concentrations within the entire model domain, with the overall minimum in 2002 and relatively high concentrations in 2010.The daily resolved MLD defines the upper integration limit for the sub-MLD mass balances, i.e., the integration depth may vary during the stratified period.
The daily MLD is defined as the vertical level of the model grid which is closest to the daily average MLD of the 2 × 2 region according to Eqs. ( 1) and ( 2).Second, we compare the O 2 mass balances of the bottom layer for regions 4 and 5 in 2002 with that of region 3 to unveil regional differences.In a last step the mass balance analysis is applied to interpret the O 2 evolution observed at North Dogger (see Fig. 2, region 2).
The O 2 concentrations and saturation concentrations shown in the different mass balances represent the average values within the analysed volume.Values of O 2 saturation concentrations were calculated according to Benson and Krause (1984) using simulated T and S. The fluxes presented are cumulative changes in the O 2 concentrations of the considered volume, i.e., the values at the end of the stratified period reflect the total net change of the O 2 concentrations due to the corresponding physical or biological process.Positive and negative values at the end of the stratification period indicate net gain and loss, respectively.The slope of each line represents the intensity of the corresponding flux at the specific moment in time, i.e., a steep positive (negative) slope implies a strong gain (loss) effect.At North Dogger, observed and simulated bottom O 2 concentrations show a steady decrease after the onset of stratification.While stratification according to Eq. (1) starts a bit earlier compared to that described by Greenwood et al. (2010), the beginning of the decrease in bottom O 2 concentrations coincides well.The simulated and observed O 2 concentrations at this time are in good agreement.Some small-scale fluctuations in the observations are not fully reproduced by the simulation, however, the general evolution is represented well by the model.The average O 2 reduction in the simulation is slightly less than in the observations, visible in the difference between the concentrations at beginning and end of the stratified period.Stratification ends a bit earlier in the simulation, with the result that simulated bottom O 2 starts to recover while the observed concentrations continue to decline.The observed O 2 concentration at the end of the stratified period is about 6.8 mg O 2 L −1 , while the simulation results in about 7.4 mg O 2 L −1 .
In 2008, we can see a similar slight overestimation of simulated O 2 concentrations, but less than in 2007.Some minor fluctuations in the observations are again not represented by the model, but the general evolution of bottom O 2 is represented well.It should be noted, that the different depths of the time series (76 m for simulation, 85 m for observation) may In 2011, persistent stratified periods derived from the simulation do not exceed 2 months at station Ems.Consequently, the temporal evolution of bottom O 2 represents mainly the temporal evolution of the O 2 saturation concentrations.Again large fluctuations of up to ±2 mg O 2 L −1 can be seen in the observations which are not fully reproduced by the model.Besides these short-term changes, the difference between simulated and observed bottom concentrations is less than 0.8 mg O 2 L −1 with higher summer values in the simulation.
The validation of bottom O 2 at the stations North Dogger and Ems shows that the HAMSOM-ECOHAM model is capable of reproducing the main features of the bottom O 2 dynamics at these two stations.The minor differences in the concentrations (< ± 0.4 mg O 2 L −1 ) at the beginning and end of the year, representing mainly the saturation concentrations, show that the general physical setting provided by the model is reasonable.The slightly slower O 2 reduction in the simulation may indicate an underestimation of the biological consumption, e.g., due to benthic remineralisation.Intra-seasonal fluctuations at both stations are not fully reproduced, due to the limited spatial resolution of the model grid.Additionally, the tides may have an effect at station Ems on the short-term.However, they are not resolved due to the daily time step of the simulated current fields.
The generally good agreement between simulation and observation is also shown by the Taylor diagram (Taylor, 2001, see Fig. 5, markers "a" for North Dogger and "b" for Ems), which presents the correlation coefficients (COR), standard deviations (SDs) and centred root-mean-square differences (RMSD) of the simulation relative to the observations.SDs and RMSDs are normalised by the SD of the corresponding observations.For analysis, the data of each data set was merged into a continuous series of data.For both stations, COR is high with values of about 0.95 and the normalised RMSD is less than 0.38.The RMSD values are mainly due to the larger range and higher (seasonal and intra-seasonal) variability in the observed bottom O 2 , which is also indicated by the normalised SDs of about 0.73 and 0.82 for Cefas North Dogger and MARNET Ems, respectively.4).Standard deviations and centred rootmean-square differences (RMSD) were normalised by the standard deviation of the corresponding observations.by low O 2 conditions and the southeastern parts, which are more vulnerable to low O 2 concentrations.In addition, the model demonstrated it is capable of capturing variations in the bottom O 2 evolution between different years.In combination with the results of the time series validation (Sect.3.1.1),this confirms that the described model setup provides reliable information on the internal physical and biological processes driving the O 2 dynamics in the North Sea.

Spatial distribution of late summer bottom O 2
The small SD of simulated bottom O 2 confirms that using the averages over a period of up to 4 weeks provides a reasonable measure for most areas.In addition, these small values imply that measurements taken late August/early September (before the breakdown of stratification) can be considered as a representative synoptic picture of the late summer bottom O 2 conditions.However, the time series validation showed that in some areas lowest concentrations of bottom O 2 may occur remarkably later in the year (see Fig. 3a).Consequently, the picture obtained from observations taken in August/September does not necessarily reflect the spatial distribution of minimum bottom O 2 concentrations, which underlines the importance of choosing the appropriate point in time for the monitoring of low O 2 conditions.
The small SD of the observations, which is a result of the data gridding, shows that in most regions vertical O 2 gradients near the bottom are negligible.The high values of 0.75 and 0.59 mg O 2 L −1 southeast of the Dogger Bank and northwest of Denmark, respectively, result from the fact that val-ues above and below the thermocline were taken into account for the averaging.
As for the time series, Fig. 5 (marker c) shows the statistical measures of the validation for the spatially resolved data.Here, COR reaches only about 0.64 which is also indicated in Fig. 4 by the variations between year 2008 and the previous years, when the simulation revealed a relative change inverse to that in the observations in the northern North Sea.The normalised RMSD of 0.77 is about twice as high as for the time series, which can be attributed to the greater regional differences in the observed bottom O 2 concentrations with higher maximum and lower minimum values.The normalised SD equals 0.67 which indicates the less strong spatial gradients in the simulation.These statistics confirm that the spatial patterns in the observed bottom O 2 concentrations are basically reproduced by the model, with only slight shortcomings with respect to the amplitude of the bottom O 2 concentrations and year-to-year variations in some regions of the North Sea.trations north of the Dogger Bank in August 2010, however, they found the minimum concentrations a bit further north around 57 • N.
The discrepancy between minimum O 2 concentrations in the 2 years and the quite similar stratification patterns demonstrates that stratification is an important prerequisite for low O 2 conditions, but other physical or biological factors do have a strong effect on the O 2 dynamics in the North Sea.

Key parameters
Table 1 shows the 2000-2012 summer (1 April to 30 September) averages of the quantities considered potentially relevant for O 2 for the regions A-D (see Fig. 2).The quantities were calculated relative to a D ref of 25 m.In addition, the stratification period (t strat ), average MLD (D mld ), average bottom depth (D bot ) and area of the regions are provided.
The mixed layer net primary production, PP mld , is strongest the coastal region A and shows decreasing values towards the central North Sea.In the SCNS region B, PP mld accounts for about 87 % of that in the highly productive coastal region A. The corresponding value for the NCNS region C and NNS region D is about 80 and 88 %, respectively.
Despite the highest PP mld in the coastal region A, the SCNS region B shows the strongest contribution of gross advection of organic matter, ADH org,in and ADH org,out .Both regions show positive net advection of organic matter, while the two northern regions C and D are characterised by negative net advection, i.e., advective loss in organic matter.The latter regions are located north of the Dogger Bank, thus, they are affected by the northern Atlantic inflow.In this region, net advection results in a loss in organic matter as the recirculated northward flowing water has higher concentrations of organic matter than the incoming Atlantic water.
The vertical export of organic matter, EXP org , below D ref consistently adds up to about 12-14 % of PP mld , indicating the clear link between these quantities.Region B yields an EXP org of 75 % of that in the coastal region A, which corresponds to a higher net advective import of organic matter of 6.7 g C m −2 in region A, compared to only 1.6 g C m −2 in Region B. EXP org in region C and D reach about 69 % and 79 % relative to region A, respectively.
The vertical mixing of O 2 , MIX O 2 , is highest within the coastal region A and adds up to 116.1 g O 2 m −2 , which is due to strong tidal mixing.The stratification period, t strat , of 151 days in region B is shorter than in region C (220 days) and does not cover the entire summer period.Thus, MIX O 2 in region B is significantly larger than in regions C and D.
The evolution of the O concentrations between the beginning and end of the summer period reveals some interesting aspects in relation to the previously mentioned parameters.The O 2 concentrations at 1 April show significant differences between the regions ranging between 9.5 mg O 2 L −1 (regions C and D) and 10.1 mg O 2 L −1 (region A).The O 2 concentrations at the 30 September yield values between 7.7 mg O 2 L −1 (region A) and 8.3 mg O 2 L −1 (region D).This implies a consistently decreasing O 2 consumption during summer from region A to D. This spatial gradient in the O 2 consumption is opposite to that in t strat , which shows a steady increase from regions A to D.
In order to give an impression of the impact of EXP org on the O 2 dynamics of the water volume below the MLD, V sub , we link the amount of exported organic matter to the amount of O 2 available within V sub assuming the organic matter is remineralised completely in the area of settlement.Based on the O 2 concentration at the beginning of April, the total amount of O 2 available is 1365 kt for region B and 4590 kt for region C.The total amount of exported organic matter is calculated as the product of EXP org and the total area of the considered region.This calculation yields 130 kt C and 115 kt C for the regions B and C, respectively.As O 2 consumption and C release occur with a molar ratio of 1 : 1 during bacterial remineralisation (Neumann, 2000), we obtain the daily O 2 consumption by dividing by the total duration of the considered 6-month period (= 183 days), yielding 0.71 and 0.63 kt O 2 d −1 for regions B and C, respectively.
The initial O 2 mass is calculated as the product of the initial O 2 concentration and V sub .Assuming the daily O 2 consumption to be constant for each region, division of this mass by the daily consumption rate calculated above provides an estimate of the amount of time required for the consumption of the entire amount of O 2 available in V sub .This calculation yields a period of about 2 years for region B, whereas the corresponding value for region C is significantly higher with almost 12 years.This great difference between the resulting periods (factor 6), compared to the relatively small difference between the daily consumption rates (factor 1.1), illustrates clearly the large influence of the sub-MLD volume, V sub , on the temporal evolution of the O 2 concentrations below the MLD.
The same calculation based on the threshold of 6 mg O 2 L −1 used by OSPAR (OSPAR-Commission, 2005) yields a consumption period of 283 days for region B, which indicates the relatively high potential for O 2 deficiency in this region.
This characteristic based on the four different North Sea regions demonstrated that the duration of stratification alone cannot explain the temporal evolution of sub-MLD O 2 concentrations.It shows the great importance of the organic matter export which drives the biological O 2 consumption.In addition, the volume below the MLD plays a key role as it governs the amount of O 2 which is available throughout the stratified period, and in combination with the organic matter export defines whether O 2 deficiency may occur or not.Thus, these three quantities can be considered as the key parameters governing the O 2 dynamics of the seasonally stratified North Sea.

The oxygen deficiency index (ODI)
The ODI resulting from the simulated stratification duration (t strat ), summer surface primary production (PP mld ) and model topography for the years 2002 and 2010 is shown in Fig. 7a and b, respectively.It can be seen that the ODI tends to be higher in 2002 than in 2010 in the region where minimum bottom O 2 is lowest in both years (see Fig. 6).North of the Doggerbank, the ODI also shows slightly higher values than in the surrounding waters which corresponds to the lowered bottom O 2 in this region.The variations of the minimum concentrations between the 2 years in this region are also well-reproduced by the ODI.Especially in 2002, the highest ODI coincides with the lowest concentrations in the entire domain.In 2010, the highest ODI is located a bit south of the minimum O 2 concentrations, which is mainly caused by the high surface production in this region.Along the northern British coast, the ODI also shows high values for both years which is in good agreement with the slightly lower minimum bottom O 2 in this area.However, ODI values tend to be too high and do not represent the slightly lower minimum O 2 concentrations off the eastern Scottish coast around 57-58 • N as the bottom depth in this area exceeds 90 m (i.e., I D = 0).Directly northwest of Denmark, the ODI also yields high values for both years with higher values in 2002.This corresponds well to the simulated bottom O 2 concentrations in this area, even though ODI values are too high, compared to ODI values in the central North Sea.
With respect to the factors selected for the calculation of the ODI, Table 1 shows that stratification alone is not sufficient to explain the North Sea O 2 dynamics.While the reduction in the O 2 concentrations is steadily decreasing from regions A to D, stratification duration is characterised by a steady increase from regions A (80 days) to D (226 days).Regarding the PP mld , the strongest O 2 reduction occurs in the regions of highest productivity A and B. In the northern regions, the higher PP mld in region D does not correspond to stronger reduction in O 2 .As t strat is also higher in this region, a further factor is needed to describe the basic O 2 dynamics.The reduced effect of surface production in the northernmost area is likely to result from the dilution effect due to the higher V sub .Considering the bottom depth D bot as a proxy for V sub , it is shown that the strongest decrease occurs in the shallower regions A and B with average depths of about 40 m, while the regions deeper than 90 m (C and D) show a weaker decrease.Net advection of organic matter, which is not taken into account in the ODI, appears to be of minor importance for subsurface O 2 relative to the local surface production as the net advective input of organic matter is significantly less than local production.The O 2 concentrations at the beginning of the stratified period were also not taken into account as they show lower values in regions with higher minimum concentrations and vice versa.
In summary, the ODI represents well the spatial and temporal variations of minimum bottom O 2 concentrations despite the small set of controlling factors.This confirms that a simple combination of only stratification duration, organic matter production and bottom depth is sufficient to reproduce the main spatial and temporal patterns of the minimum bottom O 2 in the seasonally stratified North Sea.Thus, the findings from Table 1 can be applied to most parts of the North Sea.In addition, the similarity of the ODI inside the regions analysed in Sect.3.3.1 and inside the regions selected for the mass balance analyses (see Fig. 2, regions 2-4) and their surrounding areas shows that these regions can be considered as representative, allowing for a meaningful analysis of the O 2 dynamics in these regions.

Driving mechanisms and year-to-year variability of sub-thermocline O 2 dynamics
The previous analyses showed that stratification constitutes a necessary condition for O 2 deficiency, but year-to-year variations especially in the biological factors mainly control the O 2 dynamics.For a better understanding of the processes controlling sub-thermocline O 2 a more detailed analysis is provided by the mass balances in Fig. 8.As the bottom O 2 dynamics are also influenced by the processes in the midwater, Fig. 8a and b show the O 2 mass balances for the sub-MLD volume (V sub ) in region 3 (see Fig. 2) for the years 2002 and 2010, respectively.The stratification characteristics are similar for both years with an average MLD of about 15 m and a stratification period (grey area) of 187 days in both years, only differing by a later onset (and breakdown) of stratification in 2010.The temporal evolution of the MLD (dash-dotted grey) is also similar, with deeper MLDs at the beginning and end of the stratified period and few events of enhanced mixing during the summer months (May to August).
The sub-MLD O 2 concentration (solid magenta) at the beginning of the stratified period in 2002 is about 9.81 mg O 2 L −1 , being about 0.33 mg O 2 L −1 lower than in 2010.At the end of stratification, the 2002 value of 6.85 mg O 2 L −1 is about 0.81 mg O 2 L −1 lower than in 2010.
In 2002, the clearly diverging temporal evolution of simulated O 2 and the corresponding saturation concentrations (O 2,sat ; dash-dotted magenta) reveals that the different O 2 evolution is caused not only by decreasing O 2 solubility.Hence, other factors must play an important role for the O 2 evolution below the MLD.

The influence of advection and mixing
The comparison of advection (ADV O 2 ; including horizontal and vertical components; dashed light blue) and vertical mixing (MIX O 2 ; turbulent diffusion; dashed dark blue) for the years 2002 and 2010 shows strong variations between the 2 years.ADV O 2 regularly changes its influence on the sub-MLD O 2 concentrations during stratification in both years.However, considering the temporally integrated effect, ADV O 2 causes a net gain of about 25.8 g O 2 m −2 in 2002, whereas in 2010 it results in a slight net loss in O 2 .Advection positively affects the O 2 concentrations during the last 2-3 weeks of the stratified period in both years and even causes a net O 2 increase in 2002.

Biological drivers of the sub-thermocline O 2 dynamics
In contrast to the integrated effect of the physical factors, Fig. 8a and b The integrated effect of all biological sink processes (REM pel , REM sed , RES zoo and NIT) adds up to −204.8 g O 2 m −2 in 2002 and to −136.4 g O 2 m −2 and 2010, i.e., the biological O 2 consumption in 2002 is 1.5 times higher than in 2010 and 1.2 times higher than the 2000-2012 average of 169.9 ± 21.2 g O 2 m −2 .The relative contribution of the individual processes to the biological O 2 consumption shows only minor variations during the analysed period.REM pel contributes to 53.6 ± 1.7 %, while REM sed accounts for 17.2 ± 0.9 %.For RES zoo and NIT the average contributions result in 21.6 ± 1.5 and 7.6 ± 0.8 %, respectively.The EXP org below 25 m depth (not presented; calculation analogous to Table 1) in 2002 is nearly 1.6 times larger than in 2010, which is in good agreement with the differences in the integrated effect of the biological O 2 sinks.
In late April and late June 2002 (Fig. 8a) two events of enhanced mixing reveal direct and indirect effects on the biological processes.The renewal of the nutrient pool causes short-term increases in PP around the MLD which in turn enhances RES zoo and REM pel .Consequently, only the stronger event in late June causes a net increase in O 2 .It should be noted that the shown strengthening in the different biological effects is also influenced by the change in the MLD (i.e., integration depth), however, it is also visible when considering a constant MLD (not shown).
PP shows the strongest effects on sub-MLD O 2 when the MLD is shallowest which indicates the existence of a deep chlorophyll maximum (DCM).This explains the negative influence of MIX O 2 during these periods as O 2 concentrations are highest within the DCM due to high PP.The only minor positive or even negative effect of MIX O 2 during most of the stratified period emphasises the importance of stratification for the sub-MLD O 2 dynamics as it efficiently limits the O 2 supply.
The good agreement between the variations in EXP org and the integrated effect of the biological O 2 sinks between the 2 years confirms that the supply of detrital matter to the deep layers is the driving force of sub-MLD O 2 consumption.The strong influence of pelagic remineralisation demonstrates its crucial role for the bottom O 2 concentrations as it directly affects the potential O 2 supply from the mid-water into the bottom layer.

Bottom layer dynamics of the North Sea O 2 minimum zone
Even though the dynamics in the mid-water affect the bottom O 2 levels, lowest concentrations occur in the bottom layer.In order to show which processes are the main contributors to the O 2 dynamics in this layer, Fig. 8c and d show the mass balances for the bottom layer in region 3 for 2002 and 2010.The average bottom depth in this region is 47.75 m, and the model bottom layer encompasses a volume of about 14.4 km 3 .
The O 2 concentrations at the beginning of the stratified period, 9.79 and 10.12 mg O 2 L −1 for 2002 and 2010, respectively, are similar to those in V sub .The concentrations at the end of stratification, 6.76 mg O 2 L −1 in 2002 and 7.55 mg O 2 L −1 in 2010, show larger differences to those for V sub .
The effect of the physical factors, ADV O 2 and MIX O 2 , on the bottom O 2 is different to that for V sub .While in 2002 ADV O 2 shows a similar effect on O 2 as for the sub-MLD O 2 , its effect in 2010 is opposite to that for V sub , resulting in a minor increase of about 1.1 g O 2 m −2 .During the last 3 weeks of stratification in 2002, the same positive effect of ADV O 2 as in the sub-MLD mass balance is shown, initiating the recovery of the bottom O 2 before MIX O 2 intensifies.MIX O 2 has a consistently positive effect on the bottom O 2 in both years.Its integrated effect is increased relative to V sub by the factor 1.7 and 1.4 in 2002 and 2010, respectively.
The relative contribution of the biological O 2 sinks in the bottom layer is also different to the sub-MLD volume.The 2000-2012 averages reveal that in the bottom layer REM sed accounts for 50.1 ± 1.2 % of the total biological O 2 consumption, while REM pel contributes to 32.2 ± 1.4 %.Thus, aerobic remineralisation consistently adds up to more than 80 % of the biological O 2 consumption in the bottom layer.This shift results from the different volumes considered, and the fact that REM sed only has a direct effect on the deepest pelagic layer.Average O 2 consumption due to REM sed results in values between 3.9 and 6.5 mmol O 2 m −2 d −1 .
For RES zoo , the influence on the bottom O 2 concentrations is lower than in V sub (11.3 ± 1.2 % during 2000-2012).This relates to the fact that zooplankton tends to stay in the upper part of the water column where phytoplankton concentrations are higher.NIT represents the weakest sink for bottom O 2 with an average contribution of 6.4 ± 0.6 % during these years.PP as a potential source for O 2 is negligible in the bottom layer due to light limitation.
The analysis clearly shows vertical mixing is the only efficient gain term for O 2 in the bottom layer.Benthic aerobic remineralisation constitutes the major driver for O 2 deficiency in the North Sea O 2 minimum zone, although pelagic aerobic remineralisation still has a significant effect on bottom O 2 , but accounts for a remarkably lower proportion than in the entire sub-MLD volume.The simulated benthic remineralisation rates are in the same order as those derived from observations giving a range from 7 to 25 mmol O 2 m −2 d −1 for a nearby station (station 3 in Upton et al., 1993), however, rather at the lower end of this range.The O 2 concentrations at the beginning of the stratified period in regions 4 (9.74 mg O 2 L −1 ) and 5 (9.46 mg O 2 L −1 ) are lower than in region 3.In contrast, both regions show higher concentrations at the end of stratification compared to region 3.These values reach 6.98 mg O 2 L −1 in region 4 and 7.61 mg O 2 L −1 in region 5, compared to 6.76 mg O 2 L −1 in region 3.

Spatial variability in the
In region 4, intense mixing in late June/early July, indicated by the steep increase in MIX O 2 , causes the breakdown of stratification and the complete replenishment of bottom O 2 .Integrated over the stratified period the effect of MIX O 2 is almost 1.5 times higher than in region 3.In region 5, MIX O 2 represents a significantly lower supply of O 2 , which mainly relates to the significantly greater water depth favouring more stable stratification.ADV O 2 has an opposite effect in regions 4 and 5 relative to region 3, however, showing only minor negative integrated effects.
The integrated biological O 2 consumption in region 4 is about 7 % less than in region 3. Referring to changes in O 2 concentrations, the consumption even exceeds that in region 3 by about 6 %, due to the thinner bottom layer, and corresponds to an EXP org below 25 m depth (calculation analogous to Table 1) in region 4 being 19 % higher than in region 3.In the deeper region 5, EXP org accounts for 62 % of that in region 3, while biological O 2 consumption accounts for only 35 % of that in region 3. Despite the differences in the overall biological O 2 consumption, the relative contributions of the different sink processes in region 4 are in the same order as in region 3. REM sed represents the largest contributor with about 54.6 % of the total biological consumption, while REM pel accounts for 27.2 %.Thus, the combined effect of REM sed and REM pel accounts for 81.8 % which is similar to region 3. Average daily benthic remineralisation rates are higher than in region 3 and yield 8.9 mmol O 2 m −2 d −1 .The relative contributions for RES zoo and NIT result in 13.6 and 4.6 %, respectively.
The comparison of the relative contribution of REM pel and REM sed in region 5 reveals some changes compared to regions 3 and 4. REM sed accounts for about 70 % of the total biological O 2 consumption, while REM pel contributes to only about 20 %.This relates to the generally lower amount of exported organic matter reaching the model bottom layer (63 % of that in region 3).On the one hand, this causes lower O 2 consumption due to REM pel , and on the other hand, enhances REM sed relative to REM pel as more organic matter reaches the bottom.
Considering the combined effect of stratification and biological consumption in region 4 reveals that the shorter stratification period and the strong mixing prohibit the evolution of low O 2 conditions in this region, despite the highest biological consumption.The higher benthic remineralisation rate compared to region 3 underlines the high potential for low O 2 conditions in the Oyster Grounds under persistent seasonal stratification.This is in good agreement with the findings by Greenwood et al. (2010), who observed bottom O 2 concentrations less than 6 mg O 2 L −1 in this area.As in region 3, the simulated benthic remineralisation rate lies at the lower end of the range of 5.6 to 30.6 mmol O 2 m −2 d −1 obtained from observational studies near this site (Lohse et al., 1996;Weston et al., 2008).This suggests that the model most likely underestimates benthic remineralisation rates.
In region 5, the amount of exported organic matter reaching the bottom layer (deepest pelagic model layer) is limited due to the great water depth.Thus, biological consumption in the bottom layer is low, preventing the evolution of O 2 deficiency, even though stratification lasts longer and is more stable than in the other regions.This suggests that region 5 is unlikely to be affected by low bottom O 2 concentrations.However, Queste et al. (2013) found bottom O 2 concentrations of about 6 mg O 2 L −1 near this area in 2010, which indicates that this area can be affected by O 2 deficiency.

Interpreting observed bottom O 2 at North Dogger
The O 2 observations at station North Dogger (Greenwood et al., 2010) shown in Fig. 3a and b showed similar O 2 concentrations of about 9.5 mg O 2 L −1 at the beginning of the The temporal evolution of the biological consumption processes is also similar in both years, with higher rates in 2007.The integrated effect of all biological sink processes results in −40.6 g O 2 m −2 , which corresponds to an average consumption rate of 0.18 g O 2 m −2 d −1 .For 2008, the simulation yields a 6.8 g O 2 m −2 lower biological consumption due to a 0.02 g O 2 m −2 d −1 lower consumption rate.This constitutes a relative difference of 13 % between the 2 years.For EXP org below 25 m depth during summer (calculation analogous to Table 1), the same relative difference is found.
The relative contribution of the different processes to the overall bottom O 2 consumption shows only minor changes between the 2 years.In both years REM sed accounts for about 58 %, while REM pel accounts for about 31 %.RES zoo and NIT contribute to about 3 and 8 % in both years, respectively.This shows that the variations in the O 2 dynamics between 2007 and 2008 at North Dogger are mainly driven by differences in the organic matter export.The steeper de-crease in bottom O 2 in 2007 results from the enhanced supply of organic matter and the subsequent increased degradation by bacteria.The enhanced release of ammonium due to pelagic and benthic remineralisation consequently triggers an increase in nitrification.
In contrast to the findings by Greenwood et al. (2010), who argued that the relatively strong advection at North Dogger may ventilate the bottom layers in terms of O 2 , our results suggest that advection only has a minor positive effect due to the only slightly higher O 2 concentrations in the surrounding waters.The large contribution of bacterial remineralisation (REM sed and REM pel ) accounting for almost 90 % of the overall biological consumption at station North Dogger confirms that the estimates for C remineralisation rates made by Greenwood et al. (2010) provide reasonable results.However, as NIT also accounts for about 8 % of the total biological O 2 consumption, this process should be considered to obtain more precise estimates for C remineralisation rates.

Conclusions and perspectives
The North Sea is one of the shelf regions regularly experiencing seasonal O 2 deficiency in the bottom water (Diaz and Rosenberg, 2008;Emeis et al., 2015;Rabalais et al., 2010).However, not all areas of the North Sea are similarly affected by low O 2 conditions (e.g., Queste et al., 2013) due to different characteristics with respect to stratification and O 2 consumption.Observations and model results suggest that the area between 54-57 • N and 4.5-7 • E shows the highest potential for low O 2 conditions, but also areas around the Doggerbank experience lowered bottom O 2 concentrations.The model-based analysis of different factors affecting O 2 showed that besides sufficiently long stratification (> 60 days), surface layer primary production (driving organic matter export) and sub-thermocline volume are the key parameters influencing the bottom O 2 evolution.Based on this, the North Sea can be subdivided into three different zones in terms of O 2 dynamics: (1) a highly productive, nonstratified coastal zone (region A), (2) a productive, seasonally stratified zone with a small sub-thermocline volume (region B), and (3) a productive, seasonally stratified zone with a large sub-thermocline volume (regions C and D).While the zones of types 1 and 3 are unlikely to be affected by O 2 conditions due to either continuously ongoing ventilation (type 1) or the large sub-thermocline volume diluting the effect of O 2 consumption (type 3), type 2 is highly susceptible to low O 2 conditions.This results from the specific combination of high upper layer productivity and small sub-thermocline volume, which causes a strong impact of the consumption processes on the decrease in the bottom O 2 concentrations.
The ODI demonstrates that this regional characterisation, based on only three controlling parameters, can be applied to most parts of the North Sea.The ODI is rather simple compared to the eutrophication risk index (EUTRISK; Druon et al., 2004) as it is designed to indicate regions with higher risk for O 2 deficiency.Therefore, it may also allow for an operational use as the information on stratification can be derived from operational hydrodynamical models and information on net primary production from satellite data.
The model-based mass balances showed that pelagic bacterial remineralisation constitutes the largest O 2 consuming process within the sub-thermocline volume.In the bottom layer, benthic remineralisation constitutes the major O 2 sink, which consistently contributes more than 50 % to the overall bottom O 2 consumption.Pelagic remineralisation consistently contributes to more than 20 % of the overall bottom O 2 consumption.Zooplankton respiration and nitrification are less important, however, can contribute to up to 14 and 8 %, respectively.In addition, the results suggest that the relative contribution of the different O 2 consuming processes in the bottom layer at a certain location depends on the water column depth and is independent of the overall consumption.
The mass balances also showed that differences in the surface layer primary production drive variations in the sub-thermocline and bottom O 2 evolution between different years.Increased primary production directly enhances the export of dead phytoplankton into the deeper layers.Furthermore, it enhances zooplankton growth which causes an additional increase in organic matter production and export.This enhanced zooplankton growth further increases O 2 consumption due to respiration.The overall increase in organic matter export results in stronger bacterial remineralisation which in turn triggers nitrification due to the stronger release of ammonium.
Our analysis suggests that advection usually only has a minor effect on the bottom O 2 dynamics in most North Sea regions which contradicts the interpretation by Queste et al. (2013).However, during years of especially low bottom O 2 concentrations it may play an important role for the recovery of the O 2 levels before the breakdown of stratification in autumn.In addition, we showed that during the summer period only very strong mixing results in a net increase in O 2 as the enhanced nutrient supply triggers primary production, eventually increasing the biological O 2 consumption, which balances or even exceeds the enhanced O 2 supply.

F. Große et al.: Looking beyond stratification
This study demonstrated that ecosystem models are capable of describing the key features of the O 2 dynamics as an integral part of the North Sea ecosystem.This, in combination with the provision of a spatially and temporally consistent picture is useful for the detection of regions susceptible to low O 2 conditions which therefore require enhanced management.Additionally, this is of importance for monitoring authorities as our model showed that bottom O 2 measurements taken in late summer provide a synoptic picture of the North Sea O 2 conditions.
This study provides a general characterisation and process-based analysis of the North Sea O 2 dynamics in its present state, including the actual eutrophication status in the "continental coastal region" as defined for the OSPAR assessment (Claussen et al., 2009).The question on the anthropogenic contribution to the O 2 deficiency problem in relation to elevated nutrient supply is beyond the scope of this study.However, the capability of three-dimensional models to describe the O 2 dynamics in the context of the natural variability of the ecosystem can be related to changes in anthropogenic drivers, such as increased atmospheric deposition (Troost et al., 2013) or riverine nutrient input (Lenhart et al., 2010).
Similar model studies are essential for the assessment within the Water Framework Directive (WFD), in which dissolved O 2 is used as a key parameter (Best et al., 2007).As the WFD assessment depends strongly on the description of pristine conditions, related to natural nutrient levels (Topcu et al., 2009), ecosystem models can provide a consistent picture of the North Sea O 2 dynamics under these pristine conditions and of the effects of WFD reductions (Schernewski et al., 2015).As river load reductions within the WFD regulation affect the entire North Sea ecosystem, also in terms of bottom O 2 conditions, these scenarios can also be interpreted within the frame of the Marine Framework Directive (MSFD), which involves O 2 as one of the main descriptors for the definition of the "Good Environmental Status".
Recent observational studies explicitly highlight the importance of organic nutrient loads on the O 2 dynamics in the context of nutrient reductions (Kemp et al., 2009;Topcu and Brockmann, 2015).Thus, future modelling studies on the effects of nutrient reductions on the marine environment should differentiate between the effects of organic and inorganic nutrient inputs in order to optimise measures in the catchment area with respect to cost efficiency.Diaz and Rosenberg (2008) report an exponential expansion of global O 2 deficiency (and hypoxia) since the 1960s and argue that future changes in O 2 conditions will strongly depend on the effects of climate change on stratification and riverine nutrient supply.For the North Sea, several model studies predict a rise in water temperature (e.g., Lowe et al., 2009;Meire et al., 2013;Mathis and Pohlmann, 2014) which will reduce the O 2 solubility (Weston et al., 2008).Stratification intensity may either increase (Lowe et al., 2009;Meire et al., 2013) or even decrease (Mathis and Pohlmann, 2014), implying opposed effects on bottom O 2 .Primary production could increase due to enhanced nutrient supply caused by changes in weather conditions (Rabalais et al., 2010) or due to a temperature-driven increase in metabolic rates (van der Molen et al., 2013), which could eventually aggravate the O 2 conditions (Justić et al., 2003).In contrast, Gröger et al. (2013) predicted a North Sea wide reduction in primary production by about 30 % due to reduced winter nutrient import from the Atlantic.As these potential changes in the O 2 conditions will also affect the biocoenosis of the North Sea (Emeis et al., 2015), it is important to foster the analysis of potential impacts of climate change and changes in nutrient loads on the O 2 dynamics.

Data availability
The time series data from the Cefas station North Dogger can be accessed via the Cefas Data Hub (https://www.cefas.co.uk/cefas-data-hub/) according to Greenwood et al. (2016).The time series data of MARNET station Ems can be retrieved from the Deutsches Ozeanographisches Datenzentrum (DOD Data Centre) via email query to dod@bsh.de.The spatially resolved data from the North Sea cruises in 2001 and 2005 have been released in the framework of the EU-FP6 project CARBOOCEAN.These data can be accessed via the CARBOOCEAN data portal (http://dataportal.carboocean.org/)according to Thomas et al. (2012) and Thomas and Borges (2012), respectively.The data of the North Sea cruise 2008 (R/V Pelagia 64PE294; Zemmelink, 2008) have not been published, yet, but can be requested via the CODIS data portal (http://www.nioz.nl/portals-en;registration required).The main assumption behind the rather small critical T difference of 0.05 K is, that even if the surface mixed layer is interrupted due to mixing, this does not necessarily result in a complete overturning of the water column.Thus, even a minor difference in T indicates a bottom layer unaffected from vertical mixing.Despite this significant difference to common MLD criteria (e.g Kara et al., 2000, therein Table 1), the criterion applied in this study represents the stratification conditions quite well.However, it should be noted that the end of the stratified period may be slightly overestimated.In addition, in regions with a less pronounced onset of stratification, i.e., a less distinct increase in surface T , the determined timing of the onset may be slightly too early.The use of the maximum T gradient to determine the MLD under stratified conditions yields reasonable results, and is closely related to real conditions as the thermocline is defined as the layer with the maximum T gradient.

Figure 1 .
Figure 1.Extent of observed O 2 concentrations < 6 mg O 2 L −1 in the German Bight area from 1980 to 2010.Dotted lines indicate geographical limits of the individual surveys.Light grey line marks German Maritime Area (from Topcu and Brockmann, 2015).

Figure 2 .
Figure 2. Horizontal grid and bottom topography of the HAMSOM-ECOHAM model domain.White numbers indicate depth levels.Yellow boxes A-D mark the 4 × 4 regions used for the characterisation of key features presented in Sect.3.3.Black-filled boxes (1, 2) mark the validation sites discussed in Sect.3.1.1.Redframed boxes (2-5) indicate regions used for the O 2 mass balance calculations in Sects.3.4-3.7.

Figure 3
Figure 3 shows the comparison of simulated bottom O 2 against time series data at the Cefas station North Dogger for the years 2007 (a) and 2008 (b) and the MARNET station Ems during 2010 (c) and 2011 (d).The indicated stratification period was derived from the simulated temperature fields using Eq.(1).At North Dogger, observed and simulated bottom O 2 concentrations show a steady decrease after the onset of stratification.While stratification according to Eq. (1) starts a bit earlier compared to that described byGreenwood et al. (2010), the beginning of the decrease in bottom O 2 concentrations coincides well.The simulated and observed O 2 concentrations at this time are in good agreement.Some small-scale fluctuations in the observations are not fully reproduced by the simulation, however, the general evolution is represented well by the model.The average O 2 reduction in the simulation is slightly less than in the observations, visible in the difference between the concentrations at beginning and end of the stratified period.Stratification ends a bit earlier in the simulation, with the result that simulated bottom O 2 starts to recover while the observed concentrations continue to decline.The observed O 2 concentration at the end of the stratified period is about 6.8 mg O 2 L −1 , while the simulation results in about 7.4 mg O 2 L −1 .In 2008, we can see a similar slight overestimation of simulated O 2 concentrations, but less than in 2007.Some minor fluctuations in the observations are again not represented by the model, but the general evolution of bottom O 2 is represented well.It should be noted, that the different depths of the time series (76 m for simulation, 85 m for observation) may

Figure 3 .
Figure 3. Annual time series of observed and simulated bottom O 2 concentrations at Cefas station North Dogger in (a) 2007 and (b) 2008, and at MARNET station Ems in (c) 2010 and (d) 2011.Same legend for all panels.Grey shaded areas indicate the stratification periods derived from simulated T according to Eq. (1).

Figure 4 Figure 4 .
Figure 4 shows the spatial distribution of the average simulated and observed O 2 concentrations in the model bottom layer for the years 2001 (a), 2005 (b) and 2008 (c), and the SD related to the averages in 2005 (d).The averaging period for the simulations corresponds to the complete observation period for each year, listed in the bottom right corner of each panel.In 2001, the observations show the lowest concentrations of all years with minimum values of 5.9 mg O 2 L −1 in the area 54-57 • N, 4.5-7 • E. This minimum is similarly present in the model yielding 6.96 mg O 2 L −1 .Maximum observed concentrations were found off the southern tip of Norway (9.3 mg O 2 L −1 ) and in the deepest parts of the Norwegian Trench (8.7 to 8.8 mg O 2 L −1 ).The very high observed value at the northwesternmost sampling site represents an outlier

Figure 5 .
Figure 5.Taylor diagram of simulated (×) bottom O 2 concentrations compared to observations (OBS) for time series (see Fig. 3) at (a) Cefas North Dogger and (b) MARNET Ems, and (c) spatially resolved data (see Fig. 4).Standard deviations and centred rootmean-square differences (RMSD) were normalised by the standard deviation of the corresponding observations.
Figure 6a and b show the spatial distribution of the longest persistent stratification periods (after Eq. (1) using simulated T ) for the years 2002 and 2010, respectively.Both years show similar stratification patterns with stratification periods of > 180 days in large parts of the central and northern North Sea.Comparing the corresponding minimum concentrations of bottom O 2 (Fig. 6c and d) shows significant differences.The minimum bottom O 2 concentrations in 2002 in the region from 55-56.5 • N, 4.5-7.5 • E constitute the lowest O 2 concentrations during the entire period 2000-2012 reaching values of below 5.8 mg O 2 L −1 .In contrast, the duration of stratification in this area is similar or even longer in 2010 than in 2002.The O 2 concentrations in 2002 are even below the O 2 threshold applied by OSPAR (6 mg O 2 L −1 ; OSPAR-Commission, 2005) and persist for about one month (not presented).In contrast, 2010 represents a year with relatively high minimum bottom O 2 concentrations being above 7.3 mg O 2 L −1 in the entire model domain.The areas directly north and south of the Doggerbank also reveal lower bottom O 2 concentrations in both years.The stratification periods derived from the simulation are in good agreement with the different stratification regimes described byPingree et al. (1978) andvan Leeuwen et al. (2015).The latter applied a density-based stratification criterion on model results to subdivide the North Sea into areas of different stratification characteristics, and showed that most areas of the seasonally stratified central and northern North Sea reveal stratification periods of 170 to 230 days.The increased potential for low O 2 conditions north and south of the Doggerbank corresponds well to observed bottom O 2 time series in these regions(Greenwood et al., 2010).Queste et al. (2013) also observed lower bottom O 2 concen-

Figure 8 .
Figure 8. Mass balances of simulated O 2 in region 3 (see Fig. 2) during stratification (grey shaded): (a, b) for the entire volume below the MLD (grey dash-dotted; according to Eqs. (1) and 2), (c) and (d) only for the bottom layer for the years 2002 (a, c) and 2010 (b, d).Same legend for all panels.Black y axes apply to processes, magenta y axes apply to O 2 (saturation) concentrations.Values of black y axes in (a) and (b) also apply to MLD (unit: m).Changes in O 2 due to different processes are cumulative.Text boxes list relevant stratification parameters, average volume of the analysed water body and bottom depth.Note different y axes for (a), (b) and (c), (d).
Figure 9a and b show the bottom O 2 mass balances of regions 4 (southern North Sea) and 5 (northern North Sea) in 2002 (see Fig.2for location of regions).Both regions show different stratification periods than in region 3.In region 4, the period between the first and last day of stratification accounts for only 163 days compared to 187 in region 3. Additionally, stratification is temporarily intermittent in late April and early July.In region 5, stratification lasts for 211 days without any interruptions.The water depths and bottom layer volumes also differ between the regions.Region 3 has an average water depth of 47.75 m and a bottom layer volume of about 14.4 km 3 , while region 4 is characterised by an average depth of 45 m and a bottom layer volume of 11.9 km 3 .Region 5 has an average bottom depth of 99.23 m and a bottom layer volume of 16.3 km 3 .The O 2 concentrations at the beginning of the stratified period in regions 4 (9.74 mg O 2 L −1 ) and 5 (9.46 mg O 2 L −1 ) are lower than in region 3.In contrast, both regions show higher concentrations at the end of stratification compared to region 3.These values reach 6.98 mg O 2 L −1 in region 4 and 7.61 mg O 2 L −1 in region 5, compared to 6.76 mg O 2 L −1 in region 3.In region 4, intense mixing in late June/early July, indicated by the steep increase in MIX O 2 , causes the breakdown of stratification and the complete replenishment of bottom O 2 .Integrated over the stratified period the effect of MIX O 2 is almost 1.5 times higher than in region 3.In region 5, MIX O 2 represents a significantly lower supply of O 2 , which mainly relates to the significantly greater water depth favouring more stable stratification.ADV O 2 has an opposite effect in regions 4 and 5 relative to region 3, however, showing only minor negative integrated effects.The integrated biological O 2 consumption in region 4 is about 7 % less than in region 3. Referring to changes in O 2 concentrations, the consumption even exceeds that in region 3 by about 6 %, due to the thinner bottom layer, and corresponds to an EXP org below 25 m depth (calculation analogous to Table1) in region 4 being 19 % higher than in region 3.In the deeper region 5, EXP org accounts for 62 % of

Figure 9 .
Figure 9. Mass balances of simulated bottom O 2 in regions 4 (a) and 5 (b); (see Fig. 2) during stratification (grey shaded) in 2002.Same legend for (a) and (b).Black y axes apply to processes, magenta y axes apply to O 2 (saturation) concentrations.Changes in O 2 due to different processes are cumulative.

Figure 10 .
Figure 10.Mass balances for simulated bottom O 2 at Cefas North Dogger (see Fig. 2, region 2) during stratification (grey shaded) in (a) 2007 and (b) 2010.Same legend for (a) and (b).Black y axes apply to processes, magenta y axes apply to O 2 (saturation) concentrations.Changes in concentrations due to different processes are cumulative.

Figure A1 .
Figure A1.Hovmöller diagram of simulated T and MLD according to Eqs. (1) and (2) at Cefas station North Dogger (see Fig. 2, region 2) in 2007.Depth levels represent the centre depth of model layers.

Table 1 .
Average critical quantities (2000Average critical quantities ( -2012) )characterising the O 2 dynamics in the four different 4 × 4-regions (see Fig.2, yellow boxes).Fluxes are cumulated from 1 April to 30 September and relate to a surface layer of thickness D ref = 25 m.
sidering the biological sink processes, pelagic remineralisation of organic matter (REM pel ; dashed light green) has the strongest effect on the sub-MLD O 2 , accounting for 50 % of the overall biological consumption.Benthic remineralisation (REM sed ; dashed yellow) accounts for 18.8 %, while zooplankton respiration (RES zoo ; dashed dark green) and nitrification (NIT; dashed red) contribute 22 and 8.8 %, respectively.This order in the relative importance is consistent throughout the entire period 2000-2012 (not shown).In 2002, REM pel is strongest among all years yielding −103.1 g O 2 m −2 , while 2010 represents the year of weakest REM pel .For RES zoo , 2002 yields a value of −45.1 g O 2 m −2 (1.8-fold of 2010 value).The 2002 and 2010 values constitute the highest and lowest among all years, respectively.The same applies to REM sed with a 2002 value of −28.5 g O 2 m −2 (1.2-fold of 2010 value).NIT is also strongest in 2002 resulting in −18.1 g O 2 m −2 , while in 2010 it is about 13 % below the average value of show that the biological processes cause a net loss in O 2 .The only source process for O 2 is primary production (PP; dashed light blue) which causes a gross increase of about 61.8 g O 2 m −2 in 2002 (1.4-fold of 2010 value).Conwww.biogeosciences.net/13/2511/2016/Biogeosciences, 13, 2511-2535, 2016