Sea-surface dimethylsulfide (DMS) concentration from satellite data at global and regional scales

Abstract

The marine biogenic gas dimethylsulfide (DMS) modulates climate by enhancing aerosol light scattering and seeding cloud formation. However, the lack of time- and space-resolved estimates of DMS concentration and emission hampers the assessment of its climatic effects. Here we present${\text{DMS}}_{\text{SAT}}$, a new remote sensing algorithm that relies on macroecological relationships between DMS, its phytoplanktonic precursor dimethylsulfoniopropionate (DMSPt) and plankton light exposure. In the first step, planktonic DMSPt is estimated from satellite-retrieved chlorophyll $a$ and the light penetration regime as described in a previous study (Galí et al., 2015). In the second step, DMS is estimated as a function of DMSPt and photosynthetically available radiation (PAR) at the sea surface with an equation of the form: ${\log }_{\mathrm{10}}\text{DMS}=\mathit{\alpha }+\mathit{\beta }{\log }_{\mathrm{10}}\text{DMSPt}+\mathit{\gamma }\text{PAR}$. The two-step${\text{DMS}}_{\text{SAT}}$ algorithm is computationally light and can be optimized for global and regional scales. Validation at the global scale indicates that ${\text{DMS}}_{\text{SAT}}$ has better skill than previous algorithms and reproduces the main climatological features of DMS seasonality across contrasting biomes. The main shortcomings of the global-scale optimized algorithm are related to (i) regional biases in remotely sensed chlorophyll (which cause underestimation of DMS in the Southern Ocean) and (ii) the inability to reproduce high DMS $/$ DMSPt ratios in late summer and fall in specific regions (which suggests the need to account for additional DMS drivers). Our work also highlights the shortcomings of interpolated DMS climatologies, caused by sparse and biased in situ sampling. Time series derived from MODIS-Aqua in the subpolar North Atlantic between 2003 and 2016 show wide interannual variability in the magnitude and timing of the annual DMS peak(s), demonstrating the need to move beyond the classical climatological view. By providing synoptic time series of DMS emission, ${\text{DMS}}_{\text{SAT}}$ can leverage atmospheric chemistry and climate models and advance our understanding of plankton–aerosol–cloud interactions in the context of global change.

Introduction

Ocean-emitted gases and particles control the number, size distribution and composition of aerosols in remote oceanic areas (Brooks and Thornton, 2018). These aerosols scatter sunlight and can act as cloud condensation nuclei that alter the radiative properties of clouds, both microscopic (cloud droplet number concentration and effective radius) and macroscopic (cloud abundance, albedo and lifetime). Interactions between natural aerosols and clouds are a major source of uncertainty in climate projections, confounding the calculation of natural and anthropogenic radiative forcing and the attribution of anthropogenic climate change (Carslaw et al., 2013). Therefore, there is an urgent need to better understand and model the oceanic sources of aerosols, and to better resolve their variations at relevant spatial and temporal scales, from weekly to seasonal to interannual.

The gas dimethylsulfide (DMS) is produced by marine microbial food webs in the sunlit layer of the ocean. With its emission currently estimated at 28 $\mathrm{Tg}\mathrm{S}{\mathrm{yr}}^{-\mathrm{1}}$, it contributes about 70 % of natural sulfur emissions to the global atmosphere and a major portion of the marine emission of organic volatiles (Carpenter et al., 2012; Schlesinger and Bernhardt, 2013; Simó, 2011). The cloud-seeding activity of DMS and its potential role in climate regulation were first postulated 3 decades ago (Charlson et al., 1987; Shaw, 1983). The so-called CLAW hypothesis (Charlson et al., 1987) proposed that a negative feedback could operate between marine phytoplankton, DMS emission and cloud albedo, potentially regulating the Earth's climate. Posterior research showed that the mechanisms behind the potential loop are far more complex than initially envisaged. This, and the estimated low sensitivity of each step of the feedback to changes in its forcing factors, led Quinn and Bates (2011) to refute the CLAW hypothesis.

Nevertheless, atmospheric studies powered by new analytical techniques (Kulmala et al., 2014) and modeling have shown instances where marine DMS controls ultrafine aerosol particle formation in the Arctic (Leaitch et al., 2013; Park et al., 2017), temperate North Atlantic (Sanchez et al., 2018), Antarctica (Yu and Luo, 2010) and the tropical South Pacific atmospheres (Modini et al., 2009). Moreover, Quinn et al. (2017) recently reported that DMS-derived aerosols dominate cloud condensation nuclei populations over most of the global ocean. Hence, the influence of DMS on marine stratiform cloud albedo remains in the spotlight (Brooks and Thornton, 2018), and the occurrence of a “seasonal CLAW” in remote marine atmospheres is becoming increasingly conceivable (Levasseur, 2013; Vallina and Simó, 2007a).

DMS is produced by marine microbial food webs through a complex network of biological interactions and chemical processes (Simó, 2004). Its primary source is the enzymatic cleavage of dimethylsulfoniopropionate (DMSP), a multifunctional osmolyte that accumulates at high (hundred millimolar) intracellular concentrations in some phytoplankton, especially haptophytes, dinoflagellates and some picoeukaryotes (Stefels et al., 2007). DMSP cleavage is catalyzed by a wide diversity of enzymes, called DMSP lyases, produced by some phytoplankton (Alcolombri et al., 2015) and bacteria (Curson et al., 2011). Breakage of phytoplankton cells through zooplankton grazing, viral attack and autolysis releases DMSP to the algal boundary layer and the dissolved phase and enhances DMS production (Simó, 2004; Stefels et al., 2007). Another process that contributes to DMS production is the diffusive release of DMS from phytoplankton cells, which proceeds almost instantaneously after intracellular DMSP cleavage by DMSP lyases or by photochemically produced radicals (Lavoie et al., 2015; Mopper et al., 2015). Once in seawater, DMS is removed by biotic and abiotic processes. DMS budgets in the upper mixed layer (UML) indicate that, on average, about 90 % of dissolved DMS is consumed by bacterial oxidation and UV-driven photolysis, and only 10 % is emitted to the atmosphere through turbulent diffusion (Galí and Simó, 2015).

Seawater DMS concentration controls the emission flux because the oceanic UML is largely supersaturated with respect to the atmosphere. DMS concentration in the UML is regulated by a subtle dynamic equilibrium between production and consumption processes with a typical timescale of less than 4 days (Galí and Simó, 2015). Over the seasonal cycle, DMS concentration varies mainly in response to the phenology and ecological succession of microbial species and their interplay with physical forcing factors, particularly solar exposure and nutrient supply, which are in turn regulated by vertical mixing (Galí and Simó, 2015; Lizotte et al., 2012; Simó and Pedrós-Alió, 1999). For instance, diatom blooms, typical of nutrient replete conditions at high latitudes, are characterized by low DMSP concentration per unit biomass and low DMSP-to-DMS conversion yield (Lizotte et al., 2012). The opposite is true for microbial communities typical of stratified, nutrient-depleted and highly irradiated surface waters, both at low and high latitudes (Galí and Simó, 2010; Lizotte et al., 2012). Under these conditions, two main factors act synergistically to increase DMS concentration (Galí and Simó, 2015; Vallina et al., 2008): the higher contribution of DMSP-rich species to total phytoplankton biomass (Galí et al., 2015; Stefels et al., 2007); and the higher DMSP-to-DMS conversion yield at the microbial community level, possibly caused by the effects of nutrient and irradiance stress (Galí et al., 2013; Stefels, 2000; Sunda et al., 2002, 2007; Vallina et al., 2008). As a result, similar DMS concentrations may occur in waters that differ by 1 or 2 orders of magnitude in phytoplankton biomass (Lizotte et al., 2012), and DMS tends to peak in summer across polar to tropical latitudes, lagging the annual chlorophyll peak by some months in the subtropical gyres. The mismatch between phytoplankton biomass, DMSP and DMS, termed the DMS summer paradox (Simó and Pedrós-Alió, 1999), is an essential feature that biogeochemical models strive to reproduce with mixed success (Le Clainche et al., 2010).

With nearly 50 000 DMS measurements taken between 1972 and 2010, the global sea-surface DMS database (https://saga.pmel.noaa.gov/dms/, last access: 13 April 2017) is a valuable resource for model development and validation. Gridded monthly climatologies (Kettle et al., 1999; Lana et al., 2011) calculated from this dataset are the standard DMS product used as input to atmospheric chemistry and climate models, hence emphasizing the seasonal climatological view (Mahajan et al., 2015; McCoy et al., 2015). At the other end, the climatic role of DMS is often evaluated through extreme sensitivity tests that examine the response of Earth system models to order-of-magnitude perturbations of DMS emission (e.g., Grandey and Wang, 2015). In comparison, contemporaneous decadal-scale DMS variability has received less attention. This gap can be filled using empirical remote sensing algorithms, a handful of which have been developed since the early 2000s (Tesdal et al., 2016; see also the pioneering works of Jodwalis and Benner, 1995, and Thompson et al., 1990). Interestingly, large discrepancies exist among global DMS fields estimated with interpolated climatologies, empirical algorithms or prognostic biogeochemical models (Tesdal et al., 2016). Although it is tempting to attribute these discrepancies to the poor skill of the models, they may also arise from issues in the calculation of the climatology.

Here we present ${\text{DMS}}_{\text{SAT}}$, a new remote sensing algorithm for DMS that proceeds in two steps: (i) estimation of the concentration of the phytoplanktonic DMS precursor, total dimethylsulfoniopropionate (DMSPt), from remotely sensed chlorophyll and light penetration, and from climatological mixed layer depth (MLD); (ii) estimation of DMS concentration from DMSPt and solar irradiance. This two-step empirical algorithm reflects, with a simplified formulation, the mechanistic understanding of oceanic sulfur cycling described in the previous paragraphs. The DMSPt sub-algorithm was presented by Galí et al. (2015) and is briefly described in Appendix A. Thus, here we focus on the second step, based on the nonlinear relationship between DMS, DMSPt and photosynthetically available radiation (PAR) at the sea surface. We implement our algorithm to produce a global DMS climatology, which we compare to the latest version of the interpolated DMS climatology (Lana et al., 2011; L11 hereafter) and to climatologies derived from other remote sensing algorithms that follow similar rationales (Simó and Dachs, 2002; Vallina and Simó, 2007b). Finally, we implement our algorithm using 14 years of MODIS-Aqua satellite data in the subtropical and the subpolar North Atlantic and in the northeast Pacific to illustrate and understand interannual DMS variability.

Methods

Datasets used for algorithm development and validation

In situ concentrations of DMS and DMSPt (nM) and chlorophyll $a$ (Chl, $\mathrm{mg}{\mathrm{m}}^{-\mathrm{3}}$), accompanied by ancillary data (bottom depth, temperature, salinity, wind speed), were downloaded from the global sea-surface DMS database. The latter was complemented with additional in situ datasets recently obtained by the authors' teams. After quality control, the database had 41 304, 3700 and 9182 measurements for in situ DMS, DMSPt and Chl, respectively, with 3637 DMS–DMSPt and 8141 DMS–Chl pairs.The in situ database was extended with geophysical and biogeochemical parameters, including satellite retrievals collocated in time and space (“matchups”) and gridded climatological datasets, following Galí et al. (2015; see below). Detailed information regarding data sources, quality control and processing can be found in the Supplement and in Tables S1–S3.

We performed satellite matchups using SeaWiFS (1997–2010) and MODIS-Aqua (2003–2012) retrievals of remotely sensed Chl ($\mathrm{mg}{\mathrm{m}}^{-\mathrm{3}}$), vertical attenuation coefficient at 490 $\mathrm{nm}$ (Kd490, ${\mathrm{m}}^{-\mathrm{1}}$), particulate inorganic carbon (PIC, $\mathrm{mol}{\mathrm{m}}^{-\mathrm{3}}$) and daily photosynthetically available radiation at the sea surface (PAR, $\mathrm{mol}\mathrm{photons}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{d}}^{-\mathrm{1}}$). To maximize the amount of available matchups, and after verifying the consistency between the two sensors, we produced merged variables by averaging SeaWiFS and MODIS-Aqua matchups. We employed a hierarchical search procedure whereby the matchup criteria were progressively relaxed from 1 to 8 days and from single-pixel to $\mathrm{5}\times \mathrm{5}$ pixel bins (Sect. S2 in the Supplement). These merged variables are hereafter designated with the SAT subscript (e.g., ${\text{Chl}}_{\text{SAT}}$). Daily sea-surface temperature (${\text{SST}}_{\text{SAT}}$, ${}^{\circ }$C) from the AVHRR sensors was also matched to the database.

Relationship between DMS, DMSPt and PAR across oceanic biomes and data availability. (a) DMS vs. DMSPt. (b) DMS $/$ DMSPt ratio vs. mean daily irradiance (PAR) at the sea surface. (c) Longhurst biogeochemical provinces and biomes. (d) In situ DMS data counts in $\mathrm{5}{}^{\circ }\times \mathrm{5}{}^{\circ }$ latitude–longitude bins; stippling indicates bins with measurements available for 3 or more months. In (a, b) small grey dots represent individual data points and large colored dots represent the median in a given Longhurst biogeochemical province and month and the corresponding interquartile ranges. Province–month medians are colored by biome following the map in (c), which also shows the amount of DMS–DMSPt–PAR measurements available in each biome. In (a–c), the $R^{\mathrm{2}}$ and $n$ (data counts) outside and inside parentheses correspond to non-binned data and province–month (MLongh) binned data, respectively. Regression lines in (a, b), calculated with MLongh binned data, are only illustrative.

[Figure omitted. See PDF]

The database was further extended with monthly climatological data: daily PAR from SeaWiFS (1997–2010 average); mixed layer depth (MLD, in meters) from the monthly MIMOC climatology (Schmidtko et al., 2013); bottom depth from the General Bathymetric Chart of the Oceans (GEBCO08); and sea-surface nitrate and phosphate concentrations ($\mathrm{\mu }\mathrm{M}$) from the World Ocean Atlas 2009 (WOA09). Nutricline depths were calculated from WOA09 vertical profiles as the depth where nitrate and phosphate first exceeded 1 and 0.4 $\mathrm{\mu }\mathrm{M}$, respectively. Nutricline depth estimations were robust to changes of $\pm$50 % in these concentration thresholds.

The mean daily PAR in the upper mixed layer (${\text{PAR}}_{\text{MLD}}$) was calculated as $$\begin{array}{cc}{\displaystyle }{\text{PAR}}_{\text{MLD}}=& {\displaystyle }[{\text{PAR}}_{\text{SAT}}/(\text{Kd}{\mathrm{490}}_{\text{SAT}}\text{MLD}\left)\right]\\ & {\displaystyle }& {\displaystyle }[\mathrm{1}-\exp (\text{Kd}{\mathrm{490}}_{\text{SAT}}\text{MLD}\left)\right].\end{array}$$ When satellite matchups were not available (before September 1997), we used climatological PAR from SeaWiFS (1997–2010 average) in order to increase the temporal coverage of the ${\text{PAR}}_{\text{SAT}}$ and ${\text{PAR}}_{\text{MLD}}$ variables. Statistical analyses done with climatological or matchup ${\text{PAR}}_{\text{SAT}}$ gave very similar results. This procedure was not followed with other variables (Chl, PIC, Kd490) that show wider interannual variations.

Statistical analyses and data binning schemes

All statistical analyses were conducted using (i) non-binned data; (ii) data binned by month and $\mathrm{5}{}^{\circ }\times \mathrm{5}{}^{\circ }$ latitude–longitude bins ($M\mathrm{5}\times \mathrm{5}$); and (iii) data binned by month and the 56 Longhurst biogeochemical provinces (designated MLongh; Longhurst, 2010). Data binning eliminates low-frequency variation (“noise”) below a given space or timescale. MLongh binned data were further aggregated into six biomes: two Polar biomes (Arctic and Antarctic), two mid-latitude Westerlies biomes (Northern and Southern hemispheres), one Trades biome (tropical latitudes) and one global Coastal biome (Fig. 1c). Variables with a right-skewed approximate lognormal distribution entered statistical analyses after ${\log }_{\mathrm{10}}$ transformation: DMS, DMSPt, DMS $/$ DMSPt ratio, Chl, nitrate and phosphate concentrations. To further account for non-normality, we conducted statistical explorations using both bin means and bin medians.

To develop the DMS algorithm we analyzed the relationship between DMS, the DMS $/$ DMSPt ratio and environmental variables (listed in Table 1). After an exploratory analysis based on the calculation of Pearson correlation coefficients (Table 1), we built several regression models where DMS was estimated as a function of in situ DMSPt concentration and additional variables (Tables 2 and S4). We added one variable at a time in order of decreasing data availability, and significant terms were selected using stepwise regression with entrance and removal $p$ values set at 0.001 and 0.005, respectively. The logic for adding one variable at a time, rather than building a single initial model with all the predictor variables, is that the size ($N$) of the data subset used for model fitting decreases rapidly when variables with sparse coverage are combined. Each set of initial predictors was tested across the three degrees of data binning described above and three degrees of model complexity: linear without interactions, linear with interactions and quadratic with interactions. This $\mathrm{3}\times \mathrm{3}$ nested structure provided a stringent test for the robustness of a given regression model. Improvements in model performance were assessed based on the increase in adjusted $r$ squared, $R_{\text{adj}}^{\mathrm{2}}$, and the decrease in root-mean-square error (RMSE) and the Akaike information criterion (AIC).

Correlation analysis (Pearson's $r$) for different data binning levels. Pearson's $r$ with $p$ value $<\mathrm{0.01}$ are not shown; values in italic are $\mathrm{0.0001}<p<\mathrm{0.01}$; n/a: not applicable. “Ratio” refers to ${\log }_{\mathrm{10}}$(DMS $/$ DMSPt). DMS, DMSPt, Chl, MLD, [${\mathrm{NO}}_{\mathrm{3}}$] and [${\mathrm{PO}}_{\mathrm{4}}$] were ${\log }_{\mathrm{10}}$ transformed; Pearson's $r$ calculated on ${\log }_{\mathrm{10}}$-transformed variables were higher than those calculated on the same non-transformed variables, and similar in magnitude to Spearman's rank correlations. Bottom and nutricline depth values are shown as positive values (deeper is bigger). See the text for other acronyms.

	Non-binned				Monthly $\mathrm{5}{}^{\circ }\times \mathrm{5}{}^{\circ }$				Monthly Longhurst
	data				binning ($M\mathrm{5}\times \mathrm{5}$)				binning (MLongh)
	DMS	$N$	Ratio	$N$	DMS	$N$	Ratio	$N$	DMS	$N$	Ratio	$N$
In situ data
DMSPt	0.65	3637	n/a	1442	0.58	308	n/a	157	0.46	122	n/a	87
Chl	0.45	to	$-$0.33	to	0.37	to	$-$0.28	to	0.34	to	$-$0.45	to
SST	$-$0.02	41 304	0.29	3637		1562	0.45	308		322	0.56	119
Salinity	$-$0.12		0.27								0.32
Wind speed	$-$0.12		$-$0.13		$-$0.12				$-$0.20		$-$0.27
Bottom depth	$-$0.19		0.10		$-$0.12		0.16				0.26
Day length	0.42		0.06		0.43				0.49
Climatological data
MLD	$-$0.37	35 505	$-$0.13	3433	$-$0.32	1474		298	$-$0.51	312	$-$0.24	119
$\left[{\mathrm{NO}}_{\mathrm{3}}\right]$	0.06	to	$-$0.19	to	0.16	to	$-$0.31	to		to	$-$0.45
$\left[{\mathrm{PO}}_{\mathrm{4}}\right]$	0.05	39 478	$-$0.15	3637	0.13	1535	$-$0.34	308		318	$-$0.32
N-cline	$-$0.14		0.30		$-$0.22		0.44				0.55
P-cline	$-$0.12		0.24		$-$0.14		0.45				0.37
Satellite match-up data
${\text{PAR}}_{\text{SAT}}$	0.32	16 411	0.35	1123	0.30	498	0.46	124	0.52	171	0.67	86
${\text{PAR}}_{\text{MLD}}$	0.12	to	0.37	to	0.15	to	0.49	to	0.36	to	0.66
${\text{Chl}}_{\text{SAT}}$	0.37	41 088	$-$0.42	3620	0.22	1539	$-$0.34	307	0.28	321	$-$0.39	119
${\mathrm{PIC}}_{\text{SAT}}$	0.24		$-$0.27		0.29		$-$0.30		0.33

Regression models were further optimized for global and regional domains using the bootstrap method followed by nonlinear optimization as described in Sect. S4. Selected models were then validated using an independent data subset composed of in situ DMS measurements and their satellite matchups (described in Sect. 3.1.3) and evaluated using a wide array of skill metrics (following Galí et al., 2015): $R^{\mathrm{2}}$, RMSE, mean absolute percentage error (MAPE), percentage bias and the slope of major axis (type II) linear regression between observations and model estimates (${\text{Slope}}_{\text{MA}}$). All analyses were carried out using MATLAB R2013b.

Algorithm implementation

The newly developed ${\text{DMS}}_{\text{SAT}}$ algorithm (Fig. 2) was implemented to produce (i) a monthly global DMS climatology and (ii) several regional time series with 8-day resolution for the period 2003–2016. Further details and data sources can be found in Sect. S5 and Table S2.

Scheme of the ${\text{DMS}}_{\text{SAT}}$ algorithm. The algorithm proceeds in two steps: the DMSPt sub-algorithm (described by Galí et al., 2015 see Appendix A) and the DMS sub-algorithm (this study). Dashed lines mark the PIC-based equation of the DMSPt sub-algorithm, which in practice is not used when gap-free satellite Chl fields are used as input.

[Figure omitted. See PDF]

Global ${\text{DMS}}_{\text{SAT}}$ fields were computed using ocean color data from SeaWiFS (1997–2010 monthly climatology, $\mathrm{1}/\mathrm{12}$${}^{\circ }$ grid), SST from AVHRR and the MIMOC monthly MLD climatology. We used SeaWiFS data to maximize the temporal overlap between the satellite-based ${\text{DMS}}_{\text{SAT}}$ climatology and the in situ data used to produce the L11 climatology, which span the period 1972–2010. Note, however, that ${\text{DMSPt}}_{\text{SAT}}$ climatologies derived from SeaWiFS and MODIS-Aqua are extremely similar (Galí et al., 2015). We established a reference ${\text{DMS}}_{\text{SAT}}$ run where ${\text{Chl}}_{\text{SAT}}$ was computed with a band-ratio algorithm (OC4-OCI standard NASA algorithm) and the euphotic layer depth (${\text{Zeu}}_{\text{SAT}}$) was computed as the 1 % penetration depth of 490 $\mathrm{nm}$ radiation (${\text{Zeu}}_{\text{SAT}}=\mathrm{4.6}/\text{Kd}\mathrm{490}$). The impact of this choice was evaluated with sensitivity tests where ${\text{Chl}}_{\text{SAT}}$ and ${\text{Zeu}}_{\text{SAT}}$ were calculated with the semi-analytical algorithms of Garver–Siegel–Maritorena (GSM; Maritorena et al., 2002) and Lee et al. (2007), respectively, which are more appropriate in optically complex waters. Observation gaps caused by low solar elevation at high latitudes in winter were left blank. Global monthly ${\text{DMS}}_{\text{SAT}}$ fields were averaged onto 1 and 5${}^{\circ }$ grids for mapping and comparison to other DMS climatologies: the interpolated L11 climatology (Lana et al., 2011), and the climatologies derived with the empirical algorithms of Simó and Dachs (2002; SD02) and Vallina and Simó (2007b; VS07). The procedures used to calculate these datasets are briefly described in Sect. 3.2.

Regional ${\text{DMS}}_{\text{SAT}}$ time series between 2003 and 2016 were computed using daily MODIS-Aqua data (4.64 $\mathrm{km}$) combined with the MIMOC MLD climatology. As done for the global implementation, we produced ${\text{DMS}}_{\text{SAT}}$ fields using both band-ratio and semi-analytical bio-optical products. We also performed a test comparing ${\text{DMSPt}}_{\text{SAT}}$ obtained with the MIMOC MLD climatology vs. model-derived MLD time series, showing little ${\text{DMSPt}}_{\text{SAT}}$ sensitivity (Fig. S1 in the Supplement). Since non-climatological satellite data contain gaps caused by cloudiness, we applied a binning and gap-filling procedure to obtain full coverage, such that the final regional time series had a resolution of 8 days and 27.8 $\mathrm{km}$. We produced ${\text{DMS}}_{\text{SAT}}$ time series for the Bermuda Atlantic Time-series Study (BATS) station (31${}^{\circ }$40${}^{\prime }$ N, 64${}^{\circ }$10${}^{\prime }$ W) and for the Northern Hemisphere at latitudes $>\mathrm{45}$${}^{\circ }$ N. The latter dataset was then sampled at selected North Atlantic sites and at the Ocean Station P (OSP) in the NE Pacific (50${}^{\circ }$ N, 145${}^{\circ }$ W). Satellite time series were compared to the L11 climatology and to in situ DMS and DMSPt. These in situ data, kindly provided by the BATS (Levine et al., 2016) and OSP (https://www.waterproperties.ca/linep/, last access: 29 November 2017) teams, were not used in algorithm development.

Results

Development and validation of the DMS sub-algorithm

Statistical exploration

We analyzed the correlation between potential predictor variables and ${\log }_{\mathrm{10}}\left(\text{DMS}\right)$ or ${\log }_{\mathrm{10}}(\text{DMS}/\text{DMSPt})$ (Table 1). This analysis systematically showed that (i) DMSPt was the best correlate of DMS ($r=\mathrm{0.46}$ to 0.65), and (ii) surface ${\text{PAR}}_{\text{SAT}}$ and mean PAR in the mixed layer (${\text{PAR}}_{\text{MLD}}$) were the best correlates of the DMS $/$ DMSPt ratio ($r=\mathrm{0.35}$ to 0.67). These correlation patterns remained across different binning levels, suggesting that DMS can be estimated, to first order, by the concentration of its phytoplanktonic precursor compound and by the PAR-dependent enhancement of DMSPt-to-DMS conversion. It is also noteworthy that the correlation between day length and the DMS $/$ DMSPt ratio was weak or non-significant. This supports the causal relationship between PAR and the DMS $/$ DMSPt ratio and discards other factors that might follow synchronous seasonal cycles.

Summary of fitted model coefficients and goodness-of-fit statistics. Different sets of coefficients were obtained by fitting the model ${\log }_{\mathrm{10}}\text{DMS}=\mathit{\alpha }+\mathit{\beta }$ ${\log }_{\mathrm{10}}\text{DMSPt}+\mathit{\gamma }\text{PAR}$ to observed DMS, DMSPt and ${\text{PAR}}_{\text{SAT}}$ after applying different binning schemes. Equations (2f) and (2h) were derived using a different optimization procedure, applied to the global MLongh binned dataset and to the BATS station local dataset. The models implemented to calculate a global DMS climatology are denoted in italic, regional or local time series in bold.

	Binning	Bin metric	Eq.	$\mathit{\alpha }$	$\mathit{\beta }$	$\mathit{\gamma }$	$R_{\text{adj}}^{\mathrm{2}}$	RMSE	SlopeMA	$N$
Regression,	Non-binned		(2a)	$-\mathrm{1.213}\pm \mathrm{0.028}$	$\mathrm{0.672}\pm \mathrm{0.012}$	$\mathrm{0.0136}\pm \mathrm{0.0006}$	0.50	0.35	0.62	3620
global scale	$M\mathrm{5}\times \mathrm{5}$	Mean	(2b)	$-\mathrm{1.154}\pm \mathrm{0.083}$	$\mathrm{0.669}\pm \mathrm{0.0371}$	$\mathrm{0.0130}\pm \mathrm{0.0015}$	0.55	0.28	0.67	307
		Median	(2c)	$-\mathrm{1.061}\pm \mathrm{0.084}$	$\mathrm{0.569}\pm \mathrm{0.039}$	$\mathrm{0.0130}\pm \mathrm{0.0015}$	0.46	0.28	0.58
	MLongh	Mean	(2d)	$-\mathrm{1.061}\pm \mathrm{0.115}$	$\mathrm{0.583}\pm \mathrm{0.054}$	$\mathrm{0.0155}\pm \mathrm{0.0019}$	0.57	0.24	0.70	118
		Median	(2e)	$-\mathrm{1.018}\pm \mathrm{0.100}$	$\mathrm{0.452}\pm \mathrm{0.050}$	$\mathrm{0.0163}\pm \mathrm{0.0016}$	0.57	0.21	0.69
Optimization, global scale	MLongh	Median	(2f)	$-$1.237	0.578	0.0180	0.56	0.22	0.87	118
Regression, regional scale($\mathbf{>}$ 45${}^{\circ }$ N)	$\mathbf{M}\mathbf{5}\mathbf{\times }\mathbf{5}$	Mean	(2g)	$\mathbf{-}\mathbf{1.283}\mathbf{\pm }\mathbf{0.154}$	$\mathbf{0.670}\mathbf{\pm }\mathbf{0.097}$	$\mathbf{0.0186}\mathbf{\pm }\mathbf{0.011}$	0.68	0.28	0.80	87
Optimization, local scale (BATS)	Non-binned		(2h)	$\mathbf{-}$0.898	0.316	0.0214	0.44	0.26	0.66	35

Guided by the correlation patterns, we established a base regression model expressed by the equation: $${\log }_{\mathrm{10}}\text{DMS}=\mathit{\alpha }+\mathit{\beta }{\log }_{\mathrm{10}}\text{DMSPt}+\mathit{\gamma }\text{PAR}.$$

This model explained between 50 and 57 % of ${\log }_{\mathrm{10}}\left(\text{DMS}\right)$ variance with an increasing level of data binning, and the corresponding RMSE ranged between 0.35 and 0.21 (Table 2).

We assessed whether the base model could be significantly improved by adding predictor variables and/or increasing model complexity. We started by adding a variable $X$ to a linear model without interactions of the form ${\log }_{\mathrm{10}}\text{DMS}=\mathit{\alpha }+\mathit{\beta }{\log }_{\mathrm{10}}\text{DMSPt}+\mathit{\gamma }\text{PAR}+\mathit{\delta }X$. Variable $X$ was chosen among those showing higher correlations to either DMS or the DMS $/$ DMSPt ratio: SST, nitrate concentration, nitracline depth, salinity, wind speed and ${\mathrm{PIC}}_{\text{SAT}}$ (Table 1). Although this analysis led us to discard additional predictor variables, its results are briefly described below and compiled in Table S4 for the sake of completeness. With non-binned data, all the additional variables entered regression models with significant coefficients, but only salinity, wind speed and ${\mathrm{PIC}}_{\text{SAT}}$ produced significant decreases in RMSE and AIC. With MLongh binned data, only SST and ${\mathrm{PIC}}_{\text{SAT}}$ entered with significant coefficients. Yet none of the additional variables improved simultaneously the $R^{\mathrm{2}}$ adj, RMSE and AIC skill metrics with respect to the base model. Increasing model complexity through addition of interaction and quadratic terms, or by adding a fourth variable, generally resulted in minor improvements or erratic changes in model performance (results not shown). Invariably, DMSPt and PAR were the only variables with highly significant coefficients ($p\ll {\mathrm{10}}^{-\mathrm{10}}$) regardless of the binning scheme or the inclusion of additional variables.

As a corollary, the use of ${\text{PAR}}_{\text{MLD}}$ instead of ${\text{PAR}}_{\text{SAT}}$ slightly degraded model skill (Table S4). Although ${\text{PAR}}_{\text{MLD}}$ is a priori a more realistic metric of light exposure, it is possible that the use of climatological MLD degraded the ${\text{PAR}}_{\text{MLD}}$ estimates. Another potential explanation is the episodic nature of oceanic vertical mixing, which requires the distinction between the actively mixing layer, which is defined by higher turbulence than in the ocean interior, and the mixed layer, here termed MLD and defined by the vertical homogeneity of regular temperature–salinity profiles (Brainerd and Gregg, 1995; Sutherland et al., 2014). This distinction implies that, on occasion, mean light exposure at the sea surface is better approximated by surface ${\text{PAR}}_{\text{SAT}}$ than by ${\text{PAR}}_{\text{MLD}}$. After these considerations we excluded ${\text{PAR}}_{\text{MLD}}$ from our algorithm, and focused on optimizing Eq. (2).

Finally, note that ${\text{PAR}}_{\text{SAT}}$ is used in our algorithm both as a direct driver of DMS cycling processes and as a proxy for UVR-driven processes (e.g., Archer et al., 2010; Galí et al., 2013; Royer et al., 2016). For obvious astronomical reasons, incident PAR and UVR are strongly correlated. Attenuation of PAR and UVR in the atmosphere, first, and in seawater, afterwards, also covary (Kirk, 2011), so that plankton UVR exposure is to first order well correlated to PAR.

Implications of the model structure

Here we analyze the physical and biogeochemical meaning of Eq. (2) coefficients in view of their optimization for diagnostic purposes (3.1.3).

• The intercept ($\mathit{\alpha }$) acts to adjust the magnitude of DMS concentrations by a fixed proportion everywhere; e.g., increasing $\mathit{\alpha }$ by ${\log }_{\mathrm{10}}\left(\mathrm{1.10}\right)$ would raise estimated DMS concentrations by 10 % globally.

• The ${\log }_{\mathrm{10}}\text{DMSPt}$ coefficient ($\mathit{\beta }$) expresses the nonlinear relationship between DMS and DMSPt. Fitted $\mathit{\beta }$ is smaller than 1 regardless of the binning applied. For a constant PAR, this implies that DMS increases more slowly than DMSPt (Fig. 3a) and that the DMS $/$ DMSPt ratio decreases nonlinearly and approaches an horizontal asymptote as DMSPt increases (Fig. 3a). This behavior implicitly represents the change in DMSPt-to-DMS conversion efficiency depending on the biomass and structure of the microbial plankton community.

• The PAR coefficient ($\mathit{\gamma }$) expresses the DMSPt-independent modulation of DMS concentration. Fitted $\mathit{\gamma }$ is larger than 0, meaning that the PAR sensitivity increases exponentially with PAR, regardless of DMSPt (Fig. 3b). Incident PAR is positively correlated to shallow mixing at the global scale, and to seawater transparency at low latitudes (but not at high latitudes). Thus, nonlinear PAR sensitivity possibly embodies the enhancement of sunlight exposure by shallow mixing and deeper PAR (and UVR) penetration. In summary, $\mathit{\gamma }$ represents stress-driven DMS production.

Relationship between DMS, DMSPt and PAR as represented in the ${\text{DMS}}_{\text{SAT}}$ algorithm. (a) DMS vs. DMSPt and (b) DMS $/$ DMSPt ratio vs. DMSPt, as a function of PAR. Colored circles represent the medians of in situ data binned by Longhurst province and month, and lines correspond to the model predictions. PAR levels are indicated in the color bar.

[Figure omitted. See PDF]

In biogeochemical terms, Eq. (2) implies maximal DMS $/$ DMSPt ratios when/where low DMSPt and high PAR co-occur. In biogeographic terms (Fig. 1), highest DMS $/$ DMSPt ratios are found in oligotrophic areas of the Trades biome, where low DMSPt concentrations prevail ($<\mathrm{20}$ $\mathrm{nM}$). Low DMSPt concentrations are also found in winter at high latitudes in deeply mixed waters, but the corresponding low irradiance results in DMS $/$ DMSPt $<$ 0.05. At the high DMSPt concentrations that occur at high latitudes in summer ($>\mathrm{100}$ $\mathrm{nM}$), the DMS $/$ DMSPt ratio is generally $<\mathrm{0.1}$.

As shown in Table 2, Eq. (2) coefficients change systematically as the binning spatial scale increases. To further explore the interrelationship between the model coefficients, we used the bootstrap method to produce ${\mathrm{10}}^{\mathrm{5}}$ sets of regression coefficients for the MLongh dataset. The scatter plots between $\mathit{\alpha }$, $\mathit{\beta }$ and $\mathit{\gamma }$ resulting from the ${\mathrm{10}}^{\mathrm{5}}$ bootstrapped regressions confirmed that covariation between the coefficients is non-random, such that fitted $\mathit{\beta }$ and $\mathit{\gamma }$ are negatively correlated to $\mathit{\alpha }$ (Fig. S2). These trade-offs should be kept in mind when optimizing our model for global or regional implementation.

Optimization and validation

By definition, least squares regression minimizes the RMSE, but it has been shown that regression models derived in this way do not necessarily have the best skill (Jolliff et al., 2009). Therefore, we devised an alternative nonlinear optimization procedure (Sect. S4). To obtain realistic solutions, we constrained the optimized coefficients to the confidence intervals derived from the bootstrapped regressions for the MLongh dataset (Fig. S2). The resulting optimal model (Eq. 2f) had higher DMSPt ($\mathit{\beta }$) and PAR ($\mathit{\gamma }$) coefficients and a smaller intercept ($\mathit{\alpha }$), and moved the modeled DMS concentration closer to the $\mathrm{1}:\mathrm{1}$ agreement line without degrading either RMSE or $R^{\mathrm{2}}$ (Table 2).

We validated the different versions of Eq. (2) (Table 2) by comparing ${\text{DMS}}_{\text{SAT}}$ against in situ DMS using an independent subset of the database. Since the complete DMS algorithm proceeds in two steps (Fig. 2), its validation must take into account uncertainty in variables used as input to the ${\text{DMSPt}}_{\text{SAT}}$ sub-algorithm. Galí et al. (2015) showed that, apart from the inherent algorithm uncertainty, most uncertainty in ${\text{DMSPt}}_{\text{SAT}}$ (RMSE $\le \mathrm{0.3}$ in ${\log }_{\mathrm{10}}$ space) results from error in ${\text{Chl}}_{\text{SAT}}$. Thus, the validation subset was defined according to three criteria: (i) satellite matchup data used as input to the algorithm (${\text{Chl}}_{\text{SAT}}$, $K_{d,\mathrm{490}}$, ${\text{PAR}}_{\text{SAT}}$ and ${\text{SST}}_{\text{SAT}}$) were available; (ii) in situ DMSPt was not available – thus excluding the data used for model fitting; (iii) in situ DMS and Chl were available. We used in situ Chl concentration to constrain the uncertainty in ${\text{Chl}}_{\text{SAT}}$ used as input to the algorithm. Indeed, this procedure progressively reduced the size of the validation subset as the maximum ${\text{Chl}}_{\text{SAT}}$ error tolerated decreased. Uncertainty arising from ${\text{PAR}}_{\text{SAT}}$ could not be assessed because the current database lacks in situ PAR measurements. Frouin et al. (2003) reported an error of $\pm$15 % ($<\mathrm{10}$ % for weekly and monthly periods), with negligible bias for ${\text{PAR}}_{\text{SAT}}$, suggesting it is a minor source of uncertainty.

Figure 4a summarizes the validation results for the best-performing regression model (Eq. 2e) and its optimized version (Eq. 2f). Supporting our assumption, ${\text{DMS}}_{\text{SAT}}$ skill metrics improved as the maximum tolerated ${\text{Chl}}_{\text{SAT}}$ RMSE decreased from 0.5 to 0.2 (Fig. 4). With ${\text{Chl}}_{\text{SAT}}$ RMSE smaller than 0.2, the statistics showed erratic behavior owing to reduced sample size. Other skill metrics (not shown in Fig. 4) showed comparable trends.

Algorithm validation results constrained by the uncertainty in satellite-retrieved Chl. (a) Equation (2e), derived from regular multiple regression; (b) Eq. (2f), obtained through an optimization procedure. The scatter plots compare non-binned data and model estimates, color-coded depending on the maximum tolerated error in ${\text{Chl}}_{\text{SAT}}$ with respect to Chl in situ, as shown in the $x$ axis of the center plots. The center plots show the performance of the DMS algorithm for increasing error in ${\text{Chl}}_{\text{SAT}}$, evaluated with different skill metrics: the ${\log }_{\mathrm{10}}$ space $R^{\mathrm{2}}$ and RMSE (left $y$ axis) and the linear space MAPE (right $y$ axis). $N$ increases from 86 to 1293 as the tolerated ${\text{Chl}}_{\text{SAT}}$ error increases.

[Figure omitted. See PDF]

Relative deviation (“bias”) of the global-scale optimized ${\text{DMS}}_{\text{SAT}}$ algorithm across ocean biomes when evaluated against the in situ database (a and b) and against the L11 monthly climatology (c). Comparison with the in situ database was done using either (a) satellite matchups with ${\text{Chl}}_{\text{SAT}}$ error constrained with respect to in situ Chl or (b) all ${\text{Chl}}_{\text{SAT}}$ matchups, regardless of the availability of in situ Chl. The last column shows a qualitative overall assessment of ${\text{DMS}}_{\text{SAT}}$ bias, based on columns (a–c), classifying the most likely magnitude and sign as small $\pm$ ($\pm$10 %), moderate (10 to 40 %, either $+$ or $-$) or large ($>\mathrm{40}$ %, either $+$ or $-$).

Biome	In situ database matchups					(c) L11 gridded data		${\text{DMS}}_{\text{SAT}}$ bias assessment
	(a) ${\text{Chl}}_{\text{SAT}}$ RMSE $<\mathrm{0.3}$			(b) all matchups		bias ${\text{DMS}}_{\text{SAT}}$	area${}^{\mathrm{b}}$
	bias ${\text{DMS}}_{\text{SAT}}$	bias ${\text{Chl}}_{\text{SAT}}$	$N$${}^{\mathrm{a}}$	bias ${\text{DMS}}_{\text{SAT}}$	$N$${}^{\mathrm{a}}$
Polar N	$-$28 %	$-$7 %	54	$-$39 %	2291	$-$33 %	3.0 %	Moderate $-$
Westerlies N	$-$4 %	$-$11 %	125	$-$13 %	1872	$-$32 %	10.1 %	Small $\pm$
Trades	$-$14 %	7 %	442	6 %	3489	$-$34 %	56.6 %	Small $\pm$ to moderate $-$
Westerlies S${}^{\mathrm{c}}$	$-$3 %	$-$30 %	30	$-$27 %	1585	$-$25 %	14.2 %	Small $\pm$ to moderate $-$
Polar S	$-$59 %	$-$6 %	42	$-$47 %	170	$-$74 %	5.7 %	Large $-$
Coastal	46 %	4 %	377	3 %	5673	$-$8 %	10.4 %	Small $\pm$ to large $+$
Global ocean	11 %	2 %	1053	$-$9 %	15 080	$-$33 %	100 %	Small $\pm$ to moderate $-$
Global ocean $-$ pS $-$ C${}^{\mathrm{d}}$	$-$16 %	$-$2 %	634	$-$20 %	9249	$-$30 %	83.9 %	Moderate $-$

${}^{\mathrm{a}}$ For the database comparisons, the amount of data ($N$) available for validation in a given region depends on the total amount of measurements, the proportion of data points with satellite matchups and, for the constrained case, the fraction of matchups where ${\text{Chl}}_{\text{SAT}}$ has ${\log }_{\mathrm{10}}\text{RMSE}<\mathrm{0.3}$ compared to concurrent Chl in situ. Samples with available DMSPt measurements were excluded because they were used in model fitting and optimization.${}^{\mathrm{b}}$ For the global climatology comparison we report the percentage of ocean area, excluding pixels that could not be observed by satellites (high-latitude winter).${}^{\mathrm{c}}$ Samples where DMSPt was measured (not excluded) for this biome due to the small amount of data. Removing them would leave $N=\mathrm{13}$.${}^{\mathrm{d}}$ Global ocean minus Polar S and Coastal biomes.

The optimized model coefficients (Eq. 2f) increased $R^{\mathrm{2}}$ and reduced RMSE with respect to the regression-derived coefficients, achieving a maximal $R^{\mathrm{2}}$ of 0.53 and a minimal RMSE of 0.25 for error-free ${\text{Chl}}_{\text{SAT}}$ (${\log }_{\mathrm{10}}$-space non-binned data; Table S5). Corresponding best scores in linear space were $R^{\mathrm{2}}=\mathrm{0.24}$ and RMSE $=\mathrm{3.4}$ $\mathrm{nM}$. These linear space statistics might be interpreted as a sign of poor performance, but it should be noted that they were strongly affected by a small fraction of highly biased estimates. Removing the most biased estimates (8 % of points beyond a factor of 3 from real measurements; Fig. 4) increased the linear space $R^{\mathrm{2}}$ to 0.42–0.53 and decreased the RMSE to 1.8–2.3 $\mathrm{nM}$ across the full range of ${\text{Chl}}_{\text{SAT}}$ error, with MAPE of 34–41 % and relative bias of $-$1 to 7 %. Altogether, these statistics illustrate the good performance of the ${\text{DMS}}_{\text{SAT}}$ algorithm and the greater robustness of log-space statistics. The global-scale optimized algorithm had a mean-normalized standard deviation of 1.06–1.14 (${\log }_{\mathrm{10}}$ space), meaning that the spread of ${\text{DMS}}_{\text{SAT}}$ nearly matched that of in situ DMS concentrations.

Mean area-weighted DMS concentrations for different climatological datasets and domains. Different results for the ${\text{DMS}}_{\text{SAT}}$ algorithm were obtained using alternative approaches for retrieving chlorophyll $a$ concentration (${\text{Chl}}_{\text{SAT}}$) and the euphotic layer depth (${\text{Zeu}}_{\text{SAT}}$). Results (based on $\mathrm{1}{}^{\circ }\times \mathrm{1}{}^{\circ }$ gridded data) are shown for the global domain, the Longhurst Coastal biome and, within it, areas shallower than 200 $\mathrm{m}$ (shelves); n/a: not applicable.

DMS algorithm or data product	${\text{Chl}}_{\text{SAT}}$ product	${\text{Kd}}_{\text{SAT}}$ or ${\text{Zeu}}_{\text{SAT}}$ product	Area-weighted DMS mean (nM)
			global	coastal	coastal
				biome${}^{\mathrm{a}}$	biome shelves${}^{*}$
L11 climatology (Lana et al., 2011)	n/a	n/a	2.43	2.70	2.56
SD02 (Simó and Dachs, 2002)	OC4-CI	n/a	2.12	3.07	4.27
VS07 (Vallina and Simó, 2007b)	n/a	Kd490	2.71	2.70	2.48
${\text{DMS}}_{\text{SAT}}$ (Eq. 2f and Appendix A, this study)	OC4-CI	$\text{Zeu}=\mathrm{4.6}/\text{Kd490}$	1.63	2.49	3.06
	GSM	$\text{Zeu}=\mathrm{4.6}/\text{Kd490}$	1.58	2.42	3.09
		Zeu_Lee	1.55	2.41	3.02

${}^{*}$ The Longhurst Coastal biome represents about 10.4 % of the global ocean area, of which one third is shallower than 200 $\mathrm{m}$ (3.5 % of global area).

Table 3 summarizes a detailed analysis of algorithm bias in different biomes. When ${\text{Chl}}_{\text{SAT}}$ error was constrained to RMSE $<\mathrm{0.3}$, ${\text{DMS}}_{\text{SAT}}$ showed a global mean bias of 11 % (with a global mean bias of 2 % in ${\text{Chl}}_{\text{SAT}}$ itself). ${\text{DMS}}_{\text{SAT}}$ bias was negligible in the Westerlies biomes and larger than $\pm$10 % in the other biomes, where its sign often matched that of ${\text{Chl}}_{\text{SAT}}$ bias. When ${\text{DMS}}_{\text{SAT}}$ bias was assessed using all available matchups, irrespective of ${\text{Chl}}_{\text{SAT}}$ error, we obtained a global bias of $-$9 %. In some biomes, the sign and magnitude of ${\text{DMS}}_{\text{SAT}}$ bias changed when ${\text{Chl}}_{\text{SAT}}$ error was not constrained, indicating the influence of input satellite data. Assuming that ${\text{Chl}}_{\text{SAT}}$ matchups ($N\sim \mathrm{15}\mathrm{000}$) are a random sample of the database between 1997 and 2012 ($N\sim \mathrm{24}\mathrm{000}$), this analysis suggests that global ${\text{DMS}}_{\text{SAT}}$ bias likely ranges between $-$9 and 11 %.

Global climatologies

We implemented the global-scale optimized algorithm (Eq. 2f) using the SeaWiFS climatology. As shown in Fig. 5 and 6, DMS concentrations around $\sim \mathrm{2.5}$ $\mathrm{nM}$ prevail during the astronomical spring and summer in each hemisphere, decreasing to around 1 $\mathrm{nM}$ in fall and $<\mathrm{1}$ $\mathrm{nM}$ in winter. The seasonal cycle has wider amplitude at high latitudes and is nearly flat in the tropical oceans (see also Fig. S3). Regional enhancement of DMS concentrations occurs in some coastal and shelf areas, equatorial and eastern boundary upwellings, close to the subtropical front in austral summer (40${}^{\circ }$ S) and in the subpolar North Atlantic in boreal summer (60${}^{\circ }$ N). The global mean area-weighted ${\text{DMS}}_{\text{SAT}}$ concentration is 1.63 $\mathrm{nM}$ (median and geometric mean of 1.36 $\mathrm{nM}$). This global mean decreases by less than 5 % when semi-analytical ${\text{Chl}}_{\text{SAT}}$ and ${\text{Zeu}}_{\text{SAT}}$ products are used instead of our reference products, but larger deviations occur in the Coastal biome, particularly in shallow areas, associated with optically complex waters (Table 4).

Global ${\text{DMS}}_{\text{SAT}}$ concentration fields by season. (a) December–February (DJF); (b) March–May (MAM); (c) June–August (JJA); (d) September–November (SON). Each row contains mean latitudinal profiles for the L11 climatology and ${\text{DMS}}_{\text{SAT}}$ (left); ${\text{DMS}}_{\text{SAT}}$ concentration maps (center); and maps of the percentage difference between ${\text{DMS}}_{\text{SAT}}$ and the L11 climatology (right).

[Figure omitted. See PDF]

Hovmöller diagrams comparing climatological DMS fields (monthly 5${}^{\circ }$ latitude bins). (a) In situ database means, (b) in situ database medians, (c) L11 climatology, (d) ${\text{DMS}}_{\text{SAT}}$ algorithm, (e) SD02 algorithm and (f) VS07 algorithm. Panels (a–f) share the same color scale; the black background indicates no data. Panel (g) shows the average and range of monthly $\text{mean}/\text{median}$ ratios (shown in panels a and b) in each latitude band. Geometric means of in situ binned data (b) are similar to bin medians.

[Figure omitted. See PDF]

Comparison to the L11 climatology

The L11 DMS climatology (Lana et al., 2011) was calculated using an objective interpolation procedure similar to that used in prior climatologies (Kettle et al., 1999; Kettle and Andreae, 2000). An initial template, called first-guess field, was obtained by calculating the monthly mean DMS in each Longhurst province. The gaps were filled through temporal interpolation and, in provinces with too few documented months, the seasonal cycle was extrapolated by scaling that of neighboring provinces. Objective interpolation was then applied by searching measurements within a 555 $\mathrm{km}$ radius, weighting them inversely to the distance from a given grid point, and the resulting global fields were repeatedly smoothed.

The global mean area-weighted ${\text{DMS}}_{\text{L11}}$ concentration is 2.43 $\mathrm{nM}$ (median 1.88 $\mathrm{nM}$, geometric mean 1.83 $\mathrm{nM}$). Thus, mean ${\text{DMS}}_{\text{SAT}}$ concentration is globally 33 % lower than ${\text{DMS}}_{\text{L11}}$ (Tables 3 and 4). The largest and smallest differences occur in the southern polar and the Coastal biomes, respectively (Table 3).

The disagreement between the ${\text{DMS}}_{\text{L11}}$ and ${\text{DMS}}_{\text{SAT}}$ climatologies varies depending on the regions and on the spatial–temporal scales compared. Despite the general offset, the seasonal latitudinal profiles (zonal means) of ${\text{DMS}}_{\text{L11}}$ and ${\text{DMS}}_{\text{SAT}}$ have similar shapes, and agree very well in June through August. Comparison by means of Hovmöller diagrams (Fig. 6) shows a remarkable qualitative agreement in their month-by–latitude patterns, except for the polar austral summer. Figure 6 also reveals smaller disagreements in the Arctic Ocean in winter–spring and in the equatorial band during most of the year, with lower concentrations in ${\text{DMS}}_{\text{SAT}}$ in both cases. The most striking regional disagreements appear when ${\text{DMS}}_{\text{SAT}}$ and ${\text{DMS}}_{\text{L11}}$ are compared by means of seasonal anomaly maps (Fig. 5). The sign of the ${\text{DMS}}_{\text{SAT}}-{\text{DMS}}_{\text{L11}}$ divergence changes from positive to negative in a patchy pattern, often following the boundaries of the Longhurst biogeochemical provinces.

Comparison to the SD02 climatology

The SD02 algorithm (Simó and Dachs, 2002) was designed to estimate DMS from MLD and ${\text{Chl}}_{\text{SAT}}$ using two different equations depending on the $\text{Chl}/\text{MLD}$ ratio, $$\begin{array}{}& {\displaystyle }& {\displaystyle }\text{DMS}=-\ln \left(\text{MLD}\right)+\mathrm{5.7}\text{Chl}/\text{MLD}<\mathrm{0.02},& {\displaystyle }& {\displaystyle }\text{DMS}=\mathrm{55.8}\text{Chl}/\text{MLD}+\mathrm{0.6}\text{Chl}/\text{MLD}\ge \mathrm{0.02},\end{array}$$ such that DMS increases linearly with the $\text{Chl}/\text{MLD}$ ratio in stratified productive conditions (e.g., high latitudes in summer) and inversely with MLD in typical oligotrophic conditions.

Validation of SD02 with the same dataset used for ${\text{DMS}}_{\text{SAT}}$ indicates that it explains less variance (${\log }_{\mathrm{10}}$ $R^{\mathrm{2}}$ of 0.22–0.31) but has similar RMSE and MAPE (Table S5). Globally, ${\text{DMS}}_{\text{SD}\mathrm{02}}$ is characterized by a bimodal distribution (Fig. 7), with an area-weighted mean of 2.12 $\mathrm{nM}$, 13 % lower than ${\text{DMS}}_{\text{L11}}$ (Table 4). SD02 estimates are in good agreement with the L11 climatology at tropical and temperate latitudes. In the southern Westerlies biome, however, prevailing deep vertical mixing and low Chl cause SD02 to underestimate DMS throughout the productive season (Figs. 6 and S3). Another feature of SD02 is the high DMS concentration in northern polar latitudes through late summer and fall, caused mainly by the shallow MLD due to freshwater-driven stratification. As ${\text{DMS}}_{\text{SAT}}$, ${\text{DMS}}_{\text{SD}\mathrm{02}}$ suffers from negative bias in the Antarctic biome during the productive season (November through February).

Histograms illustrating the relationship between DMS and Chl in the global ocean. (a) Global SeaWiFS 1997–2010 ${\text{Chl}}_{\text{SAT}}$ climatology and SeaWiFS ${\text{Chl}}_{\text{SAT}}$ matchups for the in situ DMS database; (b) Global L11 and ${\text{DMS}}_{\text{SAT}}$ climatologies and in situ DMS database; (c) Global SD02 and VS07 climatologies and in situ DMS database; (d) bivariate histogram of in situ DMS vs. corresponding SeaWiFS ${\text{Chl}}_{\text{SAT}}$ matchups; (e–h) bivariate histograms of the global SeaWiFS ${\text{Chl}}_{\text{SAT}}$ climatology vs. climatological DMS from (e) L11, (d) ${\text{DMS}}_{\text{SAT}}$, (e) SD02 and (f) VS07 algorithms. In (d–h), the correlation coefficient (Pearson's $r$) between DMS and Chl in each dataset is super-imprinted in white.

[Figure omitted. See PDF]

Comparison to the VS07 climatology

Vallina and Simó (2007b) reported a globally valid linear relationship between DMS concentration and the solar radiation dose (SRD) in the upper mixed layer in the global ocean, according to the equation: $$\text{DMS}=\mathrm{0.492}+\mathrm{0.019}\text{SRD}.$$

This numerical relationship was not meant to be used as an empirical algorithm but as evidence for the emerging response of ecosystem DMS production to changes in solar radiation. However, the fact that SRD explained a large part of the variance of DMS concentration across regions and seasons prompted its use in global warming projections (Vallina et al., 2007). SRD is analogous to ${\text{PAR}}_{\text{MLD}}$ (Eq. 1), but replacing ${\text{PAR}}_{\text{SAT}}$ by total shortwave irradiance (${\text{Ed}}_{\text{SW}}$; $\mathrm{W}{\mathrm{m}}^{-\mathrm{2}}$). Here we implemented VS07 with two variations: (i) we used ${\text{Kd490}}_{\text{SAT}}$ instead of a fixed Kd (note that in phytoplankton-rich and continentally influenced waters, ${\text{Kd490}}_{\text{SAT}}$ is generally higher than the fixed Kd $=\mathrm{0.06}$ ${\mathrm{m}}^{-\mathrm{1}}$ used by Vallina and Simó, 2007b); (ii) we estimated ${\text{Ed}}_{\text{SW}}$ from ${\text{PAR}}_{\text{SAT}}$ by converting the latter to units of $\mathrm{W}{\mathrm{m}}^{-\mathrm{2}}$ (Morel and Smith, 1974) and then applying a constant ${\text{Ed}}_{\text{SW}}$ $/$ ${\text{PAR}}_{\text{SAT}}$ ratio of $\mathrm{1}/\mathrm{0.43}$ (Kirk, 2011).

VS07 shows poorer performance than ${\text{DMS}}_{\text{SAT}}$ and SD02 when validated with individual measurements (Table S5). It produces rather uniform DMS fields compared to the other climatologies (Figs. 6 and 7c), with mean area-weighted concentration of 2.71 $\mathrm{nM}$ (Table 4). VS07 performs well in the Westerlies biome, especially in the Northern Hemisphere, but invariably overestimates (underestimates) DMS in the Trades (Polar) biomes (Fig. S3).

Regional ${\text{DMS}}_{\text{SAT}}$ time series

We selected different regions to test the ${\text{DMS}}_{\text{SAT}}$ algorithm based on the following criteria: the abundance of in situ data (subpolar North Atlantic), the existence of seasonal and multiannual time series (Ocean Station P and BATS station), and the challenges posed by intra- and interannual variability of phytoplankton and DMS(P) in each region.

Subpolar Atlantic and Pacific

We used MODIS-Aqua data to produce a 14-year ${\text{DMS}}_{\text{SAT}}$ time series (and the corresponding climatology) for the Northern Hemisphere at latitudes $>\mathrm{45}$${}^{\circ }$ N. In this regional implementation we used a different set of coefficients, obtained from regression of $M\mathrm{5}\times \mathrm{5}$ binned data restricted to latitudes $>\mathrm{45}$${}^{\circ }$ N (Eq. 2g, Table 2). These regional coefficients largely corrected the negative bias observed with globally optimized coefficients. In this case, further optimization did not lead to significant improvement. We then sampled the resulting time series in some representative regions: three rectangles with an area of $\sim \mathrm{200}\mathrm{000}$ ${\mathrm{km}}^{\mathrm{2}}$ each, located along the 50–56${}^{\circ }$ N band in the North Atlantic, and the Ocean Station P (OSP, 50${}^{\circ }$ N, 145${}^{\circ }$ W) in the NE Pacific.

Figure 8 shows the seasonal cycles of ${\text{DMS}}_{\text{SAT}}$, ${\text{DMSPt}}_{\text{SAT}}$ and ${\text{Chl}}_{\text{SAT}}$ in selected North Atlantic areas: (a) the deep waters of the northwest Atlantic drift, (b) the shelf break west of Ireland, and (c) the shallow southern North Sea. We observe good agreement between the 14-year ${\text{DMS}}_{\text{SAT}}$ climatology and the L11 climatology, except in the southern North Sea where ${\text{DMS}}_{\text{SAT}}$ is too high through summer and fall. ${\text{DMS}}_{\text{SAT}}$ reproduces well the east–west variation in the temporal lag between the annual peaks of DMS and Chl, with a lag of up to 4 months in the southern North Sea. The most salient result is however the wide interannual variability of the ${\text{DMS}}_{\text{SAT}}$ seasonal cycles. Estimated ${\text{DMS}}_{\text{SAT}}$ concentrations during the productive season vary by up to 3-fold between years (see variability metrics in Fig. 8), and the annual ${\text{DMS}}_{\text{SAT}}$ peak can occur within a temporal window of 2–3 months. Although years with a major peak in spring–summer are the norm, a second peak in late summer is not unusual. Additional validation supports the good performance of ${\text{DMS}}_{\text{SAT}}$ in the subpolar North Atlantic (Fig. S4), lending credit to satellite-determined variability patterns.

Interannual ${\text{DMS}}_{\text{SAT}}$ variability in the subpolar North Atlantic. (a) Northwest Atlantic drift, (b) shelf break west of Ireland, (c) southern North Sea shelf. Top panels show individual years between 2003 and 2016 determined from 8-day MODIS-Aqua data, differentiated by color, and the mean seasonal cycle according to the L11 DMS climatology (black); colored circles mark the peak of each seasonal cycle. Bottom panels show the mean annual cycles of ${\text{Chl}}_{\text{SAT}}$, ${\text{DMSPt}}_{\text{SAT}}$, ${\text{DMS}}_{\text{SAT}}$ and the L11 DMS climatology. Each variable is divided by its maximum, shown by the number in the quotient; a common scaling factor is used for ${\text{DMS}}_{\text{SAT}}$ and ${\text{DMS}}_{\text{L11}}$. Markers on the L11 line indicate the amount of in situ data on which the L11 climatology is based in a given month: no data, i.e., month filled through interpolation (no marker); 1–9 measurements in a single year (empty circles); $\ge \mathrm{10}$ measurements in a single year (crossed circles); and $\ge \mathrm{10}$ measurements distributed in different years (filled black circles). Red polygons on the map show the 3 selected areas and the larger region used for the validation scatter plots (see Fig. S4).

[Figure omitted. See PDF]

We used the same MODIS-Aqua dataset to analyze the mean seasonal cycle and the interannual variability at OSP (Fig. 9a–c), where DMS has been measured around February, June and August since 1996. ${\text{DMS}}_{\text{SAT}}$ captures in situ DMS concentrations in February and June but suffers from low bias in August. Examination of August measurements during the 2005–2016 period suggests the existence of two regimes: 8 years with in situ DMS of $\mathrm{6.6}\pm \mathrm{1.1}$ $\mathrm{nM}$, about twice as high as ${\text{DMS}}_{\text{SAT}}$, and 4 years with in situ DMS of $\mathrm{16.1}\pm \mathrm{4.8}$ $\mathrm{nM}$, about 6-fold higher than ${\text{DMS}}_{\text{SAT}}$. Local tuning of Eq. (2) using OSP data could not increase ${\text{DMS}}_{\text{SAT}}$ in August–September without degrading its performance in other months.

${\text{DMS}}_{\text{SAT}}$ vs. in situ data at long-term research stations. (a–c) Ocean Station P in the northeast Pacific (50${}^{\circ }$ N, 145${}^{\circ }$ W); (d–f) BATS station (31${}^{\circ }$40${}^{\prime }$ N, 64${}^{\circ }$10${}^{\prime }$ W). (a, d) compare ${\text{DMS}}_{\text{SAT}}$ estimates to in situ measurements; (b, e) compare monthly DMS means (i.e., climatologies) derived from ${\text{DMS}}_{\text{SAT}}$ (2003–2016 MODIS-Aqua data), in situ data (available measurements between 2003–2016) and the L11 climatology; (c, f) show the corresponding ${\text{Chl}}_{\text{SAT}}$, ${\text{DMSPt}}_{\text{SAT}}$ and ${\text{PAR}}_{\text{SAT}}$ climatologies (2003–2016) and in situ DMSPt. Lines and shaded envelopes show the mean $\pm$ range of satellite-derived data; filled circles and error bars show the mean $\pm$ range of available in situ data. DMS and DMSPt are in $\mathrm{nM}$ units, Chl in $\mathrm{\mu }\mathrm{g}{\mathrm{L}}^{-\mathrm{1}}$ and PAR in mol photons ${\mathrm{m}}^{-\mathrm{2}}{\mathrm{d}}^{-\mathrm{1}}$. In addition to public OSP data (three DMS measurements per year generally available between 2005 and 2016) we show DMSPt data reported in published studies (Asher et al., 2017; Levasseur et al., 2006; Royer et al., 2010). BATS data are comprised of the monthly 2005–2008 time series (Levine et al., 2016) and cruise data from the vicinity of BATS station in July 2004 (Bailey et al., 2008; Gabric et al., 2008).

[Figure omitted. See PDF]

Finally, note that these time series were calculated using semi-analytical bio-optical products (see section S5). Using the band-ratio Chl algorithm for MODIS (OC3) gave very similar results in deep ocean regions but 70 % higher concentrations in the shallow southern North Sea, possibly due to interference of non-algal materials on OC3 Chl retrieval (data not shown).

Bermuda Atlantic Time-series Study

Using the globally tuned coefficients (Eq. 2f), ${\text{DMS}}_{\text{SAT}}$ reproduces the shape of the mean seasonal cycle at the oligotrophic BATS station but underestimates DMS by around 2-fold between June and October (Fig. 9d and e). In August, part of this bias can be attributed to ${\text{DMSPt}}_{\text{SAT}}$ (Fig. 9f). However, replacing ${\text{DMSPt}}_{\text{SAT}}$ by in situ DMSPt raises ${\text{DMS}}_{\text{SAT}}$ by only 17 %, indicating that most of the underestimation is caused by the DMS sub-algorithm. Optimizing the coefficients using local data (Eq. 2h; Sect. S4) improves the model–data fit by decreasing the DMSPt coefficient and increasing the PAR coefficient. Indeed, different studies have shown that irradiance suffices to explain most of the DMS seasonal cycle at BATS (Galí and Simó, 2015; Toole and Siegel, 2004; Vallina and Simó, 2007b). Regarding interannual variation, it is of note that the locally tuned ${\text{DMS}}_{\text{SAT}}$ is in excellent agreement with in situ data throughout 2007 and in June through August in all years, while the underestimation persists in September and October of 2006 and 2008.

Discussion

The ${\text{DMS}}_{\text{SAT}}$ algorithm captures in situ variability (Figs. 4–9 and S3) using a small set of predictor variables (Fig. 2). Moreover, it reproduces the mismatch between DMS and Chl such that, at a given ${\text{Chl}}_{\text{SAT}}$ concentration, determined DMS can vary by up to 40-fold. This mismatch is larger than that produced by the SD02 or the VS07 algorithms, but still smaller than in the database (Fig. 7d–h). The progressive dissociation between ${\text{Chl}}_{\text{SAT}}$ and DMS imposed by the two-step algorithm, and the nonlinear relationships embodied in Eq. (2) (Fig. 3), allow ${\text{DMS}}_{\text{SAT}}$ to produce a DMS peak in the right concentration range in summer across different latitudes, despite order-of-magnitude variations in chlorophyll concentration (Le Clainche et al., 2010).

Independent validation using satellite matchups suggests that ${\text{DMS}}_{\text{SAT}}$ estimates are globally within $\pm$10 % of in situ measurements. However, this is at odds with the global mean bias of $-$33 % suggested by comparison to the L11 gridded climatology. In some regions, assessments of ${\text{DMS}}_{\text{SAT}}$ bias based on comparison to matchups and to L11 are consistent (Table 3), indicating shortcomings of the algorithm and/or in input satellite data. For instance, the globally optimized coefficients produce too low DMS in northern high latitudes, which we solved through regional tuning, building on the relative abundance of DMS(P) measurements north of 45${}^{\circ }$ N (Table 2, Fig. 8). In southern high latitudes, too-low ${\text{DMS}}_{\text{SAT}}$ primarily results from the negative ${\text{Chl}}_{\text{SAT}}$ bias (estimated at $<-$50 % by Johnson et al. (2013) and $-$56 % in our matchup dataset). Improving ${\text{DMS}}_{\text{SAT}}$ in this region would require, in the first place, the improvement of bio-optical algorithms and an important sampling effort to better document DMS(P) dynamics (Fig. 1; Jarníková and Tortell, 2016). In contrast to high latitudes, database matchups in the Westerlies and Trades biomes suggest smaller ${\text{DMS}}_{\text{SAT}}$ bias than comparison to L11 gridded data (Table 3). While some of this bias certainly arises from too-low ${\text{DMS}}_{\text{SAT}}$ in late summer at specific locations (Fig. 9), it is also plausible that divergences between ${\text{DMS}}_{\text{SAT}}$ and ${\text{DMS}}_{\text{L11}}$ arise from regional biases in the interpolated climatology.

In what follows we pinpoint the strengths and weaknesses of our novel approach focusing on two aspects. In Sect. 4.1 we examine how geo-statistical shortcomings of the in situ DMS database cause bias in the L11 interpolated climatology, highlighting the advantages of ${\text{DMS}}_{\text{SAT}}$ and paving the way towards improving gridded DMS fields. In Sect. 4.2 we speculate about potential causes behind high DMS $/$ DMSPt ratios in late summer and early fall. This exercise shows our limited capacity to account for relevant biogeochemical processes and explain their interannual variation using satellite data, and identifies knowledge gaps that need to be tackled to improve diagnostic and prognostic modeling of oceanic DMS(P). The rationale is that both kinds of issues (“geo-statistical” and “biogeochemical”) are highlighted by discrepancies among in situ data, L11, and the macroecological relationships embodied in our algorithm.

Known sources of error and bias: interpolated climatology vs. ${\text{DMS}}_{\text{SAT}}$

Global DMS fields estimated by the L11 climatology and by the ${\text{DMS}}_{\text{SAT}}$ algorithm show remarkable geographic differences (Fig. 5), and their mean concentrations differ by a factor of $\sim \mathrm{1.5}$. Particularly, changes in the sign of the ${\text{DMS}}_{\text{SAT}}$–${\text{DMS}}_{\text{L11}}$ anomaly often follow the somewhat artificial boundaries of the Longhurst biogeochemical provinces. As explained below, regional and global biases in L11 arise from the application of objective interpolation procedures to a global dataset characterized by (1) the right-skewed distribution of DMS concentrations (Kettle et al., 1999; Fig. 7), (2) the small amount of monthly data available in many biogeochemical provinces (Fig. 1), (3) the absence of repeat measurements in most oceanic regions (Halloran et al., 2010), (4) the low spatial resolution of most DMS datasets and (5) the preferential sampling of DMS-productive conditions.

A primary drawback for the calculation of interpolated climatologies is the overrepresentation of biologically productive conditions in the sea-surface DMS database. This sampling bias is clearly illustrated by comparing SeaWiFS-retrieved Chl concentration in the global ocean and in the database matchups (Fig. 7a). If we assume that SeaWiFS matchups ($N=\mathrm{11}\mathrm{600}$) represent a random sample of the DMS database ($N=\mathrm{41}\mathrm{304}$), and considering the global positive correlation between Chl and DMS (Fig. 7d), this implies a sampling bias towards high DMS concentrations. The bias is largest when the comparison is restricted to the spring–summer semester of each hemisphere, when the median ${\text{Chl}}_{\text{SAT}}$ is 0.19 in the global SeaWiFS climatology and 0.74 in the DMS database matchups. This is the period when DMS peaks and has more influence on mean annual DMS concentration.

Sampling bias is intertwined with the right-skewed distribution of sea-surface DMS concentrations (Fig. 7b) and the poor spatial resolution of most in situ DMS datasets. Spatial averaging, justified by data scarcity, is appropriate when applied over small or sufficiently homogeneous regions. However, when applied over a large and potentially heterogeneous Longhurst province, it can propagate the sampling bias over the entire province and up to the global scale. As illustrated in Fig. 6g, the mean DMS concentration in $M\mathrm{5}\times \mathrm{5}$ bins is systematically higher, by 40 % on average, than the corresponding median. Province-level averaging converts the long tail of high in situ DMS concentrations into too-large province means, such that the global mode of ${\text{DMS}}_{\text{L11}}$ is similar or even higher than that of in situ DMS (Fig. 7b).

The influence of extreme in situ DMS concentrations is maximal in productive regions, where $\text{mean}/\text{median}$ ratios of around 4 are observed in $M\mathrm{5}\times \mathrm{5}$ bins (Fig. 6g). In these regions, sharp productivity gradients and dynamic ecosystem processes complicate the task of sampling DMS through all appropriate scales (Nemcek et al., 2008), suggesting the need to apply finer-scale or dynamic regionalization (Devred et al., 2007) prior to interpolation. Emerging biogeochemical relationships, like that between net community production and DMS observed by Kameyama et al. (2013) in the northeast Pacific, might assist DMS interpolation and determination, but require validation across contrasting regions and relevant scales (Asher et al., 2017). In low-latitude oligotrophic areas where DMS shows reduced spatial variability (Royer et al., 2015), the method used to construct interpolation-based climatologies (Kettle et al., 1999; Lana et al., 2011) seems appropriate regarding spatial resolution.

In the temporal domain, the calculation of interpolated climatologies is complicated by two factors: poor interannual coverage, i.e., the scarcity of DMS measurements repeated in different years in a given region; and poor seasonal coverage, i.e., the scarcity of fully resolved seasonal cycles. Seasonal coverage is limited at the province level (see Table 1 in Lana et al., 2011) and obviously worse in 5${}^{\circ }$ bins (Fig. 1d). Regarding interannual coverage at the MLongh level, 42 % of the province–month bins contain measurements from a single year, 21 % from two years and 37 % from three or more years. Thus, data from one or two years are often assumed representative of the mean ecosystem state in interpolated climatologies, which is probably not the case in regions with wide interannual variability or long-term trends (Vantrepotte and Mélin, 2011). While this does not necessarily bias global DMS fields, it can produce artificial seasonal cycles. For example, L11 suggests the existence of early spring and fall DMS peaks in the North Atlantic drift area, which result from interpolation from neighboring regions (Fig. 8a). In contrast, ${\text{DMS}}_{\text{SAT}}$ suggests these are improbable (spring) or infrequent (fall) features. Another example is found at OSP, where ${\text{DMS}}_{\text{L11}}$, based on measurements made before 2003, is in poor agreement with measurements made between 2005 and 2016. In February and June, ${\text{DMS}}_{\text{SAT}}$ is in better accordance with in situ DMS data at OSP (Fig. 9a and b).

In summary, caution has to be taken when comparing DMS measurements, their derived climatological products, and independent model estimates that are not collocated in time. This temporal mismatch may partly explain the poor correlation between modeled DMS climatologies, on one hand, and the DMS database and ${\text{DMS}}_{\text{L11}}$ climatology, on the other (Tesdal et al., 2016). Note that the latter study compared DMS fields binned into monthly $\mathrm{5}{}^{\circ }\times \mathrm{5}{}^{\circ }$ boxes ($M\mathrm{5}\times \mathrm{5}$), such that 82 % of the bins contained measurements from a single year.

Compared to interpolated climatologies, ${\text{DMS}}_{\text{SAT}}$ provides a robust means to estimate DMS concentrations in sparsely sampled areas because it relies on satellite observations and macroecological relationships. The resulting DMS fields are in better accordance with natural gradients in plankton abundance (biogeography, phenology) and environmental forcing, as long as the models can account for the driving factors, as discussed below.

Unknown sources of error: how far can we go with remote sensing algorithms?

By testing ${\text{DMS}}_{\text{SAT}}$ in challenging biogeochemical settings (Fig. 9), we identified its main drawback: the failure to reproduce high DMS, and more specifically DMS $/$ DMSPt ratios higher than 0.3, at intermediate PAR levels (Fig. 3), as observed between midsummer and early fall at BATS and OSP in some years. This limitation can hardly be fixed without identifying the underlying biogeochemical processes, which are not necessarily the same in these contrasting biogeochemical regimes. Although this feature is probably not widespread (see Fig. 2 in Lana et al., 2011), its occurrence in emblematic time series stations warrants further discussion.

A common explanation could be the underestimation of irradiance effects, caused by the use of sea-surface ${\text{PAR}}_{\text{SAT}}$ rather than ${\text{PAR}}_{\text{MLD}}$ as DMS predictor variable. For instance, a delay of seasonal mixing, associated with deeper irradiance penetration, could enhance stress-driven DMS production well into fall. Yet examination of the BATS and OSP time series does not support this explanation. At both sites, the summer MLD is stable at about $\le \mathrm{20}$ $\mathrm{m}$ and deepens slowly in late summer (Levine et al., 2016; Steiner et al., 2012), such that ${\text{PAR}}_{\text{MLD}}$ declines faster than sea-surface ${\text{PAR}}_{\text{SAT}}$ (Eq. 1). Thus, using ${\text{PAR}}_{\text{MLD}}$ instead of surface PAR cannot delay the decline of modeled DMS $/$ DMSPt ratios through the summer, and other factors need to be invoked.

In the oligotrophic BATS station, some modeling studies proposed macronutrient limitation of bacteria (Polimene et al., 2011) and also phytoplankton (Belviso et al., 2012) as drivers of the seasonal mismatch between DMSPt and DMS, besides irradiance (Vallina et al., 2008). With this in mind, we tried to factor phosphate and nitrate limitation into our regression models using different variables: nutrient concentrations, nutricline depths (Table S4) and limitation factors estimated according to Michaelis–Menten kinetics (not shown). However, none of the tested variables improved the regression models significantly. Moreover, macronutrient availability (limitation) terms generally entered regression models with positive (negative) coefficients, even when regressions were restricted to oligotrophic low latitudes. This implies that macronutrient limitation of phytoplankton growth globally acts to decrease DMS, offsetting nutrient stress responses that increase DMS. Note also the irregular occurrence of high DMS at the BATS station in late summer in different years (Fig. 9d; Levine et al., 2016), and that a BATS-like seasonality is not observed at other sites with late summer macronutrient limitation (Archer et al., 2009; Belviso et al., 2012; Galí and Simó, 2015; Vila-Costa et al., 2008). Altogether, these observations suggest that regional macronutrient stress responses are difficult to generalize.

Analysis of the OSP time series also yields valuable information because macronutrient concentrations remain at high concentrations in late summer in this iron-limited regime (Harrison et al., 2004). While in situ DMS is accurately estimated by ${\text{DMS}}_{\text{SAT}}$ in February and June, the variable DMS peak occurring around August is strongly underestimated. Previous studies emphasized the role of iron limitation at OSP, which configures phytoplankton communities dominated by high DMSP taxa (Asher et al., 2017; Levasseur et al., 2006; Royer et al., 2010; Steiner et al., 2012). Interestingly, ${\text{DMSPt}}_{\text{SAT}}$ peaks in August at OSP, in phase with the in situ DMS peak and in good accordance with the few available DMSPt measurements (Fig. 9c), even though the ${\text{DMSPt}}_{\text{SAT}}$ sub-algorithm does not explicitly resolve phytoplankton taxonomy (Galí et al., 2015). Therefore, high DMS yields, possibly co-occurring with low DMS removal rate constants (Asher et al., 2017), are required to explain DMS $/$ DMSPt ratios observed at OSP. High DMS yields probably result from a combination of processes, including algal and bacterial DMSP metabolism (Merzouk et al., 2006; Royer et al., 2010) and microzooplankton grazing (Steiner et al., 2012). The striking late summer variability at OSP is presently not captured by biogeochemical models (Steiner et al., 2012) or empirical algorithms, and it remains unanswered whether it simply reflects too low sampling frequency, or it is caused by processes that switch on/off depending on environmental conditions in a given year, or by the variable location of oceanic fronts in response to circulation patterns.

In summary, our analysis indicates that additional factors are needed to better reproduce DMS seasonality in specific regions, where DMS $/$ DMSPt ratios are occasionally higher than the “baseline” established by Eq. (2) (Figs. 3 and 9). Biotic interactions involving phytoplankton, bacteria and microzooplankton, regionally interacting with iron and macronutrient limitation in multiple ways, are good candidates to explain strong deviations from the mean relationship between DMS, DMSPt and irradiance. However, they can hardly be represented in empirical algorithms with our current level of understanding, particularly when interannual changes are considered (Figs. 8 and 9).

From a practical standpoint, tuning the Eq. (2) coefficients is a workable alternative in certain regions (Figs. 8 and 9). If sufficient measurements were available in all oceanic areas, Eq. (2) could perhaps be generalized in a way that allows its coefficients to vary across different biogeochemical regimes, while avoiding geographic discontinuities. The inclusion of additional terms in Eq. (2) lacks strong statistical support when applied globally, at least with the current dataset (Table S4). If posterior analyses supported the addition of new satellite variables, their retrieval uncertainty and its propagation to ${\text{DMS}}_{\text{SAT}}$ should be considered. More obviously, climatological variables such as the WOA nutrient concentrations are not appropriate to produce time series, and their use in remote sensing algorithms should be minimized. The only climatological variable used in our algorithm is MLD, which enters mainly as a categorical variable (Galí et al., 2015), such that ${\text{DMSPt}}_{\text{SAT}}$ is robust to MLD uncertainties (Fig. S1).

The question of the “optimal model complexity” is a pervasive one in biogeochemistry, and the right answer may depend on the purpose of each study. The algorithms tested here showed improved qualitative and quantitative performance with increasing complexity (VS07 $<$ SD02 $<{\text{DMS}}_{\text{SAT}}$). VS07 failed to capture DMS patterns outside the subtropical band, possibly due to its inability to modulate the DMS–irradiance relationship depending on phytoplankton biomass. Inclusion of phytoplankton biomass-dependent terms in SD02, and of implicit taxonomic information through the embedded ${\text{DMSPt}}_{\text{SAT}}$ sub-algorithm in ${\text{DMS}}_{\text{SAT}}$, improved algorithm skill in productive regions, where DMS shows wider seasonal cycles and sharper spatial gradients.

More sophisticated approaches may be needed to achieve significant improvements in model skill, but they also suffer from major uncertainties. For example, neural networks were successfully used to estimate DMS in the Arctic (Humphries et al., 2012), but their robustness might be compromised by the small training datasets, the use of climatological variables and the lack of a mechanistic basis. Complex biogeochemical models with satellite data assimilation have strong potential for resolving interannual DMS variations, but reliance on several tens of poorly constrained parameters currently limits their skill (Le Clainche et al., 2010; Galí and Simó, 2015; Tesdal et al., 2016). An approach of intermediate complexity that deserves further exploration is DMS diagnosis based on a simplified steady-state budget equation (Galí and Simó, 2015). Applying empirical parameterizations for the main DMS production and removal pathways, this approach could enhance the flexibility of remote sensing algorithms across a wider range of biogeochemical settings.

Conclusions and outlook

Sensors on polar-orbiting satellites provide synoptic observations of the global ocean surface every few days, and are thus well suited to resolve spatial and temporal variations in DMS concentration. The ${\text{DMS}}_{\text{SAT}}$ algorithm presented here, based on robust macroecological relationships, reproduces the main spatial–temporal features of sea-surface DMS(P) concentrations with remarkable skill using satellite-retrieved Chl, euphotic layer depth and PAR, and climatological MLD. Other strengths of our approach are its flexibility, allowing for regional tuning, and the minimal computing cost.

When compared against the L11 interpolated DMS climatology (Lana et al., 2011), the ${\text{DMS}}_{\text{SAT}}$ climatology shows similar latitudinal profiles but disagrees in the basin-scale patterns. Examination of spatial DMS statistics highlights possible shortcomings of the L11 climatology caused by the combination of sparse and biased sampling, the right-skewed distribution of DMS and the interpolation procedures used. High-resolution measurements of DMS(P), if validated against traditional standard techniques (Royer et al., 2014), will help improve interpolated climatologies and models.

The global mean area-weighted ${\text{DMS}}_{\text{SAT}}$ concentration is 1.63 $\mathrm{nM}$, 33 % lower than ${\text{DMS}}_{\text{L11}}$ (2.43 $\mathrm{nM}$). Global-scale ${\text{DMS}}_{\text{SAT}}$ fields are insensitive to the choice of different Chl and euphotic depth satellite products, but semi-analytical products should be used in optically complex coastal waters to avoid DMS overestimation (after detailed examination at regional scale). Globally, ${\text{DMS}}_{\text{SAT}}$ suffers from negative bias for exogenous and endogenous reasons. In the Antarctic Ocean, it is affected by the negative bias in satellite-retrieved chlorophyll. At the BATS and OSP stations (at least), additional factors besides irradiance are needed to enhance DMS concentration in late summer in some years. Excluding the Coastal and Antarctic biomes, ${\text{DMS}}_{\text{SAT}}$ bias assessed using satellite matchups ranges between $-$16 and $-$20 % (Table 3). This bias is probably more realistic than the $-$33 % deduced from comparison to L11, given the evidence for DMS overestimation in the L11 climatology.

Gauging global DMS emission is critical to understand gas-to-particle conversion efficiency and the dynamics of cloud condensation nuclei (CCN) populations in the marine boundary layer. Assuming a linear relationship between global mean DMS concentration and emission (Fig. 8 in Tesdal et al., 2016), ${\text{DMS}}_{\text{SAT}}$ suggests a global emission of 16–20 $\mathrm{Tg}\mathrm{S}{\mathrm{yr}}^{-\mathrm{1}}$ (depending on the assigned bias). These emission values lie within the low range of current estimates (Lana et al., 2011; Tesdal et al., 2016).

Unlike climatologies constructed from the database, the satellite-based algorithm allows one to explore interannual change. Implementation of ${\text{DMS}}_{\text{SAT}}$ in the subpolar Atlantic between 2003 and 2016 illustrates the wide interannual variability in the timing and magnitude of the annual DMS peak(s) over large areas. This opens new avenues for studying the imprint of oceanic aerosol precursors on cloud properties using simultaneous ocean–atmosphere satellite observations (Krüger and Grabßl, 2011; McCoy et al., 2015; Meskhidze and Nenes, 2006). If coupled to atmospheric measurements and numerical models, ${\text{DMS}}_{\text{SAT}}$ enables studying the effects of contemporaneous DMS variability on atmospheric chemistry and clouds, which could lead to a better understanding of intricate aerosol–cloud interactions. Further work is warranted to analyze marine DMS emission variability patterns in regions where climate is particularly sensitive to DMS, such as the Southern Ocean and the Arctic.

The code used to perform the data analyses and produce ${\text{DMS}}_{\text{SAT}}$, SD02 and VS07 DMS fields can be provided by the authors on request.

The primary database used to develop the ${\text{DMS}}_{\text{SAT}}$ algorithm is publicly available at http://saga.pmel.noaa.gov/dms/ (last access: 13 April 2017). The database extended with satellite matchups and climatological variables and the global DMS and DMSPt climatologies derived with the ${\text{DMS}}_{\text{SAT}}$, SD02 and VS07 algorithms can be provided by the authors on request. The L11 DMS climatology and other related documents and datasets can be downloaded from https://www.bodc.ac.uk/solas_integration/implementation_products/group1/dms/ (last access: 18 July 2013).

Description of the ${\text{DMSPt}}_{\text{SAT}}$ sub-algorithm

We estimated sea-surface DMSPt concentration ($\mathrm{nmol}{\mathrm{L}}^{-\mathrm{1}}$) using the algorithm of Galí et al. (2015). This algorithm estimates DMSPt as a function of chlorophyll $a$ concentration (Chl) by switching between two different equations depending on the light penetration regime, defined by the quotient between euphotic layer depth (Zeu) and mixed layer depth (MLD), $$\begin{array}{cc}{\displaystyle }& {\displaystyle }{\log }_{\mathrm{10}}\text{DMSPt}=\mathrm{1.70}+\mathrm{1.14}{\log }_{\mathrm{10}}\text{Chl}+\mathrm{0.44}{\log }_{\mathrm{10}}{\text{Chl}}^{\mathrm{2}}\\ & {\displaystyle }& {\displaystyle }+\mathrm{0.063}\text{SST}-\mathrm{0.0024}{\text{SST}}^{\mathrm{2}},{\displaystyle }& {\displaystyle }\text{Zeu}/\text{MLD}\ge \mathrm{1},\\ {\displaystyle }& {\displaystyle }\text{stratified water column,}\\ {\displaystyle }& {\displaystyle }{\log }_{\mathrm{10}}\text{DMSPt}=\mathrm{1.74}+\mathrm{0.81}{\log }_{\mathrm{10}}\text{Chl}\\ & {\displaystyle }& {\displaystyle }+\mathrm{0.60}{\log }_{\mathrm{10}}(\text{Zeu}/\text{MLD}),{\displaystyle }& {\displaystyle }\text{Zeu}/\text{MLD}<\mathrm{1},\text{mixed water column,}\end{array}$$ with Chl in $\mathrm{mg}{\mathrm{m}}^{-\mathrm{3}}$ ($=$ $\mathrm{\mu }\mathrm{g}{\mathrm{L}}^{-\mathrm{1}}$) and sea-surface temperature (SST) in ${}^{\circ }$C. This scheme implicitly reproduces the seasonal ecological succession from low DMSP phytoplankton taxa (mainly diatoms) towards high DMSP taxa (mainly haptophytes and dinoflagellates) that are relatively more abundant in stratified conditions, when the entire sea-surface mixed layer is well illuminated. To avoid uncertain algorithm behavior at extreme Chl concentrations, ${\text{Chl}}_{\text{SAT}}$ lower (higher) than 0.04 (60) $\mathrm{mg}{\mathrm{m}}^{-\mathrm{3}}$ was set to these respective bounds. The ${\text{DMSPt}}_{\text{SAT}}$ sub-algorithm uses a third equation, based on satellite-retrieved PIC, in coccolith-laden waters where Chl algorithms fail (Galí et al., 2015). In practice, this equation is not used in climatological implementations (see Fig. 2).

${\text{DMSPt}}_{\text{SAT}}$ sub-algorithm. Dots in the background show database measurements, and lines represent fits for “stratified” and “mixed” conditions at different SST and $Z_{\text{eu}}/\text{MLD}$ ratios, respectively.

[Figure omitted. See PDF]

The supplement related to this article is available online at: https://doi.org/10.5194/bg-15-3497-2018-supplement.

MG designed the study, performed the research and wrote the paper, with input from all coauthors through the different phases. ED processed remote sensing reflectance data used as input for the Northern Hemisphere ${\text{DMS}}_{\text{SAT}}$ time series.

The authors declare that they have no conflict of interest.

Acknowledgements

We thank the NASA Ocean Biology Distributed Active Archive Center (OB.DAAC) for access to MODIS and SeaWiFS datasets, Timothy (Tim) Bates (NOAA/PMEL) for the maintenance of the GSS DMS(P) database and the DMS-GO project for the recent database update (a joint initiative of the SOLAS Integration Project and the EU projects COST Action 735 and EUR-OCEANS to Rafel Simó); Eric Rehm and Maxime Benoît-Gagné for IT support; Naomi Levine, John Dacey, Nick Bates and Scott Doney for sharing BATS data; and Marie Robert and Michael Arychuk for guidance on access to Ocean Station P public datasets. We acknowledge funding from the Canada Excellence Research Chair in Remote Sensing of Canada's New Arctic Frontier (Marcel Babin), the Canada Research Chair on Ocean Biogeochemistry and Climate, NSERC Discovery Grant Program, Northern Research Supplement Program (Maurice Levasseur), the NETCARE network (funded under the NSERC Climate Change and Atmospheric Research program) and ArcticNet (The Network of Centres of Excellence of Canada). Rafel Simó acknowledges funding from the Spanish Ministry of Economy through project BIOGAPS. Martí Galí acknowledges the receipt of a Beatriu de Pinós postdoctoral fellowship funded by AGAUR (Generalitat de Catalunya). This project is a contribution to the research program of Québec-Océan and the Takuvik Joint International Laboratory (CNRS-France and Univeristé Laval-Canada). Edited by: S. Wajih A. Naqvi Reviewed by: two anonymous referees