Biogeochemical controls of surface ocean phosphate.

Surface ocean phosphate is commonly below the standard analytical detection limits, leading to an incomplete picture of the global variation and biogeochemical role of phosphate. A global compilation of phosphate measured using high-sensitivity methods revealed several previously unrecognized low-phosphate areas and clear regional differences. Both observational climatologies and Earth system models (ESMs) systematically overestimated surface phosphate. Furthermore, ESMs misrepresented the relationships between phosphate, phytoplankton biomass, and primary productivity. Atmospheric iron input and nitrogen fixation are known important controls on surface phosphate, but model simulations showed that differences in the iron-to-macronutrient ratio in the vertical nutrient supply and surface lateral transport are additional drivers of phosphate concentrations. Our study demonstrates the importance of accurately quantifying nutrients for understanding the regulation of ocean ecosystems and biogeochemistry now and under future climate conditions.

INTRODUCTION

Dissolved inorganic phosphate (DIP) plays a central biogeochemical role in marine environments (, ), but the surface DIP concentration is commonly below the detection limit for standard methods (). We know the general pattern of high DIP in high-latitude regions, low DIP in the subtropical gyres, and intermediate levels in coastal and equatorial upwelling regions. However, we understand less about spatial DIP patterns within oligotrophic waters, although subtle DIP changes could be an important driver (–) or tracer (, ) for ocean biogeochemistry.

As a driver, DIP influences phytoplankton growth, abundance, and diversity. Phosphorus is commonly regarded as the ultimate limiting nutrient and has the longest ocean residence time of major biologically limiting nutrients (). As any N deficits may be offset by N fixation, changes in surface ocean stratification and the vertical supply of DIP are thought to be major regulators of future ocean productivity and downward carbon export (, ). DIP can also be a proximate limiting nutrient in regions such as the western North Atlantic Ocean (, ) and the Eastern Mediterranean Sea (). In addition to phytoplankton growth and productivity, subtle changes in DIP availability may influence the distribution and biogeochemical roles of plankton. Ocean heterotrophic and autotrophic plankton are adapted to relative differences in P versus N availability—often via the presence or absence of specific nutrient acquisition genes—and such adaptation can play an important role in regulating phytoplankton growth (, ). DIP availability can also be a tracer of ocean biogeochemical processes (, ). Regions with low DIP can have a very efficient biological pump (), and extreme DIP drawdown is often indicative of extensive N fixation (). Thus, mapping the detailed distribution of DIP at low concentrations is central to understanding regional differences in core ocean ecosystem and biogeochemical processes.

The underlying reason for a limited understanding in the distribution of surface DIP is that the standard methodology using an autoanalyzer combined with a spectrometer has high variance and low interlaboratory accuracy, below ~100 nM (, ). The standard methodology is routinely used for most oceanographic research programs including global hydrographic surveys such as GO-SHIP and GEOTRACES [albeit with a few exceptions (, )]. These low-sensitivity data get propagated into the World Ocean Atlas, leading to a potential misrepresentation of the geographic distribution of DIP in oligotrophic waters. Several alternative and highly sensitive methods for determining DIP in seawater have been developed and more widely used in recent years. These include the magnesium-induced coprecipitation (MAGIC), liquid waveguide capillary cell (LWCC), and solid-phase extraction (SPE) methods (–). MAGIC has a detection limit as low as ~0.5 nM, whereas LWCC and SPE have a detection limit of ~1 to 2 nM. Using these more sensitive techniques, several studies have proposed regional differences in DIP among oligotrophic environments (e.g., low DIP in the western North Atlantic Ocean) (), but DIP levels are uncertain in most oligotrophic regions. However, we now have the tools and data available to develop a global understanding of DIP variation across surface ocean environments.

Here, we conduct a global assessment of surface ocean DIP concentrations based on high-sensitivity measurements. On the basis of this assessment, we first ask whether there are systematic regional differences in DIP concentrations across oligotrophic marine environments. Second, we evaluate how well observation-based climatologies and Earth system models (ESMs) capture DIP levels in oligotrophic environments. Third, we quantify how variation in surface DIP influences important ocean ecosystem properties using shipboard measurements, satellite-based estimates, and ESM predictions. Last, we identify possible physical and biological mechanisms driving the observed DIP distribution using new simulations within the ocean component of the Community Earth system model (CESM).

Observed patterns of DIP in oligotrophic regions

To understand the biogeochemical role of P across oligotrophic regions, we compiled a global dataset of observed DIP concentrations ([DIP]obs) measured using high-sensitivity methods (Fig. 1). These data originated from 42 cruises from all major oligotrophic regions (Table 1).

Fig. 1

Global distribution of surface DIP.

(A) Observed DIP using high-sensitivity techniques ([DIP]obs). (B and C) Climatologically predicted surface DIP concentrations ([DIP]predict) from the World Ocean Atlas (WOA13) (B) and GLODAPv2 (C). (D to L) Surface DIP concentrations predicted using ESM predictions. The dataset has a total of 50,591 observations including 41,747 from the top 30 m of the water column (shown here). Climatological and ESM-predicted values are at 1° resolution and represent annual means.

Table 1

Summary of cruise data for this study.

MQ, Milli-Q water; NaOH, supernatant following alkaline precipitation of phosphate; Unp., unpublished.

Cruise no.	ID	Region	No. of samples	(<30 m)	Method	Frozen	Prefiltered	Blank	Ref.
1	Cook-Book	North East Pacific	215	65	MAGIC	Yes	No	NaOH	Unp.
2	BULA	Central Pacific	48	8	MAGIC	Yes	No	NaOH	Unp.
3	Beachbash	Central Pacific	131	42	MAGIC	Yes	No	NaOH	Unp.
4	SUPER-HI-CAT	North East Pacific	84	42	MAGIC	Yes	No	NaOH	Unp.
5	BIOSOPE	South East Pacific	288	77	MAGIC	No	Yes	NaOH	(6)
6	PROSOPE	Mediterranean Sea	137	26	MAGIC	No	Yes	NaOH	(39)
7	AMT12	Atlantic Ocean	93	45	LWCC	No	No	MQ	(10)
8	AMT13	Atlantic Ocean	60	10	LWCC	No	No	MQ	(10)
9	AMT14	Atlantic Ocean	209	59	LWCC	No	No	MQ	(10)
10	AMT15	Atlantic Ocean	220	128	LWCC	No	No	MQ	(10)
11	AMT16	Atlantic Ocean	272	117	LWCC	No	No	MQ	(10)
12	AMT17	Atlantic Ocean	276	92	LWCC	No	No	MQ	(10)
13	BOUM	Mediterranean Sea	424	167	LWCC	No	No	MQ	(40)
14	AR16	Western North Atlantic	94	45	MAGIC	Yes	Yes	NaOH	Unp.
15	COST2005	North Pacific	174	60	LWCC	No	No	MQ	(41)
16	CR 1455, CR 1487, CR 950	South China Sea	24	24	MAGIC	Yes	No	NaOH	(42)
17	POOB2	North Pacific	13	13	MAGIC	Yes	No	NaOH	(43)
18	GA03	North Atlantic	335	168	LWCC	No	No	MQ	(44)
19	KH04-5	West Pacific Ocean	13,997	13,997	LWCC	No	No	NaOH	(22)
20	KH05-2	West Pacific Ocean	6,175	6,175	LWCC	No	No	NaOH	(22)
21	SCS	South China Sea	45	12	LWCC	No	No	MQ	(19)
22	MR09-01	South Pacific	38	38	LWCC	Yes	No	NaOH	(45)
23	R/V Oceanus	North West Atlantic	171	75	MAGIC	Yes	No	MQ	(21)
24	BATS	North West Atlantic	1815	444	MAGIC	Yes	Yes	NaOH	(32)
25	AE1319	North West Atlantic	42	12	MAGIC	Yes	Yes	NaOH	(5)
26	Trophic BATS	North West Atlantic	893	183	MAGIC	Yes	Yes	NaOH	(32)
27	NH1418	Central Pacific	188	48	MAGIC	Yes	Yes	NaOH	(46)
28	BVAL	North West Atlantic	1,254	290	MAGIC	Yes	Yes	NaOH	(5)
29	AE1206/1219	North West Atlantic	110	32	MAGIC	Yes	Yes	NaOH	(5)
30	KH14-3	Central North Pacific	74	39	LWCC	Mix	No	NaOH	(47)
31	KH09-5	Central Indian Ocean	49	19	LWCC	No	No	NaOH	(48)
32	KH-06-2, MR07-01,MR07-06, KT-06-21,Nagasaki-Maru 242	West Pacific Ocean	43	43	LWCC	Mix	No	NaOH	(49)
33	Umitaka-Maru cruise	Indian Ocean	7,073	7073	LWCC	No	No	NaOH	(50)
34	HOT	North Pacific	2,545	656	MAGIC	Yes	No	NaOH	(51)
35	NH1417	North East Pacific	38	28	MAGIC	Yes	No	NaOH	(52)
36	Mixed	South China Sea	2,209	2209	SPE	No	No	MQ	(20)
37	Mixed	South China Sea	1,592	912	MAGIC	Yes	No	NaOH	(53)
38	GA06/D361	Central North Atlantic	176	176	LWCC	No	No	MQ	(54)
39	GP13	South West Pacific	777	91	LWCC	No	No	MQ	(17)
40	KT-05-24, KT-06-21, KT-07-22	North West Pacific Ocean	7,838	7,838	LWCC	No	No	NaOH	(55)
41	GEOTRACES GApr08, JC150	North Atlantic	119	119	LWCC	No	No	MQ	Unp.
42	OUTPACE	South Pacific Ocean	95	16	LWCC	No	Yes	None	(56)

Our extensive compilation revealed strong regional patterns in observed DIP concentrations (Fig. 1A). Two major Northern Hemisphere regions were ultraoligotrophic (i.e., [DIP]obs < 10 nM). As described in past studies [e.g., (3, 11)], the first region included the North Atlantic Ocean and Mediterranean Sea covering the width of the ocean basin below 40°N. The second region was detected in the western North Pacific Ocean centered at ~25°N and stretched from the Philippine Sea to the Hawaiian Islands. A low-DIP zone in the western Pacific Ocean has been previously suggested (, ), but our synthesis uncovered the broad extent of this region (Fig. 1A). In the Southern Hemisphere, gyre regions with low [DIP]obs were smaller or absent, and none of them were ultraoligotrophic. A low-[DIP]obs region centered on 28°S was seen in the southwestern Pacific Ocean. There were a few samples on the western side of the South Atlantic Ocean with low [DIP]obs, but the extent of this low-DIP region was unknown due to limited data coverage. We had incomplete data coverage in the Indian Ocean, but there appeared to be two regions with low [DIP]obs levels. The first region was the western side of the South Indian Ocean Gyre and, thus, a location mirroring other subtropical low-[DIP]obs regions. The second low-[DIP]obs region was near the Indonesian archipelago and might be linked to a combination of high iron supply () and the Indonesian throughflow supplying low P relative to N waters from the North Pacific Ocean (). In summary, the distribution of [DIP]obs suggested two consistent trends. First, there appeared to be larger low-DIP regions with more complete drawdown in the Northern Hemisphere. Second, the most intense low-DIP regions were generally found on the western side of the subtropical gyres, leading to a longitudinal gradient within basins. Thus, this analysis revealed several previously unrecognized low-DIP regions and consistent geospatial trends in the surface [DIP]obs distribution.

Bias in extant DIP climatologies and ESM predictions

We next evaluated the accuracy of current observation–based climatologies and ESM-simulated DIP concentrations ([DIP]predict). The climatological estimates were from the World Ocean Atlas (WOA13) and GLODAP (v2.2016b) (, ). The nine ESMs were from the Fifth Phase of the Coupled Model Intercomparison Project (CMIP5) using simulated values for the 1990s (table S1) (). Observational climatologies had positive phosphate biases at low [DIP]obs (fig. S1). The median bias for [DIP]predict from the World Ocean Atlas was +24× at the lowest [DIP]obs (<5 nM). The positive bias declined in regions with higher [DIP]obs and disappeared above 200 nM. The median bias for [DIP]predict from the GLODAP climatology was lower than that from the WOA climatology (+8× at the lowest [DIP]obs concentration) and disappeared above 100 nM. We found substantial positive and negative biases in [DIP]predict from ESMs. Most of the models had bias profiles resembling the climatologies, perhaps not unexpectedly, as WOA is commonly used in the development and optimization of ocean biogeochemical models. The models overestimated DIP by 10 to 100× at low [DIP]obs, and the systematic bias disappeared at higher [DIP]obs (>200 nM). The HadGEM-2ES and the two IPSL models had very different bias profiles, but both models systematically underestimated [DIP]obs. Thus, we detected biases at low DIP concentrations for both climatologies and ESMs.

The observed regional patterns generally matched WOA and GLODAP distributions, but less so for many ESMs (Fig. 1, B to K). WOA and GLODAP captured the boundaries of the low-DIP regions but failed to reveal several subtle differences across the low-DIP regions including the stronger drawdown in the Northern versus the Southern Hemisphere regions (Fig. 1, B and C). ESMs showed regional biases, and the extent of the simulated equatorial Pacific upwelling zone was generally poorly defined (Fig. 1, D to L). CESM1-BGC and MPI-ESM-LR captured many features including low concentrations in the North Atlantic and the western side of both the South and North Pacific Ocean. However, these two models had a latitudinal gradient of [DIP]predict in the Indian Ocean not seen in our observations. Furthermore, neither of them captured the large North Pacific Ocean low-DIP region and both lacked the observed longitudinal gradient in the Pacific Ocean. The GFDL models underestimated [DIP]predict in the Indian Ocean (this is also seen in MPI-ESM-MR) and lacked a clear longitudinal gradient in the Pacific and South Atlantic Ocean. The NorESM1-ME had a strong positive bias, and the two IPSL models had strong negative biases in [DIP]predict concentrations. The latter models also lacked the geographic variation seen in observations. Thus, it is clear that current formulations of ESMs struggled in simulating the geospatial distribution of surface DIP.

Phosphate as a driver of ocean ecosystem properties

We next wanted to quantify the role of DIP as a driver of key ocean ecosystem properties. Hence, we compared our high-sensitivity DIP observations to satellite-based estimates of surface chlorophyll concentration (Chlsat) and net primary production (NPPsat) (Fig. 2). We detected a similar but nonlinear relationship between [DIP]obs versus Chlsat, or [DIP]obs versus NPPsat. At low [DIP]obs, the levels of both Chlsat and NPPsat were depressed relative to a maximum observed between 20 and 40 nM. Above 40 nM, both Chlsat and NPPsat actually declined by ~30%. There is some uncertainty in the satellite estimate of NPP (), so we also compared [DIP]obs with NPPobs from a subset of cruises that included shipboard NPP measurements (fig. S2). Here, we found a similar pattern with depressed NPPobs at low [DIP]obs, with maximum values between 20 and 80 nM, and then reduced NPPobs at higher concentrations. Chlorophyll and NPP levels not only are partially linked but also include some independent regulation (e.g., lower chlorophyll per cell in stratified waters or higher phytoplankton growth rate at elevated temperature). Thus, we should not expect to see parallel patterns per se but, nevertheless, find a similar but nonlinear relationship between surface P and both ecosystem features.

Fig. 2

Relationships between DIP, chlorophyll, and NPP among observations and ESMs.

(A and B) Relationships between [DIP]obs and chlorophyll concentrations (Chlsat) (A) or NPP (NPPsat) (B). Chlorophyll and NPP are from MODIS satellite observations (n = 38,653). (C and D) Relationships between predicted DIP ([DIP]predict) and chlorophyll concentrations (Chlpredict) (C) or NPP (NPPpredict) (D) from ESMs. Chlorophyll or NPP relationships from ESMs are depicted in detail in fig. S3.

Across ESMs, we found a very different and variable relationship between DIP concentrations and ecosystem properties than observed (Fig. 2 and figs. S3 and S4). As such, none of the model relationships between [DIP]predict and Chlpredict or NPPpredict were significantly correlated. Most models had a monotonically positive relationship between [DIP]predict and Chlpredict. For some models, this relationship was well constrained, whereas others showed considerable scatter (fig. S3). Similarly, the relationship between [DIP]predict and NPPpredict varied among models as seen for chlorophyll (Fig. 2). Thus, most models had monotonically positive relationships between [DIP]predict and NPPpredict (fig. S4). This link had lowest variance for CESM1-BGC and the MPI models, but other models included a lower slope and higher variance. One exception was GFDL-ESM2G, which had a pattern with a maximum chlorophyll or NPP at intermediate [DIP]predict values and somewhat resembled the observation-based pattern. In summary, there was considerable variation in the relationships between model DIP concentrations and surface ecosystem properties, and most models showed limited agreement with the in situ relationships.

Processes regulating surface DIP in oligotrophic regions

We next evaluated the processes regulating surface DIP across oligotrophic regions. Generally speaking, surface nutrient concentrations are the balance of nutrient sources and sinks. A conceptual model commonly stated is that the external iron supply via atmospheric deposition (Featms) regulates N fixation and the relative drawdown stoichiometry of P versus N (). Thus, regional variation in Fe deposition has been proposed as the principal driver of regional DIP levels between and within subtropical gyres (, ). To test this iron control hypothesis, we performed model simulations using a newly optimized version of the CESM (v2). This optimized model captured most of the observed patterns in DIP across oligotrophic regions and represented a large improvement over CESM1 (fig. S5). First, we tested whether a uniform between-gyre Featms leads to the same degree of DIP drawdown. This was done by applying a scaling factor to the current Featms supply. Hence, the mean deposition to each southern subtropical gyre was made to match the mean deposition to the North Pacific Subtropical Gyre. This increase in the Featms supply resulted in a Southern Hemisphere stimulation of N fixation (fig. S6) and slightly lower DIP levels (Fig. 3). However, the geospatial pattern in DIP was largely conserved. Thus, the hemisphere bias in the size of low-DIP regions or the degree of DIP drawdown could not be tied to differences in Fe deposition.

Fig. 3

Changes in the global distribution of DIP as a function of dust deposition patterns.

(A) Estimated DIP in the current ocean (i.e., control scenario). (B) Estimated DIP for the scenario of increased dust deposition in the Southern Hemisphere gyre (reaching North Pacific Ocean subtropical gyre levels). (C) Estimated DIP for the scenario of even dust deposition across each subtropical gyre. (D) Differences in DIP between the “Southern Hemisphere scenario” [as shown in (B)] and current levels [as shown in (A)]. (E) Differences in DIP between “even gyre scenario” [as shown in (C)] and current levels [as shown in (A)]. Ocean biogeochemistry was simulated with the ocean component of the CESM (v2).

We next tested whether the atmospheric iron supply controlled the observed longitudinal gradients in surface DIP. We first removed all spatial intragyre gradients in deposition but retained the mean. This new deposition pattern led to N fixation shifting further east in most gyres but did not change total fixation levels notably (236 Tg N/year versus 240 Tg N/year; fig. S6). Furthermore, a uniform deposition pattern had little impact on the geospatial distribution of surface DIP levels. As a case in point, most dust and Fe are deposited on the eastern side of the North Atlantic Ocean, but the lowest DIP concentrations are observed on the western side (). Last, [DIP]obs was compared with the surface ocean vertical velocity and, thus, the physical strength of nutrient supply from deeper waters (fig. S7). Here, we observed a nonlinear relationship suggesting that changes in the vertical nutrient supply could not explain the regional differences in DIP (fig. S8). In conclusion, gradients in Fe disposition and the physical strength of nutrient supply cannot explain the inter- and intragyre bias in DIP.

We propose that regional differences in the vertical iron supply and lateral circulation are two additional factors controlling surface DIP concentrations. Thermocline waters are typically deficient in iron relative to macronutrients due to particle scavenging losses (), but the deficit is not uniform. An examination of modeled nutrient ratios revealed that nutricline waters were Fe-depleted relative to macronutrients on the eastern side of the North Pacific, South Pacific, and South Atlantic subtropical gyres (fig. S9). In contrast, both observations and our model showed that the North Atlantic Ocean has significantly higher Fe:macronutrient ratios in subsurface waters likely contributing to elevated Fe availability [fig. S8 and ()]. We do not know the full cause of the spatial differences in Fe:macronutrient ratios. However, we hypothesize that Fe-depleted water from the Southern Ocean contributes to an Fe deficit and elevated DIP concentrations in the southern hemisphere and on the eastern side of several basins (). In addition to regional differences in the Fe supply ratio, the gyre circulation transporting higher nutrient water from the eastern equatorial and coastal upwelling zones westward may contribute to a longitudinal gradient. Dissolved P relative to N is then systematically stripped from the surface ocean over time due to external N inputs (mainly N fixation). Furthermore, nutricline waters in the eastern tropical Pacific Ocean has low dissolved N:P ratios due to intense water column N loss, perhaps contributing to the observed longitudinal surface dissolved inorganic N:P and DIP gradients (fig. S9). Such a circulation/P-drawdown model has been proposed for the Bermuda Atlantic Time-series (BATS), where there is limited N fixation but still very low P levels (). Here, we extend this model as a mechanism that influences DIP concentrations in most subtropical gyres. Thus, the observed geospatial pattern of DIP suggests a suite of processes controlling surface DIP concentrations in the oligotrophic regions including the differential supply of Fe relative to other nutrients from both above and below the surface, the associated patterns of N fixation, and the stoichiometry and magnitude of lateral and vertical macronutrient inputs.

DISCUSSION

Multiple past studies as well as climate models have suggested that the availability of phosphate is an important driver of ocean biogeochemical processes. Here, it is worth considering that phosphate availability consists of two linked but not overlapping variables: the nutrient flux to the surface waters and the standing stock. The former is important for key biogeochemical exchanges between the surface and deep ocean (e.g., downward carbon export), whereas the latter (studied here) may regulate local ecosystem processes such as phytoplankton biomass accumulation and growth. The relationships between DIP availability and chlorophyll or NPP are mostly positive in models but nonlinear in the observations. Hence, we only see evidence for P limitation below ~20 nM DIP (Fig. 2), and other nutrients (probably mostly N) may be limiting at higher concentrations. A physiological threshold at this low concentration is consistent with plankton physiology studies. First, phytoplankton have expressed high affinity and uptake rates for DIP at very low concentrations (, ). Second, the cellular requirement for P appears to be more variable than for N (). Third, organic forms of P serve as alternative nutrient sources due to the efficient cycling of P in oligotrophic regions. Thus, marine phytoplankton can acclimate and adapt to very low DIP and still achieve high growth rates. This capacity for tolerating P stress implies that changes in the vertical nutrient supply will mainly limit phytoplankton biomass or primary production if diazotrophs are unable to supply sufficient N. Only in extreme cases are phytoplankton biomass and growth directly controlled by DIP as their uptake and demand are less sensitive to P limitation than currently expressed in most models (). As N fixation is tied to the Fe flux, an additional supply of iron into the subtropical gyres could provide considerable stimulation of biological processes due to the large pool of residual DIP seen in many parts of the ocean.

Our analysis demonstrated clear biases in most ESMs with important implications for climate change predictions. A future decline in the vertical P flux has been proposed as a first-order control on changes in the standing stock of phytoplankton biomass and ocean primary productivity (). However, our comparison revealed large uncertainties in how ESMs represent P cycling and regulation of surface ecosystem processes. There are likely many factors that control P cycling in the ocean, but we propose that (i) capturing complex interactions regulating the supply stoichiometry of N:P:Fe; (ii) accounting for plankton plasticity in requirements, uptake, and use of organic forms of P; and (iii) correctly simulating the global-scale patterns of N fixation and the associated drawdown of “excess” DIP in surface waters would greatly improve ESMs (). A case in point is the Mediterranean Sea, where large imbalances in the nutrient ratios are observed. The underlying cause is likely some combination of high external inputs of N and limited internal N loss (). This leads to widespread P stress with important implications for primary productivity and carbon export pathways (). Thus, the Mediterranean Sea illustrates the importance of understanding the causes and consequences of non-Redfield nutrient cycling for ocean biogeochemistry.

Our global assessment suggests that the prediction of future declines in phytoplankton biomass and NPP carries great uncertainty as models struggle to capture the impact of DIP availability on ecosystem processes in oligotrophic ocean regions. As such, we need to improve the description of phytoplankton biomass accumulation and growth based on carefully assessed nutrient concentrations. These modifications will allow us to accurately capture current and future relationships between nutrient availability, ecosystem functioning, and marine biogeochemical cycles.

MATERIALS AND METHODS

DIP observations

All abbreviations are listed in table S2. To quantify the global variation in DIP across oligotrophic environments, we compiled high-sensitivity P measurements ([DIP]obs) from 42 major cruises covering all oligotrophic regions (Table 1). The dataset comprised a total of 50,591 samples including 41,747 samples from the upper 30 m. The data are available at www.bco-dmo.org/dataset/764704. The protocols for these measurements varied including the use of LWCC, MAGIC, or SPE (Table 1). In addition, there was variation in whether the samples were prefiltered, frozen, or analyzed onboard ship, and the nature of the blank (Table 1).

We interpolated the observations by weight averaging using Ocean Data View (ODV) (). We also calculated the median and mean [DIP]obs in 1° × 1° grids from the top 30 m to create a “climatology” of 1138 unique locations for comparison with other datasets.

DIP predictions

We retrieved predicted DIP concentrations ([DIP]predict) at 1° × 1° resolution from two climatologies [World Ocean Atlas 2013 () and GLODAPv2.2016b ()] and 10-year mean outputs (1990s) from nine ESMs as part of the CMIP5 model intercomparison project (table S1). Here, the data were reformatted to fit a 1° × 1° geospatial grid.

Chlorophyll concentrations and NPP rates

We retrieved monthly MODIS satellite–estimated chlorophyll concentrations (Chlsat) and integrated NPP estimates (NPPsat) at 1/12° × 1/12° resolution. NPPsat was estimated using the Carbon-based Production Model-2 (CbPM2) (). We then matched [DIP]obs values and satellite values based on latitude, longitude, and month, resulting in 38,653 matches. Several cruises estimated in situ NPP (NPPobs) (Table 1) and were included for comparison (n = 2912). We retrieved 10-year mean outputs (1990s) of surface chlorophyll (Chlpredict) and integrated NPP (NPPpredict) from ESMs (table S1).

Vertical velocity analysis

We retrieved 10-year mean outputs (1990s) of the 1° resolution vertical velocities from CESM v1.2 and estimated the mean across the euphotic zone (surface to 105 m, top 11 layers). Then, we compared the predicted vertical velocities to [DIP]obs using Pearson correlation.

CESM analysis

A modified version of the CESM was used for the simulations here (fig. S5), which includes variable N:P in phytoplankton plus sinking particulate organic matter export, along with other modifications to the model leading up to the release of CESM v2 (www.cesm.ucar.edu) (). The model includes multiple phytoplankton functional groups and multiple potentially growth-limiting nutrients (N, P, Fe, and Si) and has been used in CESM climate simulations (). The simulations ran for 300 years, and we averaged model output from the last 20 years for our analysis. We used a 20-year repeating cycle of atmospheric forcings from the NCEP reanalysis dataset corresponding to years 1990–2009. A comparison between the old and new versions of the CESM reveals that the CESMv2 is more accurate in predicting low DIP concentrations across the Pacific Ocean (fig. S5). However, it may underestimate DIP concentrations in the Indian Ocean. As our DIPobs data coverage is low in the Indian Ocean, it is unclear what the correct level should be.

Atmospheric iron inputs in the control simulation include iron from mineral dust and combustion aerosol sources. In one experiment, we modified the deposition to the southern hemisphere gyres by applying a scaling factor to each gyre sufficient to bring the mean Fe deposition rate up to the level seen in the North Pacific subtropical gyre. The required factors to equalize atmospheric deposition were ×3.0 for the South Pacific, ×2.5 for the South Atlantic, and ×3.1 for the South Indian gyre. In a second experiment, the Fe deposition within each gyre was replaced by the mean deposition for that gyre in the control simulation. This Fe deposition manipulation removed the spatial patterns in external Fe inputs within the gyres but keeps the total Fe input the same.

Data analysis

All analyses were performed in MATLAB unless otherwise noted. The median and mean DIP values were estimated for 20 bins at 10-nM intervals (0 to 200 nM). We correlated the median relationship for DIP versus chlorophyll or NPP for observations against the same relationship in all ESMs using Spearman correlation, and no comparisons were significant (P > 0.05).