Summer Oxygen Dynamics on a Southern Arabian Gulf Coral Reef

During the summer the Arabian Gulf is the world's warmest sea, also characterized by hypersalinity and extreme annual temperature fluctuations (12–35oC), making it marginal for coral growth. Yet extensive reefs occur in all eight nations bordering the Gulf. Here we present data demonstrating recurrent summer hypoxia events [oxygen concentration (O2) <2 mg l−1] at a reef in the southern Gulf. Currently these episodes are short enough (median 3 h, max 10 h) to preclude mass mortality. Will this always be the case? Predicting future Gulf hypoxia risk for coral reef ecosystems requires diagnosing the underlying causes driving the timing and magnitude of O2 swings. To this end, we compare our data with the output of a simple coupled 1-D water column/biogeochemical model of the reef environment. This allows us to give quantitative estimates of the O2 fluxes produced by photosynthesis both in the water column and within the coral framework, by respiration processes in the benthos, and from the atmosphere. We demonstrate the role of turbulent mixing, and in particular of tides, in shaping the temporal variability of the amplitude of the diel O2 cycle. We find that, in spite of significant turbulent mixing, which maintains the temperature vertically well-mixed, the biological O2 production and consumption is dominant over the atmospheric O2 flux, and is sufficient to generate vertical differences of 1 to 5 mg l−1 between the bottom and 1.5 m above it. While estimating future trends of hypoxia frequency will require further study, the present findings single out the relevant physical and biological processes (and their interplay) which deserve further scrutiny. The Gulf today experiences temperatures expected to occur across much of the tropics by the end of the century, and the observation of recurrent hypoxia events in the Gulf suggests that similar hypoxic phenomena may represent an important, but to date underappreciated, threat to the future of global coral reefs.


INTRODUCTION
The Arabian or Persian Gulf (hereafter Gulf), despite its geopolitical importance, is one of the least studied marginal seas . In particular, there is a lack of data around its biogeochemistry, and the microbial processes governing it. Basin-wide information comes from a scant number of cruises, distributed over more than 60 years, and with about a decade in between each other (Al-Yamani and Naqvi, 2019). These few data preclude a complete picture of both a historical baseline of ecosystem functioning in the Gulf as well as decadal trends to indicate the likely impacts of ongoing climate change. Biogeochemistry as a whole concerns all relevant chemical elements; among these, oxygen (O 2 ) is necessary for aerobic respiration, critical for metazoans. When waters become deficient in O 2 , they become hypoxic (defined here as <2 mg l −1 ; Vaquer-Sunyer et al., 2008) or even anoxic (no measurable O 2 ). Sufficient exposure to these low concentrations can result in mass mortality, especially in the less mobile benthic community (Grantham et al., 2004;Levin et al., 2009). As a result, the shifting distribution of hypoxic water can significantly alter potential habitat space and community composition (Stramma et al., 2012;Wishner et al., 2020).
Globally, there are marine waters with low-O 2 concentrations, oxygen minimum zones (OMZs), which evolve and are maintained due to the complex interactions between large-scale circulation, O 2 solubility, and ecological functions such as the biological carbon pump (Resplandy et al., 2012). The prevalence of OMZs has been increasing (Breitburg et al., 2018), and is expected to continue with climate change (Bopp et al., 2013). One of these OMZs is located adjacent to the Gulf in the Arabian Sea (Lachkar et al., 2019). In recent decades, evidence suggests that the O 2 concentration in the Gulf of Oman, the body of water joining the Gulf and Arabian Sea, has been declining (Queste et al., 2018). Additionally, seasonal hypoxia has been observed in the deepest part of the central Gulf from 50 m of depth down to the bottom (60-100 m) (Al-Ansari et al., 2015;Saleh et al., 2021). Given the Gulf 's shallow bathymetry (never exceeding 100 m depth before the Strait of Hormuz) and that in its central part water masses have residence times substantially shorter than a year, the occurrence of hypoxia in the central Gulf is somewhat surprising and the mechanism responsible for producing it is unclear.
Hypoxia also occurs elsewhere in the Gulf apart from its deepest regions. Coastal hypoxia can occur in harbors such as Sulhaibikhat bay in Kuwait (Al-Mutairi et al., 2015), with the most obvious explanation being that nutrients from city pollution lead to eutrophication and excess O 2 consumption. Additionally, mangroves are well-known to withstand low O 2 concentrations within the sediment, but hypoxia in the water column and mortality can be exacerbated by pollution from various source, mostly from wastewater treatment plants, but also including oil spills in the wake of the 1991 Gulf War (Kathiresan and Bingham, 2001;Zahed et al., 2010;Etemadi et al., 2021). To date, there have been few or no studies examining the presence and impact of hypoxia in other coastal habitats within the Gulf.
Even though most of the Gulf 's seafloor is dominated by sand and mud, there are widespread regions of hard bottom that host about 6% of the global coral reefs. These ecosystems show a remarkable biodiversity, and support very valuable fisheries (Vaughan and Burt, 2016). The Gulf reefs are the most thermally tolerant in the world, resisting summer temperatures that would be unbearable for any other reef, as well as a yearly temperature cycle spanning 15 • C or more (Burt et al., 2014. In the past decade, however, extensive bleaching and reef collapse has been observed, due to increasing sea temperatures (0.6 • C/decade according to Strong et al., 2011, or 0.3 • C/decade according to Noori et al., 2019) and population-related pressures, with a complex interplay of physical and biological factors (Burt et al., 2014Paparella et al., 2019). Therefore, any additional physiological stress would represent a serious compounding effect on already vulnerable reefs.
In recent years, the crucial role of O 2 in reef ecology has been brought to the forefront. Shifts in the community structure of reef ecosystems in response to hypoxic events have been convincingly documented (Altieri et al., 2017;Bates, 2017;Hughes et al., 2020;Johnson et al., 2021). Much less attention has been given to investigating how the specific interplay of biological and physical causes of hypoxia shapes those shifts, so much that the single, broad category of "dead zone" often clumps together areas where oxygen deficit is persistent, or seasonal, or occurring on shorter time scales (Diaz and Rosenberg, 2008).
Even healthy reefs, though they generally occur in shallow, turbulent waters, and thus in places that are expected to have an oxygen content close to the atmospheric saturation level, may systematically show a diel cycle, where photosynthesis produces oxygenated waters during the day, and respiration uses up that O 2 during the night (Hughes et al., 2020). These oscillations can vary between 50% and 200% of oxygen's air saturation concentration, but it generally ranges over smaller values (see Nelson and Altieri, 2019 for a review). Some field studies have shown that shallow, isolated reefs can adapt to relatively warm, acidic, and deoxygenated waters (Camp et al., 2017), although, at the very least, these temperature conditions do not attain the extremes experienced in the Gulf.
In this paper we report on the O 2 dynamics observed during the summer and early autumn 2019 on a shallow reef located in the southern Arabian Gulf. This is an area removed from the main circulation of the Gulf (Reynolds, 1993), and too shallow to have a coherent wind-driven circulation. South of about 25 • N the bathymetry gently oscillates between 0 and 20 m of depth, on an undulated shelf with no evident northward slope. Here tides are the main sources of currents with peak speeds of about 0.2 m/s, and a residual circulation of less than 0.01 m/s , with the sole exception of the channels between coastal islands. The physical water properties are determined almost exclusively by air-sea fluxes, rather than by upwelling and lateral mixing, as elsewhere in the Gulf (Al Senafi et al., 2019;Li et al., 2020). In particular, strong evaporation produces dense, salty water, which, further North, on the continental slope located around 26 • N, slowly pours into the deepest part of the Gulf and exits through the strait of Hormuz (Swift and Bower, 2003). Lacking a coherent large-scale advection, and vigorous along-shore currents akin to those found along the Iranian and Kuwaiti coastlines, the residence time of waters in the southern shallows exceeds 600 days (Alosairi et al., 2011). Our oxygen data show very large-amplitude diel oscillations, with frequent night-time O 2 dips down to hypoxic levels. We analyze these data to determine their significance, and investigate possible relationships with temperature, tide, and wind speed. Owing to the lateral homogeneity of the physical tracers on the shallows of the southern Gulf, we feel that much can be learned by initially focusing on exploring the impact of vertical processes which may elucidate the basic mechanisms behind both FIGURE 1 | Region of study within the Gulf. Red dot indicates the location of Dhabiya reef, with nearby Gulf cities designated by blue dots. Inset shows study region (red rectangle) located within greater southwest Eurasia. This figure is generated using Blue Marble: Next Generation data provided by Reto Stöckli, NASA Earth Observatory.
the observed diel cycle and hypoxia. To this end, we interpret our data using a simple O 2 biogeochemical model representing microbial processes coupled with an underlying physical water column model.

Study Location and in situ Sampling
Data for this study come from Dhabiya Reef, located in the waters of Abu Dhabi in the southern Gulf. The reef, located ∼8 km offshore, but close to the lagoonal barrier islands, is at a depth of ∼4 m (Figure 1) and surrounded by flat sandy bottom. Hourly in situ temperature and O 2 data were collected with respectively, HOBO TidbiT V2 (±0.2 • C accuracy) and HOBO U26 Dissolved Oxygen (±0.2 mg/l accuracy) data loggers. The dissolved oxygen probe was calibrated prior to deployment, with maintenance, cleaning and battery replacement performed every 2 months. These instruments were positioned at three stations close to one another (∼50 m). Each stations consists of a 12 mm diameter vertical metal pole permanently fixed into the reef substrate with sensors placed at two depths: within the reef framework, ca. 2-5 cm above the bottom in interstices between coral colonies (these will be referred to as "bottom sensors" in the following), and 1.5 m above that (these will be referred to as "midwater sensors" in the following), for a total of 6 O 2 time series. After removing 72 h from the beginning and end of each time series to account for sensor acclimation, data from stations 1 and 3 span from June 16 until 4 November 2019, with station 2 beginning 25 August and ending the same date.
The time series data between the three stations are largely congruent with one another (Supplementary Figure 1), especially with respect to temperature, but during differing periods O 2 at either stations 2 or 3 deviated from the other two, possibly because of bio-fouling (see Supplementary Section 1).
Conversely, oxygen data from station 1 were consistently in agreement with at least one of the other sites at all times. As a result, the data from station 1 (Figure 2) have been selected for further analysis and comparison with model output.
The diel cycle was reconstructed by averaging all the oxygen concentration values corresponding to the same time of the day. Percentile values of the underlying distribution were also computed (Figure 3). The procedure was then repeated including only the days (defined as beginning and ending at noon of the local time) during which the oxygen concentration dropped below the threshold of 2 mg/l ("hypoxic event"). We used a Wilcoxon Rank Sum test (as defined by the ranksum function in Matlab, MATLAB, 2020) in order to determine during which hours the diel cycle of the subset of days with hypoxic events shows significant differences with the average diel cycle of the whole data set. The choice of the Rank Sum test over the more popular ANOVA stems from the fact that the latter requires the data to be Gaussian and homoscedastic in order to be reliable, hypotheses for which, in our case, we have little evidence, and no theoretical support (nonetheless, while not shown, ANOVA applied to these data shows results that do not contradict those found with the Rank Sum test).

Environmental Data
Surface atmospheric and tidal data were retrieved both to provide environmental context and for forcing in the physical water column model (section 2.3). Surface atmospheric data were retrieved from the European Centre for Medium-Range Weather Forecasts (ECMWF) high-resolution reanalysis data (ERA5, Hersbach et al., 2018) for all of 2019. The variables include 10 m horizontal winds, 2 m air temperature and pressure, dew point, and net shortwave radiation flux, all provided at hourly time intervals. Data were interpolated from the native grid to the reef site position using the method of Kara et al. (2007) to minimize effects of nearby terrestrial gridpoints. A comparison between the ERA5 data and in-situ observation in the Northern Gulf is available in Al Senafi et al. (2019) Tidal horizontal velocities u, v, and sea surface height (ssh) were retrieved from the Oregon State University (OSU) TPXO Arabian sea and gulf regional tidal model (Egbert and Erofeeva, 2002), available at (https://www.tpxo.net/regional, accessed 10-Nov-2020) using the TMD toolbox provided by Earth and Space Research (https://www.esr.org/research/polar-tide-models/tmd-software/, accessed 10-Nov-2020). The model provides a grid at 1/60 • resolution, incorporating 8 tidal components. Tidal model bathymetry reports the depth at the reef to be 2.8 m, resulting in unrealistically large velocities. To ameliorate this, u, v, and ssh were interpolated to a nearby position 1 km away to where the depth is 4 m, thus matching the real depth of the Dhabiya site. Comparison between the two timeseries shows a reduction in magnitude, but no significant differences between phases and timing of low and high tides. Additionally, the results of Pous et al. (2012), which used more Gulf observations, indicate that tidal magnitudes never exceed 0.2 cm s −1 . As a result, tidal velocities and height were further multiplied by the scaling factor 0.4, which brings the amplitudes computed by TPXO in line  with the findings of Pous et al. (2012). The impact of this scaling factor is discussed in the Supplementary Section 2.

Water Column Model
The in situ temperature data were matched with output from the generalized ocean turbulence model (GOTM) water column model (GOTM v.5.2, Umlauf et al. (2005), downloaded 10-Dec-2020 from https://github.com/gotm-model/code). GOTM has been previously used for reef studies in the region, and, in particular, to quantify the evaporative cooling induced by breezes and synoptic winds Paparella et al., 2019). Despite being a 1D model, its use is justified in our study because of the similarity in O 2 data between the three stations (Supplementary Section 1). GOTM depends on atmospheric surface forcing, tides, and water optical properties to simulate turbulent transport of heat. Here we use the ERA5 surface reanalysis data interpolated at the position of Dhabiya reef, and barotropic tidal current and sea surface elevation data obtained from the TPXO model with reduced tidal velocities. The optical water properties are assumed to be described by the Jerlov-3 absorption function (Stips, 2010). Lacking in-situ turbidity data, we chose this absorption function because it was found to give the best results in a previous study using a similar model . The effect of these forcings and the validation of the model are in the Supplementary Section 3. The water column has a time-averaged depth of 4m, and it is discretized into 60 equally-spaced intervals. The intregration time step is 30s.
The model is representative of the water column from the sea surface to the top of the hard-coral structure. The water enveloped within the reef framework is considered as part of the benthos. The roughness of the bottom is taken into account by GOTM's turbulence parameterizations, which we selected to be the same as in Paparella et al. (2019).

Oxygen Model
The GOTM water column model is used to estimate eddy diffusivity values as a function of depth and time. These, in turn, are used to solve the equations of a simple O 2 model, embodied by the following equations: where O 2 is the dissolved oxygen concentration and κ is the eddy diffusion coefficient of transported tracers as diagnosed by GOTM. The term aIP(I) is a crude simplification that models the phytoplankton O 2 production rate. Here, a is a constant, whose value will be fitted to match the observed data; I is the photosynthetically available radiation which depends on depth z and time t as diagnosed by GOTM; P(I) is a photo-inhibition term, for which we use the formulation of Eilers and Peeters (1988): (2) with I 1 = 500 W/m 2 and I 2 = 150 W/m 2 . This model does not have explicit phytoplankton population dynamics, a simplifying choice justified by our interest in investigating the role of physical drivers in modulating O 2 , as well as a lack of direct information on: (1) nutrient availability, (2) the dominant phytoplankton and zooplankton groups, and (3) the relative role of the microbial loop in the coastal waters of the southern Gulf. Therefore, the constant a represents the oxygen production rate per unit of solar radiation of an idealized constant phytoplankton biomass, which experiences photo-limitation at high intensities of the solar radiation. The simplifying choice of keeping the phytoplankton biomass constant explicitly rules out for the model the possibility of describing algal blooms.
Equation (1) is subject to flux boundary conditions at both ends of the water column. At the surface, the oxygen flux is modeled as in Wanninkhof et al. (2009): Here we follow the convention that upward fluxes are positive; k is the gas transfer velocity, and is modeled as a function of temperature and wind velocity: where U is the wind speed at 10m height above the sea surface and T is the sea surface temperature. The values of the positive constant k 0 and of the coefficients of a 3 rd order polynomial approximation for the dependence of the Schmidt number Sc on temperature are taken from Wanninkhof (1992). The O 2 air saturation concentration S O is a function of water temperature and salinity. We use the fitting formula of Garcia and Gordon (1992), using the in-situ measured water temperature and a climatological average salinity of 42 psu (Campos et al., 2020).
A fully realistic boundary condition for the oxygen should take into account of a number of complicated physical and biological processes that transport oxygen through the reef framework and ultimately into contact with both the individual coral polyps (Kühl et al., 1995;Larkum et al., 2003) and any other consumer/producer that may live in the benthos or in the nepheloid layer trapped in the reef framework. In our simplified model, the "benthos" begins at an imaginary surface at the top of the reef framework. Therefore our bottom boundary conditions should describe the oxygen fluxes through this imaginary surface. This allows us to summarize all the complex dynamics occurring in the reef framework with simple prescriptions for the oxygen flux at the bottom boundary of the model. We explore the effect of two distinct boundary conditions. The simplest is: where a b is a constant representing the O 2 production rate per unit area and unit solar radiation of the coral reef, I b is the photosynthetically available radiation at the bottom of the water column, and the constant F b is the oxygen consumption (respiration) rate per unit area of the benthos. A slightly more advanced version of this boundary condition introduces a temperature dependence on the respiration rate, using the Arrhenius parameterization (Kooijman, 2000), which gives: where T A is the Arrhenius temperature in Kelvin, T ref is a reference temperature in centigrades, and T is the bottom water temperature in centigrade. Finally, a more complicated boundary condition is used in order to explore the role of aerobic microbes in the benthos: which is coupled to an equation for the time evolution of the benthic microbial concentration M: Here the first term in the right-hand side of (7) represent the coral O 2 production, just as in (5). The second term assumes that the microbial O 2 consumption is modulated by O 2 availability, through the function r b . Likewise, bacterial growth is also modulated by O 2 availability, through the function r. The last term in (8) represents bacterial mortality, assumed to occur at a constant rate µ. For simplicity, we approximate the unknown functions r b and r as: where b and γ are positive constants. This microbial benthos model does not aim at being realistic, and therefore excludes many relevant factors, such as the availability of dissolved organic carbon. Its purpose is to serve as the simplest example of how a benthic microbial community may regulate the oxygen concentration of the overlying water column.
Equation (1) is discretized using centered second-order differences on the same grid used by GOTM, with the eddy diffusivity κ being evaluated on a grid staggered with respect to that of the variables O 2 and I. Time stepping is performed using a Crank-Nicholson method, without splitting the diffusion and the production term in (1). The resulting non-linear set of equations is solved by means of the Newton-Krylov method implemented in Scipy v.1.3.3 (Virtanen et al., 2020). The code is run with a time step of 5 min. We performed many test runs and found that a month is sufficient to obtain results that are independent from the choice of initial conditions. For the results discussed in the next section, the runs start on May 1 st with the water column having a depth-independent oxygen concentration of 5.9 mg l −1 , which is approximately equal to the air saturation value.

RESULTS
Two hourly timeseries of dissolved O 2 concentrations recorded at station 1 are shown in Figure 2, again produced by sensors placed within the reef framework (bottom data) and 1.5 m above the bottom (midwater data). The thick green line shows the O 2 air saturation concentration S O , computed using the hourly midwater temperature data recorded by probes co-located with the O 2 sensors.
The diel O 2 cycle is very evident at both depths, with a change in concentration from a daily minimum to the subsequent daily maximum that rarely spans less than 3 mg l −1 and can occasionally exceed 6 mg l −1 . The cycle straddles the O 2 air saturation concentration with daily maxima (minima) always above (below) the air saturation value. This is clearly shown by the enveloping thin black lines, which represent the maximum and minimum observed O 2 concentration in the previous 24 h. As expected, the maxima occur during the day, when photosynthesis dominates, and minima occur at night, when respiration dominates. On several occasions the O 2 concentration drops to hypoxic levels, that we conventionally define as 2 mg l −1 (dashed horizontal line). This occurs more often at the bottom (lower dots) than in midwater (upper dots). Within hypoxia days at the bottom, which total 26 events, duration ranges from as small as one data point (that is, 1 h), which is the most frequent case, up to 10 h, with a median of 3 h. At the bottom, the oxygen concentration drops below the conventional anoxia threshold of 0.5 mg/l twice: for 2 h on July 29th (with O 2 concentration reaching the lowest recorded level of 0.2 mg/l), and for 1 h on August 14th. This is observed at both Stations 1 and 3 (Station 2 was deployed at a later date; Station 3 also records 1 h of O 2 concentration below anoxic threshold on July 28th). Events very close to hypoxic conditions occur in September and in October, and are recorded by all three stations. The anoxia threshold is never breached at midwater. Hypoxia observations can occur as early as 9 p.m., or as late as 9 a.m. More generally, the night-time O 2 concentration is lower close to the bottom than at midwater, suggesting that a substantial fraction of the total O 2 consumption takes place in the benthos. Before October, the 1-week moving average of O 2 concentration is always lower at the bottom than at midwater, and never surpasses the air saturation value, revealing the existence of a persistent summer O 2 deficit at the base of the water column. Figure 3 shows a distribution of the diel O 2 cycle for all the observed days (solid blue lines display the median, dashed blue lines show 25/75th percentile ranges). The right panel, showing the differences between midwater and bottom O 2 , confirms that the bottom is more oxygenated (negative values) than the midwater level during the day, and less so during the night (positive values). The teal green lines in Figure 3 shows the same quantities, but limited to the days where hypoxic conditions are observed (those marked with a dot in Figure 2). The most evident difference is that in the days when hypoxic conditions are reached, the diel cycle has a larger amplitude than on average, and the night-time vertical gradient is also much larger than on average, further pointing to an important role of the benthos. In addition, the departures from the median are significant (at 99% confidence level or better, black dots in Figure 3) only during night and early morning; i.e., daytime oxygen concentration during days with an hypoxia event is not significantly different from that of the average day.
These results suggest that the processes producing hypoxia have both strong and fast onsets as well as cessations, with implications for reef recovery. The same statistical tests done in Figure 3 were undertaken for temperature (details in Supplementary Section 4). Temperatures during hypoxia events are indeed different and warmer than non-hypoxia events if one uses the entire timeseries. However, when using only data from the summer, we find no significant temperature differences between hypoxia and non-hypoxia days. This is explained by observing that hypoxia events are more common during the summer than in the autumn, but there is a significant temperature difference between the summer and the autumn. So, while temperature may indeed be a factor in favoring hypoxia events, day-to-day temperature fluctuations do not seem to play a role, and only the temperature dynamics on seasonal and longer time scales may actually be an important factor.
In searching for alternative mechanisms favoring hypoxia, environmental factors must also be considered. Figure 4 shows the moving standard deviation of the O 2 concentration measured at the bottom. The hypoxic events (black dots) generally coincide with high values of the standard deviation, that is, to days when the diel cycle has a large amplitude. This is consistent with the O 2 significance test results, showing that days with hypoxic events have similarly high daytime O 2 values while having much lower (hypoxic) nighttime values. Furthermore, the moving standard deviation has a clear fortnightly cycle, with opposite phase of the spring-neap cycle of the tides (green lines, showing the moving average of the tidal current speed, as given by the TPXO model): low-amplitude tides are associated to high values of the O 2 moving standard deviation. The dark yellow lines in Figure 4 show the moving average of the wind speed, interpolated at the reef 's position from the ERA5 reanalysis data. Large wind speeds strongly reduce the variability of the diel O 2 cycle, and when strong winds persist for several days this may even offset the effect of the tide (e.g., in June and early July). These observations suggest that the amplitude of the O 2 's diel cycle is regulated by turbulence in the water column. When turbulence is strong and mixing is high (strong tides and / or winds) the diel cycle has a small amplitude, but if both winds and tides slack, then the cycle may reach a standard deviation in excess of 2 mg l −1 . This is consistent with a scenario where O 2 consumption occurs at the bottom and a substantial production of oxygen occurs higher up in the water column: high turbulence would maintain the oxygen well-mixed in the water column, dampening the amplitude of the diel cycles observed at the bottom, while low turbulence would allow the bottom sink to build-up at night a strong vertical oxygen gradient, resulting in an oxygen-depleted bottom layer.
In order to test the turbulence-mixing scenario, we have used the model described in section 2.4. We performed several runs, with different choices for the sources of O 2 . In all cases the O 2 consumption occurs exclusively at the bottom. This oversimplifying hypothesis allows us to use only three flux parameters (a, a b , F b , cf. section 2.4) and fit their values with the following simple procedure. In all cases, we tune the parameters so as to maintain the daily moving average of the O 2 concentration slightly below the air saturation value (as in the observations), and the minimum-to-maximum range of the diel cycle close to 4 mg l −1 . This approximately matches the observations (Figure 5D, left). However, when we only impose O 2 production in the water column (Figure 5A, left, a = 3.3 · 10 −6 mg m 2 l −1 s −1 W −1 , a b = 0, F b = 3.5 · 10 −4 mg m l −1 s −1 ) the vertical O 2 gradient is always positive (Figure 5A, right). If O 2 production occurs, instead, only at the bottom (Figure 5B, left, a = 0, a b = 4.5 · 10 −6 mg m 3 l −1 s −1 W −1 , F b = 2.0 · 10 −4 mg m l −1 s −1 ) the vertical gradients are predominantly negative (Figure 5B, right). Only when we allow O 2 production both in the water column (by phytoplankton) and at the bottom (by the corals and benthic algae) (Figure 5C, left, a = 1.9 · 10 −6 mg m 2 l −1 s −1 W −1 , a b = 3.4 · 10 −6 mg m 3 l −1 s −1 W −1 , F b = 3.5 · 10 −4 mg m l −1 s −1 ) can we approximately reproduce vertical gradients that alternate between positive values during the day and negative values during the night with magnitudes that approximately match the observations (Figures 5C,D, right). Figure 6 shows the daily averaged O 2 fluxes resulting from the model of Figure 5C, using the convention that positive values correspond to O 2 added to the water column, and negative ones represent O 2 removed from it. The largest contribution comes from photosynthesis occurring in the water column [blue line, showing the vertical integral of the last term in Equation (1)], which, due to the assumed photoinhibition of phytoplankton, is insensitive to small changes in the solar irradiance, and remains roughly constant at 2 · 10 −4 (mg l −1 )(m/s). The O 2 flux due to photosynthesis occurring at the bottom is shown by the orange line. Its magnitude is smaller than the bulk production by a factor of about 2/3, and shows a slight downward trend throughout the summer months due to the declining solar irradiance after the summer solstice. The atmospheric O 2 flux (green line) has a magnitude much smaller than that of the previous two fluxes. Even though it is generally positive (O 2 flowing from the atmosphere into the water), it can occasionally attain slightly negative values (O 2 outgassing into the atmosphere), especially in June, when the solar irradiance (and thus the photosynthesis) is largest. The red line shows the constant bottom respiration flux that removes from the water column 3.5 · 10 −4 (mg l −1 )(m/s) of O 2 . The algebraic sum of these four terms cancels out to less than 2 · 10 −5 (mg l −1 )(m/s).
The inset in Figure 6 shows the instantaneous O 2 fluxes (sampled hourly) at the end of July, when several hypoxia events are observed in the data. As for the main figure, the blue, orange and green line show the contribution of, respectively, the photosynthesis occurring in the water column, photosynthesis at the bottom, and the atmospheric flux. The bottom photosynthetic flux peaks during the day at 4 · 10 −4 mg m l −1 s −1 , while the photosynthesis in the water column slightly exceeds that value for several hours during the day. Both fluxes clearly drop to zero during the night. In order to check on the realism of our model, we give an alternate estimate of the atmospheric flux (black line) by using the temperature and O 2 concentration measured by the midwater sensor in Equations (3) and (4). This procedure yields the correct O 2 atmospheric fluxes under the assumption that the temperature and O 2 concentration at the surface of the water column have a negligible difference from the values measured at midwater, and that the ERA5 wind data closely match the local winds. The agreement between the model's (green line) and the estimated (black line) atmospheric fluxes is fairly good. In both cases there is a clear diel cycle, with weak fluxes out of the water column in the late afternoons, and much stronger fluxes into the water column early in the morning, that may briefly peak at about 2 · 10 −4 (mg l −1 )(m/s). The dynamics of the atmospheric O 2 flux is strongly affected by the wind (dark yellow line): the proportionality coefficient k that appears in the expression (3) for the flux is a quadratic function of the wind speed U (cf. Equation 4). Thus, non-negligible fluxes to or from the atmosphere, in addition to an excess / deficit of O 2 in the water column, also require a sufficient wind speed. For example, from the night of July 26 th to the early morning of July 28 th , winds permanently slower than 4 m/s maintain atmospheric O 2 fluxes very close to zero, even though the O 2 dissolved in the water column at night is well below the air saturation concentration. Likewise, the early morning land breezes help replenishing the water column depleted of oxygen by night-time respiration, contributing to the morning positive maxima of atmospheric oxygen flux, while the evening sea breezes help off-gassing to the atmosphere any oxygen produced by photosynthesis in excess of the air saturation concentration, contributing to the evening negative minima of atmospheric oxygen flux.
A further comparison of the model with the observed hypoxia events at the end of July is offered in Figure 7. The agreement between the modeled and the observed O 2 is reasonably good in the days when the oxygen concentration remains above the hypoxic threshold, and less satisfactory in the days when hypoxia is observed. Furthermore, in the model, the maximum difference between the midwater and bottom O 2 concentration is generally reached in the evening or early nighttime, while the observations show persisting vertical gradients throughout most nights. This discrepancy stems from an excessive vertical mixing in the water column diagnosed by the GOTM model, in particular close to the surface, as shown by the very high values attained by the   eddy diffusion coefficient, even with weak winds and high or low tide. This suggests an inadequacy of GOTM's parameterizations (including the optical water properties) for the conditions specific to our study location, or inaccuracies in the forcing data, a topic which will be the subject of future investigations. Similarly, imperfect parameterizations might explain why the modeled O 2 concentration does not appear to be modulated by tides and wind as strongly as the observed data. A brief discussion of the model's daily heating/cooling cycle and its potential to overestimate mixing is in Supplementary Section 5.
The need for a more sophisticated treatment of biogeochemistry is highlighted in Figure 8. Figure 8A shows the daily averaged solar irradiance at the top of the water column (yellow line): from mid-June to end-October its value declines by more than 25%. Consequently, the photosynthetic O 2 fluxes decline by a similar amount. Having assumed a constant respiration flux, the model shows persistent hypoxia, and even anoxia as the fall season progresses (Figure 8B). This is clearly in contrast with the observations, which actually show a weak increasing trend of O 2 concentration ( Figure 8E). That is, in an ecosystem where O 2 fluxes appear to be predominantly produced by photosynthesis, the observed O 2 levels raise as solar irradiance decreases. While a full solution to this apparent paradox will require further in-situ measurements, we can hypothesize two distinct mechanisms (not mutually exclusive) that may explain this problem. Figure 8C shows a run with the model of Figure 8B, but where the bottom respiratory flux is corrected with an Arrhenius dependence on temperature (Equation 6, with a, a b , F b as in Figure 5C, and T A = 5, 000 • K, T ref = 34 • C): as the water temperature decreases, so does the respiratory flux, which produces an increasing trend in the O 2 concentration during October. Figure 8D shows a run where the bottom fluxes are determined by the simple benthic microbial model embodied by Equations (7)-(9), with a, a b as in Figure 5C, and b = 5 · 10 −5 m 3 (mgC) −1 s −1 , γ = 10 −5 l mg −1 s − 1, µ = 5 · 10 −5 s −1 . The inset shows the time evolution of the benthic microbial concentration. In this model, as the season progresses and light becomes less available, rather than resulting in the decrease of water column O 2 , it is the microbial population in the benthos that adjusts its abundance, pegging the O 2 to oscillate around a fixed value approximated by the ratio µ/γ = 5 mg l −1 (see details in the Supplementary Section 6). Both the models in Figures 8C,D are too simplistic, and have serious shortcomings. Introducing a temperature dependence in the respiration rates is physiologically correct, but it might not be a complete solution of the issue, because in August and early September light is already declining, while water temperatures are as high or higher than in late June and July (Figure 8A). A microbial control of the water column O 2 is an enticing hypothesis, but the simple-minded model presented here exerts too harsh a control, pegging the average O 2 to a constant, and removing nearly all variability from the diel O 2 cycle. Further work, and, most importantly, further observations, shall be necessary to give a full account of the observed O 2 dynamics on a multi-seasonal time scale.

DISCUSSION AND CONCLUSION
In this paper we report uninterrupted, high temporal resolution (1h) measurements of O 2 concentration on a coastal coral reef in the southern Gulf spanning mid-June to early November 2019. We find 26 separate occurrences of hypoxia (O 2 <2 mg l −1 ) always occurring at night and occasionally persisting through the early morning. While one hypoxia event spanned 10 h, in most cases hypoxic conditions last for about 1 h, indicating that currently this Gulf reef more often experiences near-hypoxia conditions. However, in two cases during the summer anoxic conditions have briefly been attained, and in a few other cases, even during the autumn, near-anoxic conditions have been recorded. Comparing hypoxia events to non-hypoxia conditions, we find that O 2 levels differ substantially only during the nighttime, indicating that O 2 concentration is quickly restored to normal conditions during the day, and that the processes creating hypoxia are contained within a nightly timespan. The concern that temperature is responsible for hypoxia is reasonable, but while most hypoxia indeed occurs during the summer, the events themselves do not occur during especially warm nights, indicating temperature as a contributing, but not primary factor in favoring hypoxic conditions. Even when O 2 does not drop to hypoxic levels, the diel cycle of O 2 concentration at the bottom of the water column typically spans about 4 mg l −1 from minimum to maximum (this range drops below 1 mg l −1 only once in June, but occasionally it exceeds 6 mg l −1 ). Even though the reef is only 4 m deep, is exposed to tides, and is not sheltered from the winds, we still often observe significant vertical differences in nighttime O 2 between the bottom (less oxic) and at 1.5 m above the bottom (more oxic). Differences of 1-2 mg l −1 are typical, but cases where this value exceeds 4 mg l −1 do occur. During the day smaller, reversed vertical differences are observed, with the bottom having higher O 2 than the overlying water, but this difference usually does not exceed 1 mg l −1 (Figures 2, 3).
The variability in the amplitude of the diel cycle appears to be regulated by the turbulent mixing in the water column, which, in turn, is strongly affected by the fortnightly modulation of the tidal cycle and, on shorter time scales, by the wind. Strong winds and large-amplitude tides correspond to a small-amplitude O 2 diel cycle, and vice versa (Figure 4).
A water column-model consisting in state-of-theart turbulence parameterizations, forced with reanalysis meteorological data, and equipped with a very simplified description of biogeochemical processes, suggests that the dominant contribution to the O 2 fluxes occurs in the water, through photosynthesis and respiration. Time-averaged atmospheric fluxes are about an order of magnitude smaller, although they may attain non-negligible values when the water column is both severely depleted of O 2 , and strong wind (5 m/s or more) boosts the gas transfer velocity at the air-water interface. Helped by breezes, O 2 generally flows from the atmosphere into the water column during the early morning, and from the water column into the air during the evening. The model also suggests that the O 2 produced at the bottom by the corals is about 2/3 as much that produced in the water column by the phytoplankton which, when integrated over the whole water column, amounts to about 2 · 10 −4 mg m s −1 l −1 = 2 · 10 −7 mg m −2 s −1 . This fraction is much higher than what has been reported for coral communities in the Red Sea, where the contribution of the water column is about 60 times smaller than that of the benthos (Cardini et al., 2016), and is a subject that deserves further investigation. Among the factors that may contribute to this unusual partition of the gross primary production between the benthos and the water column there is a general trend to eutrophication observed in the Gulf (Al-Yamani and Naqvi, 2019), the exponential dependence of plankton's gross primary production on temperature (López-Sandoval et al., 2019), and the extensive damage endured during the 2017 bleaching event, in the wake of which coral declined from 49 to 15% coverage, while algae coverage nearly doubled from 25 to 45% . The presence of strong vertical gradients in O 2 despite being in a turbulent water column is evidence that most O 2 -consuming processes occur in the benthos (and/or in a nepheloid layer trapped within the coral reef framework). By design, the model assumes that respiration in the water column is negligible. We are comforted in holding this assumption by the reasonable match between the fluxes at the air-sea interface computed by the model and those inferred directly from the observations (Figure 6, inset): given that a strong outflux of O 2 is required at the bottom in order to match the observed vertical gradients, significant amounts of additional respiration in the water column would require a matching increase in the photosynthetic production; this, in turn, would produce more excess O 2 at the end of each day, yielding O 2 fluxes from the water to the atmosphere larger than the observed ones.
The biogeochemical portion of the model is purposely oversimplified. Our goal is to obtain an order-of-magnitude estimate of the O 2 fluxes without introducing a large number of free parameters, whose pertinence to the specific conditions of Gulf coastal reefs is unknown. Therefore, we adopt a simple-minded model that assumes constant respiration and photosynthetic flux rates, ruling out any model sophisticated enough to represent the internal, evolving dynamics of plankton populations and the biogeochemistry of the nutrients they require. Therefore, the present model lacks the ability to represent and explain any O 2 changes due to the internal variability of biological processes.
On seasonal time scales this approach has obvious limitations. The observations demonstrate that the dissolved O 2 levels increase slightly from mid-September onward, in spite of a substantial decrease in solar radiation. This cannot be explained solely through physical processes (e.g., increased O 2 solubility with decreasing temperature) and thus there must be a biogeochemical mechanism at work. We have explored two distinct, but non-mutually exclusive hypotheses. One assumes that increasing O 2 levels with decreased sunlight are the results of temperature-dependent respiration, the other assumes that O 2 levels in the water column are regulated on seasonal time scales by the microbial population in the benthos/nepheloid layer. Both hypothesis lie on sound physiological and ecological principles, but in the form presented here neither gives fully satisfactory results. From mid-August to mid-September water temperatures remain nearly as high as in the previous month, in spite of a substantial reduction in solar insolation. If temperature-dependent respiration were the only factor, a conspicuous drop in the daily averaged O 2 levels would occur in this period (see Figure 8C), contrary to the observations. This result de-emphasizing the sole role of temperature is consistent with the statistical result that hypoxia does not occur on significantly warmer days. We demonstrate that a benthic microbial population can indeed control the O 2 level in the water column, but our model is nothing more than a proof-of-concept, and neglects important processes, such as a description of the production, sedimentation and consumption of dissolved organic carbon, which likely plays a fundamental role in determining the microbial dynamics of the benthos.
The frequency and magnitude of the hypoxic/anoxic events recorded in this study provide an insight into processes that may regulate reef communities within the southern Gulf, particularly when considered in conjunction with temperature effects. Most coral reef fishes are resilient to hypoxia, with critical oxygen concentrations varying between 13 and 34% air saturation (Nilsson and Östlund-Nilsson, 2004). This suggests that fishes can tolerate all but the most extreme hypoxia events on reefs in the southern Gulf, at which point they may be forced from their shelters within the reef framework, becoming vulnerable to predators. Yet, Nilsson et al. (2010) showed that at higher temperatures this tolerance to hypoxia is compromised, as a 3 • C increase in temperature from 29 to 32 • C resulted in an increase of 71 and 23% increase in critical O 2 concentrations in two reef fish species. Similarly, hypoxic conditions are known to reduce the upper thermal limits of fishes (Ern et al., 2017). Temperatures on southern Gulf reefs regularly exceed 34 • C and can reach 36 • C during summer months making them the hottest reefs in the world (Sheppard et al., 1992). Despite these extreme temperatures a recent study concluded that temperature alone did not explain the low diversity of fishes on southern Gulf reefs compared to reefs in the Gulf of Oman (Brandl et al., 2020). As described here fishes on southern Gulf reefs are exposed to both extremely high temperatures and repeated bouts of hypoxia. Hence, the combined action of these two stressors may exceed the tolerance limits of many fishes and may act to regulate population and species composition on these reefs.
Our findings unveil a further menace to the reef ecosystems of the southern Gulf. In this region, climate change might increase the likelihood of attaining bleaching conditions by favoring lowwind conditions that greatly reduce surface evaporative cooling during the summer . The present work shows that low-wind conditions attained during neap tides set the stage for extreme excursions of the oxygen concentration at the bottom of the water column, that can easily reach hypoxic conditions, and even approach anoxia, for a few hours during the night. Night-time hypoxia may likely be an important contributing factor to the widespread bleaching events observed in the last decade, as very recent research is unveiling (Alderdice et al., 2021). We also find that the oxygen fluxes are mainly determined by endogenous biological processes, with atmospheric fluxes playing only an ancillary role. Therefore, any factor, such as human-induced eutrophication, that increases the activity of these biological processes may greatly enhance the amplitude of oxygen concentration excursions, and hence produce longer and deeper hypoxia and anoxia events. Even if the frequency, length and intensity of hypoxia events were to increase gradually, the response of reef ecosystems is known to be characterized by sharp thresholds (Hughes et al., 2020). Coastal reefs in the southern Gulf have already lost virtually all of the Acropora corals, a key ecosystem builder . Thus, the possibility cannot be ruled out that the combination of bleaching temperatures with slightly more frequent, longer or more intense hypoxic conditions, may push these fragile ecosystem closer to total collapse. Future work aimed at investigating the outstanding questions still open on the dynamics of oxygen on the Gulf 's coastal reefs will require much additional data. It will be necessary to monitor the O 2 concentration for at least an entire year, in order to have a complete picture of the seasonal dynamics. In-situ meteorological data (wind speed, relative humidity, etc.) and vertical profiles of (at least) temperature and turbidity, as well as ADCP data will be necessary in order to assess the reliability of the turbulence representation in the GOTM water column model. Formulating a more realistic, dependable biogeochemical model will require observations of at least chlorophyll and basic nutrients (nitrate, etc.), as well as repeated surveys aimed at identifying the species composition and abundances of reef phytoplankton assemblages in the southern Gulf. Obviously, direct measurements would be required also to settle any questions on the role played by the microbial population in the benthos and/or bottom nepheloid layer in regulating the water column O 2 .
Finally, we believe that much can be learned by extending these observations beyond the reef ecosystems, and we plan to investigate the oxygen dynamics over seagrass meadows, sand/mud flats, and mangrove forests within the southern Gulf.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors upon request, without undue reservation.