By mass, dust is the largest contributor to global aerosol burden, yet long-term observational records of dust, particularly over the ocean, are limited. Here, two nearly global observational datasets of dust aerosol optical depth τd are created primarily on the basis of optical measurements of the aerosol column from 1) the Moderate Resolution Imaging Spectroradiometer (MODIS) aboard the Terra satellite spanning from 2001 to 2018 and 2) the Advanced Very High Resolution Radiometer (AVHRR) from 1981 to 2018. The quality of the new data is assessed by comparison with existing dust datasets that are spatially more limited. Between 2001 and 2018, τd decreased over Asia and increased significantly over the Sahara, Middle East, and parts of eastern Europe, with the largest increase found over the Aral Sea where emissive playa surfaces have been exposed. These daily, observational, and nearly global records of dust will allow for improvement in understanding the role of dust in climate variability.
Aeolian dust makes up the largest mass fraction of the global aerosol burden (Textor et al. 2006) and produces profound impacts on natural (Rosenfeld et al. 2001; Swap et al. 1992; Mahowald 2011; Miller et al. 2004; Ault et al. 2011; Evan et al. 2006; Strong et al. 2018) and anthropological systems (Tong et al. 2017; Griffin 2007; Shao 2008; Mani and Pillai 2010; Ai and Polenske 2008; Brown 2002). A major limiting factor in understanding the global distribution of dust, how it will change in the future, and its complex interactions with Earth’s climate, is the limited availability of measurements especially over the ocean (Prospero and Mayol-Bracero 2013). Airborne campaigns have provided measurements for select places and times (Ralph et al. 2016; Formenti et al. 2003; Chen et al. 2011; Formenti et al. 2008; Stith et al. 2009; Ryder et al. 2013; Klaver et al. 2011) and ground-based networks have provided long-term measurements for specific regions (Malm et al. 1994; Prospero and Nees 1986; Prospero 1999; Holben et al. 1998). Most satellite sensors do not isolate the dust contribution to aerosol optical depth (AOD) at the point of retrieval, so several methods have been proposed and utilized for studying spatial and temporal variations of dust using satellite measurements and models. These methods have mostly used optical properties of the aerosol column related to particle size (Ginoux et al. 2012; Kaufman 2005), color, and single scattering albedo (Ginoux et al. 2012), or have used the ultraviolet (UV) absorptive property of dust for identification (Ginoux 2003). Several studies have used satellite-based products to study dust aerosol optical depth τd over either the ocean or land in select locations. Over the Atlantic Ocean, Kaufman (2005) analyzed τd using the Moderate Imaging Spectroradiometer (MODIS) retrievals of AOD and Evan and Mukhopadhyay (2010) used a similar method to create a record of τd based on retrievals by the Advanced Very High Resolution Radiometer (AVHRR). In addition, Ginoux et al. (2012) used τd over bright land surfaces that is based on MODIS Deep Blue AOD to identify dust source regions. Monthly τd has also been used in evaluation of the dust emission schemes of global climate models (Pu and Ginoux 2018b).
More recently, studies have used a combination of satellite observations and models to investigate the distribution of dust globally and trends in dust for select regions. Ridley et al. (2016) produced a long-term mean dataset of τd between 2004 and 2008 using AOD retrievals from multiple satellite platforms and dust estimates from several global models and Chin et al. (2014) used satellite AOD from several sensors alongside the Goddard Chemistry Aerosol Radiation and Transport (GOCART) model to investigate multidecadal regional trends in aerosol species, including dust. These merged approaches have value in allowing for elucidation of aerosol sources and the ability to distinguish anthropogenic from natural dust. However, they are susceptible to model biases in dust emission and transport. Climate models display a large diversity in dust aerosol optical depth and dust emission, mostly underestimating dust emission from major sources, and do not consistently agree with observations (Evan et al. 2014; Huneeus et al. 2011; Pu and Ginoux 2018a; Kok et al. 2017; Kok 2010; Ryder et al. 2018). Therefore, it is valuable to have complimentary extensive observational records of dust for comparison.
Here, we expand on the work of Evan and Mukhopadhyay (2010) by using the extended record of satellite retrievals to produce two nearly global τd datasets developed using MODIS (2001–18) and AVHRR (1981–2017) at daily and monthly resolution, respectively, with the former dataset extended to cover both land and ocean regions. Although dust transport happens on a temporal scale of days to weeks, the previously mentioned efforts have focused on estimating monthly mean τd. This work provides the longest record of observed daily τd to date, at a temporal resolution suitable for analysis of daily transport events (daily) and 1° × 1° spatial resolution.
In the next section we describe the data and methods used in the estimation of global τd and provide an estimation of the associated uncertainty. Our results in section 3 include identification of climatological seasonal patterns in τd, a comparison of our τd estimates as derived from two different sensors, and a comparison of the MODIS-derived dataset with preexisting ground-based measurements of dust aerosol from several locations. We then identify regional trends in τd for the period of each record.
Two datasets of τd were created. The first, which hereinafter will be referred to as , is based primarily on AOD (τ) retrievals from the MODIS instrument on the Terra platform, with a companion set of estimates created for the same instrument aboard the Aqua platform, (both accessed at https://ladsweb.modaps.eosdis.nasa.gov) and extends from 2001 to 2018 for the case of Terra and 2003 to 2018 for Aqua. The includes estimates over both land and water-covered surfaces. The is based on aerosol optical thickness (AOT) retrievals from the AVHRR instruments aboard the NOAA series of satellites (Heidinger et al. 2014; Zhao et al. 2008, 2002; Zhao and NOAA CDR Program 2017) and covers the period 1981–2018. AVHRR AOT retrievals are not made over land, and thus is an overwater-only product.
a. Datasets used
MODIS level-3 (L3), Collection 6.1 (C061), daily Dark Target AOD is derived from the measured 500-m resolution reflectance from all visible MODIS bands by taking the average of the measured radiance over all scenes that are cloud free and glint free over the ocean within a 10-km grid box and fitting these values to a lookup table. Level-3 MODIS Dark Target AOD is calculated using the quality assurance confidence screened level-2 pixel-level data. MODIS 550-nm fine-mode fraction f is calculated from the ratio of the small-mode AOD to the total AOD. MODIS Collection 6 (C006) AOD over ocean has been shown to have an expected error of +(0.04 + 10%)/−(0.02 + 10%) (Levy et al. 2013). For dust-dominated regions, it has been shown that there is a bias of +5% (Kaufman 2005). There has been no systematic program for evaluation of MODIS f because its definition is somewhat ambiguous, making evaluations against other estimates difficult (Kleidman et al. 2005). At this time, there is no envelope of error for C061 f. However, Collection 5 (C005) f over ocean was found to agree with Aerosol Robotic Network (AERONET) retrieved sky radiance fine-mode fraction within approximately ±0.20 (Kleidman et al. 2005), and Bréon et al. (2011) found agreement between MODIS (τf,MODIS) and AERONET (τf,AERONET) fine-mode AOD, as evidenced by a correlation r value of 0.76, root-mean-square error (RMSE) of 0.08, bias of 0.01, and 53% of retrievals within the C005 envelope of error, δτ = 0.03 + 0.05τAERONET, which describes a confidence envelope containing 1 standard deviation σ (i.e.,68%) of τf,MODIS/τf,AERONET matchups. In our estimate of the known uncertainty in τd, we assume that the C061 f uncertainty has not changed from the C005 values. Although a quality-assured dataset does exist for τf,MODIS, we chose to use the standard dataset so as to retain additional coverage.
The equivalent set of estimates for from the Aqua platform extends from 2003 to 2018. The Aqua platform passes the equator at 1330 local time, whereas Terra passes near 1030 local time. AOD from MODIS Terra and Aqua have been shown to agree very well (Ichoku 2005). However, globally, Terra monthly mean AOD over ocean has been shown to be higher than Aqua by +0.015 (Levy et al. 2018; Remer et al. 2008). This is partially explained by diurnal cycles in cloud fraction (King et al. 2013), which may create differences in sampling between these sensors, and partially explained by issues with the total infrared radiation measured by Terra (Levy et al. 2018). This paper will focus largely on from the Terra platform, because it has a longer record. However, we will discuss decadal trends in derived from Aqua in section 4a.
Because our τd estimate is built upon AOD, it carries with it the known limitations of the sensors from which it is derived. MODIS AOD retrievals are limited to cloud-free pixels, and errors in cloud screening will impact the AOD estimate. The MODIS aerosol cloud mask uses the standard deviation of reflectance in sets of pixels to remove cloudy pixels, which increase spatial variability. Absolute reflectance at 1380 nm and the ratio of reflectances at 1380 and 1240 nm as well as several infrared tests are used for additional screening of thin cirrus. Over ocean, the brightest and darkest 25% of the remaining pixels are arbitrarily removed and the average reflectance in each channel is calculated from the remaining pixels. Over land, the brightest 50% and darkest 30% of pixels are discarded before the averaging step (Remer et al. 2012). While this usually produces a reasonable result, in the case that unscreened clouds remain before the discarding of the brightest and darkest retrievals, this can produce AOD estimates that are biased high (Kahn et al. 2007). It is known that MODIS AOD is biased very positively (high) in the Southern Hemisphere mid–high latitudes over the Southern Ocean because of extensive broken stratocumulus and cirrus cloud contamination (Toth et al. 2013). We have thus limited these datasets to extend from 60°N to 50°S.
Level-2, version 3, AOD and fine-mode fraction from select sites, as will be discussed in section 2b, retrieved from AERONET (Holben et al. 1998, 2001) sunphotometers using the Spectral Deconvolution Algorithm (SDA), version 4.1 (O’Neill 2003; https://aeronet.gsfc.nasa.gov) are used to estimate and over the ocean. In the determination of AERONET fine-mode fraction, two spectral modes of AOD are defined on the basis of the premise that the coarse-mode spectral variation is approximately neutral while the derivative of the fine-mode spectral variation is an approximate function of the fine-mode Ångstrom exponent (O’Neill 2003). From this, the ratio of fine-mode to total AOD is determined at a reference wavelength of 500 nm from a second-order polynomial fit of the natural logarithm of AOD and wavelength applied to each AOD spectrum across six bands. This definition of fine-mode fraction is different from that used in the MODIS dataset and may result in some disagreement. MODIS f has been shown to slightly overestimate the fine-mode fraction in dust or salt dominated regions by 0.1–0.2 relative to AERONET fine-mode fraction as generated from the O’Neill algorithm (Kleidman et al. 2005), which has the potential to introduce a small positive bias in our estimated τd.
MODIS C061 L3 Deep Blue products [AOD, single scattering albedo (SSA), and Ångstrom exponent (AE); https://ladsweb.modaps.eosdis.nasa.gov] were used in the estimation of over land (Hsu et al. 2013). These products use the blue channels of MODIS where surface reflectance is very low. Pixel-level retrievals are averaged over 10 km × 10 km grids, and data are then aggregated into granules. The L3 products are provided at 1° × 1° resolution. SSA at 412 nm and at 660 nm were used in addition to AE, which is inversely proportional to particle size (Angstrom 1929). Deep Blue AOD was used as the basis for estimation of over land and is reported to have an expected error of +(0.04 + 10%)/−(0.02 + 10%) (Hsu et al. 2014). The C061 Deep Blue AOD release extended coverage to all nonsnow land surfaces, whereas C005 only included bright land surfaces.
The 630-nm AVHRR AOT climate data record (CDR) (Zhao and NOAA CDR Program 2017; accessed at https://www.ncei.noaa.gov/data/) was used in the estimation of (Zhao and NOAA CDR Program 2017; Chan et al. 2013). The CDR was created using the NOAA AVHRR PATMOS-x level-2B product (Zhao et al. 2008, 2002; Zhao and NOAA CDR Program 2017). Although there is no difference between AOD and AOT, we refer to the MODIS data as AOD and the AVHRR data as AOT because this is how each is referenced in their documentation. It combines retrievals from 17 different sensors, and the overpass time and number of observations per day (4–8) varied between sensors. Large data gaps occur for years before 1985. AVHRR AOT was regridded from its native resolution of 0.1° × 0.1° to a resolution of 1° × 1° using inverse distance weighted interpolation. Aerosol retrievals are only available over ocean and require clear-sky conditions (Heidinger et al. 2014). The CDR AVHRR AOT dataset has a systematic error (positive bias) of 0.03 ± 0.006 and a random error of ±0.113 (Zhao and NOAA CDR Program 2017). The error budget of AVHRR AOT CDR is based on AERONET validation. Another product of AVHRR aerosol loading is available, derived using algorithms from modern sensors (Hsu et al. 2017; Sayer et al. 2017), but it is limited to several years. Therefore, we did not use this product.
Daily mean surface wind speed from the second Modern-Era Retrospective Analysis for Research and Applications (MERRA-2) (Gelaro et al. 2017) 1-h time-averaged surface flux assimilation product (accessed at https://disc.gsfc.nasa.gov/) was used in the estimation of marine AOD. This was regridded from its native resolution of 0.625° longitude × 0.5° latitude to a resolution of 1° × 1° using bilinear interpolation. MERRA-2 surface wind patterns have been shown to be similar to other datasets, but greater in magnitude in most regions than MERRA and ERA-Interim and weaker in magnitude than NCEP-R2 (Bosilovich et al. 2015). Similar to most reanalysis products, MERRA-2 winds are weaker than 93% of observations considered in Bosilovich et al. (2015).
b. MODIS dust optical depth over ocean
Our method for estimating over ocean surfaces is based on Kaufman (2005). Underlying this method is the assumption that the major contributions to total AOD are anthropogenic and biomass-burning aerosol, marine aerosol, and dust aerosol, and thus one can deduce τd by estimating the contributions from the other aerosol types. As anthropogenic and biomass-burning aerosols are generally dominated by submicrometer particles with AE greater than 1 (Dubovik et al. 2002; Schmeisser et al. 2017), f can be used to discriminate these types from the contributions by coarse-mode-dominant aerosol types.
The total aerosol optical depth τ is
where τm and τa are the marine and anthropogenic or biomass-burning contributions to τ, respectively. Here we have combined the biomass-burning and anthropogenic aerosols for convenience, because they both have a dominant fine mode and will thus be isolated and removed from τ simultaneously. The fine-mode aerosol optical depth is given by fτ and can be expressed as a function of the fine-mode optical depth contributions from the individual aerosol species,
Data sources for each of these variables can be found in Table 1. Characteristic values of fm, fa, and fd were determined using the average of fine-mode fraction from AERONET stations dominated by each aerosol type, weighted by the total τ (Table 2). This weighting was performed for each station according to
where f is fa, fd, or fm. For fa and fd, AERONET stations were chosen on the basis of Cazorla et al. (2013). For fm, marine aerosol-dominated AERONET island stations were selected based on Lehahn et al. (2010) and Chaubey et al. (2011) and through visual selection of additional remote island sites. MODIS retrievals were not used for determining these coefficients since f is only distinguishable over the ocean with this instrument. This makes it likely that all MODIS f estimates contain some significant contribution from sea spray and as such represents a deviation from the methods of Kaufman (2005).
Marine aerosol contribution
Kaufman (2005) parameterized marine aerosol based on NCEP surface (1000 hPa) wind speed at 2.5° × 2.5° resolution. This parameterization was based on the relationship between surface wind speed from a National Climatic Data Center meteorological station and an AERONET station on Midway Island from 14 months of data, excluding the dust season (February–May) (Smirnov et al. 2003).
Here we consider that the total marine aerosol optical depth is given by
where τm,c and τm,f are the coarse- and fine-mode contributions to marine aerosol. We make the assumption that τm,f is weakly dependent on wind speed (Satheesh et al. 2006), and that the dependence of τm on surface wind speed is primarily from τm,c (Lehahn et al. 2010). We identified a region dominated by marine aerosols (0°–25°S, 178°–130°W), and set τm,f equal to the long-term MODIS fine-mode aerosol optical depth,
where ±0.03 is the standard deviation of the daily MODIS fine-mode aerosol optical depth at 0°–25°S, 178°–130°W from 2001 to 2017 (Fig. 1c). Assuming the linear relationship,
where w denotes surface wind speed and it is implied that τm,c = 0 in the absence of surface winds, the proportionality term α was calculated via regression of MODIS Terra coarse-mode AOD onto MERRA-2 surface wind speeds, both interpolated to a 1° × 1° horizontal resolution, and over the same region used to calculate τm,f. Via the linear regression (Fig. 1b),
which is consistent with Lehahn et al. (2010), who found τm,c to be linearly related to surface wind speed with a slope of 0.009 ± 0.002 m s−1. We tested three other parameterizations for marine aerosol and found that this parameterization resulted in the lowest bias (+0.02) and lowest RMSE (0.04) when compared with AERONET AOD at Nauru, a remote island station, from 2001 to 2014 for days on which there were clear-sky measurements from both instruments. These other parameterizations are described in the online supplemental material.
Sea salt aerosol has been parameterized previously using linear (Lehahn et al. 2010; Smirnov et al. 2003; Mulcahy et al. 2008; Glantz et al. 2009; Huang et al. 2010), power-law (Mulcahy et al. 2008; Glantz et al. 2009), and exponential (Moorthy and Satheesh 2000) relationships to 10-m or surface wind speeds from in situ and remote sensing measurements. It is well known that, to first order, sea salt emission is primarily dependent on wind speed (Jaeglé et al. 2011). However, other factors including sea surface temperature (SST), organic material, vertical mixing, advection, wet and dry deposition, and relative humidity could impact sea salt emission to a lesser extent. The second most commonly used predictor, after wind speed, is SST (Jaeglé et al. 2011; Gong 2003). However, laboratory and in situ investigations of the SST dependence of marine aerosol have had contradictory results, indicating that we do not yet fully understand how temperature impacts SSA (Grythe et al. 2014; Mårtensson et al. 2003; Zábori et al. 2012). Therefore, we have chosen to retain a parameterization that is only based upon surface winds.
c. over land
The over land for the time period between 2001 and 2018 was estimated on the basis of the methods of Ginoux et al. (2012) with minor modifications. Ginoux et al. (2012) used C005 MODIS Deep Blue AOD to identify pixels over land that contained dust by their column optical characteristics including AOD of greater than 0.1, AE of less than 0, SSA at 412 nm of less than 0.95, and SSA at 412 nm that is less than or equal to SSA at 660 nm. The AE criterion selects for pixels with a dominant coarse mode to the aerosol size distribution (Ginoux et al. 2012), which should primarily be dust and marine aerosol (Dubovik et al. 2002). The criterion for SSA at 412 nm is related to absorption of solar radiation in the green channel, characteristic of dust, and is effective at removing influences of marine aerosol. The difference in SSA between 660 and 412 nm is used because an optical property of dust is a positive spectra variation of SSA with wavelength, which manifests as a sharp decrease in absorption from the red (660 nm) to blue (412 nm) channels. This criterion is effective at removing any additional influence of biomass burning, because biomass-burning particles have spectral variation with the opposite slope. If a retrieval met all of these criteria, its Deep Blue AOD contributed to the long-term mean of τd at that grid cell. To test the thresholds used in Ginoux et al. (2012), we collected these same data within three overland regions where the aerosol signal is dominated by dust, and one overland region where mineral aerosols are unlikely to be present, and binned the pixels of this data for each location by their optical properties using the time period 2006 to 2010. The boundaries of these regions were chosen visually, selecting for either desert surface boundaries or for vegetation cover (Fig. S8 in the online supplemental material). The dust-dominated locations included a region over the Sahara Desert (10°–30°E, 19°–30°N) (Fig. 2a), a region over the Arabian Desert (40°–50°E, 23°–35°N) (Fig. 2b), and a region over the Taklamakan Desert (77°–88°E, 36°–40°N) (Fig. 2c). The non-dust-dominated location chosen was over equatorial Africa (10°–32°E, 10°S–5°N) (Fig. 2d).
In all three dust-dominated regions, 97%–100% of pixels meet the criterion for SSA at 412 nm < 0.95; in the region dominated by other aerosol types, fewer than 26% of pixels meet this criterion. Similarly, in each dust-dominated region, more than 99% of pixels meet the criterion of SSA at 412 nm ≥ SSA at 660 nm, and only 40% meet this criterion for the nondust region. The criterion used in Ginoux et al. (2012) of AOD greater than 0.1 was not used in the creation of here, because we were concerned that this threshold would systematically bias our results.
Ginoux et al. (2012) used a conservative requirement of AE less than 0 to identify dust. However, the updated algorithm used to create the C006 Deep Blue AE product now limits the valid range of values to 0 ≤ AE ≤ 1.8 (Sayer et al. 2013; Hsu et al. 2013). Therefore, we chose a new threshold that is based on the distribution of AE in our three desert regions. When this threshold was replaced with a criterion of AE less than 1, 90% of pixels in the Sahara Desert region, 91% of pixels in the Arabian Desert region, and 73% of pixels in the Taklamakan Desert region met the new threshold and only 5% met it in the nondust region. The new threshold of 1 we have used here is within the values of AE reported for dust regimes in Dubovik et al. (2002). In summary, daily estimates of over land grid cells are composed of the Deep Blue AOD if that pixel met the following criteria:
MODIS C061 AE 470–670 nm, which is inversely proportional to particle size, of less than 1,
MODIS C061 SSA at 412 nm, the ratio of aerosol scattering to extinction coefficients at that wavelength, of less than 0.95, and
MODIS C061 SSA at 412 nm is less than or equal to SSA at 660 nm.
If a pixel did not meet these criteria, its τd was set equal to zero. An analysis of the accuracy of this method, through comparison with AERONET, is presented in the online supplemental materials.
The over the ocean from 1981 to 2018 was derived similarly to what was previously described for MODIS. However, AOT from AVHRR was used rather than MODIS AOD, and the climatological seasonal mean f from the MODIS period, calculated from the MODIS small-mode 660-nm AOD divided by the total 660-nm AOD, was used in Eq. (3) because there is a spectral dependence of f. Evan and Mukhopadhyay (2010) showed that over the tropical North Atlantic, a dusty region, the seasonal cycle in f is far larger than the interannual variability, such that monthly climatological mean f can be used in place of monthly-mean f without introducing significant biases. Here we also assume that the use of monthly climatological f can be used in place of daily f. We tested this assumption by calculating the using climatological f. The results were similar to calculated using daily f with a bias of 0.00 in the long-term global mean over the ocean for the period 2001–17 (RMSE: 0.03). The characteristic fine-mode fractions, fm, fd, and fan, were the same as those used in the calculation of .
Following Evan and Mukhopadhyay (2010), a time series of zonal mean stratospheric aerosol optical thickness τs as estimated in Sato et al. (1993) was removed from AVHRR AOT before calculation of for years prior to 1999 and was assumed to be zero in subsequent years. This was done to remove the impact of large volcanic eruptions, including Pinatubo (June 1991) and El Chichón (March and April 1982), on . τs was not removed in years following 1999 because there have been no volcanic eruptions of the scale of Pinatubo or El Chichón since 1999. The total column AOT peaks immediately following these major eruptions, but τs does not peak until 1–3 months later (Fig. 3). This offset results in a remaining peak in the τs-corrected AOT, which leads to a coincident peak in . In our analysis of trends in , we exclude 1982 and 1991 to avoid the impact of these events.
Our resultant global mean time series of exhibits fluctuations in years prior to 1995 that are associated with satellite drift, a known issue with the pre–K, L, and M series (KLM) NOAA satellites, which has been investigated with respect to long-term variability in cloud cover (Foster and Heidinger 2013; Norris and Evan 2015). Satellite drift, or orbital decay, results in a progressively later equatorial crossing time. Thus, biases in optical depth retrievals may be due to retrieval errors associated with larger solar zenith angles (Fig. 3), or sampling of the diurnal cycle (Norris and Evan 2015). We remove the signal associated with satellite drift via linear regression, following the methods of Norris and Evan (2015). Consequently, at each location, the resultant values are anomalies relative to the global mean.
The primary results of this work are two nearly global, daily, observational datasets of τd at a 1° horizontal resolution that can be used for studies of the global dust cycle. The long-term mean global over the ocean was found to be 0.03 ± 0.06 where ±0.06 is the ±1σ error estimate. For comparison, Ridley et al. (2016) found global dust AOD at 550 nm to be 0.03 ± 0.005 using in situ and satellite observations combined with global models. The long-term mean global over land was found to be 0.1. The remainder of this section will be used to present a comparison of and (section 3a), verification of against ground-based measurements of dust (section 3b), and a description of seasonal climatological mean τd (section 3c).
a. Comparison of and
Globally, monthly mean estimates of and agree well, as evidenced by a high correlation coefficient between the two data (r2 = 0.41), a linear least squares regression slope near unity [y = 0.90(±0.001)x], a low RMSE (0.07), and a low offset (−0.01) (Fig. 4a). Much of the disagreement between the two occurs in locations where there are few days on which both AVHRR and MODIS have retrievals (Fig. 5). This is often the case in latitudes above 30°N and below 30°S as well as in characteristically cloudy regions such as the ITCZ (0°–10°N), the southeastern Atlantic Ocean (10°–20°S, 10°W–10°E), and the southeastern boundary of the Pacific Ocean. In dust-dominated regions, and are very well correlated, but with a slope that is greater than 1. For example, in the tropical North Atlantic, the only ocean region where mean is greater than 0.3 in all seasons (Fig. 10, described in more detail below), they are linearly related with a form y = 1.21(±0.005)x − 0.03 (r2 = 0.76, RMSE: 0.08, and offset: −0.00) (Fig. 4b). Monthly mean is offset low (bias: −0.18) relative to for greater than 1. We hypothesize that the low offset in is an artifact of misclassification of optically thick dust layers as clouds in the AVHRR AOT cloud clearing algorithm (X. Zhao 2018, personal communication), although it is also contributed to by the difference wavelengths between MODIS AOD retrievals and AVHRR AOT retrievals (Zhao et al. 2008). To determine boundaries around the potential contribution to the offset by the difference in wavelength of AOD retrieval, we perform a back-of-the-envelope calculation. The largest percentage difference in AOD at each wavelength should occur when the AE is large. We use AERONET 500- and 675-nm AOD at the (dust dominated) Tamanrassett station to calculate the mean and upper and lower quartiles of AE for these wavelengths, using only days that met the criteria used for ground-truth validation of dust (see the online supplemental material). In the absence of a 550–630-nm AE from MODIS, we assume that this 500–675-nm AE from AERONET is similar enough to approximate the influence of wavelength differences on the offset between MODIS and AVHRR. Using the lower-quartile AE at Tamaransett, 0.12, we find that the wavelength difference would contribute a 1.6% offset with greater than . Using the mean AE at Tamaransett, we find a potential offset of 3.5%; using the upper quartile, it is 5.8%.
We also calculated the linear least squares best fit between the daily and daily , including only dates and pixels where both sensors retrieved AOD (Fig. 4d). We expect that if the low offset of relative to is contributed to by missing data in regions of optically thick dust, the slope of the regression between these datasets will be significantly less steep for the daily data than for the monthly data in the tropical North Atlantic (Fig. 4d). This is because for the monthly mean, MODIS would incorporate high τd days that AVHRR “missed.” Indeed, we found that the daily and daily are related by the relationship y = 0.95(±0.0002)x − 0.01 (r2 = 0.54, RMSE: 0.13, and offset: −0.01) in the tropical North Atlantic. Globally, daily and daily are related by the relationship y = 0.70(±0.0003)x + 0.01 (r2 = 0.26, RMSE: 0.09, and offset: −0.01) (Fig. 4c).
The discrepancy in retrievals over large dust plumes can be easily visualized in the case of 27 June 2014, over the tropical Atlantic (Fig. 6). In this case, indicates an optically thick dust plume, centered near 15°N, 30°W, that has been advected over the tropical Atlantic from the Sahara (Fig. 6a). This particularly strong event allowed for dust to be measured in situ as far as Colombia (Bedoya et al. 2016). In Fig. 6, grid cells that are missing estimates of τd or AOD, or grid cells over land, are masked in gray. For the same day, the (Fig. 6b) is largely missing estimates in the region of optically thick dust. This is a result of missing retrievals of AVHRR AOT (Fig. 6d), which occur due to a combination of sun glint and because of the cloud-screening algorithm classifying optically thick dust as cloud.
b. Comparison with existing datasets
It is informative to compare our τd with previously published estimates of atmospheric dust. Here, we compare with two recently published time series of dust concentrations over the western United States (Tong et al. 2017; Hand et al. 2016) and one over the Caribbean (Prospero and Nees 1986; Prospero et al. 1996). The measurements of dust over the western United States are based on surface concentrations of aerosol species. In Hand et al. (2016), iron from IMPROVE stations was used as a proxy for fine dust concentrations and compares well to these annual mean fine dust concentrations (r value = 0.61; p value < 0.01) (Fig. 7). In Tong et al. (2017) several criteria from ground-based IMPROVE data were used to confirm dust events that had been detected in satellite imagery. The criteria included high PM10 and PM2.5 concentrations, a low ratio of PM2.5 and PM10, high concentrations of crustal elements (SI, CA, K, Fe, and Ti), low concentration of anthropogenic components, and low enrichment factors of pollution elements. Annual mean is strongly correlated with estimates of annual number of dust events (r value = 0.73; p value < 0.05) from this work. Toth et al. (2014) found that globally, the average correlation between total AOD and AOD integrated from 500 m AGL to the surface was 0.61 over land. For the contiguous United States, this correlation was found to be 0.62, with 0.57 over the western time zone. Therefore, although τd is a measure of column integrated aerosol, we still expect surface concentrations to be correlated with column integrated values, especially when averaged over a large region. Interestingly, while an increasing trend in dust was noted for these studies for the period between 2001 and 2015, from these results, the extension of the time series through 2017 indicates that dust has been decreasing since 2011.
We also compared our new climatology of with long-term, independent, ground-based dust measurements of dust taken at Ragged Point, Barbados (Prospero and Nees 1986; Prospero et al. 1996; Prospero 2014). These measurements are based on ash residue from extracted filters with an adjustment factor based on average crustal abundance of Aluminum in soil dust. These data have a standard error that is approximately constant at ±0.1 μ m−3 for concentrations less than 1 μg m−3 and about ±10% for higher concentrations. We found that in addition to matching the seasonal cycle of dust at this location (Fig. 8b), monthly mean over Barbados is positively and statistically significantly correlated (r value = 0.50; p value < 0.01) with monthly mean ground-based dust aerosol measurements when the seasonal cycle is removed from both time series (Fig. 8a). Monthly mean from 1982 to 2013 is also statistically correlated (r value = 0.33) with dust concentrations at Barbados, and shows the same seasonality there as (Fig. S3 in the online supplemental material).
We also compared our dataset with a monthly time series of dust optical depth (DOD) over Syria as published by Pu and Ginoux (2016) (Fig. 9). Similar to Ginoux et al. (2012), Pu and Ginoux (2016) used MODIS Deep Blue C006 AOD, SSA, and AE to construct DOD, and thus the datasets are not independent. However, Pu and Ginoux (2016) updated their method by interpolating optical properties to a 0.1° × 0.1° grid and using a threshold of SSA at 470 nm along with an empirical relationship between AE and f to derive coarse DOD. The time series are related with a correlation coefficient of 0.87 (p value < 0.01) and show the same trend of increasing dust over Syria (Fig. 9). We are not aware of another independent estimate of against which we can make a comparison.
c. Seasonal dust AOD
Having shown that this satellite-based climatology of dust appears to reflect changes in the atmospheric concentration of dust, we next summarize some of the broader characteristics of the global dust cycle, starting with seasonality. Maxima in seasonal mean are found over the Sahara and the Arabian Peninsula in austral summer (July–August) (Fig. 10d) and over the Taklamakan and Gobi Deserts during spring (March–May) (Fig. 10c). Aeolian transport from the Sahara toward the Americas is visible during all seasons with a maximum extension westward during the summer (Figs. 10d,h), which is well known (Prospero and Mayol-Bracero 2013; Prospero et al. 1996). Seasonal mean , appears far lower than near to major dust sources in all seasons, due to limited retrievals in regions of optically thick dust as described in the previous section. The spring maximum in dust over Asia can be attributed to the seasonal cycle of the Mongolian cyclonic depression and the related high frequency of springtime cold air outbreaks (Sun et al. 2001; Ge et al. 2014). During this season, dust is transported from the Asian continent over the North Pacific Ocean toward North America at latitudes between 30° and 50°N. It is presumed that this dust is primarily from the Taklimakan Desert as the local topographical and meteorological conditions allow for it to be entrained to elevations exceeding 5 km while dust from the Gobi Desert is, for the most part, confined closer to the surface (<3 km) (Sun et al. 2001). It has also been shown that some portion of this trans-Pacific dust originates from Africa (Creamean et al. 2013).
Discrepancies between over land and over adjacent ocean are visible in some areas with significant anthropogenic aerosol sources, including the Sahel region of Africa and India during austral autumn and spring (Figs. 10a,c). This may be due to different methods and data sources used for estimating over land and ocean. The method used for over the ocean leads to a continuous range of , whereas the threshold approach used for estimation of over land may prevent these gradients when there is mixing between anthropogenic aerosols and dust aerosol. However, this may also partially be the result of the different algorithms used in the retrieval of AOD over land and ocean (Sayer et al. 2013; Hsu et al. 2013) (Fig. S4 in the online supplemental material). While these unphysical discrepancies could be avoided through spatial interpolation across these boundaries, we have chosen not to do this in order to minimize our assumptions and retain a daily dataset that is directly representative of the satellite observations. At the very least, these results suggest that direct comparison of overland and overwater τd may be problematic.
A region of high summertime mean is visible over the North Pacific Ocean near 45°N (Fig. 10d). Summertime trans-Pacific transport of dust has been observed (Yumimoto et al. 2010). However, visual inspection of MODIS Terra corrected reflectance on days when over the western North Pacific was high during summer months indicated the possibility of contamination from biomass burning over Europe and Russia. We compared summer seasonal mean and MERRA-2 organic carbon extinction optical thickness at 550 nm (online supplemental Fig. S5; https://giovanni.gsfc.nasa.gov/) and confirmed that summer seasonal mean organic carbon was highest over the western North Pacific in years when summer seasonal mean was highest. This supports the hypothesis that there may be biomass-burning contamination over the western North Pacific during these months. The fine-mode fraction f during these summer high events is approximately 0.57–0.63. In the springtime, most of the high events are associated with fine-mode fractions near 0.50, which is closer to fd.
In addition, because of extensive low cloud cover this region in the western North Pacific has few retrievals during the summer months. The limited number of retrievals is evident in the number of pixels used in the estimate of MODIS Dark Target AOD during the summer in that region. While greater than 36 pixels are incorporated on average in the equatorial Pacific estimations of MODIS Dark Target AOD, only approximately 26 on average are used in this estimation in parts of the western North Pacific. Long-term seasonal mean appears low in the western North Pacific in all seasons (Figs. 10c–f).
a. Trends in
In Fig. 11a we show regional trends in monthly mean from Terra between 2001 and 2018, where the seasonal cycle has been removed. Only trends that are significant at the 95% confidence level are shown. One region that stands out with a large upward trend in dust is the Arabian Peninsula, where values exceed 0.3 (10 yr)−1 for the MODIS period. This is consistent with previous studies that found increasing trends, of smaller magnitude, in τ in this region (Hsu et al. 2012; Alfaro-Contreras et al. 2017; Toth et al. 2016). Parts of the Saharan Desert appear to be getting dustier in from Aqua. However, coastal Northern Africa and the edges of the Sahara display decreasing trends in from Terra (Fig. 11). Our results show a statistically significant decrease in over northern China (approximately 100°E, 40°N), and over southwest Asia (approximately 180°E, 30°N). This is consistent with recent findings by Pandey et al. (2017), who found that premonsoon dust loading decreased between 2000 and 2015. This decrease has been linked to increased rainfall, and resultant wet scavenging and increased soil moisture, and, to a lesser extent, decreased wind strength. The decrease in over northeast Asia may also be the continuation of a documented decreasing trend from the late 1950s through the 1980s in the region due to reduced cyclone frequency (Qian et al. 2002).
Interestingly, the greatest increase in (1.04 units per 10 yr) is seen over the Aral Sea at the intersection of Kazakhstan and Uzbekistan (44.5°N, 59.5°E). This was once the location of one of the world’s largest lakes, spanning over 66 000 km2 (Izhitskiy et al. 2016). However, the lake has been drying up since the early 1960s, when irrigation projects diverted the rivers that fed it leaving behind a highly emissive exposed playa (Micklin 1988; Indoitu et al. 2015). The eastern basin of the lake dried up completely in 2014, but the northern parts of the sea remain filled and have largely stabilized. A dramatic increase in is also seen over Razzaza Lake in Iraq, where the shoreline is receding due to diverted irrigation water. The northern-most region of the Caspian Sea, which is also the shallowest portion of the sea, shows a smaller but significant increase in potentially due to the recession of this body of water over the past two decades (Chen et al. 2017).
τd derived from MODIS Aqua (2003–18) displays the same trends of decreasing dust over northeast China and southwest Asia and increasing trends over the Arabian Peninsula (approximately 20°N, 70°E) and central Sahara (approximately 20°N, 10°E) (Fig. 11b). There are some regional changes that are more pronounced in this record, including a statistically significant increase in τd over Southern California (approximately 34°N, 119°W). There is also a large increase in dust over Oman (0.20; approximately 20°N, 57°E) that is not present in the derived from Terra and a greater increase over Iran than is present in the derived from Terra. Differences in regional trends between these two datasets could be due the differences in the length of record, as well as the differences in overpass time for the two satellites from which the datasets are derived as it relates to the daily cycle of dust emission (Heinold et al. 2013; Kocha et al. 2013). Terra has a morning orbit while Aqua has an afternoon orbit (Ichoku 2005). However, it has been shown that although there is no significant overall trend, the offset between the Terra and Aqua C6 total AOD records has oscillated direction over time, and its variability has increased (Levy et al. 2018). This is likely to impact the differences in from Aqua and Terra.
Trends in overwater τd, from 1981 to 2018, were also calculated using monthly mean (Fig. 11c). Over these longer time periods, there is a small increasing trend in near the Arabian Peninsula [0.02 (10 yr)−1] and a decreasing trend in in the tropical North Atlantic near to the African coast [0.01 (10 yr)−1]. The decreasing trend over offshore of northern Africa is consistent with the time series of dust at Cape Verde (16.5°N, 23.04°W) created using a similar method (Evan and Mukhopadhyay 2010), and with another study that utilized satellite and coastal surface measurements (Chin et al. 2014).
We note that Zhao et al. (2013) found a slight positive trend in the CDR AVHRR AOT record that is the result of residual cloud contamination. As such, trends in the tropical North Atlantic and Arabian Peninsula could be similarly biased.
b. Comparison with a modern reanalysis dust product
To understand whether this observation-based dataset is consistent with a modern reanalysis dust product, we compared with MERRA-2 dust extinction aerosol optical thickness at 550 nm (https://giovanni.gsfc.nasa.gov/giovanni/) over two of the regions—Barbados and the southwestern United States—used previously in our comparisons with ground observations (Figs. 7 and 8). Our monthly mean , with the seasonal cycle removed, was significantly correlated at the 99% confidence level (r value = 0.73) to the monthly mean MERRA-2 dust extinction aerosol optical thickness over Barbados (Fig. 12a). However, the annual mean displays greater interannual variability and was much more similar to the ground observational record than the MERRA-2 dust extinction aerosol optical thickness over the southwestern United States (Fig. 12b). Over this location, the MERRA-2 dust record was not statistically significantly correlated (r value = 0.20) with .
These global, daily, observation-based datasets of dust aerosol optical depth are rooted in observations of optical characteristics of the aerosol column from Terra Moderate resolution Imaging Spectroradiometer and from the Advanced Very High Resolution Radiometer, but incorporate observations of characteristic fine-mode fraction in dusty, clean marine, and anthropogenically dominated regions from the Aerosol Robotic Network as well as surface wind speed from the second Modern-Era Retrospective Analysis for Research and Applications. As previously mentioned, the long-term mean global over the ocean was found to be 0.03 ± 0.06. Although there are few measurements of dust aerosol for comparison, over land compares well to ground-based measurements of dust from the Interagency Monitoring of Protected Visual Environments network in the western United States and over ocean compares well to independent ground-based measurements made at Ragged Point.
We had hoped to create a reliable dataset of extending from 1981 to 2018 using AVHRR aerosol optical thickness from the Climate Data Record AOT dataset, but we found that the AVHRR AOT had limited coverage in dusty regions. Comparison of daily AVHRR AOT and MODIS AOD over the tropical Atlantic indicated that this limited coverage is at least partially explained by optically thick dust identified as cloud in the AVHRR cloud masking process. This limited coverage results in a low offset in seasonal and monthly mean when compared with . We found that increased in the central Sahara and the Middle East and decreased over northern China and southwestern Asia between 2001 and 2018. The largest increase in τd occurred over the Aral Sea, which has receded, leaving emissive playa surfaces exposed. We found decreasing trends in near Africa in the equatorial Atlantic and increasing trends near the Arabian Peninsula from 1981 to 2018.
In 2011, the Suomi National Polar-Orbiting Partnership (SNPP) satellite was launched carrying the Visible Infrared Imaging Radiometer Suite (VIIRS), which is intended to continue the observational records that have been provided by sensors such as MODIS and AVHRR (Sayer et al. 2017). VIIRS has been shown to be capable of identifying aerosol type and partitioning between fine- and coarse-mode aerosol. As such, observational studies that require a shorter record (beginning in 2012) may be able to rely on these products. In addition, Chen et al. (2018) demonstrated that an inversion technique, utilizing the GOES-Chem model, can be applied to certain satellite aerosol optical property retrievals in order to retrieve desert dust. This was applied to one year of satellite-derived aerosol information generated by the General Retrieval of Atmosphere and Surface Properties (GRASP) algorithm from Polarization and Anisotropy of Reflectances for Atmospheric Sciences coupled with Observations from a Lidar (POLDER/PARASOL) retrievals. This technique can provide daily dust at 2° × 2° resolution with uncertainty below 25.8%, but it is computationally costly and has not been applied to longer records of remotely sensed aerosol optical properties.
The current work provides daily τd generated from a simple algorithm and observed aerosol data. Previous estimates of τd have enhanced our understanding of the global dust cycle and have been a valuable tool for evaluating the performance of global models (Pu and Ginoux 2018a; Evan and Mukhopadhyay 2010; Kaufman 2005; Ginoux et al. 2012). This longer, higher-temporal-resolution record will facilitate additional advances in our understanding of the dust cycle.
This work was funded by the California Department of Water Resources Contract 4600010378, Task Order OSCOP215, and the Army Corps of Engineers USACE (CESU) W912HZ–15–0019. We thank the AERONET principal investigators and their staff for establishing and maintaining the 39 sites used in this investigation. We are grateful for the datasets and data-archiving centers that supported this work and appreciate those who made our study possible, including the MERRA-2 team at the GMAO and staff at GSFC, and the MODIS teams at NASA, as well as Tom Zhao and the NOAA Climate Data Record team. We thank Drs. Joseph Prospero, Jennifer Hand, and Daniel Tong for their contributions of long-term datasets for comparisons made in this work, and Drs. Martin Ralph, Leah Campbell, and Nora Mascioli for helpful comments on this paper. We also thank three anonymous reviewers for helpful comments on this paper. The datasets described in this paper will be made publicly available at Pangea Open Access (https://doi.pangaea.de/10.1594/PANGAEA.909140).
Supplemental information related to this paper is available at the Journals Online website: https://doi.org/10.1175/JAMC-D-19-0194.s1.