10-year spatial and temporal trends of PM2.5 concentrations in the southeastern US estimated using high-resolution satellite data

Long-term PM2.5 exposure has been associated with various adverse health outcomes. However, most ground monitors are located in urban areas, leading to a potentially biased representation of true regional PM2.5 levels. To facilitate epidemiological studies, accurate estimates of the spatiotemporally continuous distribution of PM2.5 concentrations are important. Satellite-retrieved aerosol optical depth (AOD) has been increasingly used for PM2.5 concentration estimation due to its comprehensive spatial coverage. Nevertheless, previous studies indicated that an inherent disadvantage of many AOD products is their coarse spatial resolution. For instance, the available spatial resolutions of the Moderate Resolution Imaging Spectroradiometer (MODIS) and the Multiangle Imaging SpectroRadiometer (MISR) AOD products are 10 and 17.6 km, respectively. In this paper, a new AOD product with 1 km spatial resolution retrieved by the multi-angle implementation of atmospheric correction (MAIAC) algorithm based on MODIS measurements was used. A two-stage model was developed to account for both spatial and temporal variability in the PM2.5–AOD relationship by incorporating the MAIAC AOD, meteorological fields, and land use variables as predictors. Our study area is in the southeastern US centered at the Atlanta metro area, and data from 2001 to 2010 were collected from various sources. The model was fitted annually, and we obtained model fitting R2 ranging from 0.71 to 0.85, mean prediction error (MPE) from 1.73 to 2.50 μg m−3, and root mean squared prediction error (RMSPE) from 2.75 to 4.10 μg m−3. In addition, we found cross-validation R2 ranging from 0.62 to 0.78, MPE from 2.00 to 3.01 μgm−3, and RMSPE from 3.12 to 5.00 μgm−3, indicating a good agreement between the estimated and observed values. Spatial trends showed that high PM2.5 levels occurred in urban areas and along major highways, while low concentrations appeared in rural or mountainous areas. Our time-series analysis showed that, for the 10-year study period, the PM2.5 levels in the southeastern US have decreased by ∼20 %. The annual decrease has been relatively steady from 2001 to 2007 and from 2008 to 2010 while a significant drop occurred between 2007 and 2008. An observed increase in PM2.5 levels in year 2005 is attributed to elevated sulfate concentrations in the study area in warm months of 2005.


Introduction
Long-term exposure to PM 2.5 (particle size less than 2.5 μm in the aerodynamic diameter) is associated with various adverse health outcomes including respiratory and cardiovascular diseases (Crouse et al., 2012;Peng et al., 2009). Due to the spatiotemporally continuous nature of the distribution of fine particles, obtaining long-term and spatially resolved distribution of PM 2.5 concentrations is important to reduce exposure misclassification and facilitate epidemiological studies in the region. In addition, time-series analyses of air pollution and human health have become a common study design to compare day-to-day fluctuations of air pollution and corresponding fluctuations in health outcomes (Ito et al., 2007) and require long-term PM 2.5 concentration estimates. Previous research examined temporal trends in PM 2.5 levels. For instance, Weber et al. (2003) investigated the temporal variations in PM 2.5 concentrations at the US Environmental Protection Agency (EPA) Atlanta Supersite Experiment in August, 1999. So et al. (2007 examined long-term variation in PM 2.5 levels during two 12-month periods in Hong Kong. The EPA (2011) evaluated temporal trends of annual and 24 h mean PM 2.5 concentrations at the national level from 2001 to 2010 and reported that annual and 24 h mean PM 2.5 concentrations dropped 24 and 28 %, respectively, during these 10 years.
Accurate depiction of spatial trends of PM 2.5 levels is also important for air pollution health effects research. Stationary ambient monitors leave large areas uncovered, which makes the assessment of PM 2.5 spatial variability difficult. Using measurements from central monitors to estimate population exposure inevitably introduces measurement errors that likely have substantial implications for interpreting epidemiological studies, especially time-series analyses (Zeger et al., 2000). On the other hand, aerosol observations from satellite remote sensing could substantially improve estimates of population exposure to PM 2.5 (van Donkelaar et al., 2010). As a result, satellite-retrieved aerosol optical depth (AOD), which measures light extinction by aerosols in the atmospheric column, has been widely used to predict ground-level PM 2.5 concentrations, given its relatively low cost and large spatiotemporal coverage. A number of AOD products from sensors such as the Moderate Resolution Imaging Spectroradiometer (MODIS), the Multiangle Imaging Spectro-Radiometer (MISR), and the Geostationary Operational Environmental Satellite Aerosol/ Smoke Product (GASP) have been applied to PM 2.5 concentration prediction in previous studies (Liu et al., 2007(Liu et al., , 2009Paciorek et al., 2008;Hu et al., 2013). A limitation of those AOD products is the relatively coarse spatial resolutions. For instance, the spatial resolutions of AOD derived from MODIS and MISR are 10 and 17.6 km, respectively. Although GASP has a spatial resolution of 4 km, the AOD retrievals are less precise than those from the polar-orbiting instruments due to limited information content (one spectral band) and relatively low signal-to-noise ratio of the GOES sensor (Prados et al., 2007). Meanwhile, epidemiological studies typically have access to health data geo-coded to small geographical units (e.g., zip code and census block groups), many of which are substantially smaller than the spatial resolutions of MODIS and MISR. In addition, satellite-estimated PM 2.5 concentrations at coarse resolutions omit detailed spatial variability of PM 2.5 exposure and therefore have limited value in the investigation of spatial trends of PM 2.5 levels at urban scale (Hu et al., 2014). Hence, it is essential to use high-resolution AOD retrievals to generate high spatial resolution PM 2.5 concentration estimates. Recently, a new AOD product retrieved by the multi-angle implementation of atmospheric correction (MAIAC) algorithm based on MODIS measurements has been reported (Lyapustin et al., 2011b). MAIAC AOD has a spatial resolution of 1 km and thus has the ability to estimate PM 2.5 concentrations at that resolution. Moreover, MAIAC AOD has been demonstrated to be correlated with monitored PM 2.5 levels in the New England region (Chudnovsky et al., 2012). Hu et al. (2014) compared the performance of MAIAC with MODIS in PM 2.5 concentration prediction in the southeastern US in a case study and found that MAIAC predictions can reveal many more spatial details than MODIS. In a single 12 × 12 km 2 Community Multiscale Air Quality (CMAQ) grid cell, MODIS can only make one prediction, while MAIAC can make ∼144 predictions. As an example of the benefit gained with increased resolution, MAIAC predictions can distinctly show high concentrations along major highways, while MODIS predictions cannot.
Various statistical methods have been developed to establish the quantitative relationship between PM 2.5 and satellite-derived AOD, including linear regression (Schafer et al., 2008;Wallace et al., 2007;Gupta and Christopher, 2009). However, many of the methods do not consider day-to-day variability in the association between PM 2.5 and AOD. Lee et al. (2011) and Kloog et al. (2011) argued that the PM 2.5 -AOD relationship varies day to day, and this temporal variability needs to be accounted for in order to improve the performance of the AOD-based prediction models. As a result, both studies developed a linear mixed effects (LME) model to incorporate daily calibration of the PM 2.5 -AOD relationship and obtained predictions with high accuracy. To move one step further, Hu et al. (2014) introduced a geographically weighted regression (GWR) model as the second stage to account for possible spatial variability in the PM 2.5 -AOD relationship. This model used the MAIAC AOD as the primary predictor and meteorological fields and land use variables as secondary predictors. Hu et al. (2014) further pointed out that AOD is essential in the two-stage model framework in terms of prediction accuracy. The model can predict PM 2.5 concentrations with high accuracy and thus was adopted in this study.
The objectives of this paper were, first, to estimate spatiotemporally resolved PM 2.5 concentrations in the study domain during the period between 2001 and 2010 in the southeastern US using the two-stage model developed by Hu et al. (2014). Second, maps of annual mean PM 2.5 concentrations as well as the changes between 2001 and 2010 were generated from the daily estimates to visually illustrate the spatial trends of annual PM 2.5 levels between 2001 and 2010. Third, time-series analyses were conducted for the study domain and the Atlanta metro area specifically using the seasonal and annual mean PM 2.5 estimates to examine the 10-year temporal trends of PM 2.5 levels, and the underlying causes were discussed.

Study area
The study area is approximately 600 × 600 km 2 in the southeastern US, covering most of Georgia, Alabama, and Tennessee, and parts of North and South Carolina (Fig. 1). The domain includes several large urban centers, numerous medium-to-small cities, as well as suburban and rural areas.

PM 2.5 measurements
The 24 h average PM2.5 concentrations from 2001 to 2010 collected from the US EPA federal reference monitors (FRMs) were downloaded from the EPA's Air Quality System Technology Transfer Network (http://www.epa.gov/ttn/airs/airsaqs/). PM2.5 concentrations less than 2μgm −3 (∼ 0.2-3 % of total data records) were discarded as they are below the established limit of detection (EPA, 2008a).

Remote sensing data
MAIAC retrieves aerosol parameters over land at 1 km resolution, which was accomplished by using the time series of MODIS measurements and simultaneous processing of a group of pixels in fixed 25 × 25 km 2 blocks (Lyapustin et al., 2011a(Lyapustin et al., , b, 2012. MAIAC uses a sliding window to collect up to 16 days of MODIS radiance observations over the same area and processes them to obtain surface parameters used for aerosol retrievals. To facilitate the time-series analysis, MODIS data are initially gridded to a 1 km resolution in a selected projection. For this work, we used MODIS level 1B (calibrated and geometrically corrected) data from Collection 6 re-processing, which removed major effects of temporal calibration degradation of Terra and Aqua, a necessary prerequisite for the trend analysis.
Validation based on the Aerosol Robotic Network (AERONET) data showed that MAIAC and the operational Collection 5 MODIS Dark Target AOD have a similar accuracy over dark and vegetated surfaces, but also showed that MAIAC generally improves accuracy over brighter surfaces, including most urban areas (Lyapustin et al., 2011b). MAIAC AOD data from 2001 to 2010 were obtained from the National Aeronautics and Space Administration (NASA) Goddard Space Flight Center. Due to the lack of sufficient data records from AERONET, a comparison between MAIAC AOD and AERONET measurements in our study domain was not possible. Zhang et al. (2012) found that Terra and Aqua may provide a good estimate of the daily average of AOD. Thus, the average of the Aqua and Terra measurements can be used to predict daily PM2.5 concentrations. In this study, Aqua (overpasses at ∼ 1:30 p.m. local time) and Terra (overpasses at ∼ 10:30 a.m. local time) MAIAC AOD values were first combined to improve spatial coverage. In our study domain, the increase in spatial coverage ranged from 30.2 to 72.4 % for Aqua and from 17.2 to 26.3 % for Terra from 2001 to 2010. In a common MAIAC pixel, there might be only one MAIAC product from either Aqua or Terra, or both may be present. In the second case, when we combine them, the averaged value, as pointed out by Lee et al. (2012), is likely to better reflect daily aerosol loading, yet in the first case, AOD as an indicator of PM 2.5 abundance is biased towards the atmospheric condition either in the morning or early afternoon. To estimate the missing AOD value, Lee et al. (2011) defined a simple ratio between averaged Terra and Aqua AOD. In this study, we fitted a linear regression to define the relationship between daily mean Terra-MAIAC and Aqua-MAIAC AOD values. We used this regression to predict the missing AOD value (i.e., predict Terra-MAIAC AOD with the available Aqua-MAIAC AOD, and vice versa), and then averaged the observed and the predicted AOD values together. Finally, we set an upper bound of 2.0 for the combined AOD to reduce potential cloud contamination (∼0.05-0.1 % of total data records were excluded).

Meteorological fields
The meteorological fields provided by the North American Land Data Assimilation System (NLDAS) Phase 2 were downloaded from the NLDAS website (http://ldas.gsfc.nasa.gov/ nldas/). The spatial resolution of NLDAS meteorological data is 1/8th of a degree (∼ 13 km). Another meteorological data set used in this study is the North American Regional Reanalysis (NARR). NARR is a long-term, consistent, high-resolution climate data set for North America (Mesinger et al., 2006), with a spatial resolution of ∼ 32 km. NLDAS provides most of the meteorological fields used in this analysis, including relative humidity, U wind, and V wind, while NARR provides another critical parameter: boundary layer height. To generate daytime meteorological fields corresponding to the MODIS overpass times, 3-hourly NARR measurements and hourly NLDAS measurements from 10 a.m. to 4 p.m. standard local time were averaged.

Land use variables
Elevation data were downloaded from the national elevation data set (NED) (http:// ned.usgs.gov). NED is the seamless elevation data set covering the conterminous United States and is distributed by the US Geological Survey (USGS). The elevation data are downloaded at a spatial resolution of 1 arcsec (∼30m). The road data were obtained from ESRI StreetMap USA (Environmental Systems Research Institute, Inc., Red-land, CA). The road data at level A1 (limited access highway) were extracted. Summed length of road segments was calculated for each 1 × 1 km 2 MAIAC grid cell, and grid cells with no roads were assigned zero. The 2001 and 2006 Landsat-derived land cover maps covering the study area with a spatial resolution of 30 m were downloaded from the National Land Cover Database (NLCD) (http://www.mrlc.gov). Forest cover maps were generated by assigning one to forest pixels and zero to others. Primary PM 2.5 emissions (tons per year) were obtained from the 2002, 2005, and 2008 EPA National Emissions Inventory (NEI) facility emissions reports. Grid cells with multiple emission sources were assigned the summed value, and those with no emissions were assigned zero.

Data integration
All the data were first re-projected to the USA Contiguous Albers Equal Area Conic USGS coordinate system. For model fitting, a 1 × 1 km 2 square buffer was generated for each PM 2.5 monitoring site. Meteorological fields and AOD values were assigned to each PM 2.5 monitoring site using the nearest neighbor approach. Forest cover and elevation were averaged, while road length and point emissions were summed over the 1 × 1 km 2 square buffer by calculating the total length of road segments and total point emissions within the buffer. For PM 2.5 prediction, the same procedure was performed for each 1×1km 2 MAIAC grid cell.

Model structure
We adopted the two-stage spatiotemporal model developed by Hu et al. (2014). For the model to be valid, we assumed that particles within the boundary layer were well mixed, and the vertical distribution of particles above boundary layer was relatively smooth. The first stage is a LME model with day-specific random intercepts and slopes for AOD and meteorological fields to account for the temporally varying relationship between observed PM2.5 and AOD (Eq. 1). The model structure can be expressed as (1) where b i and b i,t (day-specific) are the fixed and random intercept and slopes, respectively.
Fixed intercepts and slopes are the same for all days and generated via conventional linear regression, while random intercepts and slopes vary independently for each individual day and are estimated via likelihood methods from the full set of observations. In this study, we generated fixed slopes for each predictor variable, but random slopes were only generated for AOD and meteorological fields, since they represent time-varying variables. The fixed slopes (e.g., b 1 , b 2 ,…, b 6 ) denote the overall relationship for all days, and the random slopes (e.g., b 1,t b 2,t ) indicate the daily relationship among PM 2.5 , AOD, and meteorological fields. PM 2.5,st is the measured ground level PM 2.5 concentration (μgm −3 ) at site s on day t; AOD st is the MAIAC AOD value (unitless) at site s on day t; Meteorological Fields st is the meteorological parameters at site s on day t and may include Relative Humidity st , Boundary Layer Height st , Wind Speed st , U Wind st , and V Wind st ; Relative Humidity st is the relative humidity (%) at site s on day t; Boundary Layer Heightst is the boundary layer height (m) at site s on day t; Wind Speed st is the 2 m wind speed (m s −1 ) at site s on day t; U Wind st is the east-west component of wind (ms −1 ) at site s on day t; V Wind st is the north-south component of wind (m s −1 ) at site s on day t; Elevation s is elevation values (m) at site s; Major Roads s is road length values (m) at site s; Forest Cover s is forest cover values at site s; Point Emissions s is point emissions (tons per year) at site s; and Ѱ is an unstructured variance-covariance matrix for the random effects. The fixed effects affect the population mean and represent the average effects on PM 2.5 concentration estimates for the entire period, while the random effects contribute to the covariance structure and account for the daily variability in associations between dependent and independent variables. Although the PM 2.5 -AOD relationship might vary by season, our first-stage LME model was able to incorporate daily variability in the relationship by generating day-specific random slopes for AOD and meteorological fields and thus should be able to capture the seasonal variability. In addition, by comparing the performances of models fitted for each season, each year, and all 10 years, we found that the models fitted for each year generally yielded higher prediction accuracy. Hence, in this study, we fitted the model for each year individually allowing the predictors used in each model to vary from year to year. Each final annual model was selected to achieve the highest prediction accuracy, and only statistically significant variables were retained. The detailed model structures can be found in the Supplement.
The second stage is a geographically weighted regression (GWR) model that can generate a continuous surface of estimates for each parameter at each location instead of a universal value for all observations. We fitted a monthly GWR model to calibrate the spatial variability within the PM 2.5 -AOD relationship, and the model can be expressed as (2) where PM 2.5 _resi st denotes the residuals from the stage one model at site s in month t, AOD st is the MAIAC AOD value (unitless) at site s in month t, and β 0,s and β 1,s denote the location-specific intercept and slope, respectively.
To assess the goodness of fit of the model, various statistical indicators such as the coefficient of determination (R 2 ), mean prediction error (MPE), and square root of the mean squared prediction errors (RMSPE) were calculated between the fitted PM 2.5 concentrations from the model and the observations. In addition, a 10-fold cross-validation (CV) technique was adopted to assess potential model over-fitting. A model that has been over-fit could perform better on the data used to fit the model than unobserved data and thus generally has poor predictive performance. The entire model-fitting data set was randomly split into 10 subsets with approximately 10 % of the total data records in each subset. In each round of cross-validation, we selected one subset (10 % of the data) as the testing samples and used the remaining nine subsets (90 % of the data) to fit the model. Predictions of the held-out subset (10 % of the data) were made from the fitted model. The process was repeated 10 times until every subset was tested. Statistical indicators such as R 2 , MPE, and RMSPE were calculated between the CV predicted concentrations and the observations. The model overfitting assessment was conducted by comparing the CV and model-fitting statistics. Crossvalidation also can provide a means to quantitatively assess prediction accuracy for areas where there are no ground observations. A relative accuracy value was also calculated for each year to make validation results comparable among different years.
The daily PM 2.5 concentrations were estimated using the final annual models for 2001 through 2010. The maps of annual mean PM 2.5 concentrations as well as the percent changes between 2001 and 2010 for the study domain and the Atlanta metro area were generated using the daily estimates to visually examine spatial trends of PM 2.5 levels between 2001 and 2010. The percent changes were calculated as follows (3) where PM 2.5 , percentchange denoted the percent changes of PM 2.5 during a study period. PM 2.5,endyear was the PM 2.5 concentrations in the end year of the study period, and PM 2.5,startyear was the PM 2.5 concentrations in the start year of the study period. Moreover, time-series analyses were conducted by year and season, respectively to quantitatively investigate the 10-year temporal trends of fine particle levels in the study domain and the Atlanta metro area.

Descriptive statistics
The descriptive statistics of variables used in fitting the models are listed in Table 1 Table 1 also shows that land use variables and meteorological fields vary from year to year within the data.

Results of model-fitting and validation
The model-fitting and CV statistics (e.g., R 2 , MPE, and RM-SPE) are listed in Table 2. The results show that R 2 ranges from 0.71 to 0.85, MPE is from 1.73 to 2.50 μgm −3 , RM-SPE ranges from 2.75 to 4.10 μg m −3 , and relative accuracy ranges from 72.9 to 80.7 %, which indicates a good fit between the predicted values from the fitted models and the observations. In addition, CV statistics results suggest that some model over-fitting is present; that is, R 2 decreases, while MPE and RMSPE increase from model fitting to crossvalidation, yet the differences are relatively small for all the years. For instance, R 2 and relative accuracy have an average decrease of 0.08 and 4.21 %, respectively, while MPE and RMSPE have an average increase of 0.39 and 0.60 μgm −3 , respectively, through all the years. Moreover, a regression of predicted values against the observations with an intercept at zero (Fig. 2) shows that, at high concentration levels, both model fitting and crossvalidation under-predicted the PM 2.5 concentrations by 3-7% (e.g., fitted/CV PM 2.5 =97% to 93% observed PM 2.5 ). Figure 3 illustrates the PM 2.5 concentration estimates at 1 km spatial resolution in the study area. The annual mean estimated concentrations are 13. 97, 13.90, 13.35, 13.31, 15.19, 13.73, 13.22, 11.34, 10.58, and 11.22 μg m −3 for year 2001 though 2010, respectively. The spatial patterns of PM 2.5 are very similar for all the years. High concentrations appear in large urban centers and along major highways, while low concentrations occur in rural and mountainous areas. In addition, high PM 2.5 levels are also seen in the southeastern part of the study domain. Hu et al. (2014) reported elevated PM 2.5 concentrations measured from monitoring sites located in this region. This area is primarily occupied by agriculture land, and high agricultural emissions may lead to elevated PM 2.5 levels. As reported by previous studies, ammonia (NH 3 ) and nitrogen oxides (NO x ) generated by agricultural activities, such as farm vehicles, domestic and farm animals, and fertilizer applications, can significantly increase the number of suspended particles (Kurvits and Marta, 1998). However, specific agricultural emissions data are needed for further validation. In addition, biomass burning also contributes to emissions of fine particles in the region, following typical seasonal variations (Zhang et al., 2010). Figure 4 shows that the pattern of ground PM 2.5 measurements from FRM monitors corresponds well with that of our estimated concentrations, and the differences between observed and estimated PM2.5 were within ±3 μg m −3 for, on average, 92 % of the monitoring sites for the 10 years, indicating a good agreement between them (Fig. 5).

Spatial trends of PM 2.5 concentrations
To take advantage of the high spatial resolution of the MAIAC data, we generated a map of PM 2.5 estimates in the Atlanta metro area for each year (Fig. 6). The annual mean estimates from 2001 to 2010 are 15. 10, 14.64, 14.00, 14.54, 15.63, 14.39, 14.14, 11.78, 10.98, and 11.65 μgm −3 , respectively. Compared to the last plot of Fig. 6, which illustrates the percentage of impervious surfaces and indicates the level of urban development, the PM 2.5 maps distinctly show that high PM 2.5 levels occur in areas with high urban land use and along major highways, while low concentrations appear in forest and recreational areas, suggesting an underlying positive relationship between air pollution levels and urban development.
As shown in Fig. 7, PM 2.5 concentrations have decreased on average ∼ 20 % for the entire domain and ∼23 % for the Atlanta metro area between 2001 and 2010. Figure 7a illustrates the spatial trend of changes in PM 2.5 levels in the study region. The results show that PM 2.5 levels in most of the areas decreased from 0 to 25 %, and large parts of the areas had decreases exceeding 25 % and as high as 50 %. Larger decreases occurred in more polluted areas such as the Atlanta metro area and along major highways, which might be due to recently enacted emission reduction programs (EPA, 2011) such as the EPA's Clean Air Interstate Rule (CAIR) issued in 2005 (http://www.epa.gov/cair/index.html), since the majority of emissions sources are located in or near urban areas and along major highways. Mitigation of fine particles has been effected by controlling direct PM 2.5 emissions from both stationary and mobile sources (e.g., through installation of scrubbers and filters and the use of alternative fuels and electric vehicles) (EPA, 2007). The mountainous area in the northeastern part of our domain with generally low PM 2.5 levels has also seen substantial decreases of PM 2.5 concentrations. PM 2.5 levels decreased from 25 to 50 % in most of the region, and some areas had decreases exceeding 50 %. By checking 2002 and 2008 point emissions data from the EPA NEI facility emissions reports (due to the lack of 2001 and 2010 data), the decreases are probably due to the dramatically reduced number of emission sources as well as the total emissions in the region during the period. Figure 7b illustrates the percent changes in PM 2.5 levels within the Atlanta metro area. Once again, the spatial trend shows that larger decreases (25 to 50 %) primarily occurred in urban built-up areas and along major highways, while smaller decreases (0 to 25 %) appeared in forest or recreational areas with generally lower pollution levels. This result is expected because the changes of emissions mostly took place in urban built-up areas and along major highways. Two pixels with unusually large changes were identified (in blue and red circles). The large decrease (> 50 %) of PM 2.5 concentration in the blue pixel was due to large emissions reduction from power plants located within that pixel during the period between 2001 and 2010. Likewise, the large increase (> 25 %) of PM 2.5 concentration in the red pixel was due to the addition of a new emission source that did not exist in 2001.
We also illustrated the percent changes between 2001 and 2007, between 2007 and 2008, and between 2008 and 2010 in Fig. 7c- Figure 7g and h illustrate the percent changes between 2008 and 2010. It revealed < 5 % increases in many areas in the eastern part of the domain with increases in some areas exceeding 5 %. On the other hand, many areas in the western part of the domain had < 5 % decreases with decreases in some parts of the mountainous region exceeding 15 %, In the Atlanta metro area, our results show decreases (< 10 %) in urban built-up areas and along major highways with some residential and suburban areas showing < 5% increases. Similarly, in the absence of 2010 emissions data, we could not examine whether these changes were associated with changes in emissions.

Temporal trends of PM 2.5 concentrations
A time-series analysis was conducted to quantitatively examine temporal trends of PM 2.5 levels in the study area as well as the Atlanta metro area during the period between 2001 and 2010 (Fig. 8). The results show our model underestimated PM 2.5 concentrations by 0.99 μg m −3 for the study domain and 1.82 μg m −3 for the Atlanta metro area. This is because satellite-estimated PM 2.5 concentrations included both urban and rural regions, while the ground measurements mostly represented urban conditions. On the other hand, our estimates over monitoring sites matched well with the ground measurements. The mean difference was 0.4μg m −3 for the Atlanta metro area and 0.41 μg m −3 for the study domain. The PM 2.5 levels in the study region as well as the Atlanta metro area had relatively small fluctuations from 2001 to 2007 and from 2008 to 2010, while there was a significant drop in 2008, which was probably due to significant emissions reduction in 2008. The results also reveal seasonal variations of PM 2.5 levels with the highest concentrations occurring in summer and the lowest appearing in winter. Between 2001 and 2010, our time-series analysis showed that the annual mean PM 2.5 concentration decreased ∼ 20 % in the study area and ∼23 % in the Atlanta metro area, which is in line with the findings documented in the US EPA report on particle pollution (EPA, 2011). Both the EPA's findings and our results illustrate a peak in PM 2.5 levels in year 2005, and this phenomenon might be attributed to the increase of sulfate concentrations emitted from electric utilities and industrial boilers during the warm months (e.g, from May to September) of 2005 (EPA, 2008b). In addition, the decrease of PM 2.5 levels after year 2005 likely is due to the emissions reduction programs that have been enacted recently (EPA, 2007, EPA, 2011 such as the EPA's CAIR issued in 2005. Figure 8 distinctly shows this sharp decrease of emissions from 2005 to 2008. In addition, the sharp decrease from 2007 to 2008 also may be partially attributed to the national financial crisis starting in late 2007. The economic slowdown had clear impacts on manufacturing productivity (e.g., the real gross domestic product (GDP) of metro Atlanta in the manufacturing sector dropped 10.2% from 2007 to 2008 (www.bea.gov)), and may also have led to decreases in PM 2.5 emissions.

Discussion
A major strength of this study is that we used high-resolution PM 2.5 estimates derived from MAIAC AOD to investigate spatiotemporal trends of PM 2.5 concentrations in the study area. PM 2.5 estimates at finer resolutions are more suitable for investigation of spatial trends than those at coarser resolutions derived from other AOD products (e.g., MODIS and MISR), because estimates at coarser scales inevitably omit local spatial details, as pointed out by Hu et al. (2014). Our results are capable of showing PM 2.5 concentrations and changes at 1 km resolution, which are very useful for air pollution studies at local scales. For instance, spatial trends of changes in PM 2.5 concentrations in the Atlanta metro area show greater PM 2.5 reduction in more polluted areas (e.g., urban built-up areas and along major highways). Some of the changes may be directly associated with the addition or removal of one or more emission sources as well as the increase or decrease of emissions from those sources. Although high-resolution PM 2.5 estimates can provide more details to examine spatial trends, difficulties lie in their validation to ground monitoring. More ground measurements at specific locations are needed to further validate the results.
Our results of temporal trends of PM 2.5 concentrations correspond well with the EPA's results (EPA, 2011). However, our results show that both ground measurements and satelliteestimated PM 2.5 concentrations over the monitoring sites were generally higher than PM 2.5 estimates over the entire study domain. This is because most of the EPA FRM monitors are located in or near urban areas with generally higher PM 2.5 levels, and therefore observed and satellite-estimated PM 2.5 levels over the monitoring sites reflect mostly urban conditions. Conversely, PM 2.5 estimates over the entire study area account for both urban and rural areas, and therefore the temporal trends of satellite-estimated PM 2.5 concentrations over the entire study domain might be more representative of the true fluctuations of regional PM 2.5 levels. Our future research will continue to explore these associations.
A limitation of our study was that only 3 years of emissions data (2002, 2005, and 2008) were available. The NEI is prepared every 3 years by the EPA primarily based on emissions estimates, emissions model inputs, and supplementary data. As a result, a more quantitative comparison between estimated PM 2.5 and emissions could not be conducted.
We realize that both meteorological fields and land use variables can be potentially incorporated in the second stage GWR model. The primary objective of this study was to investigate the spatial and temporal trends of PM 2.5 levels in the southeastern US. Hence, it was not pursued in this study and will be addressed in future research.

Conclusions
In this paper, we used a two-stage spatiotemporal model incorporating MAIAC AOD data, meteorological fields, and land use variables to estimate PM 2.5 concentrations at 1 km spatial resolution and investigated the 10-year spatial and temporal trends of PM 2.5 levels in the southeastern US. As expected, the satellite model predicted high concentrations in large urban centers and along major highways and low concentrations in rural and mountainous area with relatively high accuracy. Our time-series analysis results indicate that the PM 2.5 levels decreased ∼20 % in the study region and ∼23 % in the Atlanta metro area during the period between 2001 and 2010, and the largest drop occurred between 2007 and 2008. More polluted areas (e.g., in urban areas and along major highways) have also seen greater reduction in PM 2.5 levels, while forests and recreational areas had lower and moderate reduction. PM 2.5 estimates at high spatial resolutions can provide more details in small geographic regions and may reduce exposure misclassification in air pollution and epidemiological studies. High-resolution PM 2.5 data are also useful for air pollution monitoring in large geographic areas because they can not only greatly expand the spatial coverage of costly ground monitoring networks, but they also provide more accurate estimates of population exposure to PM 2.5 . In addition, regional transport of PM 2.5 can be better examined via comparing the spatial distributions of daily high-resolution PM 2.5 estimates, raising new questions for future research.

Supplementary Material
Refer to Web version on PubMed Central for supplementary material.    Relative accuracy (%) *