Retracing Hypoxia in Eckernförde Bight (Baltic Sea)

. In recent years, upwelling events of low-oxygenated deep water have been repeatedly observed in Eckernförde Bight (EB) situated in the Baltic Sea, Germany. Many of these events were related to massive ﬁsh-kill incidents - with negative consequences for commercial ﬁsheries and tourism. The aim of this study is to dissect underlying mechanisms and to explore the potential of existing monitoring programs to predict these events. Our main tool is an ultra-high spatially resolved general ocean circulation model which drives an elementary representation of the biogeochemical dynamics of dissolved oxygen 5 (dubbed MOMBE and EckO 2 , respectively). In addition, we integrate artiﬁcial "clocks" that measure the residence time of the water in EB along with timescales of (surface) ventilation. We present an ensemble of hind cast model simulations, covering the period from 2000 up to 2018, designed to cover a range of poorly known model parameters for vertical background mixing (diffusivity) and local oxygen consumption within EB. Our results indicate that the dynamics of low (hypoxic) oxygen concentrations in bottom waters deep inside EB is, to ﬁrst order, determined by the following antagonistic processes: (1) the 10 inﬂow of low-oxygenated water from the Kiel Bight (KB) - especially from July to October and (2) the local ventilation of bottom waters by local (within EB) subduction and vertical mixing. Biogeochemical processes that consume oxygen locally, are apparently of minor importance for the development of hypoxic events. Reverse reasoning suggests that subduction and mixing processes in EB contribute, under certain environmental conditions, to the ventilation of the KB by exporting recently-ventilated waters enriched in oxygen. A detailed analysis of the 2017 ﬁsh-kill incident highlights the interplay between westerly 15 winds importing hypoxia from KB and ventilating easterly winds which subduct oxygenated water. Finally, we explore the capabilities of - comparably computationally cheap - feed-forward artiﬁcial neuronal networks to forecast hypoxia deep in EB based on data at a monitoring site at the entrance of EB. The Boknis Time-Series Station is run by the Oceanography of the for 415 Research Kiel. The data from Boknis Eck are available from www.bokniseck.de/database-access. Author contributions. H. Dietze and U. Löptien have been equally involved in setting up and running the model conﬁgurations. Both authors contributed to the interpretation of model results, to outlining and writing of the paper in equal shares. The

to the well-ventilated surface layer) need to be known. Unfortunately, both, the respiration rates and transport timescales or subsurface residence timescales, are difficult and expensive to measure in situ. This makes this type of observations very rareeven in EB, which is renowned for its good observational data coverage. Our approach to overcome the respective limitations in this study is to integrate an ensemble of a high-resolution coupled ocean-circulation biogeochemical model configurations, that test through a range of mixing parameters (which determine residence timescales) and through respective parameter values for the biogeochemical oxygen sources and sinks. The ensemble is assessed with observations of salinity, temperature and dissolved oxygen measurements deep inside EB. The most realistic ensemble members are then analyzed in greater detail in order to dissect a mechanistic understanding of the processes involved in the dynamics of dissolved oxygen. Finally, we build an artificial neuronal network (ANN) in order to forecast dissolved oxygen concentrations deep in EB based on measurements at the entrance of EB. This approach yields a computationally cheap surrogate to the (relatively) computationally expensive 70 coupled ocean-circulation biogeochemical model to the stakeholder (cf., Figure 5). Key predictors are identified by systematic feature selection.

MOMBE (Modular Ocean Model Bight of Eckernförde) is a new configuration of a general ocean circulation model (GCM).
The GCM is coupled to a simple representation of biogeochemical processes by introducing an additional passive tracer, that is 75 advected and mixed just like the tracers temperature and salinity but, other than that, controlled by prescribed rates of oxygen production and consumption. Further, we introduce artificial tracers or "clocks" that estimate the residence times and the age (i.e. the time of last contact to the surface) of water parcels. This approach facilitates the dissection between local (i.e. inside EB) and remote (e.g., inflowing hypoxic deep water from the KB) processes that drive the oxygen dynamics. The following subsections describe the GCM, followed by a model evaluation in Section 3. The feed-forward neuronal networks designed to 80 mimic the full-fledged coupled GCM at a station deep in the Bight are described in Section 4.4.

Model Configuration
We use the Modular Ocean Model framework MOM4p1 as released by NOAA's Geophysical Fluid Dynamics Laboratory (Griffies, 2009). The model code and configurations are almost identical to those described in Dietze et al. (2014) and Dietze et al. (2020). The few exceptions are listed in the following subsections. Section 2.1.1. describes the model grid, Section 2.1.2 85 the subgrid parameterizations, and Section 2.1.3 specifies the input data (boundary conditions). Section 2.1.4 documents the representation of sea ice, Section 2.1.5 introduces the implementation of the residence time and age racers. The implementation of the oxygen module is documented in Section 2.1.6.

Grid and Bathymetry
The bathymetric data are provided by the Federal Maritime and Hydrographic Agency (BSH, https://www.geoseaportal.de/ 90 mapapps/resources/apps/bathymetrie/index.html?lang=de). We use a bilinear scheme to interpolate these onto an Arakawa 3 https://doi.org/10.5194/bg-2021-31 Preprint. Discussion started: 12 February 2021 c Author(s) 2021. CC BY 4.0 License. B model grid (Arakawa and Lamp, 1977). There are 165×103 grid boxes horizontally, each about 100 m × 100 m in size ( Figure 4). The total wet area of the model is 119 km 2 . The vertical resolution is 1 m, with a total of 31 layers. The average water depth is 11.7 m. The bathymetry was smoothed with a filter similar to the Shapiro filter (Shapiro, 1970). The filter weights are 0.25, 0.5 and 0.25. The filter essentially fills steep holes in the ocean floor which increases numerical stability of the GCM. 95 The filter was successively applied three times, as this has proven (in Dietze and Kriest, 2012;Dietze et al., 2014Dietze et al., , 2020 to be a good compromise between unnecessary smoothing on the one hand and numerical instability on the other hand.

Subgrid Parameterisations
Even a horizontal resolution as high as 100 m horizontally and 1 m vertically fails to explicitly resolve all (turbulent) processes of relevance for transport and mixing of substances in EB. Hence, effects of unresolved small-scale processes have to be 100 parameterized. We use parameterizations and setting identical to those applied by Dietze et al. (2014) in a high-resolution model configuration of the Baltic Sea. An exceptions it the parameter choice for the vertical background diffusivity: Holtermann et al.
(2012) estimates from measurements for deep water processes in the central Baltic Sea a vertical diffusivity of 10 −5 m 2 s −1 (calculated from the propagation speed of a purposely-deployed dye-like substance). Closer to coast Holtermann et al. (2012) report much higher values. Because mapping this information on conditions in EB is difficult, we decided to test a range 105 of vertical background diffusivities and to assess the respective model perfomances based on available observations. The considered diffusivities are: 5 × 10 −5 m 2 s −1 , 1 × 10 −4 m 2 s −1 und 5 × 10 −4 m 2 s −1 . This range comprises relatively low diffusivities, which are characteristic for the deep central Baltic Sea, and fairly high values, which are more representative for coastal mixing (as can be expected in the shallow Eckernförde Bight).

110
The atmospheric boundary conditions of our model are set by a reanalysis from the Swedish Meteorological and Hydrological Institute (SMHI). We use the results of the reanalysis framework as a means to interpolate (patchy) observations in time and space. The underlying atmospheric model features a horizontal resolution of 11 km. For the period 2000 to 2015 we use RCA4 (Samuelsson et al., 2015(Samuelsson et al., , 2016. RCA4 data is available only until 2015. Hence, for the period 2016 to 2018 we switched to another product: UERRA (regional reanalysis for Europe; https://cds.climate.copernicus.eu/cdsapp#!/dataset/ 115 reanalysis-uerra-europe-complete?tab=overview). UERRA is more advanced but does not include "spectral nudging" to the large-scale atmospheric circulation. This detail may allow for unrealistic shifts in the trajectories of low pressure systems. Fortunately, for the time and location under consideration here, a rough comparison with the observations from Kiel lighthouse (in position 54.3344 • N,10.1202 • E) showed a generally good agreement between reanalysis and direct observations (not shown).
Our model configuration features rigid walls in the east, where EB is connected to the KB. We mimic the respective water 120 exchange by restoring to prescribed temperature, salinity and sea surface height values at the model boundaries. For sea surface height we restore to prescribed values taken from an oceanic reanalysis carried out with MOMBA (Dietze et al., 2014). MOMBA differs from MOMBE in that it covers the entire Baltic Sea with a horizontal resolution of 1 nautical mile while MOMBE introduced here covers the EB only -albeit with much higher resolution (100 m). For the sake of consistency, 4 https://doi.org/10.5194/bg-2021-31 Preprint. Discussion started: 12 February 2021 c Author(s) 2021. CC BY 4.0 License. differs from Dietze et al., 2014). For temperature, salinity and oxygen we restore MOMBE at its horizontal boundaries with Kiel Bight to interpolated measurements from Station Boknis Eck at the entrance of EB (Lennartz et al., 2014, http://www. bokniseck.de/, http://doi.pangaea.de/10.1594/PANGAEA.855693).

Sea Ice
The focus of our investigation are ice-free seasons. We will show in Section 4.1 that the memory of the system under con-130 sideration, as given by residence times in Eckernförde Bight, is less than a month. This suggests that sea-ice dynamics are rather irrelevant to the processes and seasons examined here. Even so, for the sake of completeness, we report that our ocean component is coupled to a dynamical sea ice module, the GFDL Sea Ice Simulator (SIS). SIS uses elastic-viscous-plastic rheology adapted from Hunke and Dukowicz (1997). We use the exact same settings described in Dietze et al. (2020) (which are identical to the settings in Dietze et al. (2014), except for switching to levitating sea ice).

Artificial Clocks
In order to facilitate the dissection of local versus remote processes influencing the oceanic oxygen concentrations in EB, we introduce two artificial tracers or "clocks" to the ocean circulation model (following and approach similar to Dietze et al., 2009). Both clocks behave like dyes in that they are subject to transport processes just like like temperature, salinity and dissolved oxygen. In addition to being transported, the clocks continuously count up time in every grid box. The first clock is 140 reset to zero whenever a water parcel reaches the ocean surface. Thus, it measures the time elapsed since a water parcel had been in contact with the atmosphere. This time is also referred to as the age of the water. The second clock is reset to zero at the eastern boundaries of the model domain. Thus, it measures the time elapsed since water entered EB. This time is also referred to as the residence time of water in EB.
The ratio between residence time and age is a measure of the importance of local processes versus remote processes: if a 145 water parcel remains much longer in EB than the time has passed since the water parcel has been ventilated locally in EB, then this may be an indication of the dominance of local (inside EB) biological respiratory processes. Conversely, if the residence time is much shorter than the age, then the interplay between the inflowing water and the local ventilation of this water (by "upwelling" or mixing to the surface, where it is exposed to air-sea gas exchange) dominates.

150
Our dissolved oxygen module is dubbed EckO 2 -module (from Eckernförde O 2 ). The module is very similar to the approach of Bendtsen and Hansen (2013) dubbed OXYCON. A schematic representation of EckO 2 is given in Figure 6. Following Bendtsen and Hansen (2013), the local development over time of dissolved oxygen, ∂O2 ∂t , is defined by: where A und D denote the divergence of the three-dimensional advective and diffusive fluxes as calculated by the GCM. S 155 denotes biogeochemical oxygen sources and sinks given by the model parameters opro at the sunlit sea surface, by orewa at depth below the compensation depth zco, and by orese in the lowermost wet model grid box. These parameters determine how much oxygen is generated by primary production (opro) and how much is consumed at depth (orewa) and in the sediment (orese). The respective parameter choices are based on literature values listed in Table 1. Following Babenerd (1991) and based on AErtebjerg et al. (1981) and Jacobsen (1982) we assume that the subsurface oxygen consumption rates are rather 160 uniform throughout KB, EB and up into the Danish Straits. This assumption is necessitated by our lack of direct measurements of consumption rates in EB. EckO 2 prescribes climatological monthly mean consumption rates.

Observations
We use data from the regular monitoring program of the LLUR. Respective approx. monthly observations of salinity, temperature and oxygen covered the entire hind-cast period at the monitoring station Buoy 2a (location marked in Figure 4).

Model Evaluation
Among the challenges in simulating oxygen dynamics is that both biotic parameters (determining oxygen respiration (Section 2.1.6)), and the antagonistic abiotic parameters (that control ventilation with surface water high in oxygen such as e.g.
vertical diffusivity (Section 2.1.2)) are uncertain. Our approach is to run an ensemble of simulations encompassing a plausible range of settings. These settings are listed in Table 2. We compare low, medium and high levels of diffusivity (tagged HighMix,

170
MedMix, LowMix, respectively) and a best guess of local biotic processes versus ignoring local biotic processes altogether (tagged Rem, indicating "remote" forcing of hypoxia only). This section identifies the most realistic simulation(s). Figure 7 shows Taylor diagrams which compare simulated and observed temperature, salinity and oxygen. The simulations with high diffusivity (HiMix and HiMixRem) feature the lowest performance in reproducing the observed variability of temperature, salinity and oxygen. This is consistent with an assessment of simulated velocities by Marlow (2020). The more realistic 175 simulations LoMix and HiMix are very similar -irrespective of wether we account for local sources and sinks of oxygen or not. We conclude (from Figure 7) that the lower values for the diffusivity are more realistic and that local sources and sinks of oxygen are apparently of minor importance within EB. This suggests that hypoxic events in EB are "imported" rather than driven by local oxygen consumption.

Results
As a first step, we explore the simulated residence and ventilation timescales (Section 4.1) which provide a base for understanding the dynamics behind our hind cast, presented in Section 4.2. A complementary case study of the intense hypoxic 190 event 2017 (related to the mass fish kill incidence depicted in Figure 2) is presented in Section 4.3. Section 4.4 describes the application of artificial intelligence for feature selection and extraction of the predictive capability of monitoring data at Station Boknis Eck at the entrance of EB to forecast hypoxia within EB at the monitoring station Buoy 2a.

Residence and Ventilation Times
The estimates of residence and ventilation times are calculated with "artificial clocks", as described in Section 2.1.5. Both The distribution of ventilation times or age is similar to that of residence times in that the highest values are generally found within the Bight towards Eckernförde ( Figure 12). The horizontal gradient is more pronounced in the simulation with lower mixing, while higher prescribed vertical background mixing equalizes the effective ventilation processes horizontally. In terms 205 of vertical distribution age has, in contrast to the residence time, high values at depth and low at the surface -where it is reset to zero ( Figure 13).
In summary, we find that residence times and age are of similar magnitude. This suggests that the first order control of processes that determine oxygen concentrations in EB is an antagonistic interplay of inflowing water (generally low in oxygen) and the local aeration by vertical exchange with oxygenated surface waters. Biogeochemical processes in the interior of EB are 210 apparently of minor importance for the oxygen dynamics within EB. issue Behrenfeld, 2010;Arteaga et al, 2020;Smetacek, 1985). The biomass eventually sinks to depth where it degrades and issues oxygen consumption. Later in the season, the water column stratifies and the surface layer heats up, effectively creating a barrier to the exchange of bottom water (deprived in oxygen) and the oxygenated surface waters. As autumn approaches, the surface ocean cools again and weakens the stratified barrier to vertical mixing. This facilitates the wind-driven mixing events 220 that come along with more unstable autumn weather. In winter, convective mixing homogenizes the entire (rather shallow) water column vertically (e.g., Fennel and Testa, 2019;Petenati, 2017). Apparently the model captures this dynamic well, i.e., the ensemble mean of LoMix and MedMix features a high visual correspondence between the respective curves in Fig. 8 (see October than in August ( Figure 15). Hence, water entering EB from KB in October are more likely to "import" hypoxia.

The Typical Seasonal Cycle inside EB
We conclude: the typical oxygen deficit in late summer is imported along with water from the KB, rather than being produced locally in EB. The following Section 4.3 will elucidate the underlying succession of events by means of a detailed case study.

Hypoxic Event 2017 235
In fall 2017 a particularly pronounced hypoxic event occurred and led to a mass fish kill incidence ( Figure 2). In the following, we analyze this event in the MOMBE LoMix simulation. Figure 16 shows a sequence of snapshots of simulated hypoxia in EB, starting August 20th and ending at peak conditions on September 10th. Over the course of these several weeks, EB looses oxygen and hypoxic waters apparently enter the Bight at the bottom from the east and moves upwards. The notion of "imported" hypoxic conditions is backed by the Hovmoeller 240 Diagrams of simulated age and residence times at the monitoring station Buoy 2a in Figure 17: during the buildup of the hypoxic event in EB, the residence time features a local minimum deep inside EB. This suggests the prevalence of water masses "recently-imported" from KB (Figure 17 b). Simultaneously, the age features a maximum (Figure 17 a), indicating that the "recently-imported" hypoxic waters are well-shielded from ventilation by oxygenated surface waters. Further evidence is provided by Figure 18, showing that the oxygen decline in EB is contemporaneous with winds blowing out of the Bight. These 245 winds drive an overturning circulation, shown in Figure 19, with surface waters being pushed out of the Bight and bottom waters, for continuity reasons, being sucked into the Bight at depth. Consequently, we find in Figure 18 that the oxygen decline at the entrance of the Bight (at Station Boknis Eck) occurs earlier than the oxygen decline inside the Bight (at Station Buoy 2a) -just as expected in a system where water enters the Bight at the bottom. During the relaxation phase, that terminates the 2017 hypoxic event, the processes are reversed: Figure 20 shows that the winds are blowing consistently into the Bight for more than a week -a situation comparably uncommon in these latitudes of prevailing westerlies. Consequently, water is pushed into the Bight at the surface, having nowhere to go. Some of the welloxygenated surface water is subducted to depth and subsequently leaves EB at depth. Just as expected, the increase in oxygen at the monitoring station Buoy 2a inside the Bight occurs earlier than the corresponding oxygen increase at the entrance Station 255 Boknis Eck). The oxygen levels at Boknis Eck now lag behind Buoy 2a by approximately one week.
In summary, we identified a governing mechanism by which EB is -depending on wind direction -either: (1) impacted by imported low oxygenated waters from KB or (2) being flushed by oxygenated surface water, that is subducted to depth in the interior of EB and is exported at depth to KB -whereby EB is effectively ventilating KB.
Open question, however, remain. Of particular interest is the questions why some years are hit especially hard by hypoxia 260 and wether such events are predictable days or weeks in advance. Such predictions may, e.g., allow for netting and landing of doomed fish. The following section applies Artificial Intelligence (AI) to pursue these questions.

AI-based feature selection and time series prediction
The following section explores the statistical relations between the simulated time series at Station Buoy 2a deep in the Bight and Boknis Eck at the entrance of the Bight. The major aims are: (1) To gain further mechanistic insight.
(2) To develop a sur-265 rogate models for the stakeholder that may be implemented on off-the-shelf desktop computers, smart phones or even on very low cost (< 10,-Euros) embedded devices rather than necessitating access to a super-computing facility (as is the case with the full-fletched coupled model). This section is motivated by recent and encouraging success in emulating general circulation models (e.g Castruccio et al., 2014), ecosystem models (e.g. Fer et al., 2018), the tremendous success in machine learning and data-driven methods in fluid dynamics (as summarized e.g. by  and the sneaking suspicion that " ... the 270 most pressing scientific and engineering problems of the modern era are not amenable to empirical models or deviations of first principles ..." (Brunton et al., 2020b).
In the following, we describe the application of shallow and deep feed-forward artificial neuronal networks (ANNs) to forecast bottom oxygen concentrations deep inside EB at the monitoring station Buoy 2a two weeks in advance from the 275 atmospheric conditions and the regularly sampled monitoring station Boknis Eck at the entrance of the Bight. The forecast range is chosen as a compromise between the time needed for mitigation measures (e.g. by netting and landing of doomed fish) and forecast accuracy which typically degrades with forecasting range. During the course of this exercise we will use different combinations of predictors (or input data) and test their impact on the forecast skill -a processes also referred to as capacity estimation and feature selection (e.g., Sbalzarini et al., 2002). Note, however, that a comprehensive analysis of time 280 series forecasting, which must include traditional statistical approaches in addition to machine learning methods (Makridakis et al, 2018), is beyond the scope of this manuscript. Buoy 2a (target data). Hand-designed features are "... two edged swords" (e.g. Reichstein et al., 2019): they can be seen as an advantage because they give us control of the explanatory drivers which may be used to promote system understanding. On the other hand, hand-designed features are typically suboptimal. To this end our results here provide a lower bound on the potential of ANNs for the task at hand.
The ANN is trained using the Levenberg-Marquardt algorithm (Marquardt, 1963) applied to neuronal network training fol-290 lowing (Hagan and Menhaj, 1994) and (Hagan et al., 1996). Each training is repeated 30 times, each of which may yield (slightly) differing results because: depending on the (random) initialization of weights, the algorithm may terminate in potentially differing local optima of the cost function. As cost-function we choose mean-squared errors (calculated from MOMBE output and the ANN prediction designed to mimic the MOMBE output). Figure 21 shows respective cost as errors relative to a naive biweekly persistency forecast based on bottom oxygen concentrations at the monitoring station Boknis Eck: apparently 295 the ANN's performance converges at 45% relative to the persistency forecast. Defining this as the Pareto Frontier suggests a Pareto Optimal of 56% -which corresponds to one or two nodes. The idea of opting for a rather parsimonious two-node model that scores 80% of the Pareto Frontier rather than 100% is to reduce the risk of overfitting which may hinder generalization.
Further, parsimonious models are easier to interpret than their complex counterpart such that their robustness is easier to assess. This is especially important because we have no straightforward way to extract human semantics from the "rules" the neuronal 300 network learned during the optimization process that related our input features to the target bottom oxygen concentrations at Station Buoy 2a.
We start with a shallow (one input, one hidden and one output layer) ANN utilizing the full vertical profiles of temperature, salinity and oxygen along with a biweekly wind forecast totaling at 106 input features (given by the three 1-m resolution vertical profiles of temperature, salinity and oxygen down to 26 m depth and the 14-daily forecasts of zonal and meridional 305 winds each). This setup is based on an optimistic estimate of the number of features available to stakeholders. Specifically, we assume to have access to a correct biweekly wind forecast along with one full vertical profile of each temperature, salinity and oxygen at the monitoring station Boknis Eck located at the entrance of EB (i.e., the 106 features introduced above). Figure 21 suggests that the Pareto Frontier is at 45% corresponding to a 55% reduction in error relative to the persistence model. 80% of this yields a Pareto Optimal of 56%. This corresponds to one or two nodes. Additional tests with deeper ANN's 310 featuring up to 10 hidden layers with two nodes were unsuccessful in that respective errors were always higher than 50%. We conclude that a simple two node shallow ANN features already a reasonable performance and two input features, of the 106 tested, my suffice to capture the main effects.  Table 3). Even so, the ANN fits the training and validation data remarkably well ( Figure 22). We conclude that the ANN's biweekly forecast exploits links other than those 320 being direct consequences of the wind driven inflow versus ventilation mechanism identified in Section 4.3. Section 4.4.2 puts this exploitation to the test using independent test (model) data.

ANN Generalization
This section discusses the fidelity of the two-node ANN using bottom and surface temperature identified in Section 4.4.2 as being parsimonious but -even though -yielding reasonable results compared to more complex architectures, such as deeper 325 nets using more nodes and input data. Here, we use independent test data covering the years 2016 to 2018 of our hindcast simulation. This data has neither been used in training nor during validation so far. To rate the forecast it is compared to the "persistence model", which assumes that the oxygen concentrations at station Boknis Eck appear two weeks later at station Buoy 2a (green line in Figure 23). The first striking impression of the close-ups in Figure 23 is that all years feature a similar seasonal decline in bottom oxygen in autumn and this decline generally closely resembles the oxygen decline in Boknis Eck 330 two weeks in advance. Large interannual differences, however, occur in the onset of the trend reversal. This "return-point" in time is not captured well by the persistency model. These results are consistent with our results in Section 4.3 showing that the decline is driven by the import of low-oxygenated waters from KB. Ventilation, however, takes place in the interior of the Bight and its signal reaches Station Boknis Eck at the entrance afterwards -such that we indeed expect no predictive power of the persistency model under these circumstances. To this end, our ANN outperforms the persistency model in that 335 it predicts an earlier and more realistic recovery of oxygen values during end of summer / beginning of autumn -despite the ANN also exclusively relying on data at the entrance at Station Boknis Eck. The ANN essentially links information regarding season ("derived" from sea surface temperature) and stratification ("derived" from the temperature difference between surface and depth) at the entrance of the Bight with oxygen concentration in the interior of the Bight -without utilizing information on winds. In summary, the ANN features a remarkable (and counter intuitive) performance given that it simply relies on two 340 temperature measurement at the entrance of the Bight.

Discussion
Oxygen concentrations are controlled by the antagonistic interplay of respiration and ventilation processes. Our model-based analysis suggest that the variability in the occurrence of hypoxic conditions in EB is correlated with the a high variability in wind-driven ventilation rather than with a high variability in local respiration. This result is in agreement with AErtebjerg In our model frameworks we distinguish between two types of ventilation: for one, vertical mixing driven by isotropic turbulence and composed of a parameterization of constant background mixing complemented by a surface mixed layer model that mimics the effect of convection, shear-instability and wind-induced turbulence (more specifically we use the KPP scheme of Large et al., 1994). Vertical mixing is difficult to constrain in models because direct observations of turbulence are rare and additional complexity arises from numerical subtleties in models (e.g. Burchard et al, 2008). That said, we use the fidelity of intuition, an ANN fed with information on stratification (i.e. bottom and surface temperature whose difference is a measure of stratification) at the entrance of the Bight and season (i.e. surface temperature which is strongly correlated to season) only, performs surprisingly well without access to wind forecasts -even though the major mechanism behind the oxygen variability 375 is wind-driven. This highlights the importance of the preconditioning that has to precede a ventilating overturning event: In EB, deoxygenation continues almost monotonically until destabilizing buoyancy fluxes have eroded the stability of the water column to a point where the next shift to easterly wind can replace the denser bottom waters with lighter surface waters.
Because synoptic weather systems and associated wind directions have a lifetime of the order of a week in EB, forecasts based on state of preconditioning are, on average, accurate within a week.

380
In summary, we made an effort to explore the uncertainties that are associated with poorly constrained processes such as mixing and oxygen consumption: we tested various degrees of mixing (parameterizations) in combination with and without local sources of oxygen consumption. But caveats remain. Among those is the influence of the waste water treating facility 1mmol O 2 m −3 . This is negligible -to the extent that the assumption of instantaneous homogeneous distribution over the entire Bight holds.

Conclusions
We set out to dissect the mechanisms driving hypoxic events and associated fish-kills in EB. In order to fill data gaps, both spatially and temporally, we developed the high-resolution coupled ocean-circulation biogeochemical model MOMBE and  Smetacek (1980,1985) In-situ measurements in the western Kiel Bight with detritus traps in June (assuming negligible fraction of permanent burial) 1.6 Smetacek (1980,1985) In  . opro, orewa and orese refer to monthly (one value per month starting with the January value) oxygen production, water column oxygen respiration and oxygen consumption by the sediment, respectively (cf. Figure 6). Values for orewa and orese are derived from the published estimates listed in Table 1    Intermittent local environmental conditions of societal concern call for forecasting capabilities in order to facilitate cost-efficient mitigation and adaptation strategies of stakeholders.
A numerical prognostic model that features an explicit mechanistic representation of major processes and feedbacks is developed. A hindcast simulation retracing events of societal concern is integrated on a supercomputer in order to fill gaps in observations.
A machine learning algorithm is applied to the output of the numerical mechanistic model. The algorithm links environmental conditions with events of societal concern in an "easy-to-use" and numerically cost efficient way (i.e. no supercomputer needed).  Figure 6. Schematic of dissolved oxygen module EckO2. EckO2 calculates sinks and sources of oxygen throughout the water column for every grid box. These terms are then passed to the 3-dimensional general ocean circulation that handles the effect of advection and diffusion.
The velocity of the air-sea gas exchange is determined by the piston velocity kgt. Above the compensation depth zco, primary production produces oxygen at a rate prescribed by the model parameter opro. Below the compensation depth zco, respiration of organic matter consumes dissolved oxygen at a rate prescribed by orewa. At the bottom, prescribed oxygen fluxes orese mimic the oxygen consumption of the sediment that is fuelled by the transfer across the water-sediment boundary. Table 2