Paleoceanographic insights on recent oxygen minimum zone expansion: lessons for modern oceanography.

Climate-driven Oxygen Minimum Zone (OMZ) expansions in the geologic record provide an opportunity to characterize the spatial and temporal scales of OMZ change. Here we investigate OMZ expansion through the global-scale warming event of the most recent deglaciation (18-11 ka), an event with clear relevance to understanding modern anthropogenic climate change. Deglacial marine sediment records were compiled to quantify the vertical extent, intensity, surface area and volume impingements of hypoxic waters upon continental margins. By integrating sediment records (183-2,309 meters below sea level; mbsl) containing one or more geochemical, sedimentary or microfossil oxygenation proxies integrated with analyses of eustatic sea level rise, we reconstruct the timing, depth and intensity of seafloor hypoxia. The maximum vertical OMZ extent during the deglaciation was variable by region: Subarctic Pacific (~600-2,900 mbsl), California Current (~330-1,500 mbsl), Mexico Margin (~330-830 mbsl), and the Humboldt Current and Equatorial Pacific (~110-3,100 mbsl). The timing of OMZ expansion is regionally coherent but not globally synchronous. Subarctic Pacific and California Current continental margins exhibit tight correlation to the oscillations of Northern Hemisphere deglacial events (Termination IA, Bølling-Allerød, Younger Dryas and Termination IB). Southern regions (Mexico Margin and the Equatorial Pacific and Humboldt Current) exhibit hypoxia expansion prior to Termination IA (~14.7 ka), and no regional oxygenation oscillations. Our analyses provide new evidence for the geographically and vertically extensive expansion of OMZs, and the extreme compression of upper-ocean oxygenated ecosystems during the geologically recent deglaciation.

Introduction

Resolving records of global change through the most recent deglaciation event (18–11 ka) is one of the primary challenges to developing cohesive and robust theories regarding rapid climate change [1]. The last deglaciation was a profound event in the global climate system, wherein atmospheric [CO2] increased by 80–100 ppmv [2, 3], global average temperature rose 3–4°C [4], and sea levels rose ~110 m (Fig. 1) [5, 6]. This change was also accompanied by the pervasive loss of dissolved oxygen in the upper ocean [7], with unknown impacts on large marine ecosystems. The coupling between climate, carbon emissions and subsurface dissolved oxygen is a salient and critical element of anthropogenic climate change. The global inventory of ocean oxygen is predicted to decline between 1 and 7% by the year 2100, through stratification, ventilation reduction and decreased O2 solubility (e.g., [8]), and the hypoxic water volume in the global ocean is predicted to increase by 50% [9]. Climate model results beyond the 100-year window reveal extensive oceanic deoxygenation, on thousand-year timescales, under “business-as-usual” carbon emission scenarios, and show that oceanic deoxygenation is a fundamental and long-lasting property of anthropogenic carbon perturbation (e.g. [10]).

Figure 1

Deglacial changes in Antarctic temperature (Vostok ice core record, purple line) [146, 147], Greenland temperature (GISP2 ice core record, blue line) [106, 107], sea level (black line) [5, 6] and atmospheric pCO2 (red line) [3].

Glacial Termination IA (14.7 ka) is an event of rapid warming in the Northern Hemisphere, which initiates the warm interval of the Bølling-Allerød (B/A) from 14.7–12.9 ka. The Younger Dryas (YD), a reversal towards cool conditions from 12.9–11.7 ka, follows the B/A. The YD ends with glacial Termination IB (11.7 ka), a subsequent rapid warming event. Deglacial warming in the Southern Hemisphere begins at 18 ka.

Hypoxia substantially degrades ecosystems through mass mortality events, the alteration of food-web structures and the loss of habitat (e.g., [11]). Changes in [O2] in the ocean interior have broad consequence for global biodiversity, marine economic resources and ocean management [12]. To grasp the scale of future hypoxia disturbance, the paleoceanographic (pre-instrumental) record of natural variability provides a critical analytical and interpretive window. Ocean sediments are climatic and environmental archives, which preserve geochemical, microfaunal and sedimentary evidence that record globally relevant Earth system events, similar to other climate records such as the Greenland Ice Sheet Project [13]. Recent investigations reveal that paleoceanographic investigations hold valuable insight into modern environmental conservation and management [14, 15].

Here we synthesize published continental margin sediment core records to investigate Oxygen Minimum Zone (OMZ) changes through the last deglaciation. We build on previous syntheses of oxygenation proxy records (e.g., [7, 16, 17]), and provide a focus on regional-scale sensitivity. By integrating sediment records, sea level change, and high-resolution bathymetry, we provide geospatially analyzed paleoceanographic data that are interpretive baselines for modern oceanography and global environmental change.

The role and importance of OMZs

OMZs are tightly coupled to upwelling systems and Eastern Boundary Currents, such as the California Current, the Humboldt Current and the Benguela Current, as well as the Oman and Pakistan Margin in the Indian Ocean (Fig. 2). In these regions, respiration within the pycnocline depletes dissolved oxygen and simultaneously enriches seawater in the carbon and nitrogen byproducts of respiration [18]. Marine denitrification occurs within OMZs (e.g., [19, 20]) therefore the physical extent and intensity of OMZs is inherently coupled to the oceanic nitrogen cycle. OMZs form at shelf and upper slope depths, and are considered to be unique biological, geochemical and evolutionary environments, analogous to cold seep or deep-sea vent environments [21]. As continental margin ecosystems transition from well oxygenated surface waters to the hypoxic core of the OMZ ([O2] = 0.5–0.1 ml L-1), faunal diversity, trophic structures and physiological strategies change (e.g., [22, 23]). OMZ oxygenation gradients produce successional biological zonation and are fundamental habitat barriers for benthic and pelagic organisms [21].

Figure 2

Global Oxygen Minimum Zones, including (a) Upper depth (in meters) of intermediate water hypoxia ([O2]<1.4 ml L-1) and (b) thickness (in meters) of intermediate water hypoxia ([O2]<1.4 ml L-1).

The geospatial distributions of severely hypoxic [O2] minimums (of [O2] = 0.5 ml L-1 and [O2] = 0.2 ml L-1) are depicted on both panels as white lines. For the upper panel, regional blocks are defined by black lines to highlight where paleoxygenation reconstructions were completed. Data from World Ocean Atlas [192].

For this work, we follow the hypoxia thresholds and categories defined in [24], which synthesizes the existing hypoxia vernacular, to draw thresholds that are biologically meaningful (Table 1). Mild hypoxia begins at [O2]<2.45 ml L-1 and is the threshold where sensitive species exhibit avoidance reactions. Intermediate hypoxia, often referred to as “coastal hypoxia”, occurs at [O2]<1.4 ml L-1 and is the threshold wherein ecosystems are dominated by organisms with adaptive features. Severe hypoxia ([O2]<0.5 ml L-1) is a threshold at which mass mortality is induced for most organisms, past which only highly specialized species can survive [24].

Table 1

Hypoxia categories and associated oxygen concentrations, partial pressures and saturations for each category.

Hypoxia Category	ml O2 L-1	µmol O2 kg-1	mg O2 L-1	pO2			% O2saturation
				T = 25°C	T = 17°C	T = 12°C	T = 25°C	T = 17°C	T = 12°C
Mild	2.45	107	3.5	106	93	84	51≈50	45	40
Intermediate	1.4	61	2.0	60	53	48	29≈30	25	23
Severe	0.5	22	0.71	22	19	17	11≈10	9	8

Salinity (34 psu) and hydrostatic pressure (P = 10 bar) are assumed constant. Temperature columns indicate the temperature used for partial pressure and saturation calculations at the associated concentrations. Data table adapted from Hofmann et al., [24]. The intermediate category is described as “coastal” in Hofmann et al., [24]. We refrain from using this term here, to prevent confusion between hypoxia categories and offshore habitat locations

Oxygenation proxies in paleoceanographic records

We employ a multi-proxy approach to oxygenation reconstructions here, and include core data across sedimentary, faunal and geochemical proxies (Table 1). Fig. 3 depicts, in schemata form, how regional multi-proxy oxygenation data are interpreted, and we discuss each proxy below.

Figure 3

Schemata of a multi-proxy approach to interpreting hypoxia categories, including severe hypoxia ([O2]<0.5 ml L-1), intermediate hypoxia ([O2]>0.5–1.5 ml L-1) and mild hypoxia to oxic conditions ([O2]>1.5 ml L-1).

These hypoxia categories are detailed in Table 1, and follow Hofmann et al., [24]. Hypoxia proxies include [Re], [Mn], [U], [Cd], [Mo], δ15N, foraminiferan communities, and sedimentary laminations. Units for each proxy reflect the cited literature, which constrains the proxy to a specific oxygenation category.

Sedimentary structures. Severe hypoxic conditions preclude benthic macrofauna, preventing bioturbation and thereby allowing for the preservation of laminated sediments (Table 2) [25-27]. For example, annual sediment laminations (i.e. unmixed, fine-grained sediment displaying distinct, continuous layering) are formed in modern California margin basinal features at [O2]<0.1 ml L-1 [25, 28]. Laminations are formed under the absence of bioturbating invertebrates [29] and sufficiently high organic carbon export from the surface [16], and are one of the clearest indicators of severe hypoxia in benthic environments.

Table 2

Seafloor hypoxia proxies for paleoceanographic reconstructions, partitioned by the thresholds and capacity each proxy has to record fine-scale changes in seafloor hypoxia, as well as organic flux to the seafloor [24].

Proxy	Mild Hypoxia ([O2] = 2.45–1.4 ml L-1)	Intermediate Hypoxia ([O2] = 0.5–1.4 ml L-1)	Severe Hypoxia ([O2]<0.5 ml L-1)	Indicator of increased organic flux to seafloor?	Notes
Laminations	Not Present		Laminations only occur with extremely severe seafloor hypoxia ([O2]<0.1 ml L-1), where bioturbating benthic fauna are not present [26, 27].	Yes. Lamination formation requires high surface export [16].
Foraminifera	Above 2 ml L-1, foraminiferan composition is likely not altered by changes in bottom water oxygenation [199].	Foraminiferan communities in intermediate-severely hypoxic sediments are dense and associated with opportunistic taxa and specific morphologies [200–202]. Marker species, affiliated with a narrow oxygenation range, can be used to reconstruct seafloor oxygenation on a very fine-scale [31, 203].		Organic flux to the seafloor alters the composition and density of foraminiferan communities [204, 205], however interpretation of these community traits is not straightforward in a low oxygen setting (e.g., [206]).	Foraminifera are well adapted to the extreme chemical heterogeneity of oxygenation, methane enrichment, organic flux and sulfur-reducing environments found on continental margins [207].
δ15N	No significant denitrification, with no isotopic fractionation.		Denitrification occurs in severely hypoxic conditions ([O2]<0.23 ml L-1) [19]. The isotopic signal of denitrification reflects regional changes in biologically available N pool [35].	δ15N is an indirect record of organic flux, as photosynthesis isotopically fractures the nitrogen pool.	The δ15N of particulate nitrogen varies with the degree of surface nutrient utilization, commonly termed productivity [37, 38], and water column denitrification [128]. For continental margins, sediment denitrification contributes to the isotopic signal [19, 20]. Denitrifying bacteria are facultative anaerobes, meaning they are able to respire either oxygen or nitrogen oxides.
Molybdenum [Mo]	No accumulation increase [43].		Accumulation occurs when [O2]<0.45 ml L-1 [43]. High accumulation rates (>2 μg g-1) are related to the presence of anoxia [36].	Yes, non-lithogenic particulate Mo, associated with sinking particles, contributes ≤15% of authigenic Mo accumulation [43].	Mo is a biologically-essential, and is conservative in behavior. The pathway for authigenic Mo accumulation in sediment is unresolved, however it is clear that sulfide concentrations must be >2 μM for active Mo formation [43] Concurrent enrichment of Mo and Re is an indicator of an anoxic depositional environment [41].
Cadmium, [Cd]	No accumulation increase	Some accumulation at lower levels in “coastal hypoxia” [208].	Significant Cd enrichment (>2 ppm) is associated with laminated depositional records, a suboxic to anoxic seafloor and the presence of trace presence of H2S [42].	Composition of Cd in sinking particles is similar to plankton material; therefore changes in productivity may affect levels of enrichment [41].	Cd is involved with nutrient cycling and scavenging. It is released in sediment during the diagenesis of organic matter, forming insoluble sulfides when trace amounts of H2S are present [209]. Cd enrichment occurs in surface sediments [41].
Uranium, [U]	No accumulation increase	U is enriched in “suboxic” and anoxic sediments as a result of reduction and precipitation (>5μg g-1) [210].	U enrichment occurs in bioturbated sediments (>0.1 ml L-1), indicating enrichments at an unclear “suboxic” threshold [42].	Unclear, as the composition of U in plankton is not known.	U is conservative in behavior, and the accumulation of U may be kinetic in nature. U accumulation may be promoted by kinetic effect, in areas of low sedimentation rates [17, 211, 212]. U enrichment occurs in subsurface sediments [41].
Rhenium, [Re]	No enrichment	Re enrichment occurs in bioturbated sediments ([O2]>0.1 ml L-1), indicating enrichment at an unclear “suboxic” threshold [42].	Re is highly enriched in “suboxic” and anoxic sediments [42, 210].	Unclear, as the composition of Re in plankton is not known.	Re is conservative in behavior, diffuses across the sediment-water interface, and the direct mechanism for accumulation is unknown. Free H2S is not critical for the accumulation of Re [210]. Concurrent enrichment of Mo and Re is an indicator of an anoxic depositional environment [41].

A diagram for how multiproxy <span class="Chemical">records a<span class="Chemical">re interpreted into paleoxygenation categories is shown in Fig. 3.

Microfossils. Oxygen is a physiochemical parameter that acts as a limiting factor, thereby determining the distribution of where organisms can live based on their physiological requirements. Within the optimal range of a limiting factor, a species will reach maximum competitiveness and maximum abundance [30], resulting in tight faunal affiliations with specific oxygenation ranges. Benthic foraminifera have species-specific oxygenation thresholds and therefore marker taxa function as oxygenation proxies [31]. Metrics of foraminiferan community structure, such as diversity and density, also record changes in seafloor oxygenation [22, 23, 32]. Due to their opportunistic responsiveness to environmental change, Foraminifera are ideally suited for high-resolution oxygenation reconstructions [32-34].

Geochemistry. Isotopic and trace element geochemical proxies of oxygenation are complex, as the cycling, preservation and attribution of these elements can be controlled by processes such as of surface productivity and flux, water column processes, sediment-water interface flux and diffusion, and sub-surface sedimentary processes. Denitrification in continental margin regions occurs under hypoxic conditions in both the water column and sediment ([O2]<0.23 ml L-1) [19, 20], and changes in the nitrogen isotopic signal are a product of regional changes in the biological available N pool [35]. Although intense denitrification is often indicative of severe hypoxia, more subtle changes in δ15N can be attributed to a range of related localized factors, such as source water and nutrient availability [19, 36–39]. Redox sensitive trace elements such as cadmium (Cd), uranium (U), chromium (Cr) and rhenium (Re) occur in bulk sediment and carbonate fossils, and can be used to characterize the redox state of the seafloor as well as functioning as proxies of paleoenvironmental redox changes [40-43]. Paleoproductivity proxies are useful supportive proxies in OMZ reconstructions, due to the mechanistic coupling between surface export and subsurface respiration [44]. Additionally, the δ13C record of planktonic carbonate reflects the surface ocean isotopic pool, and can be used to reconstruct surface water productivity and carbon export [45-47].

Methods

Sedimentary archives

We include cores in this review that met criteria for further geospatial analysis (Table 3). These criteria dictated sediment records that:

Table 3

Deglacial sediment core metadata and references for sites used in geospatial analyses.

Core ID	Modern water depth (m)	Latitude	Longitude	Sedimentation rate (cm kyr-1) or Mass Accumulation Rate (MAR; g cm-3 kyr-1)	14 C dating	δ 18 O stratigraphy	Laminations	Microfossils	Geochemistry [δ15N]	Geochemistry [Mo]	Geochemistry [Cd]	Geochemistry [U]	Geochemistry [Re]	References
MV0811–15JC	418	34° 22.2′ N	120° 7.8′ W	42.2–100	x	x		x						[103]
ODP 893A	576.5	34° 17.251′ N	120° 2.2′ W	140	x	x	x	x	x	x	x	x	x	[27, 29, 83, 85, 86, 88, 90, 100, 213]
MD02–2503	570	34° 16.48′ N	120° 2.24′ W	100		x	x	x						[89, 102, 103]
MD02–2504	481	34° 13.998′ N	119° 52.122′ W	80–100		x	x	x						[89, 102, 103]
ODP 1017E	955	34° 32.1′ N	121° 6.43′ W	23	x	x		x	x	x	x	x	x	[17, 94, 214]
F2–92–P3	786	35° 37.39′ N	121° 36.28′ W	~17	x		x			x				[40, 95, 96, 105]
F2–92–P34	610	35° 1.85′ N	121° 13.54′ W	~10	x		x			x				[96, 105]
F2–92–P40	760	35° 25.09′ N	121° 24.95′ W	~42	x		x			x				[40, 96, 105]
ODP 1019E	980	41° 40.9′ N	124° 55.8′ W	41	x	x	x		x	x	x	x	x	[42, 105, 215, 216]
F–8–90–G21	1650	37° 13.4′ N	123° 14.6′ W	4–23	x	x		x			x			[95–97]
Core JT96–09	920	48° 54.067′ N	126° 53.034′ W	5–169	x			x		x	x	x	x	[92, 217]
EW0408–85JC	682	59° 33.32′ N	144° 9.21′ W	9–500	x		x	x	x	x	x	x		[58, 67, 72, 74]
EW0408–11JC	183	55° 37.6′ N	133° 30.69′ W	37–280	x					x		x		[72]
CH84–14	978	41° 44′ N	142° 33′ E	55–208	x		x	x		x		x		[68, 179]
RAMA 44PC	2980	53° 0′ N	164° 39′ E	4–29	x			x				x		[68, 70, 71, 74]
GGC-55/JPC-56	818	27° 28.16′ N	112° 6.26′ W	83	x	x	x	x	x					[16, 128, 182, 218]
MV99-PC08	705	23° 28.2′ N	111° 36′ W	30	x		x		x	x	x			[16, 105, 121, 128, 182, 218–220]
MV99-PC14	542	25° 12′ N	112° 43.2′ W	108	x		x							[16, 105, 121, 219, 220]
NH15P	420	22° 4′ N	106° 28.8′ W	4–17	x	x	x	x	x	x	x	x	x	[16, 119, 120, 122, 219, 221, 222]
NH8P	1018	22° 23.3′ N	107° 4.5′ W	11.5–16	x	x	x	x	x					[16, 119, 120, 122, 126, 127, 221, 222]
MD02–2508	606	23° 27.91′ N	111° 35.74′ W	15–43	x		x		x		x			[16, 17, 120, 126, 127, 129, 221]
DSDP Site 480	655	27° 54′ N	111° 39′ W	77–114	x		x	x						[16, 17, 26, 129]
ODP Site 1242	1364	7° 51.352′ N	83° 36.418′ W	~10	x	x			x					[16, 26, 124, 181, 223]
W7706–37K	370	-13° 38′ S	76° 51′ W	~100	x		x		x					[124, 145, 151, 181, 223, 224]
W7706–40K	186	-11° 15′ S	77° 58′ W	~55	x		x		x					[145, 151, 224–226]
W7706–41K	410	-11° 21′ S	78° 7′ W	N/A	x				x					[145, 151, 224–226]
ME0005A-27JC	2203	1° 51.2′ N	82° 47.2′ W	~3 MAR	x	x			x			x		[149, 151, 227, 228]
ME0005A-24JC	2941	0° 1.2′ N	86° 27.6′ W	~4.5–12 MAR	x				x			x		[124, 149, 181, 227–229]
Y69–71P	2740	0° 6′ N	86° 28.8′ W	~3–6 MAR	x	x			x			x		[124, 149, 181, 229, 230]
P7	3085	2° 36′ N	83° 59.4′ W	~1.5–2.5 MAR	x				x			x		[148, 149, 230, 231]
TR163–31P	3209	-3° 36′ S	83° 57′ W	~3–9 MAR	x	x			x			x		[148, 149, 227, 228, 231–235]
TR163–19P	2348	2° 15.6′ N	90° 57′ W	~2–4 MAR		x			x			x		[149, 227, 228, 232–235]
GeoB7139–2	3269	-30° 12′ S	71° 59′ W	~13	x				x					[123, 149, 150, 181, 228, 232, 235, 236]
V19–30	3091	-3° 22.8′ S	83° 30′ W	~7		x			x			x		[123, 144, 150, 181, 236, 237]
RC13–140	2246	2° 52.2′ N	87° 45′ W	>5	x	x	x		x			x		[144, 237, 238]
RC11–238	2573	-1° 31.2′ S	85° 49.2′ W	>5	x	x	x		x			x		[144, 238]

exist in a geographic <span class="Chemical">region that includes a modern <span class="Chemical">OMZ,

range in age from approximately 20–0 ka,

have well-constrained chronologies (radio<span class="Chemical">carbon and/or isotope stratigraphy),

include one or mo<span class="Chemical">re primary proxies for <span class="Disease">hypoxia, as described above.

Associated data provided here include modern water depth of core extraction (m), latitude and longitude, sedimentation (cm kyr-1) or mass accumulation rate (g cm-2 kyr-1), the presence of 14C dating or δ 18O stratigraphy, and the presence of published hypoxia proxy records (including laminations, microfossils, or geochemistry).

We rely on published chronologies for the deglacial and post-glacial core material, and do not reinterpret any chronologies. Regional archives are assembled on a unifying age axis, which is determined only by previously published chronologies. Sedimentation rates for sites <1000 mbsl span a minimum of 9 cm kyr-1 to a maximum of 500 cm kyr-1, with the majority sedimentation rate ranging from 20–100 cm kyr-1. Sedimentation rates for sites >1000 mbsl span a minimum of 4 cm kyr-1 to a maximum of 16 cm kyr-1. Mass Accumulation Rates (MARs) are published for deep sites (>2200 mbsl) and they range from 2–9 g cm-3 kyr-1.

Paleoxygenation was assessed based on proxy evidence available for each for core. If proxy evidence indicated hypoxia, we partitioned that signal into intermediate or severe hypoxia groups, based on the biologically-relevant classification scheme identified by Hofmann et al., [24]. To ensure hypoxia designations are as conservative as possible, we refrained from any hypoxic designation without evidence. In practice, this meant we only designated hypoxia where the proxy evidence was explicit. In the absence of that evidence, we did not designate any hypoxia (i.e., in the absence of laminations, we did not designate a potential range from oxic-intermediately hypoxic). Importantly, we found no instance where multiple proxies within one archive produced conflicting hypoxia reconstructions.

Sedimentary archives were then assigned paleodepths based on the modern water depth for the core site, the age model of that core, and deglacial eustatic sea level change. Paleodepths were calculated using estimated global eustatic sea level fluctuations constructed from multiple deglacial sea level records (Fig. 1) [5, 6]. We acknowledge eustatic sea level is a simplification of local sea level trends through the deglaciation, which are also impacted by glacio-hydro-isostatic effects [48, 49].

Bathymetric and Geospatial Analyses

Data visualizations applied to bathymetry, or submarine topography, provide unique perspectives into the distribution of seafloor ecosystems. Here we apply bathymetric masks associated with reconstructed paleoxygenation, wherein mask depths were chosen based on the extent of regional hypoxia across two or more cores at key temporal deglacial intervals. If only a single hypoxia record was available for a deglacial event, a hypoxic mask ±50 m from the core depth was applied. Oxic and hypoxic benthic surface area (km2) and water volume (km3) were quantified for temporal intervals through the deglaciation. To produce geospatial maps of paleoxygenation, we synthesized regional scale patterns in oxygenation, took into account eustatic sea level change, and analyzed the deglacial changes in both hypoxic seafloor surface area and hypoxic water volume from the SRTM30_PLUS global bathymetry dataset [50] using ArcGIS geospatial software [51]. Analyses were limited to the continental margin within a 400 nautical mile buffer offshore of the continental coastline. One exception to this is the inclusion of the seafloor within 400 nautical miles around the Galapagos Islands, in order to capture the associated cluster of deep core sites nearby. Within each region, the analysis was limited to the seafloor and water column above a specified isobath, selected below the deepest regional core depth.

We include four regions for which geospatial analyses were conducted: the Subarctic Pacific (SP; 0–3,200 mbsl; southern latitude limit at 38° 30′ N in the Western Pacific and 49° 30′ N in the Eastern Pacific), the California Current (CC; 0–2,400 mbsl; 31° 40′-49° 30′ N), Mexico Margin (MM; 0–1,200 mbsl; 20°-30° N) and the Humboldt Current and Equatorial Pacific (HC; 0–3,300 mbsl; 10° 30′ N-32° S) (Fig. 4). The Benguela Current (BC) and the Oman and Pakistan Margin (OPM) are discussed, however sediment records from these regions did not meet the criteria established to warrant geospatial analysis. We limit our statements of specific volume and surface area changes to Tables 4 and 5. These analyses represent the state of paleoceanographic records of deglacial OMZ expansion, including caveats associated with proxy records and limitations of the spatial resolution of core sites; our work may also be used to highlight existing knowledge gaps.

Figure 4

Deglacial sediment core locations for four Eastern Pacific regions, including (a) the Subarctic Pacific, (b) the Mexico Margin, (c) the California Current, and (d) the Equatorial Pacific and Humboldt Current.

Table 4

Surface area calculations of deglacial oxygenation changes for four Eastern Pacific continental margins: Subarctic Pacific (SP), California Current (CC), Mexico Margin (MM), and the Humboldt Current and Equatorial Pacific (HC).

Region	Ka (Calendar Age)	Regional margin surface area (km2), corrected for sea level rise	Total vertical extent of hypoxia (in meters)	Oxic upper ocean, from surface to upper subsurface hypoxic boundary (km2)	Total Oxic (km2)	Oxic (%)	Severe Hypoxia (km2)	Severe Hypoxia (%)	Intermediate Hypoxia (km2)	Intermediate Hypoxia (%)
SP	14	1,177,000	2,298	426,000	603,000	51	574,000	49	0	0
	10	1,259,000	0	N/A	1,259,000	100	0	0	0	0
	4	1,325,000	100	335,000	1,167,000	88	158,000	12	0	0
CC	18	274,000	0	N/A	274,000	100	0	0	0	0
	14	285,000	1,232	60,000	153,000	54	50,000	17	83,000	29
	12	292,000	0	N/A	292,000	100	0	0	0	0
	10	297,000	1,232	80,000	165,000	56	10,000	3	122,000	41
	4	309,000	528	99,000	254,000	82	10,000	3	45,000	15
MM	18	158,000	100	109,000	148,000	94	10,000	6	0	0
	14	181,000	603	79,000	120,000	66	61,000	34	0	0
	10	216,000	518	111,000	154,000	71	61,000	28	0	0
	4	264,000	398	151,000	221,000	84	43,000	16	0	0
HC	18	2,648,000	100	716,000	1,370,000	52	0	0	1,278,000	48
	13	2,572,000	3,022	98,000	621,000	24	54,000	2	1,897,000	74
	10	2,507,000	2,749	185,000	421,000	17	18,000	1	2,068,000	82
	4	2,428,000	100	428,000	2,388,000	98	0	0	40,000	2

Table 5

Volumetric calculations of deglacial oxygenation changes for four Eastern Pacific continental margins: Subarctic Pacific (SP), California Current (CC), Mexico Margin (MM), and the Humboldt Current and Equatorial Pacific (HC).

Region	Ka (Calendar Age)	Depth of seafloor analyzed (m), corrected for sea level rise	Regional margin volume (km3), corrected for sea level rise	Oxic upper ocean, from surface to upper subsurface hypoxic boundary (km3)	Total Oxic (km3)	Percent of margin oxic (%)	Severe hypoxia (km3)	Percent of margin severely hypoxic (%)	Intermediate hypoxia (km3)	Percent of margin intermediately hypoxic (%)
SP	14	3,214	1,524,000	456,000	479,000	32%	1,045,000	68%	0	0%
	10	3,255	1,500,000	N/A	1,500,000	100%	0	0%	0	0%
	4	3,299	1,480,000	144,000	1,397,000	94%	83,000	6%	0	0%
CC	18	2,280	316,000	N/A	316,000	100%	0	0.00%	0	0%
	14	2,314	302,000	75,000	119,000	39%	86,000	28%	98,000	32%
	12	2,336	295,000	N/A	295,000	100%	0	0%	0	0%
	10	2,355	288,000	83,000	118,000	41%	18,000	6%	152,000	53%
	4	2,399	277,000	91,000	189,000	68%	17,000	6%	70,000	25%
MM	18	1,080	65,000	60,000	62,000	96%	2,000	4%	0	0%
	14	1,114	66,000	35,000	38,000	58%	28,000	42%	0	0%
	10	1,155	69,000	41,000	43,000	62%	26,000	38%	0	0%
	4	1,199	74,000	49,000	54,000	73%	19,000	26%	0	0%
HC	18	3,180	6,155,000	4,731,000	4,794,000	78%	0	0%	1,360,000	22%
	13	3,222	5,885,000	265,000	306,000	5%	531,000	9%	5,047,000	86%
	10	3,255	5,651,000	731,000	740,000	13%	225,000	4%	4,686,000	83%
	4	3,299	5,245,000	2,330,000	5,151,000	96%	0	0%	193,000	4%

Continental margin surface area is corrected for eustatic sea level rise [5, 6] and calculated using the SRTM30_PLUS global bathymetry dataset [50]. Total vertical extent of hypoxia is stated, including intermediate and severe. Hypoxic and oxic margin surface area (km2) and percent (%) are calculated for the seafloor above a region-specific isobaths: SP (0–3,300 mbsl), CC (0–2,400 mbsl), MM (1–1,200 mbsl), and HC (0–3,300 mbsl). Upper ocean oxic surface area (km2), from the surface ocean to the upper subsurface hypoxic boundary, is included.

Continental margin volume is corrected for eustatic sea level rise [5, 6] and calculated using the SRTM30_PLUS global bathymetry dataset [50]. Hypoxic and oxic margin volume (km3) and percent (%) are calculated for the water column above a region-specific isobath: SP (0–3,300 mbsl), CC (0–2,400 mbsl), BM (1–1,200 mbsl), and HC (1–3,300 mbsl). Upper ocean oxic volume (km3), from the surface ocean to the upper subsurface hypoxic boundary, is included.

Results

Subarctic Pacific

The SP is dominated at the surface by the eastward flowing North Pacific Current, which branches into the Alaskan Gyre and the southward-flowing CC (e.g., [52]). North Pacific Intermediate Water (NPIW) forms in the Sea of Okhotsk and mixes laterally into the Pacific subtropical gyre [53]. Deep water in the North Pacific is poorly ventilated, nutrient-rich and oxygen-depleted [18], and plays a critical role in maintaining extensive intermediate water OMZs [54, 55]. The SP has a modern, seasonal OMZ in the Gulf of Alaska (Fig. 2), attenuated in the summer, with the core of the OMZ from 670–1060 m in the water column [56]. Hypoxic waters spread westward, south of the Aleutian ridge and the Russian Kamchatka Peninsula, during winter months and form a seasonal OMZ that spans the Subarctic Pacific [56].

The SP is a region with inherent complexities, due to the influence of the Cordilleran Ice Sheet [57], the variability of sea ice formation [58], and circulation sensitivity to the episodic closures of the Bering Straight [59]. Previous investigations have characterized the Last Glacial Maximum (LGM) (Fig. 1) as a time of low surface ocean productivity, high ice-rafted debris sediment flux, and cold surface temperatures [60-65]. Archives track the deglacial oscillations of the warm Bølling-Allerød (B/A) and cool Younger Dryas (YD), consistent with the expected atmospheric teleconnections between the North Pacific and the North Atlantic [66]. The B/A and Holocene are warmer and more productive across the region, from the Gulf of Alaska [67] to the Western margin of Japan [68]. In step with the changes in surface productivity, seafloor hypoxia developed during the warm, productive intervals of the deglaciation across intermediate (600 mbsl, site RAMA 44PC) and deep (2,900 mbsl, site EW0408–85JC) water depths [58, 67–72]. The productivity and oxygenation oscillations of margin and shelf environments to Termination IA, the B/A and the YD reveal the substantial changes this contiguous high latitude environment recently underwent.

Paleoxygenation reconstructions for the Subarctic Pacific. Four sediment cores meet the criteria for deglacial reconstructions: EW0408–85JC, EW0408–11JC, RAMA 44PC, and CH84–14 (Fig. 4; Table 3). Cores EW0408–85JC and EW0408–11JC were both collected from the Gulf of Alaska, RAMA 44PC was collected east of Kamchatka Peninsula, and CH84–14 was collected east of Hokkaidō Island, Japan. Shifts in seafloor oxygenation throughout the SP exhibit a cohesive, though limited, hydrographic picture (Figs. 5 and 6). We limit our SP reconstruction to mid-way through the deglaciation (14–4 ka), wherein the regional margin is extensively hypoxic, followed by a contraction to shallow, upper intermediate water hypoxia in the Holocene. The Bering Sea and Sea of Okhotsk are not included in the analysis, due to their unique, regional-scale oceanography.

Figure 5

Subarctic Pacific (SP) deglacial core data synthesized into hypoxia categories.

Changing deglacial core depths reflect global eustatic sea level change. The encircled number adjacent to each core label corresponds to the number of available oxygenation proxies, which are enumerate in Table 3. Vertical grey bars correlate to temporal intervals in OMZ geospatial reconstructions for this region.

Figure 6

Subarctic Pacific bathymetric seafloor masks and surface area (km2) histograms of deglacial hypoxia impingement for (a) 14 ka, (b) 10 ka, and (c) 4 ka.

The Bering Sea and Sea of Okhotsk are not included in the analysis, due to their unique, regional-scale oceanography. Seafloor is selected between 0–3,200 mbsl, with a southern latitude limit at 38° 30′ N in the Western Pacific and 49° 30′ N in the Eastern Pacific with a northern limit along the Aleutian Arc. Analyses were limited to the continental margin within a 400 nautical mile buffer offshore of the continental coastline and the Aleutian Arc. The changing gray shoreline through the panels depicts the paleo-shoreline. At 14 ka, severe hypoxia ranged from 596–2,894 mbsl. At 10 ka, the water column was oxygenated, and at 4 ka severe hypoxia was found between 132–232 mbsl.

Records from ODP Site 882 (Detroit Seamount, water depth 3244 mbsl) indicate that the deep glacial North Pacific was “suboxic” (based upon [U]) [73]. Although lacking a high enough sedimentation rate necessary to be included in this analysis, ODP Site 882 indicates that the deep SP (>3000 m) gained oxygen through the deglaciation. At Termination IA, the three deepest sites record the development of severe hypoxia, associated with laminations [67] and high authigenic [U] and [Mo] [74] in EW0408–85JC, and high sedimentary δ15N values in CH84–14 and RAMA 44PC [68]. These proxy records indicate that at 14 ka severe hypoxia ranged from 596–2,894 mbsl (Fig. 5). Regional hypoxia was followed by subsequent oxygenation in intermediate and deep waters through the YD, Termination 1B and the early Holocene. In the mid-Holocene, the shallowest site recorded the development of severe hypoxia in upper intermediate waters ([Mo]) [74].

Paleoceanographic reconstructions from the SP indicate extensive and severe hypoxia developed during the warming of the deglaciation, and suggest severe hypoxia as shallow as ~600 mbsl and across ~2,298 m of contiguously hypoxic water column (Fig. 6). The absence of hypoxia during the YD cooling of the Northern Hemisphere indicates that the subsurface SP is sensitive to rapid climatic oscillations. Coherent oxygenation oscillations, from the Gulf of Alaska to Eastern Japan, indicate that geographically widespread changes occur in response to global climate change.

California Current

The CC is the eastern limb of the North Pacific Subtropical Gyre (e.g., [75-77]), and is driven onshore in the Northern Hemisphere spring [54]. The CC is characterized by seasonal upwelling cycles, a highly productive continental margin, and an acute and extensive subsurface OMZ [54, 55, 78]. Additionally, a strong nearshore, poleward Counter Current, know as the Davidson Counter Current, has a surface and subsurface expression [79] and brings warmer, saltier, less oxygenated, and equatorially-influenced water up the continental margin. The CC OMZ exhibits an onshore-offshore gradient, wherein the lowest oxygen concentration values are found within 200 km from the coastline [54]. The upper boundary of the CC OMZ (as defined by [O2 = 1.4 ml L-1) is relatively deep (~600 mbsl), and the OMZ is thick (~1200–1500 m) (Fig. 1). Oxygen delivery to the CC subsurface is balanced by contributions from well-oxygenated, northern-sourced NPIW and poorly oxygenated, southern-sourced intermediate waters [53].

Early investigations into CC continental margin sediments described the presence of laminations during the late Pleistocene to Holocene transition [80-82]. ODP Site 893A (576.5 mbsl) stands out amongst the region’s extensive suite of cores [27, 83–88] as a globally recognized site with high-resolution synchroneity to Northern and Southern Hemisphere climate events [40, 88, 89]. The expansion and contraction of the regional CC OMZ occurred with remarkable fidelity to the glacial terminations, warming intervals, and cooling oscillations of the Northern Hemisphere, with clear regional-scale deoxygenation associated with Termination IA, the B/A, Termination IB and the start of the Holocene (e.g., [27, 29, 40, 42, 90, 91]). Sediment cores from Vancouver Island record OMZ expansion from 13.5 to 12.6 ka, slightly delayed from more southern CC sites [92].

Paleoxygenation Reconstructions for the California Current. Eleven deglacial sediment cores were included in our geospatial analysis, from Santa Barbara Basin (ODP Site 893A, MD02–2504, MD02–2504, MV0811–15JC), Pt. Conception (ODP Site 1017E), Central California (F-8–90-G21, F2–92-P3, F2–92-P40, F2–92-P34), the Oregon-California border (ODP Site 1019E), and Vancouver Island (JT96–09) (Fig. 4, Table 3). These data indicate the extensive and shallow influence of OMZ waters on the CC continental margin during the deglaciation, and the regional-scale sensitivity of abrupt hydrographic change to Northern Hemisphere rapid warming and cooling events.

Intermediate waters (>1,000 mbsl) in the CC were oxygenated during the LGM, up until Termination IA at 14.7 ka [27, 29, 40, 88, 93, 94]. The three deepest sites exhibit evidence of hypoxia in the LGM, including ODP Site 1017E and ODP Site 1019E, and F-8–90-G21 (Fig. 7). The ODP cores have elevated [Mo] values in the LGM [42], and hypoxia-associated foraminiferan marker species are preserved in the deepest core [95-97].

Figure 7

California Current (CC) deglacial core data synthesized into hypoxia categories.

Changing deglacial core depths reflect global eustatic sea level change. The encircled number adjacent to each core label corresponds to the number of available oxygenation proxies, which are enumerate in Table 3. Vertical grey bars correlate to temporal intervals in OMZ geospatial reconstructions for this region.

Geochemical ([Mo], %Corg, and δ15N) data for all ODP Sites in the CC region, including Santa Barbara [42, 98–100], Pt. Conception [94], and the California/Oregon border [42, 101] collectively reconstruct the presence of denitrifying and severely hypoxic waters along the California margin at Termination 1A and through the Bølling-Allerød (Fig. 7). The Santa Barbara Basin sites record the OMZ intensification as laminations, severely hypoxic foraminiferan communities, and elevated redox metal concentrations [29, 40, 89, 93, 102, 103]. Archives from central California margin exhibit preserved laminations and elevated concentrations of [Mo] and [Cd] [40, 95, 96, 104, 105]. Core F2–92-P34 is anomalous to the broad CC pattern, with no enrichment of [Mo] at Termination IA and only a slight increase in %Corg [105]. At 14 ka, the records in this region indicate the seafloor was severely hypoxic from 395–869 mbsl and was bracketed above (332–395 mbsl) and below (869–1,564 mbsl) by intermediate hypoxia (Fig. 8).

Figure 8

California Current bathymetric seafloor masks and surface area (km2) histograms of deglacial hypoxia impingement for (a) 18 ka, (b) 14 ka, (c) 12 ka, (d) 10 ka, and (e) 4 ka.

Seafloor is selected between 0–2,400 mbsl and latitudinally constrained between 31° 40′-49° 30′ N. Analyses were limited to the continental margin within a 400 nautical mile buffer offshore of the continental coastline. The changing gray shoreline through the panels depicts the paleo-shoreline. The seafloor was oxic at 18 ka. At 14 ka the seafloor was severely hypoxic (395–869 mbsl), and bracketed by intermediate hypoxia above (332–395 mbsl) and below (869–1,564 mbsl). The seafloor returned to an oxic state at 12 ka. At 10 ka severe hypoxia (436–525 mbsl) was bracketed by intermediate hypoxia (373–436 mbsl, 525–1,605 mbsl). At 4 ka severe hypoxia (525–625 mbsl) was bracketed by intermediate hypoxia (417–525 mbsl,625–954 mbsl).

Intermediate waters returned to oxygenated conditions midway through the deglaciation (12.9–11.7 ka), synchronous with the ephemeral YD cooling observed in other Northern Hemisphere climate records (Fig. 8) [106, 107]. This margin-wide oxygenation is evidenced by the unanimous absence of laminations across all depths, foraminiferal assemblages that reflect increased oxygen concentrations, and low concentrations of redox metals in all regional cores. Termination IB (11.7 ka) initiated a return to regional hypoxia, detected in eight of the cores. Sites within the Santa Barbara Basin preserved laminations at Termination IB [29] and were dominated by hypoxic, OMZ-associated foraminiferal communities [93, 102, 103]. Termination IB was a secondary expansion of the CC OMZ, and is associated with a narrow band of severely hypoxic water at 436–525 mbsl, bracketed by intermediate hypoxia from 417–525 mbsl and 625–954 mbsl (Fig. 8).

Paleoceanographic reconstructions from the CC reveal the extraordinary shallowness (severe hypoxia <300 mbsl) and extensity (1,233 m of contiguously hypoxic water column) of the regional OMZ during the recent deglaciation. Oxygenated upper ocean ecosystems are dramatically compressed at both Termination events to <300 m from the ocean surface. The YD cooling mid-way through the deglaciation is a remarkably ephemeral event of regional oxygenation. The analysis presented here reveals how geographically pervasive and temporally responsive subsurface oxygen concentrations in the CC system are to global climate change.

Mexico Margin

The MM OMZ is remarkable in thickness and intensity (Fig. 2), and is a product of high surface production, a sharp pycnocline that inhibits local ventilation of subsurface waters, and the “sluggish and convoluted deep circulation” of regional subpycnocline waters [108]. The eastern Pacific sea surface temperature warm pool is centered off of southern Mexico and Guatemala, a product of the large seasonal net heat flux and weak wind mixing of the region [108, 109]. The CC and NPIW travel equatorward along Baja California and at 25° N turn westward with gyral circulation [110-112], maintaining an oceanographic transition zone between 20°-30° N. In the north, the OMZ sits between 500–1,000 mbsl [56] and in the south shoals upwards (100–700 mbsl), due to the influence of southern sourced intermediate water (including AAIW and “Equatorial Waters”) [113]. The oxygen concentration of these waters is low (~0.2 ml L-1), and the upper hypoxic boundary is extremely shallow. The term “Equatorial Waters” refers to a composite of Subtropical Underwater found below the equator [114, 115], north and south subtropical surface water [108], 13°C waters [116, 117] and subtropical mode water [118].

Sediment records in this region indicate major reductions in bottom water oxygenation through the last deglaciation, however the deglacial intensification of the OMZ does not appear to be synchronous across multiple locations [16]. Primary production and carbon export is thought to have been low during the LGM [119], potentially due to reduced upwelling [17, 120]. Changes in the depth of the equatorial thermocline and nutricline have been hypothesized to drive deglacial oscillations in oxygenation [121]. Unlike SP and CC records, oxygenation oscillations are not a comprehensive feature of the margin record, and only appear in records from the Sea of Cortez, which is isolated from the open margin and considered a more continental record of surface and atmospheric conditions [16, 26]. The equatorial Pacific is a major global denitrification zone. δ15N records from MM sites exhibit glacial-interglacial variability [122] and synchroneity to sites along the Peru-Chile margin [123-125].

Paleoxygenation Reconstructions for the Mexico Margin. Seven paleoceanographic records were selected for deglacial reconstructions, including two cores within the Sea of Cortez (GGC-55/JPC-56 and DSDP Site 480), three cores along the southwestern margin of Baja California Sur (MV99-PC14, MV99-PC08, MD02–2508 and GC31/PC08) and two cores along the Mexican Margin (NH8P, NH15P) (Fig. 4; Table 4). The MM exhibits reduced regional coherency, with asynchronous hypoxia and oxygenation oscillations (Fig. 9). However, this region exhibits a clear trend of higher oxygenation in the LGM, transitioning to severe hypoxia in intermediate waters at the start of the Holocene.

Figure 9

Mexico Margin (MM) deglacial core data synthesized into hypoxia categories.

Changing deglacial core depths reflect global eustatic sea level change. The encircled number adjacent to each core label corresponds to the number of available oxygenation proxies, which are enumerate in Table 3. Vertical grey bars correlate to temporal intervals in OMZ geospatial reconstructions for this region.

The two deepest, southernmost sites from the Mexican margin suggest the seafloor was hypoxic during the LGM, as indicated by the presence of laminations (Site NH8P) [126, 127] and high [Mo] concentrations (Site NH22P) [119]. At 18 ka, hypoxia was detected at only one site in intermediate waters (898 mbsl). Hypoxic waters spatially expanded along the Mexican margin and in the Sea of Cortez at 14 ka, and ranged from 334–932 mbsl. Site NH8P records the expansion of severe hypoxia along the margin, with laminations preserved between 15–7.5 ka [126, 127]. Site GGC-55/JPC-56 preserved laminations from 14.8–12.5 ka, concomitant with high δ15N values (>14‰), and high opal accumulations rates [128]. Lamination preservation occurred from 13–11.2 ka in DSDP Site 480, subsequently followed by bioturbation from 11.2–10.5 ka [26].

Five sites, including MV99-PC14, MD02–2508, DSDP 480, GC31/PC08, and GGC-55/JPC-56 (all from Sea of Cortez and Western Baja), recorded a synchronous northward expansion of hypoxic intermediate waters at ~11–10 ka (Fig. 10), with severe hypoxia extending from 375–973 mbsl. The regional shift was recorded in laminations at DSDP Site 480 [128]. At the southern tip of Baja (GC31/PC08 and MD02–2508), lamination preservation and steep increases in redox metals occurred [16, 17, 105, 121, 129]. Core MV99-PC14, extracted from Soledad Basin (which is oxygenated by sill depth waters at 290 mbsl) [16], became strongly laminated at 10 ka, indicating that extremely shallow upper intermediate waters deoxygenated at the start of the Holocene. The Holocene remained severely hypoxic at 4 ka from 419–817 mbsl.

Figure 10

Mexico Margin bathymetric seafloor masks and surface area (km2) histograms of deglacial hypoxia impingement for (a) 18 ka, (b) 14 ka, (c) 10 ka, and (d) 4 ka.

Seafloor is selected between 0–1,200 mbsl and latitudinally constrained between 20°-30° N. Analyses were limited to the continental margin within a 400 nautical mile buffer offshore of the continental coastline. The changing gray shoreline through the panels depicts the paleo-shoreline. At 18 ka, severe hypoxia was limited from 937–1,037 mbsl. Severe hypoxia ranged from 334–932 mbsl at 14 ka. At 10 ka, severe hypoxia was found from 375–973 mbsl, and at 4 ka severe hypoxia was contracted slightly to 419–817 mbsl.

Paleoceanographic reconstructions for MM provide evidence for the regional intensification of subsurface hypoxia through the recent deglaciation. The MM exhibited reduced sensitivity (Fig. 10) to the first Northern Hemisphere glacial termination event (14.7 ka), as compared to the NP (Fig. 6) and CC sites (Fig. 8). A regional-scale expansion of hypoxia was recorded at ~10–11 ka, wherein ~600 m of the water column became severely hypoxic (Figs. 9 and 10).

Eastern Equatorial Pacific and Humboldt Current

The HC, also known as the Peru Current, is the eastern limb of the South Pacific subtropical gyre, characterized by upwelling and extreme biological productivity [130]. A thick (~500 m), intense ([O2]<0.2 ml L-1) and shallow (upper boundary ~50 mbsl) OMZ characterizes the HC (Fig. 2)[131]. Cold HC surface waters move north along the South American continental margin, and then become a part of the equatorial cold tongue. The southward-flowing Peru-Chile Undercurrent is found 75–500 mbsl, in association with the OMZ [132]. This equatorially-derived water is transported as far as 48° S [133], and is derived from Equatorial 13°C Waters [116, 117], Subtropical Underwater [115], and Eastern South Pacific Intermediate Water [134]. The spatial distribution of OMZ thickness is correlated with upwelling conditions [135]. An intense upwelling system is located off of Peru [136, 137], which brings nutrient-rich AAIW to the surface and stimulates high production in the equatorial cold tongue [138-140]. A more seasonal upwelling cell exists southward off of Chile [141]. A functional break in deep water properties exists at 15–25° S, where tropical and subtropical deep waters are dynamically separated, as deep waters are directly connected in the Western Pacific [142]. Abyssal waters (>4,000 m) are a mixture of water from the Weddell Sea and the North Atlantic, and are relatively high in oxygen concentrations [143].

The deglaciation is thought to have co-occurred with extensive ventilation of the deep-sea, renewing oxygen concentrations in the deep ocean interior (e.g., [2, 63, 73, 144]). Recent syntheses document global increases in deep ocean [O2] through the deglaciation, reflecting the transfer of respired carbon from the deep ocean to the atmospheric and surface ocean carbon pool [7]. Numerous high quality records of deep-sea oxygenation are present in the Eastern Equatorial Pacific, and are included in the HC analysis. Intermediate water records from the HC are limited; however, it is clear that intermediate water (~300–1,400 mbsl) deoxygenated at these depths prior to the Northern Hemisphere glacial Terminations [124, 145].

Paleoxygenation reconstructions for the Equatorial Pacific and Humboldt Current. Thirteen paleoceanographic cores met the criteria for deglacial geospatial reconstructions, including cores from the Panamanian margin (ODP Site 1242), the Galapagos (TR1–163–19P, RC13–140, Y69–71P, ME0005A-24JC, RC11–238), Ecuador (P7, ME0005A-27JC), Peru (V19–30, TR163–31P, W7706–41K, W7706–40K, W7706–37K), and the Chilean margin (GeoB7139–2) (Fig. 4; Table 4). Sites exhibit regional synchroneity from deep to upper intermediate water depths, lack coherent timing with Northern Hemisphere climate records (Fig. 11), and exhibit a deglacial temporal signature similar to the Southern Hemisphere (e.g., [146, 147]). Broad deglacial patterns include increasing deep ocean oxygenation, severe and shallow hypoxia shoaling at ~17 ka in upper intermediate water, and a ~6 kya temporal overlap between deep and intermediate hypoxia (Fig. 12).

Figure 11

Equatorial Pacific and Humboldt Current (HC) deglacial core data synthesized into hypoxia categories.

Changing deglacial core depths reflect global eustatic sea level change. The encircled number adjacent to each core label corresponds to the number of available oxygenation proxies, which are enumerate in Table 3. Vertical grey bars correlate to temporal intervals in OMZ geospatial reconstructions for this region.

Figure 12

Equatorial Pacific and Humboldt Current bathymetric seafloor masks and surface area (km2) histograms of deglacial hypoxia impingement for (a) 18 ka, (b) 13 ka, (c) 10 ka, and (d) 4 ka.

Seafloor is selected between 0–3,300 mbsl and latitudinally constrained between 10° 30′ N-32° S. Analyses were limited to the continental margin within a 400 nautical mile buffer offshore of the continental coastline and the Galapagos Islands. The changing gray shoreline through the panels depicts the paleo-shoreline. At 18 ka, severe hypoxia was limited from 937–1,037 mbsl. At 18 ka, intermediate hypoxia was found in deep water (2,082–3,088 mbsl). At 13 ka, severe hypoxia was found between 108–331 mbsl, and intermediate hypoxia was deeper (332–3,130 mbsl). At 10 ka, severe hypoxia was again found in shallow waters (315–415 mbsl) and intermediate hypoxia was deeper (415–3,164 mbsl). At 4 ka, intermediate hypoxia retracted to between 1,313–1,413 mbsl.

For the eight equatorial deep waters sites, from 2,203–3,091 mbsl, the LGM is associated with high concentrations of [U] [144, 148, 149]. These records reconstruct the presence of intermediate hypoxia at 18 ka from 2,082–3,088 mbsl. The deepest, most southerly record (GeoB7139–2) stands out as anomalous to the deep equatorial cores, wherein δ15N values rapidly increased (~5–6‰) at ~17 ka, suggesting a shift towards denitrification and an associated decrease in oxygen concentration [123, 150]. Deep-sea hypoxia continues to be a regional feature through the deglaciation and into the mid-Holocene, although the timing of post-glacial deep ocean oxygenation is debated. Concentrations of [U] remain high into the mid-Holocene in the equatorial cores [144], beyond 14.7 ka, which differs from the timing of oxygenation in the North Pacific [73]. This temporal difference may be an artifact of increases in productivity and organic matter flux, or it may be a true signal of the temporal differences in oxygenation across the Pacific. For our purposes here, we follow existing interpretations of high [U] in the post-glacial Equatorial Pacific as a signal of intermediate hypoxia (Fig. 11).

Deoxygenation in intermediate waters occurred at ~17 ka, as exhibited by laminations, a ~5‰ δ15N increase along the Chilean Margin (W7706–41K, W7706–40K, W7706–37K) [145, 151], and high δ15N values adjacent to the Panamanian Margin (ODP Site 1242; Fig. 11) [124]. Sedimentation is discontinuous through the deglaciation in W7706–37K and may be due to sediment disturbance rather than a signal of oxygenation reversal [145]. Together, these cores indicate severe hypoxia at 13 ka ranged from 108–331 mbsl and intermediate hypoxia ranged from 332–3,130 mbsl. At 10 ka, severe hypoxia is recorded in upper intermediate waters (365±50 mbsl), and is bordered below by intermediate hypoxia from 415–3,164 mbsl. At 4 ka, sediment at the deep equatorial sites shows dramatic reductions in [U] concentrations. Concurrently, the shallow central Chilean sites were oxygenated, and the remaining location with a clear hypoxia signal is limited to Panamanian ODP Site 1242. As such, hypoxia at 4 ka is substantially attenuated to 1,363 mbsl.

The HC provides a unique record of synoptic changes occurring in both deep and intermediate water through the events of the deglaciation. The HC exhibits a coherent regional signal of oxygenation, while lacking synchroneity to the Northern Hemisphere. Striking features of the record include the regional deoxygenation of upper intermediate waters at 17 ka, the extreme vertical expansion of hypoxic water at 13 ka, and the ~6 kyr overlap of deep and intermediate hypoxia creating ~3,000 m of contiguous hypoxic water column.

Benguela Current

The BC system, also referred to as the Angola-Benguela Current, is the equatorward flowing Eastern Boundary Current of the South Atlantic subtropical Gyre (Fig. 2), associated with high productivity, organic-rich sediments, and seasonal upwelling [152]. Upwelling events, and the export of surface productivity, are linked directly to modern subsurface OMZ structure [153, 154]. The HC OMZ is significantly shallower and less geographically extensive than other systems reviewed here (Fig. 2). The core of the BC OMZ is between 300–400 mbsl, and hypoxic waters extend to as shallow as 50 m [155].

Sediment records from the BC reveal mixed deglacial productivity signals, wherein a few single sites indicate decreased productivity during the LGM [156], while more recent and contradictory work indicates a decrease in surface productivity in the Holocene as compared to the LGM [157-160]. Contributing to this mixed signal, it appears that the productivity center may have moved offshore since the LGM [160]. From the evidence available, it does not appear that the OMZ associated with the BC follows a glacial/interglacial cycle like that which dominates the Pacific, but is both more heterogeneous and directly linked to regional cycles of productivity and upwelling [160].

Oman and Pakistan Margin

The OPM, within the Indian Ocean, has the most globally extensive OMZ in terms of water column vertical extent (>1000 m) [161] (Fig. 2). Subsurface dissolved oxygen content of the region is controlled by nutrient additions, changes in intermediate and deep water ventilation, and convective mixing of deep waters [162]. Oxygen depleted waters sourced from small inlet seas contribute to the maintenance of severe hypoxia in the region [163]. Monsoonal seasonality of the Arabian Sea drives OMZ oscillations, wherein strong onshore winds force the upwelling of nutrient-rich waters [164, 165]. OPM sediments are more organic rich than other low-latitude OMZ regions, attributed to high primary productivity and reduced remineralization at depth [166].

The OPM is a sensitive recorder of climate variability and an archive of the connectivity between low latitude monsoonal dynamics and high latitude temperature variability [167, 168]. Paleoceanographic investigations show that the strength of the OPM monsoon has high-resolution, sub-millennial cycles that exhibit synchrony to Greenland climate variability [168]. Regional surface water productivity has co-varied with OMZ intensity, such that OMZ intensity is especially weak during cold stadial events when summer monsoon effects are reduced [167, 169–172], while OMZ intensity is enhanced during warm interstadials (e.g., the B/A) [171]. Deglacial cores demonstrate rapid (102 year timescales) shifts from intermediate to severe hypoxia due to climate forcing [173]. Weakening of the OMZ occurred during Heinrich 1 (~18–14 ka) and the YD [174, 175]. Monsoon intensity peaked with a corresponding increase in surface productivity [172, 176] and ocean hypoxia from 9.5 to 5.5 ka [177, 178]. Geospatial reconstructions of deglacial hypoxia were not conducted for the OPM, as the limited high-resolution deglacial records were not sufficient to warrant further analyses.

Discussion

Mechanisms and implications of OMZ variation

Large-scale physical and biogeochemical processes in the ocean drive OMZ formation, and have high-latitude and low-latitude derived features. As such, the deglacial changes in dissolved oxygen may be explained by proximal and distal mechanisms. High latitude sources for deglacial OMZ variability have been hypothesized to include changes in intermediate water production/ventilation [17, 27, 40, 68, 71, 84, 93, 150, 179] and changes in deep ocean circulation [70, 123, 170]. Deglacial micronutrient [Fe] enrichment has been hypothesized to intensify hypoxia, particularly in the Subarctic Pacific [58, 67, 151]. Additionally, oxygen consumption at the site of NPIW subduction has been hypothesized as a high latitude mechanism [68]. Low latitude mechanisms for deglacial OMZ variability include the relative importance of equatorial counter currents [87, 122], the intensity and location of upwelling systems [127, 180] changes in the oceanic preformed nutrient inventory [73, 124, 126, 128, 181], surface ocean productivity [119, 121], atmospheric structural changes [17, 85, 182], and the strength of monsoonal systems [167, 168].

In some cases, independence between mechanisms that force hypoxia has been demonstrated, for example where local productivity and the high-latitude ventilation of intermediate waters are shown to decouple during the deglaciation [90, 183]. However, many of the physical and biological processes that drive the development of subsurface hypoxia are ultimately linked. Indeed, circulation models from the latter half of the 20th century have revealed that changes in the rates of surface productivity can ultimately be driven by physical perturbations in circulation, through interactions between thermocline depth, nutrient flux, and particulate export into hypoxic waters [184].

Major changes in the distribution of [O2] in Eastern Boundary Currents occurred during the recent deglaciation. The substantial oxygen changes arose concomitantly with increasing atmospheric carbon concentration, surface warming and rising sea levels [2-6]. The analyses presented here reveal that rapid oscillations in oxygen distribution are an inherent feature of large ocean regions, as subsurface [O2] exhibits extreme sensitivity to hemisphere-specific warming and cooling. Mid-way through the deglaciation, hypoxic waters shoaled vertically, shifting the upper OMZ boundary towards the ocean surface, with shoaled depths ranging from 596 mbsl in the SP (Fig. 6), 332 mbsl in the CC (Fig. 8), 334 in the MM (Fig. 10), to 108 mbsl in the HC (Fig. 12). Subsurface hypoxia expanded and intensified, creating OMZ environments which ranged in thickness from 2,298 m in the SP, 1,233 m in the CC, 598 m in the MM, to 3,022 m in the HC (Table 4). These data provide evidence of the capacity for OMZs to exhibit extreme shallowness and water column extensity to states that have no analogue in the modern ocean.

The deglaciation offers an extremely informative case study of the sensitivity and coupling between OMZs and global climate oscillations. Our analyses show that during global-scale warming events, vast expanses of the upper global ocean deoxygenate, resulting in the vertical compression of oxygenated ocean ecosystems. These analyses of past changes in dissolved oxygen are fundamentally relevant to the loss of [O2] observed in the modern ocean [185-187], and provide bounds and reasonable expectations for expanding OMZs in the future.

Changes in the distribution of oxygen translate directly to the structure of marine ecosystems (e.g., [11, 12]) and OMZ expansions have implications for modern oceanography, biodiversity conservation, ocean management and sustainable fisheries. Modern oceanographers can anticipate that Eastern Boundary Current OMZs have the capacity to expand to hydrographic structures that have never been instrumentally observed. These analyses underscore the continued need for high-quality instrumental hypoxic time series data and predictive models of modern oceanographic change. The challenge to resource managers and conservationists is how to plan for the great deal of uncertainty introduced when ecosystems undergo changes that only have precedent in the geologic record.

Modern subsurface oxygen variability

In the modern ocean, the vertical expansion and intensification of OMZ regions has been detected in the California Current [185, 188], equatorial waters [186], the Subarctic [189, 190] and Subtropical Pacific [191], the North Atlantic [192], the Indian Ocean [193] and the Southern Ocean [194]. The loss of subsurface dissolved oxygen is an acute perturbation to coastal ecosystems and both benthic and pelagic communities [187]. These data indicate that deoxygenation in the 20th and 21st Centuries is a feature of every global ocean basin. However, unique physical and biological processes within oceanographic provinces provide a more nuanced and complicated view of climate-forced ocean deoxygenation.

OMZs exhibit high-frequency variability on diurnal and semi-diurnal timescales [195], on the intra-annual timescales of upwelling and relaxation cycles [196], and in tight association in the California Current with La Niña events [197]. These OMZ oscillations illustrate the potential for undescribed scales of variability in the coastal ocean. Additionally, as longer instrumental time series are developed, novel questions regarding intrinsic high-frequency oceanographic variability arise. For example, conflicting interpretations of California Current oxygen trends exist, such that the data can be analyzed to reveal a long-term deoxygenation trend [185] or a 20–25 year undescribed oscillation [188]. This interpretive conflict highlights the importance of understanding high frequency variability of oxygen concentrations, which may be overprinted by anthropogenic climate forcing. Complicating the picture even more so, recent analyses of historical sediments from the Eastern Equatorial Pacific indicate that the regional OMZ has contracted, and that this contraction is correlated with the slackening of trade winds in the tropical Pacific [198]. These regional trends in dissolved oxygen illustrate the potential for OMZ oceanographic provinces to respond differently to anthropogenic climate forcing, however 20th century oxygenation in the Equatorial Pacific is thought to be anomalous in the context of a broader global picture. On centennial time scales, processes associated with a warm surface ocean, including gas solubility reduction and physical stratification, are predicted to substantially reduce dissolved oxygen concentrations in the ocean interior [8, 184].

Conclusions

We integrate existing deglacial geochemical, sedimentary, and microfossil oxygenation proxies to reconstruct the timing, depth and intensity of seafloor hypoxia in Eastern Boundary Currents, principally in the Eastern Pacific. These analyses illustrate the high degree of coupling between the global climate system and OMZ environments and provide the most comprehensive window, to date, into the spatial capacity of OMZ ecosystems to expand and contract due to climate change. The recent deglaciation was accompanied by the dramatic shoaling of the upper hypoxic boundary toward the ocean surface, the compression of upper ocean oxygenated habitat, and the expansion of the subsurface hypoxic water column. Subarctic Pacific and California Current continental margins exhibit tight correlation to the oscillations of Northern Hemisphere deglacial events, whereas the Mexico Margin and the Equatorial Pacific and Humboldt Current exhibit hypoxia expansion prior to Termination IA (14.7 ka), and no regional oxygenation oscillations. Oxygenation changes occurred in synchrony across ocean basins, revealing the extensive sensitivity of upper ocean systems to changes in global climate.