Evaluating MODIS Dust-Detection Indices over the Arabian Peninsula

Sand and dust storm events (SDEs), which result from strong surface winds in arid and semi-arid areas, exhibiting loose dry soil surfaces are detrimental to human health, agricultural land, infrastructure, and transport. The accurate detection of near-surface dust is crucial for quantifying the spatial and temporal occurrence of SDEs globally. The Arabian Peninsula is an important source region for global dust due to the presence of extensive deserts. This paper evaluates the suitability of five different MODIS-based methods for detecting airborne dust over the Arabian Peninsula: (a) Normalized Difference Dust Index (NDDI); (b) Brightness Temperature Difference (BTD) (31–32); (c) BTD (20–31); (d) Middle East Dust Index (MEDI) and (e) Reflective Solar Band (RSB). We derive detection thresholds for each index by comparing observed values for ‘dust-present’ versus ‘dust-free’ conditions, taking into account various land cover settings and analyzing associated temporal trends. Our results suggest that the BTD (31–32) method and the RSB index are the most suitable indices for detecting dust storms over different land-cover types across the Arabian Peninsula. The NDDI and BTD (20–31) methods have limitations in identifying dust over multiple land-cover types. Furthermore, the MEDI has been found to be unsuitable for detecting dust in the study area across all land-cover types.


Introduction
Dust storms and sandstorms are defined as events during which constituent particles of dust and sand are raised to altitudes of up to 3000 m by strong winds [1][2][3].Dust storms have been shown to have direct detrimental effects on human health, the economy and the environment [4][5][6].An example of their negative impact on health is the associated increase in instances of respiratory diseases [7].Dust storms may also cause an assortment of other issues, such as reduced visibility, which limits different activities and raises the risk of fatal transport crashes.
The Arabian Peninsula, classified as an arid and semi-arid climate [8,9], is characterized by vast sand and gravel deserts located at high-elevation plateaus.These factors mean that sand and dust storms constitute a significant natural hazard to communities across the Arabian Peninsula.One of the most damaging dust storms in the last decade swept over Riyadh city on 10 March 2009, resulting in dense airborne dust, causing zero visibility and a shutdown of the airport.In addition, buildings, vehicles, electricity poles and trees were damaged during this event [10][11][12].During another event, on 18 March 2012, a dust storm swept over the Arabian Peninsula, closing schools across Saudi Arabia and sending many people to hospitals with breathing problems [13].
Despite the large detrimental effects these airborne dust events have on the region, a comprehensive study of the spatio-temporal incidence of dust storms required to understand their key drivers has yet to be undertaken.Therefore, the main aim of this paper is to evaluate the suitability and efficiency of extant satellite dust methods over the Arabian Peninsula.In addition, we aim to quantify how well these different tests perform considering a range of land-use type.

Previous Work
Traditionally, dust storms have been detected using observations made at fixed monitoring stations [14].However, the discrete locations of these in-situ observations result in sparsely sampled records, especially because large areas of the Arabian Peninsula such as the Rub'al Khali desert, which covers approximately 560,000 km 2 [9], do not contain monitoring stations.Where these stations do exist, records are made at high temporal resolution, typically every three hours and, in recent years, every hour, but, in practice, the archives tend to be incomplete due to technical failures and/or when the observing system is being updated or maintained.Hence, there is keen interest in using satellite images to identify dust storms, because such images have the potential for providing regularly repeated, consistent observations of airborne dust events over the entire Arabian Peninsula.As a consequence, satellite observations can provide observations of the instantaneous extent of a dust storm, which in-situ measurements cannot.
Satellite observations are a powerful tool for examining the characteristics of large-scale dust storms.More precisely, the application of indices and algorithms for detecting airborne dust at the time the satellite passes over provides instantaneous information about the location and extent of individual dust storms.Previous studies have used a range of satellite observations to remotely detect airborne dust such as the TOMS aerosol index to monitor dust storms [15][16][17][18][19].Other studies have used the SEVIRI instrument to detect and identify dust sources [20][21][22].While this instrument provides observations at very high temporal resolution (15 min), the high-altitude geostationary orbit and fixed acquisition geometry result in insufficient spatial resolution required to understand spatial variations in airborne dust and its drivers.Due to its lower orbit, wide swath and broad spectral bandwidth, the Moderate Resolution Imaging Spectroradiometer (MODIS) instrument, on board NASA's Terra and Aqua satellites, is particularly well-suited for detecting and monitoring airborne dust.The Moderate Resolution Imaging Spectroradiometer (MODIS) instrument is carried on board the Earth Observing System (EOS) Terra and Aqua satellites, which were launched in December 1999 and May 2002, respectively.The satellites follow polar orbit with a repeat period of 16 days, and have north-to-south equator crossing times of about 10:30 a.m. and 1:30 p.m. respectively.While MODIS gives almost full worldwide coverage each day, they do not provide the same level of high temporal sampling as e.g., the SEVIRI instrument on board the Geostationary orbiting Meteosat.This means MODIS observations of individual dust storm events are discontinuous, limiting the ability to observe how individual events evolve over time.However, this is compensated by the much higher spatial and spectral resolution of MODIS, which has 36 spectral channels and generates imagery at spatial resolutions 250 m, 500 m and 1 km, depending on the channel [23].Hence, several previous studies have developed and described dust-detection tests based on MODIS data [24][25][26][27][28].Consequently, the aim of this study is to evaluate which of the existing MODIS dust detection methods are most effective over the Arabian Peninsula.

Study Area
Figure 1 presents the location and the elevation map of the study area.Most of the Arabian Peninsula lies within a wide band of deserts; some of these dry lands are covered with drifting sand dunes [29].The high-elevation plateaus of the Arabian Peninsula are intersected by numerous shallow valleys, while the southern part of the peninsula is dominated by the Rub' al Khali (Empty Quarter) desert, the largest sand desert in the world, covering an area of approximately 560,000 km 2 and characterized by an average total annual precipitation of less than 100 mm [30].In addition, the large An Nafud desert is situated in the northern part of the peninsula with an area of 72,000 km 2 [9], while, in the east of Saudi Arabia, there is the Ad Dahna desert, which is about 1450 km 2 [31] (Figure 1).The dominance of sand-covered areas makes the Arabian Peninsula a challenging region to detect airborne dust.
Figure 1 presents the location and the elevation map of the study area.Most of the Arabian Peninsula lies within a wide band of deserts; some of these dry lands are covered with drifting sand dunes [29].The high-elevation plateaus of the Arabian Peninsula are intersected by numerous shallow valleys, while the southern part of the peninsula is dominated by the Rub' al Khali (Empty Quarter) desert, the largest sand desert in the world, covering an area of approximately 560,000 km 2 and characterized by an average total annual precipitation of less than 100 mm [30].In addition, the large An Nafud desert is situated in the northern part of the peninsula with an area of 72,000 km 2 [9], while, in the east of Saudi Arabia, there is the Ad Dahna desert, which is about 1450 km 2 [31] (Figure 1).The dominance of sand-covered areas makes the Arabian Peninsula a challenging region to detect airborne dust.

Figure 1.
Overview map of the Arabian Peninsula with annotation of the major deserts, cities and vegetated areas with a topographic map (grey scaled) as background.

Materials and Methods
In this paper, we analyzed 17 MODIS images capturing major dust storm events that occurred over the Arabian Peninsula between 2000 and 2015.More precisely, Saudi Arabia and the surrounding area covered by the bordering states located on the Arabian Peninsula were considered (as indicated in Figure 2).The a priori extent of these events were defined using visual interpretation of the MODIS true color (RGB) images (Figure 2), which was then validated using in-situ observations METAR (Meteorological Terminal Aviation Routine) observations issued at hourly or half-hourly intervals [32].Moreover, the same method was applied to the dust-free observations we used in-situ measurements to validate dust-free images.We compared each 'dust storm' image with a 'dust-free' MODIS image from a nearby point in time and in the same location with the exact extent of the dust storm as shown in Figure 2 in order to evaluate the various MODIS dust-detection

Materials and Methods
In this paper, we analyzed 17 MODIS images capturing major dust storm events that occurred over the Arabian Peninsula between 2000 and 2015.More precisely, Saudi Arabia and the surrounding area covered by the bordering states located on the Arabian Peninsula were considered (as indicated in Figure 2).The a priori extent of these events were defined using visual interpretation of the MODIS true color (RGB) images (Figure 2), which was then validated using in-situ observations METAR (Meteorological Terminal Aviation Routine) observations issued at hourly or half-hourly intervals [32].Moreover, the same method was applied to the dust-free observations we used in-situ measurements to validate dust-free images.We compared each 'dust storm' image with a 'dust-free' MODIS image from a nearby point in time and in the same location with the exact extent of the dust storm as shown in Figure 2 in order to evaluate the various MODIS dust-detection approaches.Moreover, we evaluated the suitability of dust-detection indices for each land-cover type.We used MODIS Level 1B Calibrated Radiances 1 km (MOD021KM) observations for our analysis.Once downloaded, we converted the pixel values using the Top of Atmosphere (TOA) reflectance calibration for the reflective bands, and the TOA brightness temperature calibration for the thermal emissive bands.Hence, each MODIS image was corrected after the atmospheric and solar angle effects were removed accurately (but not atmospheric distortions or dust).These MODIS-calibrated reflectance data were computed using the MCTK: the MODIS Conversion Toolkit plugin for ENVI [33].
The present study aims to define the suitability and threshold for each of the considered dust detection indices.In this context, we also aim to reveal the possible existence of temporal trends.Moreover, we compared the results with in-situ observations of airborne dust from Saudi Arabia meteorological observation stations as reported in METAR aviation weather reports issued at hourly or half-hourly intervals.The main aim of this validation exercise is to obtain a general understanding of the effectiveness and suitability of the studied dust-detection indices, and therefore, it has only been applied to one major dust storm event as an example.
approaches.Moreover, we evaluated the suitability of dust-detection indices for each land-cover type.We used MODIS Level 1B Calibrated Radiances 1 km (MOD021KM) observations for our analysis.Once downloaded, we converted the pixel values using the Top of Atmosphere (TOA) reflectance calibration for the reflective bands, and the TOA brightness temperature calibration for the thermal emissive bands.Hence, each MODIS image was corrected after the atmospheric and solar angle effects were removed accurately (but not atmospheric distortions or dust).These MODIS-calibrated reflectance data were computed using the MCTK: the MODIS Conversion Toolkit plugin for ENVI [33].
The present study aims to define the suitability and threshold for each of the considered dust detection indices.In this context, we also aim to reveal the possible existence of temporal trends.Moreover, we compared the results with in-situ observations of airborne dust from Saudi Arabia meteorological observation stations as reported in METAR aviation weather reports issued at hourly or half-hourly intervals.The main aim of this validation exercise is to obtain a general understanding of the effectiveness and suitability of the studied dust-detection indices, and therefore, it has only been applied to one major dust storm event as an example.

Moderate Resolution Imaging Spectroradiometer (MODIS) Images
MODIS satellite observations have been used to detect airborne dust in previous studies.In 1997, Ackerman introduced the combination of the Brightness Temperature (BT) method for multiple channels.As this method considers the emissive properties of dust, it is sensitive to variations in chemical structure and size distribution [34,35].Multiple studies have attempted to further develop this index by introducing more robust versions based on the Brightness Temperature Difference (BTD) between two wavelength channels (e.g., [36,37]).However, during recent years, researchers have considered more advanced techniques focusing on dust storm detection through the use of various indices, such as the Normalized Difference Dust Index (NDDI) [26,38] and the Middle East Dust Index (MEDI) [24,39].

MODIS Cloud Mask
As clouds affect the successful detection of dust storms, we used the MODIS cloud mask Level 2 product (MOD35) to indicate the presence of cloudy pixels.The reason for using the MODIS cloud mask is to restrict the analysis to areas free of clouds [40].However, the cloud mask tends to misclassify thick dust storms as clouds.This might affect dust detection [41].

MODIS Land Cover
Most previous work has only distinguished between bright land cover and dark land cover (Barren-Vegetation); however, the performance of the observation in reflective bands are different over various land cover types [27].In this study, we go further by evaluating how the effectiveness of the various dust detection tests varies between urban, barren and vegetated areas within each region.We used the MODIS Land Cover Type product, MCD12Q1 at 500 m resolution [42] in order to evaluate which of the considered indices is the most effective test for detecting airborne dust over these different land-cover types.

MODIS
We evaluated the following existing MODIS dust detection approaches: • The normalised difference dust index: The NDDI, used in several previous studies [26,27,38], is calculated using reflective solar bands 3 and 7 of MODIS (Equation ( 1)).The spatial resolution of these bands is 500 m and is characterized by a spectral reflectance at a wavelength range of 459-479 nm for band 3 and 2105-2155 nm for band 7.These bands are often used for land surface studies, for detection of aerosols and to obtain information on clouds' optical thickness, phase and effective radius.
• Brightness Temperature Difference (BTD) indices: BTD methods utilize differences between multiple thermal emissive bands and have been used successfully across a wide range of regions [24,27,43,44].The most commonly used thermal emissive bands are BT20, BT31 and BT 32, which are characterized by a wavelength range of 3.66-3.84µm, 10.78-11.28µm and 11.77-12.27µm, respectively, and a spatial resolution of 1000 m.The equations for the existing BTD-related indices considered in this study are given by Equations ( 2) and (3).The unit of these indices is Kelvin (K).
BTD (31)(32) • The Middle East Dust Index: As Karimi [24] found that bands 31 and 32 are having issues in distinguishing airborne dust from the desert surface, they presented a new method called the Middle East Dust Index (MEDI) [24,39,45].This index is based on Ackerman's dust detection technique but includes also band 29.Hence, the MEDI is calculated using thermal emissive bands BT29, BT31 and BT32 (Equation ( 2)) with wavelength ranges of respectively 8.40-8.70µm, 10.78-11.28µm, 11.77-12.27µm at a spatial resolution of 1000 m, and has the specific purpose of differentiating airborne dust and desert surfaces.
• The Reflective Solar Band-derived index: Finally, we also consider an index derived from the MODIS reflective solar bands B2 (841-876 nm) at a spatial resolution of 250 m and B18 (931-941 nm) at a spatial resolution of 1000 m (Equation ( 5)) as developed by Samadi [27].This index has been used over the western part of Iran to detect dust over bright and dark surfaces [27].

Statistical Analysis
We have analyzed differences in index values between pairs of dust and associated non-dust images closest in time via the Wilcoxon signed rank test for paired samples.More precisely, as the Kolmogorov-Smirnov test indicated that most of the considered datasets were non-normally distributed, we chose the Wilcoxon signed rank test for paired samples to conduct this analysis, as it is the non-parametric alternative of the more commonly used paired t-test.
To obtain a more profound understanding of the variability in index values under both dust and dust-free conditions, we also investigated the potential existence of temporal trends across the 17 observations between 2000 and 2015 by conducting regression analyses.This regression analyses allowed us to (i) detect any significant systematic changes in the considered index over time and (ii) assess the difference in the index value between dust and dust-free conditions in a temporal context.Model uncertainty, as well as the uncertainty on the slope parameter of the linear regression model, were used to evaluate the goodness of fit of the regression lines and to draw conclusions regarding the suitability of the indices in a temporal context.The threshold value for dust-storm detection determined whether the dust and non-dust observations are significantly different, enabling an evaluation of the suitability of the considered method as a function of time.More precisely, the threshold value was defined as either the upper or the lower 95% confidence interval line of the model error bound of the trend line fitted through the non-dust observations, whichever is closest to the trend line fitted through the dust storm observations.We applied the Wilcoxon tests and temporal regression analysis across the entire selected dust storm-affected areas considering the pixels regardless of their land-cover classification as well as separately stratified by land-cover type (i.e., urban, barren and vegetation).

Dust Detection Validation
We used in-situ observations of a dust storm over the Arabian Peninsula on 2 February 2012 to validate the MODIS dust-detection results.This intense dust storm initially occurred in the northeast of the peninsula, and moved later further towards the southwest (Figure 3).The in-situ observations from 10 sites over the Arabian Peninsula were used to validate the results obtained by applying the various MODIS dust-detection indices.In essence, by comparing the dust detection results with the observation data the effectiveness of the threshold in identifying dust was tested.Note that the cloud mask was used to eliminate the pixels for which the presence of cloud might affect the detection result.Figure 3 shows a MODIS true-color image illustrating the total extent of this particular dust storm event.Table 1 presents the outcome of the dust detection by using the threshold values determined in this paper.We compared the dust detection results with the observation data to test the effectiveness of the threshold in identifying the dust.Figure 4 presents the area within the total extent of the dust storm for which a given index indicated the presence of dust.Sub-panel (a) shows the result of the land cover independent dust detection index, whereas sub-panel (b) shows the outcome of the land cover specific dust detection indices.Within the latter, only the land use pixels of which the considered index is able to detect dust are displayed.
determined in this paper.We compared the dust detection results with the observation data to test the effectiveness of the threshold in identifying the dust.Figure 4 presents the area within the total extent of the dust storm for which a given index indicated the presence of dust.Sub-panel (a) shows the result of the land cover independent dust detection index, whereas sub-panel (b) shows the outcome of the land cover specific dust detection indices.Within the latter, only the land use pixels of which the considered index is able to detect dust are displayed.

Land Cover Independent Analysis
When considering all the land uses together, the results show that the NDDI values under dust and dust-free conditions are significantly different (Wilcoxon p value = 0.000, Table 2), suggesting that NDDI is a useful index for dust-detection.Nevertheless, Figure 5 shows a clear increasing trend in mean NDDI values over time under non-dust circumstances, as indicated by significant slope

Land Cover Independent Analysis
When considering all the land uses together, the results show that the NDDI values under dust and dust-free conditions are significantly different (Wilcoxon p value = 0.000, Table 2), suggesting that NDDI is a useful index for dust-detection.Nevertheless, Figure 5 shows a clear increasing trend in mean NDDI values over time under non-dust circumstances, as indicated by significant slope parameter values (p < 0.01, Table 3) representing yearly index increases of 0.0109 ± 0.0035 year −1 under non-dust conditions (Table 3).Figure 5 shows that the trend lines fitted through the dust and non-dust situations for the BTD indices are close to each other in the early years (2000)(2001)(2002)(2003)(2004), in which the 95% uncertainty bounds overlap.However, the Wilcoxon signed-rank test indicates a significant result with p values of 0.002 for BTD (31)(32) and for BTD (20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30)(31).Furthermore, the associated linear regression slope parameter for the BTD (31)(32) index has been found to be significant only at the 90% but not at the 95% level of confidence for the dust value (i.e., −0.0174 ± 0.0100 year −1 , Table 3), whereas, under dust-free conditions, this was found not to be significant.However, the threshold value to detect dust in the atmosphere for the BTD (31)(32) and BTD (20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30)(31) varies over time (between 2000 and 2016) from below ca.0.45 K to 0.7 K and above ca.−8 K to −9 K, respectively (Figure 5).Our results show that differences in MEDI values between dust and non-dust values are not significant when considering all land uses (i.e., a Wilcoxon signed-rank test with a p value 0.758; Table 2).In addition, the analysis did not show a significant linear regression slope parameter (Table 3) underlining the absence of a temporal trend in MEDI value as shown in Figure 5.The Wilcoxon test showed that the RSB index provides significantly different values for dust versus non-dust values detected (p = 0.005; Table 2).The trend line in the RSB index demonstrates clear discrimination between the dust and non-dust data-series (see Figure 5).

Land Cover-Specific Analysis
We investigated the inter-annual variability of dust and dust-free values over different land cover (Figure 5).This approach was chosen to evaluate which of the dust indices is the most useful test for detecting airborne dust over the different major land-cover types present in Saudi Arabia (i.e., urban, barren and vegetation), building on the findings of previous studies that found differences in dust-detection results depending on the land-cover settings [27].
Considering the NDDI index for urban areas, the trend line and the data overlap in recent years.This convergence of the lines seems to be the consequence of remarkably different slope parameter values of the linear regressions, which has been found to be significant at p < 0.1 (but not at p < 0.05) under dust conditions (i.e., 0.0097 ± 0.0058 year −1 , Table 3), but not significant under non-dust conditions (i.e., −0.0019 ± 0.0056 year −1 , Table 3).In addition, the Wilcoxon signed-rank test shows that the result is not significant at p ≤ 0.05 in NDDI values under dust versus non-dust circumstances (p value = 0.381; Table 2).
As shown in Figure 5, the high inter-annual variability in the index values in MEDI does not allow the determination of a clear threshold value for dust detection.Additionally, the p value is 0.210, which means it is not significant at p ≤ 0.05 for dust and non-dust values (Table 2).
Finally, the RSB index result shows a significant difference under dust versus non-dust conditions with a p value = 0.002), especially in recent years, as the difference between dust and non-dust image increases considerably over time in particular between 2012 and 2015.Moreover, the threshold value remains approximately constant at 0.175 (Figure 5).Also, our analysis shows a significant difference in the linear regression slope parameter only over dust conditions at p < 0.1 (but not p < 0.05) (i.e., 0.0051± 0.0030 year −1 ) (Table 3).
Figure 5 also presents the indices for vegetated areas.The associated differences in the NDDI index values under dust and the non-dust condition is significant (p value = 0.049; Table 2).In addition, the associated linear regression slope parameter values for this index are found to be not significant under dust free and dust conditions (i.e., 0.0057 ± 0.0086 year −1 and 0.0063 ± 0.0070 year −1 , respectively) (Table 3).
Our results suggest that the RSB index is effective at detecting dust when considering vegetated areas (with significantly different index values under dust versus non-dust conditions at a p value = 0.001; see Table 2).Moreover, the associated linear regression slope parameter has been found to be non-significant (see Table 3, with a threshold value of ca.0.16-0.18)(Figure 5).

Validation
From the validation results, as shown in Table 1, we can see that the NDDI index has limitations in distinguishing airborne dust from the desert areas.The BTD method reveals an ability to detect dust and distinguish dust from different lands; however, the BTD (31-32) did not detect dust over two locations (Table 1).Furthermore, the BTD index failed to detect dust because the barren area was falsely detected as a dust storm [46].In addition, there is some false detection over the urban areas (Figure 4).Moreover, the RSB method highlighted false detection over the vegetation area; however, it did detect 8% of the dust storm area (Table 1).Also, it is interesting to note that the BTD (20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30)(31) and RSB detect dust in the urban area (Bisha), which is not in the extent of the dust storm as shown in Figure 4. Our results show that each dust detection method has some limitations in detecting dust over different land covers.

NDDI
The threshold value to detect dust in the atmosphere (defined in this case as the lower 95% confidence interval line of the model error bound of the trend line fitted through the non-dust observations) increases over time from below ca.0.25 in 2000 to 0.4 in 2015 (Figure 5).Hence, the strong significant temporal trend in NDDI values illustrates that the NDDI can only be considered as a successful method to detect dust storms when a dust-free image from a period in time nearby (preferentially <1 year apart) is available.In addition, the dust in the NDDI index will be detected whenever the index value exceeds the identified threshold.

MEDI
As few studies have been published using the MEDI test (e.g., [24,39]), in our analysis, we found that band 29 has an error, which might have affected the results in some of these cases.However, the results clearly indicate that the use of the MEDI is not advisable for dust detection in Saudi Arabia due to the lack of a clear distinction between dust and not-dust values as indicated in Figure 5.It is worthwhile noting that this poor performance of this index was not expected because it is specifically designed to detect airborne dust in arid areas.

RSB
For this test, our results reveal a significant linear increase in index values under dust conditions (p < 0.1, 0.0062 ± 0.0033 year −1 ) and a significant increase in non-dust values (i.e., p < 0.05, 0.0054 ± 0.0026 year −1 ), which highlights the need for a contemporaneous dust-free image in order to be able to detect the presence of airborne dust.Therefore, the threshold value also increased over time from ca. 0.21 to 0.29 between 2000 and 2016 (Figure 5).Our analysis indicates that the RSB index can detect dust whenever the index value is over the fixed threshold.
The findings of the four dust detection indices (BTD methods, NDDI and RSB) confirm the temporal trends in the observed values.We suggest that these may indicate changes in prevailing atmospheric conditions and/or changing land use, although more research is required to identify the drivers behind these temporal trends.In addition, our analysis emphasizes that all indices require a cloud mask to avoid misclassification between the dust and clouds.Our analyses suggest that all MODIS indices except MEDI are useful in determining dust detection thresholds when considering all land cover types.However, during the period 2000-2003, the NDDI and BTD tests require a contemporaneous dust-free image to identify the optimal threshold.This is also the case for the RSB index after 2014 (Figure 5).

Urban
Our results demonstrate that the NDDI is not useful index in detecting dust over urban areas (Figure 5).On the contrary, the results seem to indicate that the BTD (31)(32) and BTD (20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30)(31) are efficient indices for detecting dust over urban areas, with significantly different index values under dust versus non-dust conditions.However, it is important to note that, in early years (2000)(2001)(2002)(2003), this index seems less suitable as the curve's 95% confidence intervals are overlapping.Furthermore, the results illustrate that the RSB index is effective at detecting dust over the urban area.Importantly, our findings suggest that no comparison with a dust-free image from a comparable period of time is required in order to detect dust.

Barren
The results of the barren area are similar to the results obtained when considering the overall land cover results because the majority of the Arabian Peninsula is barren (Figure 5).Hence, the same conclusions can be made as presented in the previous section, which considered the land cover independent analysis.Temporal trends in the NDDI index values under dust and non-dust storm conditions, as well as an associated threshold values, have been identified for the indices.The NDDI result shows that the associated linear regression slope parameter has been found to be significant for the non-dust images (i.e., 0.0101 ± 0.0031 year −1 ); therefore, the threshold value does not remain constant but varies from ca. 0.25 to 0.4 (Figure 5).Our findings agree with previous studies indicating that NDDI can be used to detect dust over barren areas (e.g., [38,48]).Our results reveal the existence of a temporal trend in this particular index value, which makes the NDDI only suitable when a non-dust image from a nearby moment in time is available between 2000-2003.However, in more recent MODIS images, there is no need to compare a dust-storm image with a non-dust image when attempting to detect airborne dust over the Arabian Peninsula using the NDDI.
When considering the BTD (31-32) output, our findings are in line with Samadi [27], showing a detection threshold increasing over time from ca. 0.6 K to 0.7 K (Figure 5).However, the associated linear regression slope parameter has been found to be significant for the dust values (p < 0.1, −0.0198 ± 0.01120 year −1 ) (Table 3).In addition, for the BTD (20-31) index, the associated slope parameter value for dust and dust-free circumstance has been found to be non-significant (Table 3), but the threshold values tend to decline over time from ca. −8 K to −9 K.
As mentioned in the section above, the MEDI index is unable to make a distinction between dust and non-dust circumstances for barren areas, as indicated by the Wilcoxon signed-rank test p values of 0.868 (Table 2).Hence, this analysis confirms the unsuitability of this index under this specific land cover.For the RSB index, the Wilcoxon test reveals significant different values between dust and non-dust values (p = 0.000; Table 2) However, the associated slope parameter values have been found to be significant (0.05 < p < 0.1) under both dust circumstances (i.e., 0.0060 ± 0.0033 year −1 ) and non-dust circumstances (0.0050 ± 0.0027), resulting in an increasing threshold value over time from ca. 0.2 to 0.3 (Figure 5), underlining the need for a non-dust image from a nearby moment in time in order to make this index suitable for dust detection over barren areas.

Vegetation
Our result shows that the NDDI values representing dust and non-dust conditions almost overlap in more recent years (especially from 2008 to 2015) (Figure 5).Hence, the converging character of these curves indicates that the suitability of this index is restricted to the early years only (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008).On the contrary, the BTD (20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30)(31) and MEDI methods have failed to distinguish dust values from non-dust values over the vegetated areas.Despite the significantly different index values under dust versus non-dust conditions (i.e., p value = 0.007; Table 2), the pattern of inter-annual variability in the index values indicates that the determination of a clear threshold value for dust detection is impossible.This is demonstrated by the fact that the uncertainty bounds are overlapping.However, our findings confirm previous studies (e.g., [27]) which demonstrate that RSB index can detect airborne dust over dark surfaces such as vegetation.
Unexpectedly, our results clearly show temporal trends for dust and non-dust values.Moreover, as can be seen in Figure 5, the RSB trends are mostly consistent across the different land cover in dust and dust-free values.In contrast, the other indices show that there is a difference in the trends observed for dust and non-dust conditions.Previous studies have noted that MODIS sensor degradation can introduce temporal artefacts into time series of MODIS C5 observations [49,50], we suggest this is unlikely in our study as we are using C6 data, which includes an improved calibration to take this into account [51,52].However, the observed temporal trends may indicate long-term changes in atmospheric conditions such as the increased temperature and the decrease in the precipitation over the Arabian Peninsula, which may alter the reflectance properties of the surface over time [53][54][55].This is supported by studies which have sown increased desertification in the study area [56][57][58].It is also possible that land use changes in the region have changed the surface reflectance properties, such as the expansion of the urban areas in Saudi Arabia [59,60].In addition, the various trends of the dust storm values may reflect variations in atmospheric conditions (i.e., temperature, wind speed, and precipitation), the human activities of urbanization and the changes in the feature of the land cover [61,62].While we note that these suggested temporal trends are very interesting and warrant further investigation, an in-depth analysis of them is beyond the scope of this study.

Conclusions
We have evaluated the suitability of multiple MODIS dust-detection indices across the entire Arabian Peninsula.The main finding of this study is that the BTD (31)(32) test and the RSB test are the most useful MODIS indices to detect airborne dust over different land cover settings.Furthermore, the NDDI is a suitable method to distinguish dust from non-dust circumstances over barren areas.Our results show unambiguously that the detection of airborne dust above urban areas is not possible when using the NDDI index.Also, the BTD (20-31) is a useful method, except for vegetation areas.Our findings suggest that the determination of a suitable dust-detection threshold, to be used with these indices, requires the use of dust-free 'reference images' from the same period.Further, the dust-validation analysis reveals that the BTD (31-32) method is a successful index for detecting dust storms over different land cover.However, more work is needed to improve the accuracy of detection thresholds as a function of time.Finally, we identified a clear need for further research to investigate the possible drivers of the temporal trends we observe in our results.
Author Contributions: S.P. and J.M. designed the research and revised this manuscript.S.A. analyzed the data and wrote the manuscript.S.P., T.W. and J.M. provided comments and suggestions for the manuscript.

Figure 1 .
Figure 1.Overview map of the Arabian Peninsula with annotation of the major deserts, cities and vegetated areas with a topographic map (grey scaled) as background.

Figure 2 .
Figure 2. The extent of the 17 dust storms over the Arabian Peninsula as indicated by light yellow-colored polygons with a topographic map (grey scaled) as background.

Figure 2 .
Figure 2. The extent of the 17 dust storms over the Arabian Peninsula as indicated by light yellow-colored polygons with a topographic map (grey scaled) as background.
situ dust recorded, -No in-situ dust records, √ The MODIS dust detection index indicates the same outcome as the in-situ observations, X The dust detection index indicates a different outcome to the in-situ observations, + Test not suitable.

Figure 3 .
Figure 3.An RGB true colour image of MODIS observation indicating the extent of the dust storm over the Arabian Peninsula on 2 February 2012.

Figure 3 .Figure 4 .
Figure 3.An RGB true colour image of MODIS observation indicating the extent of the dust storm over the Arabian Peninsula on 2 February 2012.

Figure 4 .
Figure 4. Validation results within subpanels (a) the detected dust-storm area irrespective of land cover types, and in subpanels (b) the detected dust-storm area using different threshold values for each land cover type (the dashed polygon indicates the extent of the dust-storm as detected on the RGB true color image).

Figure 5 .
Figure 5. Dust indices generalized for all land-cover types and land cover-specific analysis (urban, barren and vegetation).The scatter plots showing for all indices the temporal trends fitted through the mean index value per event as a result of a linear regression, with associated model error bounds (at 95% confidence interval) considering all pixels.The error bars represent the standard deviation for each individual observation.The orange dot represents the dust data and the blue non-dust.The units for all indices are dimensionless, except for BTD where this is in Kelvin.

Figure 5 .
Figure 5. Dust indices generalized for all land-cover types and land cover-specific analysis (urban, barren and vegetation).The scatter plots showing for all indices the temporal trends fitted through the mean index value per event as a result of a linear regression, with associated model error bounds (at 95% confidence interval) considering all pixels.The error bars represent the standard deviation for each individual observation.The orange dot represents the dust data and the blue non-dust.The units for all indices are dimensionless, except for BTD where this is in Kelvin.

Table 1 .
The Validation of Index Performance over Different Land Cover (Dust Event: 2 February 2012).

Table 1 .
The Validation of Index Performance over Different Land Cover (Dust Event: 2 February 2012).

Table 2 .
Wilcoxon Signed Ranks Test p Values for Each Index Considering All Pixels and Land Cover.

Table 3 .
Slope Coefficient Values of the Temporal Linear Regression Analysis with Associated Uncertainty Value (and Level of Confidence).

Table 3 .
Slope Coefficient Values of the Temporal Linear Regression Analysis with Associated Uncertainty Value (and Level of Confidence).