New Research In
Physical Sciences
Social Sciences
Featured Portals
Articles by Topic
Biological Sciences
Featured Portals
Articles by Topic
- Agricultural Sciences
- Anthropology
- Applied Biological Sciences
- Biochemistry
- Biophysics and Computational Biology
- Cell Biology
- Developmental Biology
- Ecology
- Environmental Sciences
- Evolution
- Genetics
- Immunology and Inflammation
- Medical Sciences
- Microbiology
- Neuroscience
- Pharmacology
- Physiology
- Plant Biology
- Population Biology
- Psychological and Cognitive Sciences
- Sustainability Science
- Systems Biology
Deep ocean nutrients imply large latitudinal variation in particle transfer efficiency
Edited by David M. Karl, University of Hawaii, Honolulu, HI, and approved June 14, 2016 (received for review March 16, 2016)

Significance
Organic particles that sink to the deep ocean release carbon dioxide in waters that remain out of contact with the atmosphere on long timescales. This paper reconstructs particle flux profiles from large-scale ocean nutrient distributions, revealing systematic regional variations that have proved difficult to detect from direct observations alone. We show that the “transfer efficiency” of organic matter to depth is greatest in high-latitude regions and lowest in the subtropics, and is well explained by variations in phytoplankton community structure. These results suggest a greater sensitivity of atmospheric CO2 to high-latitude carbon export, and predict reduced ocean carbon storage as subtropical communities expand in response to climate warming.
Abstract
The “transfer efficiency” of sinking organic particles through the mesopelagic zone and into the deep ocean is a critical determinant of the atmosphere−ocean partition of carbon dioxide (CO2). Our ability to detect large-scale spatial variations in transfer efficiency is limited by the scarcity and uncertainties of particle flux data. Here we reconstruct deep ocean particle fluxes by diagnosing the rate of nutrient accumulation along transport pathways in a data-constrained ocean circulation model. Combined with estimates of organic matter export from the surface, these diagnosed fluxes reveal a global pattern of transfer efficiency to 1,000 m that is high (∼25%) at high latitudes and low (∼5%) in subtropical gyres, with intermediate values in the tropics. This pattern is well correlated with spatial variations in phytoplankton community structure and the export of ballast minerals, which control the size and density of sinking particles. These findings accentuate the importance of high-latitude oceans in sequestering carbon over long timescales, and highlight potential impacts on remineralization depth as phytoplankton communities respond to a warming climate.
Sinking organic particles deliver carbon from the surface euphotic zone (upper ∼100 m) of the ocean into deeper layers that do not exchange with the atmosphere (1). The longevity of oceanic carbon storage by this “biological pump” depends on the depth at which particulate organic matter (POM) decays and releases CO2 into seawater (2). Most POM is consumed within the mesopelagic zone (100−1,000 m) of the water column, which recirculates rapidly to the surface, leaving only a small fraction to remineralize in the deep ocean where carbon can be sequestered on centennial and longer timescales (3). The “transfer efficiency” of particulate carbon from the euphotic zone to depth is therefore a critical determinant of atmospheric pCO2 (4), but its underlying controls are poorly understood and crudely represented in Earth system models used to project global carbon cycling and climate.
The depth scale over which POM fluxes attenuate hinges on both the sinking speed of particles and their rate of decomposition (5), each governed by a range of factors. Decomposition rates are thought to depend on the abundance of heterotrophic microbes (6) and the temperature sensitivity of their metabolism (7, 8), as well as the palatability of the organic matter itself (9, 10). Particle sinking speeds depend on their size and density, which may be ultimately dictated by the plankton community structure and trophic web of the euphotic zone where they are produced (11, 12). Because these factors exhibit distinct regional variations, their relative importance might be discerned by detecting large-scale patterns of transfer efficiency in the ocean.
Arrays of neutrally buoyant sediment traps deployed at multiple depths provide the most direct estimate of particle fluxes from the euphotic zone and their attenuation in the upper ocean at local scales (13). This approach has recently revealed a pattern of rapid flux attenuation to ∼500 m in the subtropics, and slower attenuation in subarctic regions, that is consistent between the Atlantic and Pacific Oceans (8, 14). Transfer efficiency to greater depth has been estimated by combining particle fluxes captured in bottom-moored sediment traps and thorium (Th)-based measures of surface export production, yielding an opposite pattern with highest values in the low latitudes (15). Direct observations are still too sparse to constrain particle fluxes through the mesopelagic zone on the scales required to discriminate between these patterns, especially given the well-known limitations of most sediment trap systems (13).
Diagnostic Method
We reconstructed particle fluxes to the deep ocean by exploiting the large-scale tracer signatures of POM decomposition, which consumes oxygen (O2) and converts organic nitrogen (N) and phosphorous (P) to inorganic nutrients. Remineralization can therefore be identified by the accumulation of apparent oxygen utilization (AOU), nitrate (NO3−), and phosphate (PO43−) along subsurface circulation pathways (16). Each tracer integrates the remineralization signal over sufficient timescales to determine time mean rates, and has a global distribution that is characterized by orders of magnitude more observations than the particle fluxes themselves (17).
Two factors confound this tracer-based approach. First, in poorly oxygenated waters, microbial respiration is fueled by NO3− rather than O2, so remineralization leaves no AOU signature and acts as a sink, rather than source, of NO3− (18). Second, remineralization tracers reflect the decomposition of both POM and dissolved organic matter (DOM), so isolating the particulate signal requires an independent estimate of DOM decomposition. Global DOM distributions are not well sampled, but multiple observations indicate that dissolved organic phosphorous (DOP) cycles rapidly to PO4− and does not accumulate in the subsurface like dissolved organic carbon and nitrogen (19⇓–21). For these reasons, we focus on PO43− as the most promising tracer of particle remineralization, and restrict our analysis to the ocean interior below 300 m, where DOP concentrations are negligible throughout the stratified low latitudes (19). In high-latitude regions, we further excluded seasonal mixed layers from our analysis, which entrain DOP from the surface (20) and where the nutrient concentrations vary considerably over time (17).
We estimated particulate organic P (POP) remineralization from the global PO43− distribution using a data-constrained ocean general circulation model (OGCM) that accurately captures the isopycnal structure of the ocean, and the mixing and aging of water masses (22). The observed PO43− distribution in the ocean interior was first divided into two components,
Pcirc represents the PO43− carried by water masses as they pass through the 300-m horizon or the base of the wintertime mixed layer, whichever is deeper. Pcirc was determined using the model circulation to propagate observed PO43− concentrations at the boundary into the interior region without alteration, except by mixing. The difference between this transported value and the observed concentration at each location, Prem(POP), represents the PO43− produced within water masses by POP remineralization in the ocean interior (Fig. 1A). The remineralization rate, Jrem, was then diagnosed from the accumulation rate of Prem(POP) over time, as water masses spread and age along isopycnals away from the boundary in the circulation model (SI Materials and Methods). Finally, profiles of the vertical POP flux (FPOP) were reconstructed by integrating Jrem from each depth level to the seafloor, assuming that benthic POP burial is negligible. We used a Monte Carlo approach to quantify the uncertainty in these fluxes, repeating our calculation 10,000 times using eight different versions of the tracer-constrained OGCM and propagating uncertainty in the observed PO43− distribution.
Tracer method to reconstruct large-scale particle fluxes. (A) Zonal mean Prem(>500m) in the Pacific, identified by subtracting PO43− that is circulated through the 500-m boundary from the observed distribution below 500 m. The accumulation rate of this tracer over time (white contours are years since passing 500 m) as water masses spread along isopycnals (black contours) is used to diagnose the particle remineralization rate in an OGCM. (B) Averaging regions for POP fluxes, which capture differences in factors that control particle remineralization and sinking, including mean upper-ocean temperate (colormap), and phytoplankton community composition. Diagonal hatching shows communities comprising >50% picoplankton, and cross-hatching indicates regions with ballast mineral export ratios exceeding 5 grams per gram of POC (SI Materials and Methods).
To distinguish broad spatial patterns from gridpoint noise, diagnosed POP fluxes were averaged over regions that capture large-scale variations in those factors hypothesized to influence Teff (Fig. 1B). The Atlantic and Pacific Oceans were divided into warm subtropics [subtropical Atlantic (STA) and subtropical Pacific (STP)] with picoplankton-dominated communities and rapid nutrient recycling; cold northern subarctic regions [North Atlantic (NA) and North Pacific (NP)] with seasonal blooms of large microplankton (23); and productive eastern tropical upwelling zones [eastern tropical Atlantic (ETA) and eastern tropical Pacific (ETP)]. The Southern Ocean was divided between the productive subantarctic zone (SAZ), with deep winter mixed layers, and the upwelling Antarctic zone (AAZ), with short growing seasons dominated by heavily silicified diatom blooms (24).
Results
Diagnosed profiles of FPOP (Fig. 2) generally follow the same asymptotic shapes that are characteristic of observed particle fluxes (3), attenuating quickly in the upper ocean and more gradually at depth. The diagnosed fluxes are broadly consistent with the range of rates measured by bottom-moored sediment traps in the deep ocean (Fig. S1) (25), but also exhibit systematic variations between regions that are not apparent in the sparse observations. In both the Atlantic and Pacific Oceans, the subtropics stand out as regions of particularly low POP flux, which decreases rapidly from 4−5 mmol⋅m−2⋅y−1 at 300 m to values indistinguishable from zero at the base of the mesopelagic zone, defined here as the 1,000-m depth horizon (Fig. 2). In all other regions, FPOP is significantly greater than zero down to at least 2,500 m [90% confidence interval (C.I.)]. The largest particle fluxes to the deep ocean are found in the subarctic Atlantic and Pacific, where FPOP at 1,000 m is 6.4 ± 1.16 mmol⋅m−2⋅y−1 and 7.3 ± 0.63 mmol⋅m−2⋅y−1, respectively, and slightly lower fluxes (4.4–5.6 mmol⋅m−2⋅y−1) are found in the tropics and Southern Ocean.
Reconstructed profiles of POP flux (FPOP). In each panel, the central red line is the regional mean flux profile for the listed ocean province (labels match Fig. 1B). Dark-red and light-red envelopes represent the 1σ interval and 90% confidence interval, respectively, of the regional mean fluxes, quantified using a Monte Carlo method that samples uncertainty in the PO43− distribution and model circulation (SI Materials and Methods).
Comparison of diagnosed and observed particle fluxes. Central red lines are the regional mean diagnosed flux profiles, as in Fig. 2. Red envelopes represent the regional variability of diagnosed fluxes within regions (1σ), as well as uncertainties associated with the PO43− distribution and model circulation. Blue dots (with error bars) are the mean and SD of particle flux observations in each region, binned to the model’s vertical grid. The model–data comparison is discussed in SI Materials and Methods, Comparison with Observed Particle Fluxes.
Regional patterns in particle flux attenuation between the mesopelagic zone and the deep ocean can be identified from the ratio between FPOP at 1,500 m and 500 m in each profile (Fig. 3A). These depth horizons straddle the putative base of the mesopelagic zone, and 500 m is the shallowest depth at which FPOP estimates are available for all regions (Fig. 2). The fastest attenuation is found in the subtropics, where only ∼10% of the 500 m fluxes reaches 1,500 m, followed by the two tropical zones, where FPOP at 1,500 m is reduced to ∼20% of its 500-m value. Slower attenuation is observed in the high latitudes: Both subarctic regions retain ∼35% of the 500-m flux down to 1,500 m, and the Antarctic region retains >40% (Fig. 3A).
Particle transfer efficiency over depth. (A) Attenuation of reconstructed FPOP profiles, quantified as the ratio of 1,500-m and 500-m fluxes. (B) Transfer efficiency of POP from zeu to 1,000 m, combining diagnosed FPOP at 1,000 m (Fig. 2) with satellite and CMIP5 estimates of regional mean export at zeu (Fig. S2A). Large dots use algorithms and models that are in good agreement with tracer-based export estimates within the given region (SI Materials and Methods). Vertical bars propagate uncertainty of FPOP into Teff. (C) Transfer efficiency estimated using composite export estimates (Fig. S2C) that weight satellite and CMIP5 ensemble members according to their skill at reproducing observed rates. Error bars propagate the uncertainty in both diagnosed FPOP (Fig. 2) and export (Fig. S2). Red and blue dots are the Teff estimates derived from neutrally buoyant sediment trap observations (8, 14). Thick dashed line in B and C is the prediction of the Martin curve relationship: FPOP(z) = FPOP(zeu) (z/zeu)−0.86, given z = 1,000 m and zeu = 100 m.
Because our diagnostic method cannot be extended shallower in the water column, the reconstructed flux profiles alone cannot provide estimates of transfer efficiency (Teff) from the base of the euphotic zone (zeu) to the deep ocean (1,000 m). Instead, we require an independent estimate of POM export at zeu, which can only be provided by remote sensing approaches and global biogeochemical/ecosystem models on the large spatial scales considered here. To account for the full range of uncertainty in these methods, we compiled 14 different estimates of POP export, nine from satellite algorithms and five from state-of-the-art Coupled Model Intercomparison Project 5 (CMIP5) model predictions (Fig. S2A and Table S1). Satellite estimates of particulate organic carbon (POC) export were generated using all combinations of three net primary production (NPP) climatologies and three algorithms for the export efficiency (e ratio) of organic matter (SI Materials and Methods), then converted to POP using a spatially varying relationship for particulate P:C ratios (26).
Export estimates and data constraints. (A) Regionally averaged P export from nine satellite algorithms and five CMIP5 models (SI Materials and Methods, Export Estimates). These export values were used for the Teff estimates in Fig. 3B. (B) Annual mean C export predicted by satellites and models at locations within each region (SI Materials and Methods, Export Estimates), compared with the tracer-based observational constraints compiled by refs. 27 and 28. (large black dots). This comparison was used to derive the weighting factors listed in Table S2. (C) Composite P export estimates for satellite and CMIP5 ensembles, derived by taking weighting averages and SDs of the individual ensemble members (in A) using the weighting factors listed in Table S2. These composite values were used for the Teff estimates in Fig. 3C.
CMIP5 models used in this study
Export estimates were averaged over regions and combined with the reconstructed POP flux to determine Teff from zeu→1,000 m (Fig. 3B). One robust pattern to emerge from this analysis is that Teff is consistently low in the subtropics relative to the other regions, and lower than expected from the canonical “Martin curve” relationship often used to predict flux attenuation over depth (3). In general, Teff estimates cluster close to the Martin curve prediction of ∼15% in the tropics, but span a wide range between 10–50% in the high latitudes, where export estimates vary considerably (Fig. S2A). Robust differences between these regions therefore cannot be established without additional constraints on export.
We evaluated each satellite algorithm and CMIP5 model against tracer-based observations of export production at locations within seven of the eight ocean regions considered here (27, 28) (Fig. S2B). Composite export estimates for each region were then derived by weighting members of the satellite and CMIP5 ensembles according to their skill at reproducing the empirical estimates at these locations (Table S2). This approach reduces the uncertainty in our estimated POP export at zeu (Fig. S2C), and better clarifies the systematic regional differences in Teff (Fig. 3C). In addition to confirming low subtropical Teff (∼5%), this result distinguishes tropical upwelling zones with intermediate Teff (∼15%) from high-latitude oceans with high values (∼25%) that exceed the Martin curve prediction (Fig. 3C). The same pattern is observed in transfer efficiency from zeu→2,000 m (Fig. S3) and is consistent with the attenuation estimates derived from diagnosed FPOP profiles alone (Fig. 3A). Taken together, these results establish a coherent geographical pattern of particle flux attenuation through the mesopelagic zone and beyond: The subtropics deliver the smallest fraction of the surface export flux to depth, and high latitudes deliver the largest, with the tropics falling in between.
Regional weighting factors for satellite and model export
Transfer efficiency to 2,000 m. As in Fig. 3C, but using diagnosed fluxes at 2,000 m rather than 1,000 m. Red and blue dots are the predictions of power-law fits to the neutrally buoyant sediment trap data of refs. 8 and 14, respectively. Transfer efficiency to 2,000 m shows the same ordering (high latitude > tropics > subtropics) as 1,000 m (Fig. 3C), demonstrating a consistent pattern of flux attenuation over depth.
Our findings are at odds with previous attempts to constrain deep ocean transfer efficiency (15), which estimated that >20% of the surface export flux reaches 2,000 m throughout the low latitudes, compared with <5% in the Southern Ocean. The large-scale distribution of remineralized PO43− does not support this pattern, which might be compromised by the collection biases of bottom-moored sediment traps (13) (Fig. S1) and/or analytical errors associated with the Th-based export observations used (29). On the other hand, our results are highly consistent with the few available POC flux profiles to 500 m compiled from neutrally buoyant sediment traps (8, 14) (Fig. 3C), and extend the observed regional difference between subtropics and high latitudes over a much larger depth range. Furthermore, the close quantitative agreement between these approaches (Fig. 3C) indicates that large-scale patterns detected in our reconstructed POP profiles also apply to organic carbon fluxes.
Discussion
To investigate the mechanisms that control remineralization depth, we compared our tracer-constrained Teff estimates (Fig. 3C) to large-scale patterns in the factors hypothesized to influence particle sinking and decomposition rates (Table 1 and SI Materials and Methods). The upper-ocean temperature (T100–500m) exhibits the opposite trend from Teff (subtropics > tropics > high latitude, Fig. 1B), consistent with the effect of temperature on microbial respiration rates. This relationship explains 78% or 67% of the regional variations in Teff (using satellite and CMIP5 export, respectively, Table 1), but underpredicts the difference between the tropics and subtropics while overpredicting the difference between subarctic regions and the AAZ (Fig. S4). Preferential decomposition of “fresh” POM, exported rapidly during blooms and not subjected to intensive recycling in the surface, has also been hypothesized as a dominant control on remineralization depth (10). The seasonal ranges in NPP and export ratio of organic matter are therefore expected to correlate negatively with Teff (9, 30) but we find the opposite trend, suggesting that POM lability does not drive the large-scale pattern diagnosed here.
Statistical predictors of Teff
Statistical analysis of Teff. As in Fig. 4 A and B but correlating Teff against (A and B) upper-ocean temperature, averaged between 100 m and 500 m; (C and D) seasonality of primary production; (E and F) carbon export:NPP ratio; and (G and H) mass ratio of ballast minerals (calcite+silicate) to POC export. The region symbols are the same as in Fig. 4, and the derivation of regional values for each variable is discussed in SI Materials and Methods, Statistical Analysis of Teff. A, C, E, and G use Teff derived from satellite C export; B, D, F, and H use Teff derived from CMIP5 C export.
Of all of the variables we considered, the size structure of the surface phytoplankton community was the best predictor of particle transfer efficiency. The fraction of photosynthetic biomass contributed by picoplankton (Fpico), which varies between >55% in the subtropics and <30% in subarctic regions (Fig. 1B and ref 23), explains up to 93% of the variance in Teff between regions (Table 1 and Fig. 4 A and B). Mechanistically, this relationship might arise from the size spectrum of exported particles, which skews toward small, slowly sinking particles when picoplankton dominate the community (12, 31). Particle models with size-dependent sinking rates (32) predict that the observed difference in size spectra between regions with high and low Fpico (31) can fully account for the range of Teff that we found. It is therefore possible that the strong relationship between T100–500m and Teff identified here (Table 1) and by others (8) emerges only from the spatial covariation of upper-ocean temperature and phytoplankton community structure, or that both factors combine to explain the diagnosed pattern.
Statistical model of transfer efficiency. (A and B) Scatter plots of Fpico vs. Teff estimated using (A) satellite and (B) CMIP5 export, with best-fit regression line (thick black line) and 90% confidence intervals (dotted lines). The Antarctic zone of the Southern Ocean is an outlier to the strong correlation in both cases (lies outside 90% C.I.). (C) A multiple linear regression (MLR) model that combines Fpico with ballast (Eq. 2) significantly improves the prediction of Teff (P < 0.05), and no further significant improvement can be gained by adding additional predictors.
The Antarctic zone of the Southern Ocean is an outlier (90% C.I.) to the relationship between Fpico and both sets of Teff estimates (Fig. 4 A and B). To reconcile this anomaly, we pursued a multiple linear regression model by combining Fpico with additional predictors, in turn. The most significant improvement (P < 0.05) was achieved by including the “ballast ratio” (RB) of mineral (silicate+calcite) to POC export in the statistical model (Fig. 4C),
Ballast minerals alone are a poor predictor of Teff (Table 1), but the elevated Si:POC export ratio in AAZ (Fig. 1B) may account for the anomalously high Teff in that region, raising the total explained variance to 97–98%. Mineral content likely contributes to Teff variations through its ballasting effect on particle density and sinking rate (11), as observed in AAZ (24), but might also protect soft organic matter from microbial degradation (33).
Applying our statistical description of Teff (Eq. 2) to global distributions of Fpico and RB (Fig. 5A and SI Materials and Methods) highlights the biases of the commonly used Martin curve relationship: It dramatically underpredicts POM transfer to the deep ocean throughout the high latitudes and overpredicts it in the low latitudes, with the exception of narrow coastal and equatorial upwelling zones.
Biogeochemical implications of variable transfer efficiency. (A) Transfer efficiency is mapped by applying a multiple linear regression model (Eq. 2) to global distribution of Fpico and RB. High latitudes and coastal/equatorial upwelling zones have Teff greater than Martin curve prediction (15%, thin black contour). (B) Fraction of POC exported from euphotic zone that will be sequestered for at least 100 y, before carbon is recirculated to the surface. (C) Mean areal C sequestration rates compared with cases with uniform Martin curve flux profile (Fig. S5B), averaged globally and over high- and low-latitude regions (poleward and equatorward, respectively, of thick black lines in B).
To quantify the importance of spatially varying transfer efficiency to the ocean carbon cycle, we estimated the sequestration of remineralized carbon in waters that remain out of contact with the atmosphere for at least 100 y. We used our OGCM to map the minimum depth at which the transport timescale to the surface exceeds 100 y (Fig. S5A). The fraction of exported POM that reaches this “sequestration depth” was then computed based on our reconstructed flux profiles (Fig. 5B and SI Materials and Methods), revealing the tropics and high latitudes as regions of efficient carbon sequestration (>50%). The global mean sequestration rate predicted in this scenario (∼0.9 mol C⋅m−2⋅y−1) is very similar to the prediction of a globally uniform Martin curve flux profile (Fig. 5C and Fig. S5 B−D). However, the scenario with variable Teff sequesters ∼30% more carbon in high-latitude regions (>40°N or <40°S) compared with the uniform case. Models with a constant remineralization scale likely underpredict the sequestration of carbon exported in these regions, suggesting an even greater sensitivity of atmospheric pCO2 to high-latitude nutrient drawdown than previously recognized (34).
Calculation of carbon sequestration. (A) The minimum depth at which remineralized carbon will be sequestered from the surface mixed layer for >100 y, calculated from the first passage time (51) in our OGCM (SI Materials and Methods, Carbon Sequestration). (B) The fraction of sinking organic matter that reaches the sequestration depth, given a uniform Martin curve flux profile. The equivalent map for spatially varying flux profiles is shown in Fig. 5B. (C and D) The rate of carbon sequestration (particle flux through the sequestration depth) in scenarios with (C) spatially varying and (D) uniform Martin curve flux profiles, given a satellite composite estimate of surface export (SI Materials and Methods, Export Estimates and SI Materials and Methods, Carbon Sequestration). Thick black contours separate high- and low-latitude regions used for the averages in Fig. 5C.
Our findings also have important implications for the global carbon cycle in a warming climate. As the upper ocean warms and stratifies, subtropical biomes dominated by picoplankton that form small, slow-sinking particles are expected to expand. Combined with direct effects of warming on microbial respiration, this might shoal the mean depth of organic particle remineralization and recirculate carbon to the surface mixed layer on shorter timescales. Quantifying and predicting these climate feedbacks will first require a deeper mechanistic understanding of the factors responsible for varying transfer efficiency in the modern ocean than our statistical analysis provides. This might be achieved using realistic models of particle size spectra, ballasting, and remineralization (11, 32), constrained by our diagnostic flux estimates and further guided by novel in situ respiration measurements (35) and size-resolved particle counts provided by visual profilers (31).
SI Materials and Methods
Reconstructing Deep Ocean Particle Fluxes.
Circulation model.
We reconstructed organic particle fluxes by analyzing nutrient distributions within the Ocean Circulation Inverse Model. A complete description of the model can be found in refs. 22 and 36; here we briefly summarize the main features of the model.
The underlying dynamical model is based on the linearized Navier–Stokes equations under a steady-state assumption, with baroclinic pressure gradients prescribed from climatological mean density fields, and with barotropic pressure gradients computed from modeled sea surface height fields. DeVries and Primeau (22) used observed temperature, salinity, radiocarbon, sea surface height, and sea surface heat and freshwater fluxes to invert for error terms in the horizontal momentum balance. These error terms account for all errors in the horizontal momentum balance that may result from errors in the diagnosed pressure gradients or wind stress fields, as well as errors resulting from neglected nonlinear terms and non-steady-state flows, and errors resulting from discretization of the momentum equations. This inversion procedure results in a close agreement to observed tracer distributions while maintaining a momentum balance that is close to geostrophic balance throughout most of the ocean. Because of the good agreement of the model with water mass and age tracers, the model faithfully reproduces the isopycnal structure and ventilation age distribution of the ocean interior; this is critical for our study, because diagnosing slow rates of nutrient accumulation along transport pathways requires that ocean circulation rates are accurately represented.
Here we use the most recent version of this model (36), which has a horizontal resolution of 2° × 2° and 24 vertical layers, and includes standard parameterizations for subgrid-scale processes including a rotated eddy diffusion tensor and a surface mixed layer scheme. This version of the model also includes CFC-11 as an additional tracer constraint. Observed tracer fields (temperature, salinity, radiocarbon, and CFC-11) are matched within observational uncertainty, generally with an R2 value of at least 0.95. Details of the model fit to observations can be found in ref. 36. Steady-state advective–diffusive flow fields of this model are stored in matrix form, facilitating efficient tracer simulations with the Transport Matrix Method (37).
Eight different configurations of this observationally constrained circulation model (36) were used to perform our calculations, propagating uncertainty in large-scale circulation rates into the reconstructed flux profiles. These configurations vary in their horizontal and vertical diffusivities and their weighting of observational and dynamical constraints, and propagate uncertainty in the distributions of observed tracers. The circulation models can be obtained from www.geog.ucsb.edu/∼devries/Home/Models.html.
Regionalization.
The global model domain was divided into two regions—an “interior” region, where particle fluxes were reconstructed from the observed phosphate (PO43−) distribution, and a “surface” region, where no reconstruction was attempted due to confounding factors. In regions where DOP is actively cycling, it is not possible to distinguish the geochemical signatures of particle remineralization from DOP remineralization. The interior and surface regions were therefore divided along the 300-m-depth horizon throughout the stratified low-latitude ocean, as DOP concentrations are observed to decline to negligible levels below this depth (19). In high-latitude regions, where seasonal mixed layer depths can exceed 300 m, we extended the surface region to the base of the wintertime mixed layer for two reasons: First, the mixed layer entrains DOP from the surface (20); second, nutrient concentrations in the seasonal mixed layer vary substantially over time (17), confounding the diagnosis of their accumulation rates in our annual mean circulation model. At each latitude/longitude, the boundary between the interior and surface regions was therefore defined as the maximum of 300 m or the deepest mixed layer depth from a monthly climatology (38).
Having assigned each model grid cell to either the surface or interior region, we divided the transport matrix into two components, following the procedure detailed in ref. 37. Matrix B (the “boundary propagator”) propagates observed nutrient concentrations across the boundary from the surface region into the interior. Matrix A circulates nutrients within the interior region and back across the boundary.
Diagnosing remineralization rates.
We divided the PO43− distribution in the interior region ([PO43−]interior) into two components: (i) a “circulated” component (Pcirc), which represents the PO43− in the surface region ([PO43−]surface) that is transported by water masses across the boundary and then mixed conservatively within the interior region; and (ii) a “remineralized” component [Prem(POP)] that represents the accumulated product of POP remineralization within the ocean interior, assuming that DOP remineralization is negligible within this region (see SI Materials and Methods, Regionalization). The two components follow the continuity equations
Regional flux profile reconstruction.
The solution to Eqs. S2 and S4 yields an estimate of POP remineralization rate (Jrem) at each gridpoint within the interior region of our ocean model. The downward flux of POP (FPOP) at each grid point was then reconstructed assuming it balances all remineralization occurring in grid cells directly below it in the water column. At depth horizon z′, the flux is therefore given by
Reconstructed fluxes were averaged over the large provinces illustrated in Fig. 1B, yielding a regional mean flux profile for each province. In Fig. 2, the upper limit of the reported flux profiles is set by the deepest depth of the boundary separating the ocean surface and interior (see SI Materials and Methods, Regionalization) in each province. They therefore extend up to 300 m in the low latitudes, and to the base of the deepest wintertime mixed layers in other regions.
Uncertainties and Error Propagation.
Monte Carlo approach.
We quantified the error associated with our reconstructed regional mean POP fluxes by propagating various sources of uncertainty into our calculation using a Monte Carlo approach. Each Monte Carlo calculation proceeded as follows:
i) A “realization” of the global PO43− distribution was generated (see SI Materials and Methods, PO43− Realizations).
ii) A circulation model was randomly chosen from a set of eight configurations of the Ocean Circulation Inverse Model (SI Materials and Methods, Reconstructing Deep Ocean Particle Fluxes).
iii) A remineralization restoring timescale (τ) was randomly chosen from the following set: 1 mo, 3 mo, 6 mo, 1 y, 2 y, 5 y.
iv) Regional mean POP fluxes were reconstructed following the methods outlined in SI Materials and Methods, Reconstructing Deep Ocean Particle Fluxes.
This procedure was repeated 10,000 times to ensure that all variables were well sampled, but the probability distributions of FPOP converged after just a few thousand iterations. For each province (Fig. 1B), this results in 10,000 independent estimates of the regional mean flux profile. The SD and 90% confidence intervals across these 10,000 estimates were used to place uncertainty bounds on the regional mean fluxes (Fig. 2).
PO43− realizations.
Each Monte Carlo calculation used a new synthetic global PO43− dataset, which adds smoothed error fields to the objectively analyzed WOA09 monthly climatology. Error fields were generated following the procedures outlined in ref. 22, which are briefly reviewed here.
For each month, we mapped the “total uncertainty” in PO43− on the WOA09 grid, which combines a contribution from the raw data itself (data uncertainty) and a contribution from the objective mapping procedure (mapping uncertainty). In grid cells that contained raw data, the data uncertainty was taken as the variance of all observations within that cell. Cells that lacked raw data were assigned the mean data uncertainty for their depth level. The mapping error was defined as the mean difference between the objectively analyzed and raw data for each depth level.
To construct a synthetic dataset, we first selected a random month from the WOA09 climatology, and generated an error field by taking a random draw at each grid cell from a normal distribution with zero mean, and with variance defined by the total uncertainty map for that month. The error field was then smoothed over a correlation length scale of 200 km and added to the objectively analyzed PO43− distribution for the selected month. Finally, the synthetic dataset was interpolated onto the circulation model’s grid for analysis.
Comparison with Observed Particle Fluxes.
Our new particle flux reconstructions were designed to overcome the shortcomings of bottom-moored sediment trap data traditionally used to quantify the flux of organic matter to the deep ocean: They are subject to well-known collection biases (13) and are too sparsely deployed to accurately capture time mean, large-scale patterns. Nevertheless, we compared our results to a large database of bottom-moored sediment trap measurements (25) as a reality check on the magnitude our reconstructed rates (Fig. S1).
To facilitate comparison with our results, measured POC fluxes were converted to POP using the empirical stoichiometric algorithm of ref. 26 in which organic P:C is related to the ambient [PO43−] where the organic matter formed. A P:C conversion factor was determined for each sediment trap observation using the closest surface [PO43−] (WOA09) to the latitude/longitude of the trap deployment. After conversion, the observations were binned to the model’s vertical grid and into the eight ocean provinces considered here (Fig. 1B). Because the sediment traps measure local rates and are too sparse to capture regional mean fluxes, we expanded the uncertainty estimates on our reconstructed flux profiles for this comparison (Fig. S1) to reflect spatial variability in diagnosed FPOP within each region (envelopes in Fig. 2 represent uncertainty in regional mean fluxes only).
For the most part, observed particle fluxes fall well within the range of model-diagnosed rates in each region but tend to cluster toward the lower end of that range (Fig. S1). This might be explained by a number of factors. First, sediment traps may have been deployed at locations or times with particle fluxes lower than the large-scale time mean rates, although it seems unlikely that this would occur systematically across regions. Second, diagnosed flux rates in the deep ocean may be systematically too large, which would occur if the model’s overturning circulation were unrealistically fast. We find this unlikely given the observational constraints imposed on the model and its ability to accurately reproduce water mass age tracers (2). Finally, the discrepancy might arise from the incomplete collection efficiency of bottom-moored sediment traps. Given that such a systematic collection bias has been observed for these instruments and discussed previously in the literature (13, 41), we find this explanation most likely.
Export Estimates.
Global export ensembles.
To derive estimates of particle transfer efficiency through the mesopelagic zone, we combined our reconstructed particle fluxes with independent estimates of annual mean POP export from the euphotic zone. We compiled two ensembles of global export estimates, one derived from remote sensing algorithms (satellite ensemble) and one from global models in the CMIP5 archive (CMIP5 ensemble).
The satellite ensemble contains nine members, comprising all permutations of three NPP algorithms and three estimates of the export/NPP ratio. The NPP algorithms used were: (i) the chlorophyll-based Vertically Generalized Production Model (VGPM) (42); (ii) the VGPM–Eppley model (43), which modifies the temperature dependence of production in the standard VGPM; and (iii) the Carbon-based Production Model (CbPM), which uses backscatter-derived particulate carbon as the metric for plankton biomass, rather than chlorophyll (44). We downloaded monthly data from www.science.oregonstate.edu/ocean.productivity/index.php for 1997–2009, and averaged the results into monthly climatologies. For the export/NPP ratio, we used the empirical relationships of refs. 45⇓–47, each of which takes NPP and temperature as inputs. Monthly export ratios were computed and combined with NPP to yield nine monthly climatologies of carbon export. These were averaged over the seasonal cycle and finally converted to POP using the spatially varying P:C relationship of ref. 26.
The CMIP5 export ensemble comprises five members, derived from the Earth System Models listed in Table S1. For each model, we downloaded particulate organic carbon sinking fluxes at 100 m (the putative base of the euphotic zone) from “historical” simulations in the CMIP5 Data Portal, and averaged across the time interval 1951–2000. Each of these models contains a fully prognostic ecosystem, in which primary production is driven by temperature, light, and multiple limiting nutrients, and export is modulated by higher trophic levels. Differences in time mean export between the models stems from their biological parameterizations and underlying physical models, which control the nutrient supply from depth. POC was converted to POP using the stoichiometric assumptions of each model [“Redfield” stoichiometry in most; variable stoichiometry in GFDL-ESM2M (Table S1)].
All satellite and CMIP5 POP export estimates were interpolated to a common grid and averaged regionally (Fig. S2A) using the same provinces as the reconstructed particle fluxes (Fig. 1B).
Data constraints.
We derived composite regional export estimates by evaluating each satellite algorithm and CMIP5 model against tracer-based constraints. In the Atlantic and Pacific Oceans, we used the compilation of ref. 28, which combined a range of tracer-based estimates to define likely ranges of annual carbon export at five locations, spanning five of our averaging regions (Fig. S2B): (i) Ocean Station Papa (OSP), within the NP region; (ii) Station ALOHA, within the STP region; (iii) Bermuda Atlantic Time Series (BATS), within the STA region; (iv) Subarctic Atlantic (40°N–65°N, 10°W–60°W), within the NA region; and (v) Equatorial Pacific (170°W–95°W), within the ETP region.
In the Southern Ocean, we used carbon export rates estimated from O2/Ar ratios by ref. 27. We used only the observations between 120°E and 150°W and poleward of 40°S, where data density was highest, and averaged them within our AAZ and SAZ regions for each month. The monthly values were then integrated over the seasonal cycle assuming zero export during unsampled wintertime months, consistent with the methods used in other regions (28).
To compare against these constraints, annual mean satellite and CMIP5 carbon export rates were averaged in 2° × 2° bins centered around time series sites (OSP, ALOHA, BATS), and over the latitude/longitude ranges listed above for other regions (Fig. S2B). For the purpose of Fig. 3B, satellite and model export was defined to be in good agreement with the tracer constraint if it fell within 1σ of the observational constraint.
Composite export estimates were generated by assuming that satellite algorithms and CMIP5 models that best reproduce observational constraints at locations within a given region also provide the most realistic rates for the region as a whole. For each location, we quantified the probability of drawing each satellite/model prediction from a normal distribution with mean and SD defined by the tracer-based export estimates (Fig. S2B). The relative magnitude of these probabilities was then used to define weighting factors for each member of the satellite and CMIP5 ensembles independently for each region. These weighting factors are listed in Table S2 and were used to derive regional weighted means and SDs of POP export (Fig. S2C). In the Tropical Atlantic region (ETA), where data constraints were lacking, we used the same weighting factors as the Tropical Pacific region (ETP). Given the relatively narrow spread of satellite and CMIP5 export in the tropics (Fig. S2A), a similar result is obtained by weighting each ensemble member equally in ETA.
Statistical Analysis of Teff.
Transfer efficiency estimates were correlated against five environmental and ecological factors that have been hypothesized to control particle sinking and remineralization rates (Table 1, Fig. 4, and Fig. S4). We obtained global estimates of these factors and averaged regionally over the provinces in Fig. 1B. The previous hypotheses are reviewed below, and the global estimates we compiled are discussed.
Upper-ocean temperature (T100–500m).
Heterotrophic bacteria decompose organic matter faster in warmer waters that stimulate their metabolism (7). The temperature of waters through the mesopelagic zone has therefore been hypothesized as a dominant control on particle transfer efficiency to the deep ocean (8). Following ref. 8, we averaged the annual mean temperatures in the WOA09 climatology between 100 m and 500 m, where the majority of particle remineralization occurs.
NPP seasonality.
Large seasonal phytoplankton blooms rapidly export organic matter without significant recycling and processing in the surface ocean (10). The lability of this fresh organic matter should result in a negative correlation between the seasonality of NPP and transfer efficiency (9, 30). We quantified a dimensionless metric of seasonality in each of the three NPP climatologies (VGPM, VGPM–Eppley, CbPM) by normalizing the monthly range (maximum–minimum) to the annual mean rate. These metrics were then averaged across the regions in Fig. 1B, and the average of the three values for each region was correlated against Teff.
Export ratio.
If the degree of recycling within the euphotic zone strongly influences the lability (and remineralization rate) of exported organic matter, the ratio of carbon export to NPP (e ratio) should be negatively correlated with Teff (15, 30). We calculated regional averages of the annual mean export ratio using the three export ratio climatologies discussed in SI Materials and Methods, Export Estimates. Because all of the empirical relationships yield similar large-scale patterns of the e ratio, they result in a similar correlation to Teff, and only one is reported in Table 1 and considered in the MLR model.
Size structure of phytoplankton community (Fpico).
Recent observations have revealed a strong relationship between the size structure of phytoplankton communities and the size distribution of organic particles sinking from the euphotic zone (31). Communities with a larger contribution of picoplankton to photosynthetic biomass (Fpico) produce a “steeper” spectrum that is more skewed to toward small (<100 μm) particles. Because particle sinking speeds are strongly controlled by their size (48, 49), remineralization should occur at a shallower depth when these small particles are more abundant, resulting in a negative correlation between Fpico and Teff. We created a monthly climatology of Fpico using the chlorophyll-based algorithm of ref. 23, which was then averaged annually and over regions (Fig. 1B) for comparison with Teff.
Ballast ratio (RB).
The mineral (SiO2+CaCO3) content of particles has been hypothesized to control remineralization depth, either by setting their density and sinking speed (11) or by protecting “soft” organic matter from microbial degradation (33). We obtained independent global estimates of the SiO2:POC and CaCO3:POC mass ratios of exported organic matter. The former was derived by diagnosing the Si:P export ratio from the global distributions of PO43− and Si(OH)4 in our OGCM following the procedure of ref. 50 and was converted to Si:POC following ref. 26. Because the alkalinity distribution (tracer of CaCO3 export) is less well sampled, we obtained the CaCO3:POC export ratio from the predictive ecosystem model in GFDL-ESM2M (Table S1), where it is governed by the plankton community composition. Both mineral:carbon ratios were averaged regionally and summed to give the total ballast ratio (RB), which was correlated against Teff. Because the two minerals are thought to differ in their ballasting potential (11), we also correlated SiO2:POC and CaCO3:POC against Teff individually, but neither showed a significant improvement over RB.
MLR model.
Of the individual variables considered, only T100–500m and Fpico were significantly (P < 0.05) correlated with Teff in the expected direction. T100–500m underpredicts the difference in Teff between tropics and subtropics and overpredicts the difference between subarctic regions and the Antarctic zone (AAZ), whereas Fpico leaves one significant outlier in AAZ. We took an MLR approach to determine if other variables could explain these residuals. MLR models were used to combine T100–500m and Fpico with the other variables described above, in turn, and f tests were used to determine whether the statistical fit improved significantly (P < 0.05). No additional variable could significantly improve the Teff prediction when combined with T100–500m, but RB resulted in a significant improvement when combined with Fpico.
Carbon Sequestration.
We investigated the implications of our results for ocean carbon storage, by computing the amount of remineralized CO2 that remains sequestered within the ocean interior for at least 100 y. First, we used our OGCM to determine the mean transport timescale to the surface mixed layer from each interior grid cell (the “first passage time”) using the methods described in ref. 51, and found the minimum depth at which this timescale exceeds 100 y (Fig. S5A). This was defined as the sequestration depth. We then mapped the fraction of organic matter exported from the euphotic zone that would reach this depth (the “sequestered fraction”). This requires a higher-resolution estimate of spatially varying particle flux profiles than is afforded by our regionally averaged reconstructions (Fig. 2). We obtained this by defining a power-law flux profile at each latitude/longitude,
Finally, to convert to actual carbon sequestration rates, we combined our estimates of sequestered fraction with a map of carbon export rates from the euphotic zone. The export map was derived by averaging nine satellite estimates (SI Materials and Methods, Export Estimates) according to the regional weighting factors listed in Table S2 (each algorithm was weighted equally outside of these regions, e.g., in the Indian Ocean). The resulting carbon sequestration rates (Fig. S5 C and D) were averaged over low latitudes (equatorward of 40°) and high (poleward of 40°) latitudes for a simple comparison between the Martin curve and variable Teff cases (Fig. 5C).
Acknowledgments
We acknowledge the World Climate Programme’s Working Group on Coupled Modeling, which is responsible for the CMIP5 archive, and the modeling groups listed in Table S1 for producing the model output we used. This work was supported by the Marine Microbes Initiative of the Gordon and Betty Moore Foundation (GBMF 3775).
Footnotes
↵1Present address: Department of Earth and Environmental Sciences, University of Rochester, Rochester, NY 14627.
- ↵2To whom correspondence should be addressed. Email: tsweber{at}uw.edu.
Author contributions: T.W., T.D., and C.D. designed research; T.W. performed research; T.W., J.A.C., and S.W.L. analyzed data; and T.W., J.A.C., S.W.L., and C.D. wrote the paper.
The authors declare no conflict of interest.
This article is a PNAS Direct Submission.
This article contains supporting information online at www.pnas.org/lookup/suppl/doi:10.1073/pnas.1604414113/-/DCSupplemental.
Freely available online through the PNAS open access option.
References
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵.
- Kiørboe T,
- Grossart H-P,
- Ploug H,
- Tang K
- ↵
- ↵.
- Marsay CM,
- Sanders RJ,
- Henson S,
- Pabortsava K,
- Achterberg EP
- ↵
- ↵.
- Lam PJ,
- Doney SC,
- Bishop JKB
- ↵
- ↵
- ↵
- ↵.
- Buesseler KO, et al.
- ↵
- ↵
- ↵.
- Garcia HE, et al.
- ↵
- ↵
- ↵.
- Torres-Valdés S, et al.
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵.
- Galbraith ED,
- Martiny AC
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵.
- Armstrong R,
- Lee C,
- Hedges JI,
- Honjo S,
- Wakeham SG
- ↵
- ↵
- ↵
- ↵.
- Khatiwala S
- ↵.
- de Boyer Montégut C,
- Madec G,
- Fischer AS,
- Lazar A,
- Iudicone D
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵.
- Laws EA,
- D’Sa E,
- Naik P
- ↵
- ↵
- ↵.
- Jin X,
- Gruber N,
- Dune JP,
- Sarmiento JL,
- Armstrong R
- ↵
Citation Manager Formats
More Articles of This Classification
Physical Sciences
Related Content
- No related articles found.
Cited by...
- No citing articles found.