An Improved Aerosol Optical Depth Map Based on Machine-Learning and MODIS Data : Development and Application in South America

In zones where aerosol properties have been poorly characterized, satellite-based (MODIS) and ground-based (AERONET) aerosol optical depth (AOD) values typically differ. In this work, we use machine-learning based methods (artificial neural networks and support vector machines) to obtain corrected AOD values taken from MODIS in regions that are positioned far from AERONET stations. The method has been validated using several approaches. The area suitable for improvement covers 62% of the South American continent, and the degree of improvement compared to MODIS values, expressed in terms of the fraction of data within the MODIS error, was found to be 38% and 86% for the Terra and Aqua satellites, respectively. The results show absolute monthly average differences between the MODIS and the proposed method of up to ± 0.6 AOD units. The MODIS AOD distribution for the analyzed period shows a mode of –0.04, while that for the method presented here is 0.08.


INTRODUCTION
Tropospheric aerosols play a key role in the Earth radiation budget, directly through the scattering and absorption of solar radiation and indirectly by modifying cloud microphysical and radiative properties.When located in the lower troposphere, aerosols cause poor air quality, reduce visibility, and pose public health hazards.
As aerosol lifetimes last on the order of days, aerosol distributions are inhomogeneous in time and space with much higher concentrations found near the emission source.Aerosols also vary in size by several orders of magnitude, and their properties vary as they age in the atmosphere.It is thus very difficult to characterize aerosols on a global scale.
The total light extinction effect caused by aerosols over a vertical column in the atmosphere is known as the aerosol optical depth (AOD).Considerable efforts have been made to observe AODs from space and from ground-based instruments due to the effects of aerosols on global climate change and environmental research.This parameter is required in several atmospheric applications and models (e.g., the monitoring of sources and sinks of aerosols, volcanic eruptions, biomass-burning, and radiative transfer model calculations).
Numerous studies have shown that AOD is proportional to PM 2.5 and it can be used as a surrogate dataset to address spatial and temporal gaps of ground-based monitoring networks (e.g., Weber et al., 2010).Therefore, due to its global coverage and the wide range of instruments available, satellite remote sensing may be used to fill a number of these gaps.A number of efforts have been made to develop satellite-derived indicators for several issues such as air pollution and biomass-burning (i.e., de Sherbinin et al., 2014;Zeeshan and Kim Oanh, 2014).A series of satellite observational schemes was applied at the end of last century to obtain the global distribution of aerosol optical characteristics and seasonal and annual variations in the direct and indirect radiative forcing of aerosols (Kaufman et al., 1997).Among these, the Moderate Resolution Imaging Spectrometer (MODIS) is an important sensor installed onto the Earth Observing System satellites Terra and Aqua, with 36 channels ranging from visible to infrared bands.The scanning information provided by each MODIS instrument covers the entire earth once each day, offering abundant data that could be used to retrieve products related to land, clouds, aerosols, water vapor, ozone, sea color, phytoplankton, biogeochemistry features, etc.
Satellite products are limited in that they are based on important assumptions about aerosols and surface properties (Levy et al., 2007).Recently, improved methods for the treatment of surface reflectance have been proposed (He et al., 2012), or techniques to adjust AOD based on terrain and canopy (Liu et al., 2001).On the other hand, surface measurements do not present such constraints and therefore provide highly accurate aerosol optical properties.This is why aerosol optical properties retrieved from satellites are usually validated against surface measurements (e.g., Estellés et al., 2012).
On a global scale, several studies have compared MODISbased AOD measurements with ground-based measurements retrieved through the Aerosol Robotic Network (AERONET) (Holben et al., 1998) for different areas and have generally revealed strong correlations for several regions (e.g., Remer et al., 2008;More et al., 2013).However, studies on correlations between satellite and ground-based measurements for South American stations have been limited and have focused mostly on the biomass burning season in the Amazonas region.It has been noted that MODIS products present some limitations when covering South America and other regions where aerosols have been poorly characterized.Hoelzemann et al. (2009) presented a multiyear comparison between AERONET and MODIS observations for South America.Key under and over estimations or biases were observed for some South America stations when performing a one-to-one correlation (MODIS-AERONET).The relevance of properly-calibrated satellite observations in South America is crucial, as only nine AERONET stations exists in the continent (at the time that we started this study) with at least three years of high-quality data.
Determinations of bias in AOD are rather problematic due to unknown sources of errors, dataset nonlinearity and the presence of non-normal distributions.Therefore, commonly used linear methods are not always adequate for estimating aerosol characteristics.A powerful approach for managing non-ideal cases is the non-parametric approach, i.e., machine learning techniques.Machine learning is a subfield of artificial intelligence that involves developing algorithms that can empirically learn from the behaviors or properties of a given dataset.When employing these methods, a training dataset is separated into "inputs" and "outputs" and an algorithm is designed to find a connection between them.Outputs are the variables to be predicted and inputs are the variables that will feed the algorithm (Vapnik, 1995).
In this regard, we have recently presented a method for correcting MODIS AOD retrievals at 550 nm (Lanzaco et al., 2016).The method is based on machine learning techniques (Artificial Neural Networks, and Support Vector Machines) and significantly improves AOD values retrieved from satellites, obtaining stronger correlations for AOD from the AERONET than from the MODIS.A similar methodology was used by Palancar et al. (2016) to obtain an AOD value at 340 nm from the MODIS and to in turn calculate the aerosol radiative forcing for the UV-B region.The obtained values were found to be in good agreement with results obtained using AERONET AOD as input.
In the present work, we used MODIS collection 6 to develop an empirical methodology that allows us to define areas surrounding stations where the above listed method for correcting MODIS AOD at 550 nm remains valid.Thus, the final goal here was to use models trained at one station to obtain better AOD values for distant regions.

AERONET
AERONET (Holben et al., 1998) is a global ground-based network of sun-sky photometers established and maintained by NASA and PHOTONS.It covers several hundreds of locations around the world, and it provides a long-term public database of aerosol optical, microphysical and radiative properties.Due to the high quality of the measurements, they are used to validate any other AOD measurement method.
The photometers make direct sun irradiance measurements at eight spectral bands over a sequence of three measurements taken 30 seconds apart to create a triplet observation for each wavelength.AOD is calculated from the spectral extinction of direct beam radiation for each wavelength based on the Beer-Bouguer Law after removing the contributions of other gaseous absorbing pollutants.The method reported uncertainty is ± 0.01 (Duvobnik et al., 2002).
In this work, AERONET Level 2 data were used to validate the models results.AERONET does not provide an AOD measurement at 550 nm, and so it was obtained by interpolation using a log-log quadratic regression (Eck et al., 1999).

MODIS
The MODerate-resolution Imaging Spectroradiometers (MODIS) started measuring aboard the Terra and Aqua satellites in 2000 and 2002, respectively.The satellites complete sun-synchronous orbits around the Earth in such a way that the instruments can generate at least one data point per day for the Earth's surface.
The instruments measure over a broad wavelength range from 0.4 µm to 14.4 µm.Due to low aerosol absorption levels in the mid-IR regions, the surface reflectance can be determined, thus allowing for AOD retrievals.Two different algorithms can process the signal and compute the AOD depending on the surface where a measurement is made: dark target (Kaufman et al., 1997) retrieve the AOD over dark, vegetated areas and deep blue (Hsu et al., 2004) retrieve the AOD over bright areas.In this work, we present only results for dark target algorithm because it is known (Sayer et al., 2014) that deep blue algorithm provides around 20% less matchups than dark target in Amazonas region.Also, the previous work by Lanzaco et al. (2016) showed a better correlation of dark target than deep blue against AERONET for all the sites included in this study.
MODIS Collection 6 Level 2 data were obtained from the "Level 1 and Atmosphere Archive and Distribution System (LAADS)" website (https://earthdata.nasa.gov/about/daacs/daac-laads).We used Optical Depth Land and Ocean variables, which present AOD values at 550 nm based on the best quality data (QA Confidence Flag = 3) measured over a 10 km × 10 km pixel area (at nadir) and processed via the dark target algorithm.All available data up to the year 2014 were used.

METHODOLOGY
In describing the methodology applied, we first describe the method used to improve the MODIS/AERONET AOD correlation.We used the method employed for stations other than the one used for its development, and improvements from the model were measured.Once the improved stations were identified, a common characteristic among them was identified to define the limits of regions where the models can be applied, thus generating an improved result.A flow diagram of the procedure employed is presented in Fig. 1.

Local AOD Correction
To obtain corrected AOD values (AOD C ) based on MODIS AOD observations (AOD M ), we developed a method based on machine-learning algorithms.The method was trained to obtain AOD C values that are as similar as possible to ground truth data, i.e., AOD values provided by AERONET (AOD A ). Machine-learning algorithms are designed to learn from the behaviors or properties of a given dataset.They are trained providing truthful outputs (in this case, AOD A ) so that algorithms can correct themselves until they achieve a minimum error function, thus predicting the desired answer.To determine the potential of the trained algorithm, a portion of the dataset was separated beforehand and it was never used during the training process.Then, this unseen portion of the dataset was used to test the previously trained algorithm.If a method correctly predicted an unseen dataset, it was deemed capable of generalizing.We used Artificial Neural Networks (ANN) (Basheer and Hajmeer, 2000) Matlab Neural Network Toolbox version 2011Rb implementation, and Support Vector Machines (SVM) (Vapnik, 1995;Melgani and Bruzzone, 2004), LIBSVM implementation (Chang and Lin, 2011), as machine learning algorithms to train the models; in both cases the cost function to be minimized is the MSE function.
To train the algorithms, AOD M and AOD A were spatiotemporally collocated following the methodology presented by Petrenko et al. (2012).In addition to AOD M , meteorological variables recorded at the airport closest to the AERONET station were also included as input as well as a variable related to the day of the year.Meteorological variables evaluated for each station included temperature, relative humidity, wind speed and direction.All the possible inputs combinations were used for training both ANN and SVM algorithms and the set of variables that better reproduced AOD A were selected for that station.
This process was carried out separately at nine South American AERONET stations.These stations are located in Brazil (Alta_Floresta, Campo_Grande_SONDA, CUIABA-MIRANDA, Rio_Branco and Sao_Paulo), Argentina (Cordoba-CETT and CEILAP-BA in Buenos Aires), Chile (Arica) and Bolivia (SANTA_CRUZ).For each site, the best resulting algorithm trained with satellite data, station meteorology and ground AOD values will be herein after referred to as model.Thus, there is one model for Alta_ Floresta/Terra data, one model for Alta_Floresta/Aqua data, etc.
Further details on the methodological results, training methodology, criteria used to select the best machine learning algorithm and combination of inputs for each station are presented in Lanzaco et al. (2016).As a summary on the performance of the models, it can be concluded that based on all of the South American AERONET stations where models were trained, the percentage of AOD values within MODIS error (defined as ± (0.05 + 0.15 × AOD A ) by Levy et al. (2013)) increased from 57% for AOD M to 91% for AOD C .The R 2 of the linear regression against AOD A also showed an improvement for all of the stations.
However, the ultimate purpose of these models was not only to obtain AOD C values at the AERONET stations where they were trained, but to also extend their range of applicability to distant regions.Therefore, an AOD C value for other locations can be calculated using a previously trained model but by replacing the input variable AOD M measured at the training station with the AOD M measured at the pixel corresponding to the new location.Thus, to obtain an AOD C value for a given pixel using an already trained model, the following variables are required:  The averaged AOD M for every valid measurement inside a circle with a diameter of 55 km centered on the desired pixel, an approximation which is commonly used to spatio-temporally collocate MODIS and AERONET measurements (Petrenko et al., 2012).Here, it is very important to note that this value must not exceed the highest or lowest AOD values used to train the model, as this would imply an extrapolation from the range where the algorithm was trained. Temporally collocated meteorological data from the station used to train the model and the variable related to the day of the year.Meteorological data at the remote location (named as remote meteorology) were not used as model inputs due to the uncertainties and low resolution of the global meteorological databases in South America.
In addition, we tried to keep the method as simple as possible.Further discussion about this decision is presented in the Sensitivity analysis section.

Improvement Evaluation
To assess the extent of the applicability area of each model trained at one AERONET station, the models were separately used to calculate AOD C values for all of the other stations in South America.Both AOD C and the original AOD M were correlated against AOD A values over the entire measurement period.This was done to compare the correlations and to then determine whether the model improved or not satellite values at that station.AOD M values recorded from Terra and Aqua were treated separately.
The three parameters used to evaluate the correlations and to then determine whether the AOD C was better than the AOD M were: the slope of the linear regression against AOD A , the R 2 of the linear regression against AOD A and the percentage of data within the MODIS error.
These parameters were computed for the whole dataset at every station and also separately for a subset formed only by cases where AOD M values were equal to 0.3 or larger.This different subset was used to test higher AOD values, as its examination is important from an air quality point of view.However, as they usually represent only a small fraction of the entire dataset, information on higher AOD values can be missed when an entire dataset is analyzed, e.g., when analyzing the percentage of data within MODIS error.
In summary, the results derived from each model were evaluated over each station by considering the three parameters of the two datasets (six parameters in total).The output AOD was considered to be improved (i.e., AOD C is closer to AOD A than AOD M ) when at least four out of the six parameters were improved and when none of the remaining parameters presented a decrease of more than 10%.When any of the parameters decreased by more than 10%, the model's application to the given station implied a worsening of raw satellite data.If none of the mentioned conditions were met, the model's application had no significant effect (neither worsening nor improving AOD M ), and such cases were defined as neutral.

Applicability Area Construction
Once all of the stations with results improved by the models were identified, we identified characteristics that they had in common based on similarities in statistical AOD values.In doing so, we considered as many improved stations as possible and none of the worsened ones.This was done to extend the model's applicability by determining which MODIS pixels over South America could be improved by each model.
The data analysis was carried out by dividing the South American area into cells of 0.25° × 0.25° and by then computing several statistical parameters for all of the AOD M pixels inside each cell.We computed mean, standard deviation, 25 th percentile, median, 75 th percentile, maximum and minimum values using all data for the entire measurement period of each MODIS instrument.The statistical parameters were compared with those for the cell where the model was trained by computing the difference between them.Such comparisons were carried out to determine the common statistical requirements that describe cases with improved results while excluding those presenting worsened results.These requirements are the conditions to be met by a MODIS pixel to be suitable for correction by a model.Further information on the identified results is presented in section 4.2.

RESULTS AND DISCUSSION
To summarize the methodology, constructing the applicability area of each model first involved validating the results against the ground truth data, i.e., AERONET, to determine which sites could be improved (valid) and which could not (invalid).Once valid and invalid stations were identified for each model, a common statistical feature among the valid stations was identified to differentiate them from the invalid stations.Finally, the MODIS pixels sharing these validated characteristics were deemed suitable for improvement with a given model.In this way, applicability areas were defined for each model.
All AOD M values falling within the applicability area of a model were corrected to obtain an AOD C value, and these results were analyzed to assess their validity when there were overlapping applicability areas.Then, an AERONET station not used in the study was used to compare the methods.Finally, differences with the original AOD M values were also calculated.The results obtained from this process are presented and discussed in the following sections.

Performance of the Models for Other Stations
From the nine models developed for South America and presented in Lanzaco et al. (2016), six of them were chosen to build their applicability area.The chosen models were those trained at Alta_Floresta, Campo_Grande_SONDA, CUIABA-MIRANDA, Cordoba-CETT, SANTA_CRUZ, and Rio_Branco while Arica, CEILAP-BA (Buenos Aires) and Sao_Paulo were set aside.Arica was not used in this study because of its known local pollution features (Mélin et al., 2010), which make it ineligible to define an area with similar aerosol characteristics or for comparisons with other AERONET stations in South America.Although the CEILAP-BA model results show a considerable improvement over AOD M values, the model was not applied to other locations due to the particular AOD M characteristics of the site, where a strong overestimation of AOD values and a poor correlation with AOD A was observed (Lanzaco et al., 2016).For all of the other stations presented here, underestimations of most data were often found, and CEILAP-BA models would exaggerate such underestimations even more because they were trained to correct overestimations.In Sao_Paulo, the improvement accomplished by the models was not significant enough to justify the construction of their applicability area.
The six chosen models were tested for validation at all of the other available AERONET stations with the exception of Arica but including CEILAP-BA and Sao_Paulo.Both the Buenos Aires and Sao Paulo stations show local and regional aerosol effects, thus making them useful for validating the models.Therefore, as each model was not tested in the same place where it was trained, seven validation stations were used for each model.
As stated in section Improvement evaluation, the correction was evaluated based on AOD M and AOD C correlations against AOD A values for each station.For the 42 cases (6 models, with each applied to the other 7 stations), the result worsened for 23, 10 presented no significant changes, and for 9 the AOD was improved for the Terra satellite.When Aqua data were used, these values were recorded as 21, 5, and 16, respectively.More detailed results are presented in Fig. 2.
The Cordoba-CETT/Aqua model results are presented in detail in Fig. 3 as an example of how decisions on the performance of the model in each site were made.AOD values at Alta_Floresta, CUIABA-MIRANDA, Campo_ Grande_SONDA, and Sao_Paulo were considered after correcting the model as better than the unmodified AOD M values, as at least four parameters were improved and none of the remaining parameters worsened by more than 10%.The Rio_Branco results did not worsen after the correction was applied (at a threshold of lower than 10%), but they did not present significant enhancements either (only three out of six parameters were improved), so it was concluded that the model's application did not significantly affect the results.Finally, even when high CEILAP-BA AOD values were clearly improved when considering the R 2 and the slope of the linear regression, the percentage of data within MODIS error decreased considerably when analyzing the entire dataset.Hence, AOD M values were labeled as worsened by the model.This result was expected, as the Cordoba-CETT models were trained to correct underestimations and they are not capable of coping with overestimations found in CEILAP-BA.
The results derived from applying the Cordoba-CETT models in the Amazonian stations were surprising, as it was not expected that a model trained in a completely different environment could improve the results of that region.Moreover, this result was unexpected due to the large distance that separates the station were the model was trained and the Amazonian stations; the distance between Cordoba-CETT and Campo_Grande_SONDA (the closest Amazonian station) is approximately 1600 km.Even when the Amazonian stations present one of the highest R 2 when correlated with AERONET (Hyer et al., 2011), these values result from compensation between the underestimation of lower AOD values and the overestimation of high values (Lanzaco et al., 2016).As is shown in Table 1, the Cordoba-CETT/Aqua model was trained with a maximum AOD M value of 0.737 as input, and it is only applicable within that range.Thus, when narrowing the AOD range for Amazonian stations to fall within the valid range for this model, a lack of compensation by high AOD values generated R 2 values lower than 0.6 in the correlation with AOD A , and these were the R 2 values that were improved.Both Aqua and Terra models developed at Cordoba-CETT were trained to correct underestimations of low AOD values and they can thus correct this deviation in any site that presents it.
The performance of the other models (other than Cordoba-CETT) was mostly correlated with the distances to the stations that they were able to improve.This was only untrue for the improvement found in CEILAP-BA when using the Rio_Branco model, the only model to  improve AOD in this city.In Rio_Branco, we found one of the highest overestimations among all of the analyzed stations and it is because of this condition that the models trained here can correct systematic overestimations found in CEILAP-BA in spite of the great distance that separates them (roughly 2,900 km).This was also the case for the Cordoba-CETT models that improved the Amazonian stations.

Applicability Areas Construction and Validation
At first, it was thought that the geographical area for which the models could improve AOD M would be coincident with the anisotropic representativeness of AERONET stations defined by Hoelzemann et al. (2009) during the biomassburning season.These authors defined representative areas for each AERONET station in South America according to the R 2 of the linear correlation between AOD A and a MODIS pixel located in a given position.When the R 2 was higher than 0.5, they concluded that the AERONET station could represent the AOD in the location of that MODIS pixel.Li et al. (2014) found similar results for the Alta_Floresta station using the correlation between anomalies in the time series of AOD M and AOD A measurements.
However, there was no coincidence between the representative areas described in the above listed works and the results presented here.The stations that could be improved by a model trained at another station were occasionally located outside of their representative areas, and stations within these representative areas were sometimes not improved.A high R 2 value in the linear correlation between the AOD A (where the model was trained) and an AOD M pixel (where the model was used) was not coincident with the results presented in section Performance of the models for other stations, and thus this approach was discarded.Thus, a different parameter for defining which pixels are suitable for model improvement was investigated.
According to the results presented in section Performance of the models for other stations, the capacity for a model to improve AOD at another station seemed to be correlated with some AOD M statistical parameters.Therefore, several statistical variables for AOD M data were compared between each cell for South America and the cell where the model was trained (see section Improvement evaluation).Each case was classified as better, similar or worse according to the results of the previous section (Fig. 2), observing that in 9|16 (Terra|Aqua) cases the results improved, in 22|21 worsened, and in 11|5 they were neutral.Then, a correlation between this classification and the computed statistical parameters was identified.
We found that all of the stations worsened by the models presented differences in the 75 th percentile of more than 0.1 and differences in the standard deviation of more than 0.075.Thus, if criteria for the improvement of a cell are defined as a difference in the 75 th percentile of less than 0.1 and the difference in the standard deviation of less than 0.075, 10 out of 16 cases with improved results meet them and results were not worsened for any of the 21 cases.Some improved stations present greater differences than the above mentioned values, but as a priority, when defining the criteria, they were required to exclude all of the worsened stations to the detriment of losing areas that are suitable for improvement.All of these tasks were completed using models, data, and classification results for Aqua data.When the same criteria were applied to the Terra results, we found that 4 of the 9 improved stations and 0 of the 22 worsened stations met the criteria, verifying that the criteria are also valid for Terra data.Once the criteria were defined, the standard deviation and 75 th percentile measures were computed for all of the available data in each cell for South America to classify them as suitable or not for improvement through a given model.In this way, applicability areas were constructed for each model by determining which cells met the aforementioned criteria.
The applicability areas are presented in Fig. 4 for the models trained with data from the Terra satellite, and Fig. 5 presents applicability areas trained with data from the Aqua satellite.The same information presented in Fig. 2 is also plotted to compare the cross validation data with the areas covered by each model.It must be noted that none of these areas include worsened stations, even though the plot resolution suggests the opposite in some cases.
We did not find significant differences between areas obtained from the Terra and Aqua data, indicating that the aerosol loading does not change significantly between the morning and afternoon, at least over a large scale.Some regions of the applicability areas change depending on the satellite considered, but their relative sizes are small.In almost all cases, the anisotropy observed is very noticeable.The sizes of the applicability areas are presented in Table 1 in addition to the maximum AOD M value that can be used to determine AOD C .When all of the areas are overlaid, the resulting area reaches 12/14 millions of km 2 (Aqua/Terra), meaning that 62% of the South American surface is suitable for improvement using only six models.
To further validate the methodology, MODIS pixels presenting several models that can predict AOD C values were analyzed (Fig. 6(a)).Overlapping applicability areas allow one to test model accuracy levels by determining whether all of the independent models provide the same AOD C value.For each of the 117,157 cases presenting temporal and spatial overlap, the average and standard deviation was calculated for the AOD C obtained from the different models at a given latitude, longitude, and time.Fig. 6(b) shows a histogram of the standard deviation values, and in 89.3% of the cases, the standard deviation of the AOD C values was less than 0.05 (i.e., the MODIS algorithm accuracy level).This result confirms the robustness of the model.

Sensitivity Analysis
To assess the sensitivity of the method, three different meteorological datasets (other than the meteorological values recorded where the model has been trained, or model meteorology), were used as inputs, in addition to AOD M .The first dataset corresponds to values obtained from an airport at the same place of AOD M pixel (remote meteorology); the second corresponds to a grid point from the Global Data Assimilation System dataset (GDAS) (Kanamitsu, 1989)  with a 1° × 1° grid, and the third dataset was built with randomly generated numbers replacing the meteorological values.All three sets have been used to obtain AOD C , following the same steps previously described, to assess the closeness to AOD A as stated in section Improvement evaluation by means of statistical parameters (data within error, R 2 , and slope).The difference between AOD C obtained with meteorology from the three datasets and AOD C obtained with model meteorology was evaluated and presented as a histogram in Fig. 7.The histogram includes the 17 cases where a station was located inside the applicability area from other site.Eq. ( 1) is presented as an example of the criteria of evaluation, where RM and MM stand for remote meteorology and model meteorology, respectively.
Regarding to Fig. 7(a), the difference in all the statistical parameters using model meteorology and remote meteorology was less than ± 10% for all the stations.Therefore, and considering that the use of remote meteorology is not always an option due to the lack of meteorological data in remote places, its use has been discarded.
The use of a global meteorological database like GDAS as input is an option that can overcome the lack of reliable meteorological data in remote sites.However, in this case the AOD C is worse (for up to 40%) than the AOD C generated with model meteorology (Fig. 7(b)), therefore its use was not considered as a suitable option.Finally, when using random values as meteorological data for the input (Fig. 7(c)) the improvement was reduced up to 25% for all the cases considered, pointing out that although the importance of the meteorology in the model is not large, it is still relevant.
In conclusion, despite the similar results obtained by using remote meteorology at some sites, the meteorology from the site where the model was trained was chosen as the best set of inputs to the models.

Applicability Areas Analysis
Overall, when considering all of the cases where a station was located within other applicability station areas (thus, AOD C and AOD M can be compared to AOD A ), we found that the fraction of data within error values increased by 38%from 0.59 to 0.81 (4,475 points) for the Terra satellite.For the Aqua satellite, an 86% increase from 0.43 to 0.80 (4,525 points) was found.
The areas constructed for the six models can be roughly grouped in pairs according to their shape and/or location (Figs. 4 and 5).These groups are: Alta_Floresta/Rio_Branco, CUIABA-MIRANDA/SANTA_CRUZ and Campo_Grande_ SONDA/Cordoba-CETT.For all three pairs, both stations are located within each applicability area.It should also be noted that the three different group areas almost never overlap.
An analysis of the shape of these groups suggests a correlation with the Enhanced Vegetation Index (EVI).EVI is a MODIS product that was developed to obtain a vegetation signal with improved sensitivity in high biomass regions through a de-coupling of the canopy background signal and a reduction in atmospheric influence (Huete et al., 2002).To illustrate this correlation, Fig. 8 shows the EVI value for South America and for the applicability area of each group for the average of all MODIS/Aqua measurements (2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014).The EVI is quite homogeneous within each area, though differentiating pixels within and outside of the applicability area (as defined by the AOD M statistical parameters) based strictly on EVI values was not straightforward.An analysis of the areas and associated EVI values is given in the following paragraphs (Table 2 presents a description of EVI statistics for all of the stations).
Alta_Floresta and Rio_Branco, two stations located close to the Amazonas rainforest, showed very similar applicability areas, covering a large region in the southern area of Amazonas (0.978 and 1.24 million of km 2 , respectively) in Brazil and Bolivia.The median EVI for the pixels included within the area are 0.767 and 0.765, respectively.These nearly identical values and their similarities within 25 th and 75 th percentiles suggest a strong correlation between applicability areas and EVI.The region can be characterized as a transitional zone between the Amazonas rainforest and the savanna area found in the Cerrado region of Brazil.
CUIABA-MIRANDA and SANTA_CRUZ shows an applicability area surrounding the Alta_Floresta/Rio_Branco area to the south, including Paraguay, central Brazil, parts of Bolivia and Peru, and areas of Ecuador and Colombia in the case of the SANTA_CRUZ area.The median EVI for the area is 0.657 and 0.666 for both station areas, respectively.Vegetation in the regions can be mainly characterized as savanna.
Finally, the Campo_Grande_SONDA area is fully included within the Cordoba-CETT area (roughly three times larger) that surrounds the CUIABA-MIRANDA/SANTA_CRUZ areas to the south, east and west.Cordoba-CETT covers the largest total area of all six stations (11.1 million km 2 ).The large geographical extent of the Cordoba-CETT model obeys to the same reason that makes the model suitable for improving the Amazonian stations as described in section 4.1.The model covers the western region of Brazil, Uruguay, and much of Argentina while excluding arid areas of Patagonia and the Andes mountain range.Areas of Brazil of up to 4,000 km are also included within this applicability area.

Applicability Areas Use Petrolina Station
As an example of the method's application to a different station, we analyzed the Petrolina_SONDA AERONET station (9.38°S, 40.50°W, 370.0 m.a.s.l.), which at the start of this study had less than 1,000 days of Level 2 AERONET data and which was thus not included in the analysis.More days of measurements were included later, thus rendering it a reliable station that has not been previously used for model training.As the station is located within the applicability area of the Cordoba-CETT station, we compared the AOD M and AOD C (obtained from the Cordoba-CETT model) against AOD A for 298 measurements taken from the Aqua satellite and for 258 measurements taken from the Terra satellite.Fig. 9 shows a box plot of all of the cases, from which it can be observed that the average values for all of the points recorded at the Petrolina_SONDA station, measured via AERONET, were (0.091 ± 0.060) and (0.078 ± 0.057) when the Terra or Aqua satellites passed over the station, respectively.The MODIS reported average AOD values of (0.046 ± 0.087) and (0.020 ± 0.069), whereas the corrected method generated average AOD values of (0.080 ± 0.063) and (0.072 ± 0.051).Clearly, the model generated more reliable results than the MODIS instrument not only in regards to the average values but also in terms of the representation of cases with high AOD values as denoted by the 75 th percentile plot.

Changes in the South America AOD Map
As an example of the real change in the AOD distribution over South America after applying the method, four plots are presented in Fig. 10: one showing the average Aqua AOD M values for September of 2007 (Fig. 10(a)) and other showing the difference between average AOD C and AOD M values for the applicability areas (Fig. 10(b)) for the same dates.If the AOD C of a pixel can be estimated using two or more models, then their average value is used.As expected, the Amazonas region shows the largest AOD values in South America with values of up to 2 due to the biomassburning processes.In Fig. 10   In order to quantify the differences between AOD M and those obtained from the corrected method, in terms of AOD distributions in South America, Fig. 11 shows a histogram of every available AOD value obtained from both methods for pixels within all applicability areas for all of the available dates (9.3 × 10 7 points).AOD C presents a mode of 0.08, whereas AOD M presents a mode of -0.04, a consequence of the unrealistic negative AOD values measured by MODIS, which are permitted in order to avoid the appearance of an artificial bias in long term statistics.After excluding negative values, the AOD M mode is 0.02, which is still lower than that of AOD C .These distributions reinforce the fact that the MODIS routinely underestimates lower AOD values and that a method for obtaining values closer to AERONET values is very useful.

SUMMARY AND CONCLUSIONS
Despite the great potential of satellite data products, MODIS and AERONET values often differ, and especially for regions where aerosols are poorly described.This difference has been partially corrected through the latest MODIS collection 6, but it is still noticeable.
In this work, we presented a means of correcting AOD values obtained from MODIS instruments.The method is based on machine learning algorithms, which are based on AOD data from the MODIS together with meteorological data for a given site.The model was applied and validated for other AERONET stations to define an applicability area, where the model is still suitable.The validity of the assumptions was verified using data for surface stations not included in the model's development and through a statistical analysis of AOD values for areas where two or more models made predictions.The shape of these areas was correlated with the MODIS Enhanced Vegetation Index, suggesting that vegetation types and degrees of coverage have a strong effect on the aerosol loading in South America.
The proposed methodology can be applied to other regions and may prove useful when applied to zones where aerosols are poorly characterized (e.g., Africa, Asia and Australia according to Shi et al., 2011).Its computational requirements are minimal, and data required are AERONET and MODIS AOD values and simple meteorological station data.Criteria used for determining applicability areas have been proven to be independent of sites and satellites.
When all applicability areas overlap, MODIS measurements for 62% of the South American area can be improved using this method.When an AERONET measurement was available for comparison, an improvement of 38% for Terra satellite values and of 86% for Aqua satellite values was obtained, measured as the fraction of data falling within the MODIS error.The difference between the AOD corrected using the method and MODIS is very noticeable according to the histogram of AOD values, which routinely shows monthly differences of ± 0.1 and of up to ± 0.6 for some places.
The differences found are of great relevance, as many climatic and radiative transfer models depend heavily on AOD values to obtain reliable results.They thus require consistent data to propose scenarios for the climate change framework, and even small variations in the AOD can lead to changes in the calculated variables.As an example, computed changes in AOD values after the publication of MODIS Collection 5 showed weakened radiative forcing compared to Collection 4 of -0.8 to -0.65 W m -2 (Myhre, 2009).Weihs et al. (2008) also reported that an increase in AOD at 368 nm of 0.1 produces a decrease in the UV index of approximately 2.5%, and (Kim et al., 2013) showed that a 1% increase in AOD forces a decrease of (0.29 ± 0.06)% in surface erythemal ultraviolet irradiance.
This methodology may be used to make estimations of PM 2.5 and PM 10 based on satellite measurements.By using local measurements of particulate matter and applying a similar methodology to that described in this work, it would be possible to obtain better values of PM 2.5 and PM 10 from satellite data in the surface sites (likely using other variables as well and not only meteorological ones).Then, this methodology could be extended to remote regions.

Fig. 1 .
Fig. 1.Flow diagram of the procedure for obtaining corrected AOD values.

Fig. 2 .
Fig. 2. Assessment of the improvement in the AOD C values compared to AOD M values for a given station (vertical axis) using a model trained at another station (horizontal axis).Solid symbols represent MODIS/Terra models and unfilled symbols represent MODIS/Aqua.

Fig. 3 .
Fig. 3. Evaluation of the Cordoba_CETT/Aqua model for all of the other stations in terms of R 2 and the slope of the linear regression against AODA, and the fraction of data within the MODIS error.

Fig. 4 .
Fig. 4. Applicability areas (in gray) constructed using models trained at different (a-f) AERONET stations (blue circles) based on data from the Terra satellite.Circles mark the station where the model was developed, green triangles mark stations where the model improved the MODIS results, and red triangles mark where the model worsened the MODIS results.Black diamonds mark stations where no significant changes were observed.

Fig. 5 .
Fig. 5. Applicability areas (in gray) constructed using models trained at different (a-f) AERONET stations (blue circles) based on data from the Aqua satellite.Circles mark the station where the model was developed, green triangles mark stations where the model improved the MODIS results, and red triangles mark where the model worsened the MODIS results.Black diamonds mark stations where no significant changes were observed.

Fig. 6 .
Fig. 6. a) Map showing the number of models that can be used in a given cell.Circles represent the AERONET stations where models were trained.White zones represent regions outside of the applicability areas where no model can improve the MODIS AOD values.b) Histogram of the AOD C standard deviations when at least two models can correct AOD values at the same time and location.

Fig. 7 .
Fig. 7. Histogram of the difference between AOD C obtained with meteorology of the remote site (a-b) or random values (c) and AOD C obtained with meteorology in the site where the model was trained.Data within error, R 2 , and slope refer to the comparison of each AOD C against AOD A .The area shaded in grey indicates the cases where AOD C obtained with model meteorology was worse than the one obtained with the meteorology indicated in the figure.
(b), an almost systematic

Fig. 10 .
Fig. 10.a) Average Aqua AOD M values for September of 2007.b) Average AOD C minus average AOD M for the same period.Black areas represent the zones where differences fall within the MODIS algorithm accuracy level (± 0.05).c) and d) are the equivalent plots for July 2007.

Fig. 11 .
Fig. 11.Histogram of the AOD C and AOD M values for the Aqua satellite for all of the dates and locations where an AOD C value was calculated.

Table 1 .
Maximum MODIS AOD values for the models corresponding to each AERONET station, together with the size of the applicability areas for both Terra and Aqua satellites.

Table 2 .
25 th percentile, median and 75 th percentile of the EVI values (2002-2014) located inside the applicability areas of each model.