Hypoxic volume is more responsive than hypoxic area to nutrient load reductions in the northern Gulf of Mexico—and it matters to fish and fisheries

While impacts of low oxygen on marine organisms have been reviewed from physiological and ecological perspectives, relating broad population- and ecosystem-level effects to the areal extent of hypoxia (dissolved oxygen concentration below 64 μM, or 2 mg l−1) has proven difficult. We suggest that hypoxic volume is a more appropriate metric compared to hypoxic area because volume better integrates the effects of hypoxia on ecological processes relevant to many marine taxa. In this paper, we compare the volume-based load responses from a simple biophysical model with results from an independent three-dimensional hydrodynamic-biogeochemical model, and discuss the implications with respect to potentially more ecologically-relevant hypoxia management goals. We also show that hypoxic volume appears more sensitive than hypoxic area to nutrient load reductions. Model simulations indicate that even under a modest 25% nitrogen load reduction, the thickness of the hypoxic layer in the northern Gulf of Mexico decreases markedly, and hypoxia remains localized to a relatively thin layer near the bottom that most fish and other mobile organisms can more effectively avoid. This finding should be considered when reviewing and potentially setting hypoxia management goals.


Introduction
Nutrient inputs to marine ecosystems have increased globally, leading to enhanced primary production, eutrophication, harmful algal blooms, and increases in the number and sizes of low oxygen zones. Coastal hypoxic areas (with dissolved oxygen (DO) concentrations below 64 μM or 2 mg l −1 ), in particular, have increased over the past few decades, with large hypoxic regions in the Baltic Sea, Black Sea, northern Gulf of Mexico (nGOM) and East China Sea (Diaz and Rosenberg 2008, Breitburg et al 2018. Substantial hypoxic and anoxic areas also occur in US estuaries (Bricker et al 2008) and the Laurentian Great Lakes (Scavia et al 2014). Many of these systems also harbor large commercial and recreational fisheries because moderate levels of external nutrient inputs typically enhance fishery production (Caddy 1993, Nixon and Buckley 2002). While the potential impacts of these low oxygen zones have been reviewed from physiological and ecological perspectives (Vaquer-Sunyer andDuarte 2008, Breitburg et al 2018), relating fish and fisheries effects to the areal extent of hypoxia has been difficult (Breitburg et al 2009a(Breitburg et al , 2009b. One reason might be that hypoxic volume, which is rarely used to develop relationships for systems experiencing hypoxia, is a more appropriate metric compared to hypoxic area. From a biogeochemical perspective, hypoxic volume better accounts for the integrated ecosystem metabolism of the hypoxic zone. For example, while sediment respiration appears to have a controlling influence on the areal extent of hypoxia in the nGOM (Obenour et al 2015, Fennel et al 2016, water column respiration, which is often a dominant component of the overall oxygen uptake (e.g. Murrell and Lehrter 2011), can exert a stronger volumetric influence. Li et al (2016) estimated that water column respiration to be 2 to 4 times as large as the sediment oxygen demand in the Chesapeake Bay.
From an ecological perspective, hypoxic volume may be a better integrator of overall habitat quality for highly mobile organisms such as fishes and some crustaceans, compared to hypoxic area. Many species use both bottom waters and pelagic habitats and so are influenced by both the areal and vertical extent of hypoxia, though to differing degrees depending on species (Hazen et al 2009, Roberts et al 2009, Zhang et al 2009, Craig et al 2010, Craig 2012, Reeves et al 2018, LaBone et al 2018. One of the key uncertainties in quantifying the effects of hypoxia on mobile species is how variation in the areal and vertical extent of hypoxia influences spatial overlap between predators and prey, and the attendant effects on trophic interactions, fish growth rates, fishery interactions, and ultimately, fishery production. In the nGOM, for example, the vertical extent of hypoxia is highly variable and can encompass from 10% to 80% of the water column . Growth rate potential of Gulf Menhaden is lowest in hypoxic areas, but overall hypoxia effects on growth are relatively small compared to other factors (Zhang et al 2014). In contrast, hypoxia extends over a larger portion of the bottom water column in Chesapeake Bay, and has much stronger effects on Atlantic Menhaden growth (Brandt and Mason 2003). These results suggest the thickness (i.e. vertical extent) of the hypoxic layer relative to water column depth is a key factor mediating the effects of hypoxia on processes (i.e. growth) underlying fish production (Zhang et al 2014). Likewise, a 'temperature-oxygen squeeze,' whereby fish avoiding the bottom low oxygen layer are forced into surface waters where temperatures are sub-optimal for growth has been demonstrated in lakes, estuaries, and coastal seas (Coutant 1985, DeStasio et al 1996, Kraus et al 2015a, 2015b. The severity of this effect on pelagic fish production will depend on the thickness of the bottom hypoxic layer that determines the degree of vertical displacement in the water column. The thickness of hypoxia may also determine the extent to which hypoxia represents a physical barrier to migration (Zimmerman and Nance 2001), as well as susceptibility to exploitation by mobile fisheries (Goodyear et al 2008, Langseth et al 2016, Froehlich et al 2017, Purcell et al 2017. Hence, a better understanding of variation in the vertical (in addition to areal) extent of low oxygen is necessary to fully evaluate the consequences of hypoxia for processes underlying fish production and the associated consequences for fisheries and fisheries management.
However, hypoxic area is the most commonly measured, reported, and used management metric.
For example, in the case study developed here, a federal-state-tribal Action Plan for the nGOM (Task Force 2001 set a goal of reducing the midsummer hypoxic extent to 5000 km 2 . For over three decades, the relative areal extent of the hypoxic region on the Louisiana-Texas continental shelf has drawn scientific and policy attention , Rabalais et al 2002, O'Connor and Whitall, 2007, Langseth et al 2014. During that time, models ranging from regressions (Turner et al 2006, 2012, Greene et al 2009, Forrest et al 2011, Obenour et al 2012 to relatively simple (Scavia et al 2003, Obenour et al 2015 and more complex biogeochemical models (Hetland and DiMarco 2008, Wang and Justić 2009, Fennel et al 2011, Justić and Wang 2014 have been used to explore hypoxia dynamics and to provide management guidance relating the size of the hypoxic zone to key drivers. However, to date, only one model has provided load-response curves for hypoxic volume , and that analysis suggested that for the Gulf of Mexico hypoxic volume is more responsive than hypoxic area to load reductions. That model used hypoxic volume estimates determined by Obenour et al (2013) through a geostatistical reanalysis of shelf-wide hypoxia cruises conducted between 1985 and 2011. In this paper, we expand on that earlier analysis by comparing the volume-based response curves of Scavia et al (2013) with results from an independent three-dimensional hydrodynamic-biogeochemical model (Justić and Wang 2014), and discuss the implications with respect to a more ecologically-relevant hypoxia management goal.

Hypoxic extent
Hypoxic area, thickness, and volume estimates used here are from Obenour et al (2013), and are based on a rigorous geostatistical modeling framework for analyzing data collected during midsummer, quasisynoptic monitoring cruises from 1985 through 2011 (Rabalais 2011). Their simulation-based approach, which yields both areal and volumetric extent estimates with quantified uncertainty, also makes use of covariate information (i.e. trend variables such as depth and spatial position) to reduce estimation uncertainty. In their analysis, adjustments accounted for observational bias resulting from use of different sampling instruments, which could be lowered to within different proximities of the sea floor in different years. For our analyses, we used east and west shelf areal extent estimates, as described by Obenour et al (2015), and east and west shelf volumetric estimates are unique to this study. The west shelf section extends from 342.5 to 672.5 km Universal Trans Mercator (UTM) Zone 15 easting, and the east shelf section extends from 672.5 to 837.5 km UTM easting.

The models
The Scavia model used here is the same as that developed originally to relate nGOM hypoxic area to loads from the Mississippi and Atchafalaya rivers (Scavia et al 2003). It has been compared to other models (Scavia et al 2004) and used to explore nitrogen versus phosphorus control and climate change effects (Scavia andDonnelly 2007, Donner andScavia 2007). It provided guidance for the 2001 and 2008 Action Plans (Task Force 2001 and combined with three other models to provide ensemble forecasts and policy scenarios (Scavia et al 2017). The model formulation is an adaptation of the Streeter-Phelps (1925) river model that simulates DO concentration downstream from point sources of organic matter based on mass balance equations for oxygen-consumption. The model produces a subpycnocline DO concentration profile stretching from the mouth of the Mississippi River toward the Louisiana-Texas border. From that profile, total hypoxic length is determined for which DO< 3 mg l −1 because that average DO concentration roughly corresponds to a bottom water DO concentration of 2 mg l −1 and hypoxic conditions. Hypoxic length is converted to area based on historical measurements of hypoxic extent: A=57.8×length . Hypoxic volume, in turn, is determined in the model as V=A×τ 1 +A 2 ×τ 2 , where τ 1 and τ 2 are empirical coefficients that are estimated along with other coefficients through Bayesian inference , considering the fit of the model to geostatistical estimates of both hypoxic area and volume .
The model explained 69% of the variability in the Gulf of Mexico's hypoxic area and 60% of the variability in hypoxic volume over 1985-2011. This model has been used to provide annual forecasts of hypoxic area, which have compared well with measurements, especially for years without storms or high winds (R 2 =70%). Its load-response relationship indicates that reaching the hypoxia target goal set by the action plan (Task Force 2015) of 5000 km 2 requires a 62% (95% credible interval: 54%-71%) decrease in May total nitrogen load relative to the 2007-2011 average. That 62% load reduction leads to a predicted a 70% reduction in hypoxic area from the 2007−2011 mean and an 84% reduction in hypoxic volume.
The Justić and Wang model is a coupled, threedimensional hydrodynamic-biogeochemical model of the nGOM (Justić and Wang 2014). The physical underpinnings of the system are simulated using the high-resolution, three-dimensional, unstructuredgrid hydrodynamic Finite Volume Coastal Ocean Model (FVCOM). This model accurately describes the offshore circulation generated by the westerly winds during summer months, as well as the prevalent westward flow along the coast caused by the easterly winds during the rest of the study period (Wang and Justić 2009). The seasonal stratification cycle is also well represented, and simulations support the conclusion that local wind forcing and buoyancy flux resulting from riverine freshwater discharges are the dominant mechanisms affecting the circulation and stratification. DO dynamics were modeled using an expanded and revised version of the Water Analysis Simulation Program (WASP) that was fully coupled with FVCOM.
The coupled FVCOM-WASP model (Justić and Wang 2014) was driven by surface wind forcing, tidal forcing, offshore remote forcing, heat fluxes, oxygen exchanges at the air-sea interface, solar radiation, and freshwater and nutrient (nitrogen and phosphorus) fluxes from the Mississippi and Atchafalaya Rivers. Model skill was assessed for 2002 using multiple sets of observations that included time series of DO concentrations from a station within the core of the nGOM hypoxic zone, DO measurements collected during the mid-summer shelfwide cruise, and vertical DO profiles throughout the year. The agreement between the observed and simulated bottom DO levels was very high (Willmott index=0.91; R 2 =0.72; p<0.01), and the simulated hypoxic area was within 15% of the observed extent (Justić and Wang 2014).
The model results indicate that hypoxia originates in bottom waters on the mid-continental shelf, where isolated pockets of hypoxic water develop during early spring and later coalesce into a larger continuous hypoxic zone, and that the short-term variability in hypoxic extent is significant. The dynamics of bottom water hypoxia are clearly influenced by the bathymetric features of the shelf, namely the presence of both shallow shoals in the Atchafalaya Bay region and deeper shoals in the northwestern section of the region.

Responses to load reductions
The models were run independently and their estimates of changes in hypoxic area and volume resulting from 25% and 50% nitrogen load reductions from a baseline of the 2002 load were compared. The 2002 baseline was used because that was the year used to calibrate the Justić and Wang model, and a year in which the Scavia model matched both hypoxic area and volume estimates well. Response curves were generated for both models by varying the nitrogen load, keeping all other inputs constant, and reported results are the percent change from the modeled 2002 values. We find there are also differences in hypoxic thickness between the two shelf segments, with the east segment generally thicker than the west segment. In addition, there are systematic, nonlinear relationships between hypoxic area (A, km 2 ) and volume (V, km 3 ) ( figure 1(a)), and they differ between the east and west shelf segments, with the east shelf having a steeper slope ( figure 1(b)). These relationships appear well represented by a quadratic trend, consistent with the finding that thickness is positively correlated with area ( figure 2(a)). In addition, the increasing trend in hypoxic thickness over time  appears to apply equally to both shelf segments ( figure 2(b)).

Results and discussion
Meteorological conditions can affect the extent and geometry of the hypoxic region. Scavia et al (2013) used quantitative criteria to identify unusual weather years for purposes of model calibration. They identified 'storm' for 4 years with tropical storms that tend to reaerate hypoxic waters prior to the shelf wide cruise (1988,1989,1997,2003), and 'wind' for 2 years with unusually strong westerly winds (1998,2009) that tend to shoal water to the east such that, while hypoxic volume may be the same, area would be unexpectedly low. They reported that 'wind' years resulted in a larger estimated value of τ 1 , the coefficient relating volume to area (the τ 2 coefficient did not change), such that 'wind' years have relatively large volumes relative to area ( figure 1(a)). Reported 'storm' years, on the other hand, had a higher reaeration coefficient, reflecting increased mixing. However, because 'wind' and 'storm' years tend to compensate for each other in this relationship ( figure 1(a)), removing them did not significantly affect the overall relationship developed here ( figure 1(a)).
The nonlinear relationship between volume and area suggests that volume should respond more strongly to changes in load than area. Results from the two models here confirm that expectation. The Scavia model shows that reducing the 2002 load by 25% results in a 16% reduction in hypoxic area and a 27% reduction in volume (table 1). A 50% load reduction results in 46% and 60% reductions in area and  (table 1).
Because the models have very different time scales, structures, and process details, it is not surprising that their absolute values differ. One model simulates two bulk properties (annual maximum hypoxic area and volume as determined from a single shelfwide cruise) driven by nutrient load, organic matter decomposition, DO rearation, and advection. In contrast, the other model simulates the three-dimensional, daily dynamics of DO, algal production and settling, and nitrogen and phosphorus biogeochemical cycles across the shelf, driven by winds, heat flux, tides, boundary fluxes, and river discharges. Daily hypoxic area and volume are then estimated from the three-dimensional outputs. Despite the absolute differences, it is clear that in all cases, the percent volume reductions are larger than the percent area reductions. The ratio of the percent change in volume to the percent change in area (for a 25% load reduction) are all greater than one, ranging from 1.06 in the Justić and Wang model to 1.5 in the Scavia model (table 1).
Understanding variability in hypoxic thickness and volume is relevant to assessing the effects of hypoxia on fish and fisheries. Enhanced susceptibility to fishing due to increased aggregation of fish avoiding low oxygen waters has been demonstrated in a number of diverse systems (Stramma et al 2012, Kraus et al 2015a, 2015b, Froehlich et al 2017. This suggests accurately capturing the characteristics of hypoxia (e.g. area, thickness, volume, patchiness) relevant to  Table 1. Midsummer (21-26 July) hypoxic areas (km 2 ) and hypoxic volumes (km 3 ) in the northern Gulf of Mexico hypoxic zone simulated by the Justić and Wang and Scavia models for present day condition (2002), for a 25% N reduction scenario and under a 50% N reduction scenario. For a given nutrient reduction scenario, the % change with respect to 2002 is shown in the parentheses. The δV/δA ratio denotes the relative change in hypoxic volume to hypoxic area ratio. marine organisms and the fishing process is important for understanding hypoxia effects on harvested resources. For example, the Gulf Menhaden purse seine fishery, the largest volume fishery in the Gulf of Mexico, peaks in the summer months and is concentrated on the nearshore Louisiana shelf (<30 m depth) concurrent with the Gulf hypoxic zone. Menhaden primarily occupy surface waters above and inshore of the hypoxic zone to avoid low oxygen If hypoxic volume has more ecological importance than previously thought, and if hypoxic layer thickness tends to be increasing over time across the shelf (figure 2), it is worth considering volume in the context of management goals. An ensemble of four different models (Scavia et al 2017) previously indicated that reaching the 5000 km 2 action plan goal (Task Force 2015) requires a 59% reduction in nitrogen load. However, that goal and the large load reduction required to meet it may need to be re-evaluated. Because hypoxic thickness decreases with decreasing hypoxic area ( figure 2(a)) and, therefore, hypoxic volume is more sensitive than area to N reductions (table 1), volume should be considered in any future goal setting. For example, the Scavia model indicates that a 25% load reduction would reduce the hypoxic area to 14 000 km 2 (a 27% reduction), similar to the area reduction from the Justić and Wang model ( figure 3, left panel), but a 57% reduction in volume. In addition, the Justić and Wang simulations indicate that the 25% load reduction decreases volume markedly and hypoxia remains localized to a relatively thin layer near the bottom (figure 3 right panel) that most fishes and other mobile organisms can avoid more effectively.
The most recent action plan (Task Force 2015) sets an interim target of a 20% load reduction. While the model ensemble (Scavia et al 2018) suggested that a 20% reduction would have minimal impact on hypoxic area, it is in the right range to reduce volume significantly. The original 5000 km 2 goal was established through Task Force negotiations based on reducing the hypoxic extent to areas hindcast to the 1970s (Scavia et al 2003). Therefore, it may be time to reconsider the area goal, relax the current 45% load reduction target, and consider setting a more ecologicallyrelevant hypoxic volume goal.
While such a change might be welcome news because after US Farm Bill conservation programs spent more than $28 billion in the 20 Mississippi Basin states since 1995 (EWG 2017) river nitrate concentrations have not declined (Sprague et al 2011, Murphy et al 2013, and the current 5-y running average nitrate load to the Gulf is not significantly different from the 1980-1996 baseline (Task Force 2015). Clearly, to even reach a 25% load reduction target will take significantly more effort, or a new approach.
While our work focused on the nGOM, the results have broad implications given the numerous models of other systems that are also capable of simulating hypoxic volume. For example, the Scavia model described here was applied to Chesapeake Bay (Scavia et al 2006) where it was used to evaluate hypoxia controls (Stow and Scavia 2009, Evans and Scavia 2010, Liu and Scavia 2010, Liu et al 2011, produce annual forecasts of hypoxic volume (e.g. http://scavia.seas. umich.edu/hypoxia-forecasts/), and advance public understanding and participation (Testa et al 2017). While that particular model is limited to systems that are strongly one-dimensional longitudinally, other regression-based and simple biophysical models originally calibrated to hypoxic area (e.g. Turner et al 2006, 2012, Greene et al 2009, Forrest et al 2011, Obenour et al 2015 Scavia D 2016 Simulating the effect of nutrient reduction on hypoxia in a large lake (Lake Erie, USA-Canada) with a three-dimensional lake model J.