Evaluation and application of multi-decadal visibility data for trend analysis of atmospheric haze

There are few multi-decadal observations of atmospheric aerosols worldwide. This study applies global hourly visibility (Vis) observations at more than 3000 stations to investigate historical trends in atmospheric haze over 1945–1996 for the US, and over 1973–2013 for Europe and eastern Asia. A comprehensive data screening and processing framework is developed and applied to minimize uncertainties and construct monthly statistics of inverse visibility (1/Vis). This data processing includes removal of relatively clean cases with high uncertainty, and change point detection to identify and separate methodological discontinuities such as the introduction of instrumentation. Although the relation between 1/Vis and atmospheric extinction coefficient (bext) varies across different stations, spatially coherent trends of the screened 1/Vis data exhibit consistency with the temporal evolution of collocated aerosol measurements, including the bext trend of −2.4 %yr −1 (95 % CI: −3.7,−1.1 %yr) vs. 1/Vis trend of−1.6 %yr (95 % CI: −2.4, −0.8 %yr) over the US for 1989–1996, and the fine aerosol mass (PM2.5) trend of −5.8 %yr −1 (95 % CI: −7.8, −4.2 %yr) vs. 1/Vis trend of−3.4 %yr (95 % CI:−4.4, −2.4 %yr) over Europe for 2006–2013. Regional 1/Vis and Emissions Database for Global Atmospheric Research (EDGAR) sulfur dioxide (SO2) emissions are significantly correlated over the eastern US for 1970–1995 (r = 0.73), over Europe for 1973–2008 (r ∼ 0.9) and over China for 1973–2008 (r ∼ 0.9). Consistent “reversal points” from increasing to decreasing in SO2 emission data are also captured by the regional 1/Vis time series (e.g., late 1970s for the eastern US, early 1980s for western Europe, late 1980s for eastern Europe, and mid 2000s for China). The consistency of 1/Vis trends with other in situ measurements and emission data demonstrates promise in applying these quality assured 1/Vis data for historical air quality studies.


Introduction
Atmospheric aerosols have broad implications for air quality and climate change.The Global Burden of Disease (GBD) assessment attributed ambient exposure to aerosol particles with an aerodynamic diameter below 2.5 µm (PM 2.5 ) as the sixth largest overall risk factor for premature mortality with 3.2 million premature deaths per year (Lim et al., 2012).Aerosols are also considered to be the most uncertain component for global radiative forcing (IPCC, 2013).Aerosols are formed from a variety of emission sources and chemical processes with a short tropospheric lifetime against different removal mechanisms, yielding a highly variable spatiotemporal distribution that is not well understood (Fuzzi et al., 2015).Information on long-term aerosol temporal evolution is crucially needed across a range of disciplines.Historical PM 2.5 exposure and its trends are needed to understand changes in Global Burden of Disease (Brauer et al., 2012), and to guide mitigation actions (Apte et al., 2015;Wong et al., 2004).Observations are needed to evaluate historical emission inventories that are crucial to accurately represent the changes in aerosol sources and its consequent feedbacks on climate (Lu et al., 2011; S. J. Smith et al., 2011a;Xu et al., 2013).Aerosol trend analysis is also fundamental to assessing radiative forcing, evaluating model processes, and projecting future changes (Chin et al., 2014;Leibensperger et Published by Copernicus Publications on behalf of the European Geosciences Union. C. Li et al.: Evaluation and trend analysis of multi-decadal visibility data al., 2012;Li et al., 2014).Various studies have been carried out to investigate aerosol trends using in situ measurements (Collaud Coen et al., 2013;Hand et al., 2012a;Murphy et al., 2011), satellite/ground remote sensing (Hsu et al., 2012;Li et al., 2014;Zhang and Reid, 2010), and analysis of measurements with models (Boys et al., 2014;Chin et al., 2014;Pozzer et al., 2015;Turnock et al., 2015).However these studies are mostly limited to the recent 2 decades, since few satellite or in situ aerosol observations exist over land prior to the 1990s.Long-term observations of aerosols at the global scale are needed to place current knowledge of their spatial distribution and temporal evolution in a historical context for all these applications.
Visibility observations offer an alternative information source for the investigation of historical aerosol trends.Horizontal visibility (Vis) from worldwide meteorological stations and airports is mainly determined by the optical extinction (b ext ) of the atmospheric boundary layer, and has been recognized as a proxy of the atmospheric aerosol burden/loading (Husar et al., 2000).Historical Vis data from more than 3000 stations have been applied to characterize decadal trends in global aerosol optical depth (AOD) from 1973 to 2007 (Wang et al., 2009).Regional trend studies of Vis were also conducted for populated areas e.g., the US (Husar et al., 1981;Schichtel et al., 2001), Europe (Vautard et al., 2009) and China (Che et al., 2007;Chen and Wang, 2015;Lin et al., 2014;Wu et al., 2012Wu et al., , 2014)), and the inferred trends were usually attributed to changes in anthropogenic emission.Another study employing Vis over desert regions (Mahowald et al., 2007) found an association of Vis with meteorology factors such as drought index (based on precipitation and temperature) and surface wind speeds.Trends in Vis data interpreted with other data sets also supported studies of several aerosol related climate trends such as the western Pacific subtropical high (Qu et al., 2013) and precipitation (Rosenfeld et al., 2007;Stjern et al., 2011).
Despite the abundance of the above-mentioned studies, the interpretation of Vis data and their trends might be limited by insufficient data processing or poor data quality.Multidecadal Vis data might contain possible variation or even reversal in haze trends as expected from historical emission and surface solar radiation (SSR) data (Lu et al., 2010;Stern, 2006;Streets et al., 2006;Wild et al., 2005).It is of particular interest how these changes would associate with the trends of air quality, and would be captured by the Vis data.Detailed variation in global Vis trends are rarely reported in these previous studies.On the other hand, Vis data are inherently uncertain because most Vis are recorded through human observations with variable protocols.For example, an increase in inverse visibility (1/Vis) has been reported over the US during 1993-2010 (Wang et al., 2012) that is opposite in sign with the significant decline (> 10 % decade −1 ) of observed PM 2.5 , sulfate and b ext (Attwood et al., 2014;Hand et al., 2012aHand et al., , 2014;;US EPA, 2012), and raises questions about the quality of Vis observations.This study revisits the Vis observations to characterize historical trends of atmospheric haze by asserting two major efforts: a more comprehensive data quality assurance processing and a more detailed trend analysis for separate periods.This analysis provides multi-decadal information about air quality evolution and its connections to emission trends over major industrialized regions.To facilitate interpretation, the theoretical relationship between Vis and atmospheric extinction is reviewed in the following section.Section 3 describes the data and processing methods, followed by an evaluation of the screened monthly 1/Vis and its trends using in situ measurements in Sect. 4. Section 5 provides an extensive discussion of the resultant spatial distribution and temporal variation of the derived 1/Vis trends for three highly populated regions (i.e., the US, Europe and eastern Asia), and comparative analysis of these trends with sulfur dioxide (SO 2 ) emission data.The final section summarizes this work and its implications.

Relationship between Vis and b ext
Visibility is a measure of the transparency of the atmosphere, and is defined as the greatest distance at which a black object can be recognized against the horizon sky ( WMO, 2008).The visibility of a particular object (i.e., visibility marker) is determined by the contrast C between the radiation intensity I of the background b and of the object o reaching an observer at distance x from the object: Under assumptions of a plane-parallel atmosphere and homogeneous background intensity (i.e., constant sky brightness), C exhibits an exponential decay based on Beer's law, where b ext is the extinction of the atmosphere (including extinction of aerosols and molecules).Since Vis represents the furthest distance corresponding to a minimum critical contrast C crit below which the observer cannot discern the object, we have Rearranging to solve for b ext yields where K = − ln C crit C 0 .This is the Koschmieder equation (Griffing, 1980), representing a linear relationship between 1/Vis and b ext .The slope K of this relationship is mainly determined by two factors: the inherent contrast at the object's position C 0 and the critical contrast of the observer's eye C crit .This equation is only valid for a plane-parallel and homogeneous atmosphere.For situations with high gradients of b ext (e.g., smoke plumes), this could readily break down.Even for ideal conditions, this relationship could vary due to the variation of C 0 (change of markers or observing conditions) and/or C crit (change of observer or protocol).It is sometimes assumed that the object is perfectly black (C 0 = 1) so that K is only determined by C crit .Nevertheless, K still varies from 1.5 to 3.9 (e.g., Husar and Wilson, 1993;Schichtel et al., 2001;Wang et al., 2009) because of different C crit values or different observing conditions.Below we similarly find that even where 1/Vis is highly correlated with b ext data, K still varies significantly for different stations.

Visibility data
We begin with raw Vis data from synoptic observations over 1929-2013 in the Integrated Surface Database (ISD, https://catalog.data.gov/dataset/integrated-surface-global-hourly-data) archived at the NOAA's National Centers for Environmental Information (NCEI).ISD data are generated through merging hundreds of data sources (A.Smith et al., 2011).The data from different networks have different report frequencies (e.g., hourly, 3-hourly, 6-hourly, etc.).We reject the daily averaged data called "global summary of the day" (GSOD) since an arithmetic mean could bias the daily and monthly statistics because of threshold and discreteness issues, as discussed in Sect.3.1.2.Each processing step is described below.

Conventional screening
We begin with "conventional screening" using algorithms adapted from prior studies.We eliminate effects on Vis of weather conditions such as fog, precipitation, low cloud and high relative humidity (RH > 90 %, estimated from temperature and dew point) following the description in Husar et al. (2000).A sensitivity test that limited conditions to RH < 80 % reduced data density but yielded similar trend results without changing the main findings in this study.Potential human errors are reduced by statistical checks of daily spikes and non-repeating values following Lin et al. (2014).Duplicate stations with different names are combined, and stations lacking geolocation information are removed following Willett et al. (2013).After this screening step, 21 703 stations remain from the 30 895 original ISD sites.

Threshold filtering
We develop a filter to address spatial and temporal variation in the threshold of reported Vis.The "threshold" is the maximum reported Vis at a station that often depends on the furthest employed Vis marker.Vis above this threshold is not resolved.Thus the threshold acts as an artificial detection limit.The ability of Vis data to capture the variation of b ext is weak when the air is clean and/or the adopted threshold Vis at the station is low.We identify the 99th percentile of reported Vis in each year as the threshold for each station, and reject months with ≤ 50 % of the data below the threshold.This approach differs from eliminating stations with low thresholds (e.g., Husar et al., 2000).Observations could still be meaningful at heavily polluted stations even if the threshold is low, while for clean stations with high thresholds most of the reported Vis could remain unresolved.To further ensure data representativeness and variability, data are removed for any month with less than 4 different days of data or with nearly identical percentile values (i.e., the ratio of 50th and 25th percentile Vis is less than 1.07 or the ratio of the 25th to 10th percentile Vis is less than 1.1) following Husar et al. (2000).This data screening step further reduces the number of qualified station to 10 446.
We describe the monthly Vis level with nonparametric statistics rather than arithmetic mean for a few reasons.First, an arithmetic mean would have biased monthly statistics due to the variable fraction (50-100 % after the threshold filtering) of Vis reported under the threshold in 1 month.Second, Vis is recorded as discrete values with coarse and uneven increments, and is not normally distributed (Schichtel et al., 2001).The protocol of reporting Vis varies across stations, depending on local regulations and available Vis markers.Both issues would affect the GSOD data or the monthly mean 1/Vis so we work with the raw data.We follow the convention to adopt the 75th percentile 1/Vis as the monthly representation of haziness (Husar et al., 2000;Qu et al., 2013).Other statistics, such as 50th and 90th percentile 1/Vis lead to similar trends and do not alter the conclusion of this study.However, the 50th percentile is closer to and more vulnerable to the detection limit, while the 90th percentile tends to be more susceptible to extreme events.Husar and Patterson (1987) assessed the effects of different choices of statistics.Below we commonly refer to the 75th percentile as "monthly 1/Vis" unless stated otherwise.

Completeness check
Completeness criteria are applied for further screening.A year of data is removed if less than 6 months in this year is available to guarantee annual representativeness.Short-term time series covering less than 7 years are also removed since they offer little information on trends.A total of 6466 stations comply with these standards and remain in the data archive.

Change point detection
Sudden discontinuities in characteristics of the derived monthly time series of 1/Vis are frequently found even after the comprehensive filtering.Any change of the Vis marker (i.e., change of C 0 ) or observing standard (i.e., change of Figure 1.An example of change point detection and determination based on the time series of 50th (red) and 75th (black) percentiles of monthly 1/Vis.Automatically detected change points are represented by vertical lines.Text in the inset lists the dates of automatically detected points.In this example, five significant change points are identified, in which February 1988 is determined as the separation point for further analysis, while other reported breaks are considered to be false detections.
C crit ) could alter the relationship (K) between b ext and 1/Vis, introducing inconsistency in the time series unrelated to actual b ext change.For example, instrumentation (e.g., telephotometers, transmissometers and scatterometers) has replaced human observers at many sites in the US (Kessner et al., 2013) and to a lesser extent in Europe (Vautard et al., 2009), but there is a lack of documentation recording when and at which stations this switch occurred.Such artificial changes could seriously bias the inferred trends if not addressed.Various methods have been proposed to detect abrupt "change points" (Costa and Soares, 2009;Reeves et al., 2007).For example, the RHtest software package developed for multiple change point detection is based on penalized maximal t and F test (Wang, 2008a;Wang et al., 2007) embedded in a recursive testing algorithm (Wang, 2008b).We adopt the FindU function in the RHtest (version 4, available at http://etccdi.pacificclimate.org/)software to detect "type-1" change points (without reference time series).We manually examine all reported change points for possible false detections.By visually inspecting each remaining station from Sect.3.1.3,we retain only obvious structural discontinuities in the time series of 50th or 75th monthly percentiles from the candidate change points provided by the RHtest results.Figure 1 shows an example of change point detection based on the time series of 50th and 75th percentiles of monthly 1/Vis at one ISD station.The change points are reported in three different types (95 % confidence): significant change, possibly significant (undetermined) change and insignificant change.In this example, although four significant changes for the 50th percentiles 1/Vis and two significant change points for the 75th percentiles 1/Vis are reported, only one candidate (February 1988) indicated by both time series is considered to be an obvious discontinuity and chosen as the actual change point.
The candidate change points provided by RHtest allow greater efficiency than pure manual detection, which is prohibitive for thousands of stations.Any gap of more than 4 years in a time series is also considered to be a change point.Such a large gap could obscure protocol changes and introduce uncertainties in the derived trends without separation.We analyze separately the 1/Vis time series before and after the determined change points.Finally, we eliminate any year of data with annual 1/Vis (average of monthly 1/Vis) less than 40 Mm −1 to address the poor data variation and representativeness of clean environments, as will be discussed in Sect.4.1.
We acknowledge that, although guided by RHtest results and a synthetic analysis based on the time series of 50th and 75th percentiles, this is still a subjective method.A small fraction of determined change points could be extreme events, while a few undetected change points missed by this subjective judgement might remain in the analysis.Several time series with irregular temporal variation are also removed during the visual examination.In summary, only 1/Vis time series considered to be consistent and continuous are analyzed here.
A total of 3930 stations (5320 time series) remain after this processing step, in which 856 sites (22 %) are diagnosed as containing change points and thus separated.This small fraction of structural discontinuities generally has minor impacts on the large-scale trend features and regional trends in Sect. 5 according to our sensitivity test using data without separation.But the separated data reduce spatial incoherency in the derived trend maps, and are more reliable for studies over small areas or independent stations, as shown in Fig. 1.
The threshold filtering (Sect.3.1.2)and change point detection (Sect. 3.1.4)are designed to ensure basic representativeness and continuity of the derived monthly 1/Vis time series, and are the main differences of this processing from prior investigations.Figure 2 (bottom) shows that most US stations are screened after the mid 1990s.This is because more than 90 % of the ISD stations gradually switched to employ a low Vis threshold of 10 miles (∼ 16 km) after the mid 1990s (Fig. A1), probably due to the introduction of unified instrumentation (Kessner et al., 2013).A maximum Vis of 16 km may be sufficient for airport navigation and weather reports, but this threshold Vis under clear sky conditions represents a moderate pollution level, and clean cases are not resolved.Thus most of the US stations with such low thresholds are rejected during the threshold screening.In contrast, screened stations remain densely distributed with long-term data over Europe and eastern Asia after the mid 1990s because the adopted thresholds are generally higher and more consistent (Fig. A1).

Complimentary in situ data
We adopt complimentary data to evaluate and interpret the constructed monthly 1/Vis time series and trends.The measured and calculated aerosol optical data from the Interagency Monitoring of PROtected Visual Environments (IMPROVE) program (http://vista.cira.colostate.edu/improve/Data/data.htm)are employed to evaluate the screened 1/Vis data and its trends after 1988.IMPROVE applies empirical mass extinction and RH growth factors to measured mass of aerosol components to calculate and report ambient b ext in a 3-4-day frequency (Pitchford et al., 2007), and for several stations concurrent measurements of aerosol scattering coefficient (b sp ) are also made at hourly frequency using nephelometers.We generate monthly mean total b ext (including aerosol extinction and site-specific Rayleigh scattering) and b sp from data with RH < 90 % and status flags as "V0" (valid).Any month with less than 4 available days for averaging is abandoned.Pitchford et al. (2007) demonstrated that the estimated b ext is consistent with measured b sp .We also find high correlation (r = 0.90, N = 3439) between monthly b ext and b sp across IMPROVE stations (Fig. A2).
The measurement of b ext or b sp is sparse outside the US.Therefore we obtain long-term measurements of fine particulate matter mass (PM 2.5 ) from the European Monitoring and Evaluation Programme (EMEP, http://ebas.nilu.no)for comparison over Europe (Tørseth et al., 2012).Forty-five stations of data collected by filter-based ambient samplers are used.Similarly, these daily PM 2.5 data are averaged monthly provided at least four valid measurements are available.

SO 2 emission data
We apply bottom-up total anthropogenic SO 2 emission inventories to interpret historical 1/Vis trends.This approach exploits the close relation of sulfate aerosol concentration with SO 2 emission due to the short time scale of SO 2 oxidation (Chin et al., 1996(Chin et al., , 2014;;Daum et al., 1984;Hand et al., 2012a), the major PM 2.5 contribution from sulfate aerosols over land for most populated areas (Chin et al., 2014;Philip et al., 2014), and the dominance of sulfate for light extinction due to its hygroscopicity (Hand et al., 2014).We employ three different SO 2 emission data sets, including country-level data for 1850-2005 (S.J. Smith et al., 2011a, b), gridded data from EDGAR (Emissions Database for Global Atmospheric Research) version 4.2 (EC-JRC/PBL, 2011) at 0.1 • resolution for 1970-2008 (http: //edgar.jrc.ec.europa.eu/),and data from Lu et al. (2011) at 0.5 • resolution for 1996-2010 over China.The data from S. J. Smith et al. (2011a) are referred to as "Smith emissions" below.The data from Lu et al. (2011) are referred to as "Lu emissions".

Trend analysis
In this study, we separately calculate trends for several periods of 8-10 years to allow possible trend reversal, and to include stations with short-term data.The choice of study periods is mainly based on the historical SO 2 emission data.Figure A3 shows the Smith emission data for several representative countries.SO 2 emission trends in the US changed direction at ∼ 1944, ∼ 1954, and again at ∼ 1973.Also, for most eastern European countries, there is a sharp reduction of SO 2 emission starting from ∼ 1989 after the breakdown of the communist system, while the 1997 Asian financial crisis affected the SO 2 emission trend in Korea.It is of particular interest to examine how Vis is affected by these emission changes.Data for most ISD stations outside the US start from the year 1973, and representative coverage of Vis stations over the US starts from the year 1945, although the earliest records after screening start from 1929.Based on these transition points of SO 2 trends and Vis data availability, eight periods (1945-1953, 1954-1963, 1964-1972, 1973-1980, 1981-1988, 1989-1996, 1997-2005, 2006-2013) are chosen to be analyzed in detail over the US, while the latter five periods are studied for Europe and eastern Asia.We also briefly examine two short periods before 1945 (1929-1934 and 1935-1944) over the US where stations are less spatially representative (not included in regional quantitative analysis) but still show prominent trend information in 1/Vis.
We assess the linear trend and its significance (p value, two-tail test) in the deseasonalized monthly anomalies using Sen's slope (Sen, 1968) and the Mann-Kendall (MK) test (Kendall, 1975;Mann, 1945).All monthly data are deseasonalized by removing multi-year monthly means of each period before trend estimation.Pre-whitening is introduced to reduce the effect of lag-1 autocorrelation (Yue et al., 2002), and 95 % confidence interval (CI) of the slope is calculated (Li et al., 2014).This nonparametric trend estimation method is insensitive to missing values and outliers in the time series, and does not require a normal distribution, thus it has been widely adopted to study aerosol trends in previous studies (Collaud Coen et al., 2013;Papadimas et al., 2008).Least square trends (Weatherhead et al., 1998) are also calculated, and are found to be consistent with the MK-Sen trends.For all the 8027 calculated slopes in 1/Vis, 88 % are unanimously diagnosed as significant (90 % confidence, p < 0.1) or insignificant by both methods.For the significant trends 76 % of their differences are within 20 %.Relative trends are cal-culated by normalizing the absolute MK-Sen slopes to the multi-year mean of monthly 1/Vis in the corresponding period to facilitate the comparison and interpretation with other in situ data.
Short-term trends of 8-10 years are expected to be less statistically robust and more sensitive to extremes.For each period, a time series is required to contain at least half of the total months and two-thirds of the total years (e.g., at least 60 monthly data in at least 7 years for a 10-year period) for the calculated trend to be representative.This step only reduces the number of stations at which trends are reported, but does not further screen the data.
The meaning and observing methods of daytime and nighttime data differ.According to WMO (2008), Vis at night, as determined using illuminated objects, also depends on the light source intensity, the adaptation of the observer's eyes to darkness and the observer's illuminance threshold.We compare the relative trends calculated using daytime and nighttime data to the combined trends adopted in this paper, over all remaining sites and the eight periods.The 5183 daytime trends have a correlation of 0.85 with the combined trends, in which 84 % of the differences between significant trends (p < 0.1) are within 50 %.For the comparison between 4109 nighttime and combined trends, the correlation is 0.80 and 78 % of the differences between significant trends are within 50 %.Therefore, after representing the data as a monthly resolution and normalizing the changes in 1/Vis to relative trends, the daytime and nighttime data show generally consistent trends in haze level compared to the combined data, and do not meaningfully alter our results and conclusions.
We calculate composite trends based on monthly 1/Vis averaged from an ensemble of stations (e.g., for the time series of collocated stations in Sect. 4 or defined regions in Sect.5).To ensure temporal representativeness, a station is considered in the average only if two-thirds of the total months of data are available for the study period.Qualified stations are gridded to 1 • resolution before averaging to avoid biased averaging towards more densely distributed areas.To ensure spatial representativeness, only monthly data derived from at least 75 % of the total grids (i.e., number of unique grids covered by all the monthly data) for each study period are used in the composite trend estimation.This strategy reduces sampling difference within each periods; however, the composite 1/Vis for different periods might be averaged from a different distribution of stations.We expect the uncertainty from spatially variant K and data quality to be random, and to be reduced by spatial averaging and by normalizing the slopes into relative trends.Over these regions, we also calculate several time series and trends for longer merged periods with consistent station coverage and similar trends, to assess the consistency of the short-term trends.We compare the monthly IMPROVE b ext data with the quality controlled monthly 1/Vis from Sect.3.1.Collocations are considered between IMPROVE and ISD time series over 1988-2013 within the distance of less than 1 • and altitude difference of less than 500 m.One IMPROVE station could pair with more than one ISD station and vice versa.A total of 59 collocations (each with at least 20 paired monthly values) are made.We expect a maximum correlation of 0.9 given the relation between measured b sp and calculated b ext (Fig. A2).Similarly, we create collocations between ISD 1/Vis and EMEP PM 2.5 on a monthly basis, and expect a weaker correlation due to variation of aerosol water and mass extinction efficiency.
Figure 3 shows the comparison results between collocated 1/Vis and b ext over the US.This evaluation highlights the following major findings.2. As shown in the top right panel, correlation coefficients of monthly values vary from ∼ 0 to 0.85.About half of the collocations (29 out of 59) have r < 0.5, while 10 collocated ISD stations have r > 0.7.The overall moderate correlation is not unexpected, as is similarly found in previous studies (Mahowald et al., 2007;Wang et al., 2012).Correlations are expected to differ from station to station, due to the inherent difference in observing conditions, protocols, and residual uncertainties.This preliminary evaluation suggests that Vis data at individual stations can be unreliable, and in the following discussion we focus on interpreting regionally coherent observations.
3. Correlations generally exceed 0.5 in the eastern US, where the mean b ext is higher due to higher aerosol concentration (Hand et al., 2012b;van Donkelaar et al., 2015) and to a larger fraction of hygroscopic sulfate aerosols (Hand et al., 2012b).4. The slope of fitted linear relationship (bottom left) between 1/Vis and b ext varies from ∼ 0.8 to ∼ 2 even over the eastern US where correlations are higher.This supports the expectation that this slope (K) would differ spatially with observing conditions (Griffing, 1980;Husar et al., 2000;Schichtel et al., 2001), as discussed in Sect. 2. Thus in the later analysis we focus on the relative trend of 1/Vis which is independent of K.
Figure 3 (bottom right) also shows the correlation between monthly 1/Vis and PM 2.5 over Europe.Although the relation of 1/Vis with PM 2.5 is expected to be more uncertain than with b ext , we find more stations with high correlation (r > 0.5) over Europe (93 out of 129, 72 %) than over the US (51 %).Wang et al. (2012) similarly found higher correlation of 1/Vis with PM 10 over Europe and China than over the US and Canada.The higher thresholds and higher concentration of fine aerosol over Europe (van Donkelaar et al., 2015) allow 1/Vis to better resolve PM 2.5 variation there.These findings suggest more reliability of Vis observations at areas with both higher aerosol loading and sufficiently high thresholds to resolve b ext variation, e.g., the three populated regions investigated in this study.

Trend evaluation
Figure 4 shows the spatial distribution of relative trends in 1/Vis, in IMPROVE estimated b ext and in measured b sp over the US for 1989-2013.Overall, the trend maps of 1/Vis, b ext and b sp show a dominant trend of decreasing haziness over the continental US after 1988, which reflects reduction of aerosol sources (Hand et al., 2014;Leibensperger et al., 2012).The overall decrease across the US is consistent with recent trend studies employing IMPROVE b ext (Hand et al., 2014) and b sp (Collaud Coen et al., 2013) data, and is determined by the reduction of both aerosol mass and hygroscopicity (Attwood et al., 2014).For the last two periods   (1997-2013), the number of available ISD stations for trend analysis is dramatically reduced by their detection limit and improved air quality.Although the remaining sparse ISD stations still show overall consistency in trends with nearby b ext and b sp , they cannot provide spatially coherent and aggregated trend information.We thus suggest that the ISD Vis data over the US are not appropriate for studying haze trends after the mid 1990s, and limit our analysis to data before 1996 for this region.Over 1989-1996, the 1/Vis trends still reproduce the b ext trends, with decreasing tendencies in the eastern and western US.For this period, 15 ISD stations and 9 IMPROVE stations with significant trends are collocated and labeled.Thus the apparent discrepancy in sign of trends in 1/Vis (Wang et al., 2012) with trends in other aerosol measurements (Attwood et al., 2014;Hand et al., 2012aHand et al., , 2014; US EPA, 2012) is resolved by more comprehensive data processing and screening.
Figure 5 shows the spatial distribution of relative trends in 1/Vis and PM 2.5 over Europe for 2006-2013.There is a tendency of greater reductions in 1/Vis over western Europe than over eastern Europe as examined further in Sect.5.2.The dominant decreasing trends of PM 2.5 are adequately captured by the 1/Vis trends, especially at the 19 ISD and 10 EMEP collocated sites with significant trends, as discussed further below.
Figure 6 (top) shows the composite time series of the collocated 1/Vis and b ext stations over the US for 1989-1996.The seasonal variation of the averaged b ext is well reproduced by that of collocated 1/Vis, with a correlation of 0.77 between these two time series.Both composite 1/Vis and b ext show a peak in summer months, due mostly to increased aerosol con-centration in warm months because of increased photochemical activity and biogenic emission (Chen et al., 2012;Hand et al., 2012b).The trend of collocated 1/Vis (−1.6 % yr −1 ; 95 % CI: −2.4,−0.8 % yr −1 ) is within the confidence intervals of the decrease of b ext (−2.4 % yr −1 ; 95 % CI: −3.7, −1.1 % yr −1 ).The slight underestimation may reflect the weak sensitivity of discrete 1/Vis data to the continuous decrease of b ext in clean environments due to the threshold and discreteness issues.
Figure 6 (bottom) shows composite time series of PM 2.5 and 1/Vis of these collocated 1/Vis and PM 2.5 stations over Europe for 2006-2013.High correlation (0.80) between these two time series indicates consistent seasonal variation.The winter maximum in the composite 1/Vis over Europe well represents the PM 2.5 seasonality at most collocated EMEP sites, which could be attributable to near surface inversion and low surface winds (Yttri et al., 2012), to greater nitrate aerosol formation (Aas et al., 2012;Yttri et al., 2012), and to higher carbonaceous aerosol emission from residential wood combustion (Denier van der Gon et al., 2015).The CI of the 1/Vis trend (−3.4 % yr −1 , 95 % CI: −4.4,−2.4 % yr −1 ) overlaps with that of the PM 2.5 trend (−5.8 % yr −1 , 95 % CI: −7.8, −4.2 % yr −1 ), but underestimates the relative decrease of PM 2.5 .In addition to the weak sensitivity of discrete 1/Vis to resolve aerosol variation under clean environment (the collocated EMEP stations are mostly in the cleaner western Europe), the inclusion of Rayleigh scattering in 1/Vis and the non-linear association between ambient 1/Vis and dry PM 2.5 (fixed at 50 % RH) also contribute to this bias.In summary, 1/Vis exhibits spatially variant K (i.e., relationship with b ext ) and data quality that suggests uncertainty in the information of one station especially at clean locations.However the aggregated 1/Vis time series successfully capture the seasonal variation and trends of collocated in situ data.The high correlation between composite time series and the overall consistency of composite trends suggest that the interpretation value of 1/Vis data benefits from averaging over multiple stations.

United States
Figure 7 presents the calculated relative trend of 1/Vis of all qualified stations over the US for 1945for -1988 (Fig. 4 (Fig. 4 contains 1/Vis trends over 1989-2013).Figure 8 shows the regionally averaged time series and trends of 1/Vis over the eastern US for 1945-1996, superimposed with the evolution of SO 2 emission data.Historically, 1/Vis in the eastern US experienced a pronounced decrease (−2.8 % yr −1 , p < 0.001) after World War II until the mid 1950s, a consistent upward trend afterwards (0.9-1.8 % yr −1 , p < 0.001) during the following 2 periods until the early 1970s, variable tendencies during 1973-1980, and a significant decreasing trend (−1.1 to −2.0 % yr −1 , p < 0.005) from the early 1980s until 1996.Over 1954Over -1973, the long-term trend of 1/Vis is 1.2 % yr −1 (p < 0.001), lying between the separated short-term trends.This 1/Vis trend evolution resembles the SO 2 emission trend.Industrial activity gradually decreased after World War II un-til the mid 1950s, followed by economic growth until the early 1970s with the emergence of both the oil crisis and the Clean Air Act (Greenstone, 2001).The emission of SO 2 starts to consistently decrease after 1973 for the Smith inventory, and after 1977 for the EDGAR inventory.For the period 1973-1980 the regional 1/Vis is generally consistent with these two inventories except for an anomalous peak of annual 1/Vis in 1977-1979.The NOAA Climate Extremes Index (http://www.ncdc.noaa.gov/extremes/cei/)describes the winters of 1977-1979 as the coldest during 1945-1996 across the US.Increased emissions from domestic heating, as well as stagnant weather may contribute to the 1/Vis peak.After 1978, the three annual time series uniformly exhibit a downward tendency.
Table 1 contains the correlation of annual 1/Vis with SO 2 emissions.Annual 1/Vis over the eastern US exhibits a correlation of 0.66 with the Smith SO 2 emissions over the entire US (1946US ( -1995)), and of 0.73 with the EDGAR SO 2 emissions over the eastern US (1970US ( -1995)).The 1/Vis trends over the western US (where SO 2 emissions are much lower than in the eastern US, organic aerosols dominate in PM 2.5 and forest fires are more prevalent) are less consistent than over the eastern US with the SO 2 emission data, given the influence of other sources.In summary, the 1/Vis time series successfully capture large-scale haze evolution over the eastern US from 1945 to 1996, which is consistent with changes in SO 2 emissions as well as previous investigations on 1/Vis for this region (Husar and Wilson, 1993;Schichtel et al., 2001).
Figure A4 shows the calculated 1/Vis trends over the US for two short periods prior to 1945.Although the stations  of the year, where all annual values are plotted.Light green dots represent the average monthly 1/Vis of all qualified stations (error bars showing the 25th and 75th percentile) in the defined region.Red dots show the number of grid cells for averaging, and data gaps indicate months with less than 75 % of the total grids for each period.Blue lines and text represent the 1/Vis trends calculated using the monthly anomalies for each period.Trends in parentheses are the 95 % confidence intervals.Black lines are the annual 1/Vis averaged from at least 8 monthly values.SO 2 emissions for the entire US from S. J. Smith et al. (2011a) are in orange.Purple indicates EDGAR SO 2 emissions for the entire US (dashed) and for the defined region (solid) in Fig. 7. are sparsely distributed, the nearly uniform trends in 1/Vis strongly suggest a prominent decrease over 1929-1934, and then a rapid increase over 1935-1944.This evolution reflects the significant drop in industrial activity following the 1929 Great Depression, and the economic recovery after ∼ 1933 during the New Deal programs and World War II.The Smith SO 2 emissions of the US (Fig. A3) also reflect these socioeconomic events.

Europe
Figure 9 presents the spatial distribution and temporal evolution of haze trends over Europe as derived from the 1/Vis data for 1973-2005.The historical trend pattern of 1/Vis is quite different between western and eastern Europe.The large-scale 1/Vis trend over western Europe is consistently decreasing for the four periods after 1981 (also in Fig. 5).Some countries such as the UK and France begin decreasing prior to 1981, consistent with the SO 2 emission decrease over these countries (Fig. A3).Prior analysis also indicated Vis improvements after ∼ 1973 for most sites over the UK (Doyle and Dorling, 2002).Meanwhile stations over eastern Europe have significantly increased 1/Vis for 1973-1980, a mostly decreasing trend in its western part for 1981-1988, and then a decrease-dominant trend after 1989.Red dots show the number of grid cells for averaging, and data gaps indicate months with less than 75 % of the total grids for each period.Blue lines and text represent the 1/Vis trends calculated using the monthly anomalies for each period.Trends in parentheses are the 95 % confidence intervals.Black lines are the annual 1/Vis averaged from at least 8 monthly values.The Smith SO 2 emissions in orange are the total emission of all countries listed in Table 2 for each region.The EDGAR SO 2 emissions in purple are summed from all pixels inside the defined region (Fig. 9).
Figure 10 shows the regionally composite time series of 1/Vis as well as SO 2 emissions over western and eastern Europe for 1973-2013.Table 2 lists the specific country names included in the Smith emissions for the two regions.The evolution of 1/Vis over western and eastern Europe is broadly consistent with the SO 2 emissions, and reflects the lag of emission reduction in eastern vs. western Europe.Stjern et al. (2011) similarly reported later improvement in Vis over eastern vs. western Europe.The SO 2 emission reduction extends from the 1980s to the end of the data record for western Europe, and primarily over 1989-2000 for eastern Europe.The composite 1/Vis time series successfully capture the significant reduction of haze over western Europe (−1.1 to −1.7 % yr −1 , p < 0.08).Long term 1/Vis trend over western Europe for 1981-2011 (insufficient qualified stations after 2011) is −1.8 % yr −1 (p < 0.001), consistent with the separate short-term trends.For eastern Europe the decrease of 1/Vis is stronger before 1997 (−2.0 % yr −1 , p < 0.001) than after 2006 (−1.1 % yr −1 , p = 0.03), and the calculated trend over 1997-2005 is insignificant, consistent with the SO 2 emission evolution.There is an obvious peak in 1/Vis from October 1995 to March 1996 especially over eastern Europe, which is consistent with the peak sulfate concentration that Stjern et al. (2011) attributed to the anomalously cold winter of 1996 with stagnant air.
Table 1 shows that the annual 1/Vis time series exhibit a correlation of 0.91 (0.92) with the Smith Emissions for 1973-2005, and of 0.92 (0.92) with the EDGAR emissions for 1973-2008 over western (eastern) Europe, respectively.Such high correlations suggest a major role of SO 2 emissions to determine the decadal trends of haze over Europe.Table 2. List of countries included to calculate regional SO 2 emission from the country-level emission data (countries with most parts inside the defined region) of S. J. Smith et al. (2011a

Eastern Asia
Figure 11 shows the calculated relative trends of 1/Vis over eastern Asia after 1973.A persistent increasing trend of 1/Vis dominates over eastern China for more than 30 years.A prominent feature in the trends over China is more hetero-geneity in the spatial distribution compared to the trend maps over the US and Europe.This could be a result of asynchronous economic development, as several studies reported "lagging" of Vis impairment in rural sites (from ∼ 1990s) compared to urban sites (from ∼ 1960s) in China (Quan et al., 2011;Wu et al., 2012).The overall increasing trend in 1/Vis reverses in the last period of 2006-2013, when most stations in southern China and many in northern China show a statistically significant decreasing trend of 1/Vis.This is consistent with the implementation of fuel-gas desulfurization facilities in power plants after ∼ 2007.This recent reduction was also supported by satellite observations of SO 2 (Li et al., 2010;Lu et al., 2010Lu et al., , 2011;;S. Wang et al., 2015;Zhao et al., 2013).
Figure 11 also shows a consistent increase of 1/Vis over Korea from 1973 to 1996.After 1997 when the SO 2 emission transits to decrease (Fig. A3), the increase in 1/Vis levels off and reverses.The aerosols over China also affect areas downwind through long-range transport (Aikawa et al., 2010).For the 1997-2005 period, most eastern stations of Korea show a downward trend, in contrast with the increasing 1/Vis Figure 12.Regional time series analysis of 1/Vis and SO 2 emission over southern and northern China.The long ticks on the horizontal axis indicate January of the year, where all annual values are plotted.Light green dots represent the average monthly 1/Vis of all qualified stations (error bars showing the 25th and 75th percentile) in the defined region.Red dots show the number of grid cells for averaging, and data gaps indicate months with less than 75 % of the total grids for each period.Blue lines and text represent the 1/Vis trends calculated using the monthly anomalies for each period.Trends in parentheses are the 95 % confidence intervals.Black lines are the annual 1/Vis averaged from at least 8 monthly values.The SO 2 emission in Lu et al. (2011) in orange and the EDGAR SO 2 emission in purple are summed from all pixels inside the defined region (Fig. 11).over the west, which is more strongly influenced by pollutant transport from China. Lee et al. (2015) also discovered insignificant improvement of Vis over urban areas of Korea after late 1990s despite the national emission reduction policy launched in early 2000s, which was attributed to the regional transport from upwind continental areas.Longterm aerosol measurement over Gosan Island, Korea showed rapid increase of sulfate and nitrate concentrations from early 2000s to ∼ 2006, which were closely related with the trends of China's emission (Kim et al., 2011).Similarly, stations over the western and coastal areas of Japan consistently exhibit an upward 1/Vis trend before 2006, despite the continuous decrease of local SO 2 emission and concentration since 1970 (Wakamatsu et al., 2013).Aikawa et al. (2010) found a zonal gradient in terms of both the magnitude and trend of measured SO 2 and sulfate concentrations over Japan, and in the modeled contribution from China to the sulfate concentration in Japan.Lu et al. (2010) reported that most EANET (Acid Deposition Monitoring Network in East Asia) stations over Japan and Korea have increasing trends in SO 2 and sulfate aerosols from 2001 to 2007.For the last period 2006-2013, 1/Vis shows a dominant decreasing trend over Japan and Korea that may reflect in part China's SO 2 emission controls.Itahashi et al. (2012) reported a trend reversal of MODIS (Moderate Resolution Imaging Spectroradiometer) fine aerosol optical depth (AOD) over the Sea of Japan from increasing to decreasing at ∼ 2006 that is more consistent with China's SO 2 emission than the local emission.This analysis highlights the sensitivity of 1/Vis to long range transport, and the value of international collaboration for air quality improvement over eastern Asia.
Figure 12 presents a regional analysis of averaged 1/Vis time series over northern and southern China, and the evolution of SO 2 emissions from two inventories.The overall Vis impairment trend in China for 1973-2005 reflects the consistent SO 2 emission increase.Both the north and south show C. Li et al.: Evaluation and trend analysis of multi-decadal visibility data a steady and significant (p < 0.001) increase of haziness for the 1973-1980 period, and southern China shows an even faster impairment (2.9 % yr −1 ) than the north (1.2 % yr −1 ).For the next 2 decades  the 1/Vis increase slows down in both the south and the north, in accordance with other investigations using Vis and SSR data (Chen and Wang, 2015;Luo et al., 2001;Wu et al., 2014).The south exhibits a slower (0.2 % yr −1 ) and less significant (p > 0.3) increase than the north (0.5-0.6 % yr −1 ).The long-term trend over 1981-1996 for northern China (0.5 % yr −1 , p < 0.001) also exceeds that for southern China (0.2 % yr −1 , p = 0.04).This difference is determined not only by the slower increase of SO 2 emissions in the south (Lu et al., 2010), but also by more precipitation and ventilation in the south that favors the removal of aerosols and their precursors (Xu, 2001;Ye et al., 2013).The decline of SO 2 emissions from 1996 to 2000 reflects both the 1997 Asian financial crisis, and a decline in coal use and sulfur content (Lu et al., 2011).Both regions show a leveling off or even reversal of 1/Vis increase during this short period, which is again more significant in the south.The period 2000-2006 exhibits significant growth (> 1 % yr −1 ) of 1/Vis in both the north and south, resembling the steady growth in SO 2 emissions.The recent reduction of SO 2 emissions is reflected in the Lu emissions while not in the EDGAR emissions.After 2006, significant (p < 0.05) decreasing trends in 1/Vis are apparent (−0.9 to −1.6 % yr −1 ) for both northern and southern China, which is more consistent with the Lu emissions.As shown in Table 1, the annual 1/Vis time series exhibit a high correlation of 0.78 (0.87) with the Lu emissions (1996Lu emissions ( -2010)), and of 0.91 (0.88) with the EDGAR emissions  over northern (southern) China, respectively.

Connections to SSR and AOD trends
Long-term records of surface solar radiation (SSR) and columnar aerosol optical depth (AOD) serve as complimentary data resources to study and interpret changes in air pollution during the last few decades, especially for regions with fewer ground-based aerosol measurements.SSR is determined by the total columnar extinction of aerosols and clouds while 1/Vis represents the extinction level at the surface.Moreover, the direct scattering and absorption of solar radiation by aerosols could be amplified in less polluted regions or dampened over highly polluted stations, due to aerosol-cloud interaction (Fuzzi et al., 2015;Wild, 2009).Despite these uncertainties, the observed reversals of SSR from "dimming" to "brightening" in 1980-1990 over the US and Europe (Streets et al., 2006;Turnock et al., 2015;Wild, 2012) generally agree with the reversals around the 1980s of 1/Vis trends in this study.Over China, the recently reported decadal SSR variation shows dimming before the 1990s and no significant trend afterwards (Tang et al., 2011;K. Wang et al., 2015).The latter phenomenon may reflect compensation of more aerosol extinction by less cloud cover (Norris and Wild, 2009).
Reliable AOD data over land are limited to the recent 2 decades, but exhibit even greater consistency with 1/Vis trends.The recent decrease in 1/Vis after late-1990s over the US and western Europe in this study is consistent with previous studies on AOD trends based on both ground based (e.g., Li et al., 2014;Yoon et al., 2012) and satellite (e.g., Chin et al., 2014;Hsu et al., 2012;Pozzer et al., 2015) observations.Over China, several studies on AOD trends in the 2000s showed notable increasing tendency (e.g., Hsu et al., 2012;Pozzer et al., 2015;Yoon et al., 2012), while some recent studies also discovered that separating AOD time series could reflect the plateauing and reversal of trends in recent years due to emission control strategies (Che et al., 2015;He et al., 2016;Lu et al., 2011).PM 2.5 trends derived from satellite AOD over 1998-2012 have decreasing tendencies over North America and Europe, and increasing tendencies over eastern Asia (Boys et al., 2014;van Donkelaar et al., 2015), similar to the 1/Vis trends found here.

Conclusion
This study examines Vis observations as a trend indicator of haziness and air quality over the US , Europe (1973-2013), and eastern Asia (1973-2013).We comprehensively process the raw data from over 20 000 stations considering effects from meteorological factors, protocol design, and human errors.We develop filters to exclude relatively clean cases (i.e., months with ≤ 50 % records below the threshold Vis, or years with annual 1/Vis ≤ 40 Mm −1 ) with weaker sensitivity to b ext variation, and apply change point detection and separation to largely reduce the intrinsic discontinuities.Nearly 4000 stations remain after the processing with 753 over the US, 1625 over Europe, and 791 over eastern Asia.The composite time series of 1/Vis over the US for 1989-1996 generally agrees with the collocated IMPROVE b ext in terms of both seasonal variation (r = 0.77) and trends (−1.6 % yr −1 , 95 % CI: −2.4,−0.8 % yr −1 ) in 1/Vis vs. b ext (−2.4 % yr −1 , 95 % CI: −3.7, −1.1 % yr −1 ).Similarly, for 2006-2013 over Europe, the seasonal variation (r = 0.80) and significant decrease (−5.8 % yr −1 , 95 % CI: −7.8, −4.2 % yr −1 ) in PM 2.5 are captured by collocated 1/Vis (−3.4 % yr −1 , 95 % CI: −4.4,−2.4 % yr −1 ).This consistency highlights the benefits of thorough data screening to reduce uncertainties brought by the inherent issues in Vis observations such as threshold choices, discreteness and discontinuities.As discussed in Sect.3.1, the inclusion of unresolved values in the mean 1/Vis and the contaminants of discontinuities could dampen the ability of 1/Vis to correctly resolve aerosol trends.Admittedly, the derived 1/Vis trends are still subject to several uncertainties, e.g., the spatially variant K and data quality, the less robust short-term trends, sampling differences and direct averaging in composite time se-ries.Nevertheless, the interpretation value of 1/Vis data is shown to be enhanced by the comprehensive screening and spatial averaging.Therefore we focus on the trend results that are regionally coherent and aggregated, and avoid drawing strong conclusions based solely on the 1/Vis trends.Although at individual stations the 1/Vis changes might be affected by these above-stated artificial factors, regionally coherent trend signals suggest these derived 1/Vis trends represent actual changes in b ext .Our filtered monthly 1/Vis data are freely available as a public good (http://fizz.phys.dal.ca/~atmos/martin/?page_id=2527).
Analysis of the 1/Vis trends for several short periods reveals haze trend evolution and reversals.These historical 1/Vis trends and their evolution also exhibit compelling consistency with SO 2 emissions and SSR studies.For example, 1/Vis shows statistically significant decreasing trends from the late 1970s to the mid 1990s over the eastern US (−1.1 to −2.0 % yr −1 ), from the early 1980s to 2013 over western Europe (−1.1 to −1.7 % yr −1 ), in the early 1990s (−2.0 % yr −1 ) and after the mid 2000s (−1.1 % yr −1 ) over eastern Europe, and after the mid 2000s over China (−0.9 to −1.6 % yr −1 ).These recent decreases in 1/Vis are attributable to emission changes in these populated areas.Reversal points of 1/Vis trends also consistently reflect several historical socioeconomic events e.g., the New Deal programs (from decrease to increase at ∼ 1934), the end of World War II (from increase to decrease at ∼ 1945) and the Clean Air Act (from increase to decrease at ∼ 1979) in the US, the collapse of communism in eastern Europe (from increase to decrease at ∼ 1989), and the 1997 Asian financial crisis.
Therefore, the constructed 1/Vis data are applicable to resolve historical aerosol trends on a regional and annual basis, and provide complementary information about the historical changes in air quality.For instance, the annual 1/Vis time series exhibit high correlations (0.7-0.9) with SO 2 emissions for five large domains (Table 1).Apart from verifying the historical 1/Vis trends, this consistency also provides an evaluation of emission inventories.For example, after ∼ 2006 1/Vis trends agree better with Lu et al. (2011) than the EDGAR emissions in capturing the SO 2 emission controls over China.Emission inventories differ significantly (S.J. Smith et al., 2011a), and 1/Vis data offer constraints on these inventories.
However, SO 2 emission inventories cannot fully explain the trends in ambient haze due to the influence of other emissions and meteorological factors.Notable reductions in emissions of nitrogen oxides and black carbon have been reported over North America and western Europe (Bond et al., 2007;Lu et al., 2015;US EPA, 2012;Vestreng et al., 2009), while steady increase in emissions of nitrogen oxides, organic carbon and black carbon were identified over China (Lu et al., 2011;Zhao et al., 2013).Observed (Leibensperger et al., 2012;Murphy et al., 2011) and simulated (Lin et al., 2010;Wang et al., 2013) changes in various aerosol chemical species suggest increasing importance of emissions other than SO 2 on air quality trends in recent years.We have also shown that occasional cold winters in the US and Europe, and the long-range transport of China's pollutants into Korea and Japan could affect the association between 1/Vis and local emission.Future work includes applying a chemical transport model to further interpret the observed 1/Vis (b ext ) trends, as well as the contribution from meteorology and emissions.

Figure 2 (
Figure 2 (top) shows the ISD stations and the number of years with available data for 1929-2013 before and after data processing.Most of the remaining stations are located in the US (753), Europe (1625) and eastern Asia (791).More than 6000 removed stations have less than 7 years of data as indicated in the left panel.Many other removed stations have small population density or harsh observing environment (e.g., islands and polar regions), which might correspond to poor observing conditions or maintenance.

Figure 2 .
Figure 2. Distribution of Integrated Surface Database (ISD) stations before (left) and after (right) data screening.Colors indicate the number of years with available visibility data for (upper) 1929-2013 and (lower) 1995-2013.
1.The mean b ext level of collocated IMPROVE stations after 1990 is below 50 Mm −1 for the western US, and below 120 Mm −1 for the eastern US (top left).As discussed in Sect.3.1, the low threshold Vis of ∼ 16 km (equivalent to b ext ∼ 100-240 Mm −1 depending on K) recently adopted by most US stations fails to resolve actual b ext variation under this relatively clean environment.Thus many stations are rejected by the threshold filtering.

Figure 3 .
Figure 3. Spatial distribution of (top left) average of the collocated b ext of IMPROVE stations, (top right) Pearson correlation coefficients between collocated pairs of monthly ISD 1/Vis and IMPROVE b ext , (bottom left) slope of monthly b ext against monthly 1/Vis after linear fitting through the origin point using the reduced major-axis linear regression (Ayers, 2001) and (bottom right) Pearson correlation coefficients between collocated pairs of monthly ISD 1/Vis and EMEP PM 2.5 .

Figure 4 .
Figure 4. Spatial distribution of relative trends in 1/Vis (top row), IMPROVE b ext (middle row), and IMPROVE b sp (bottom row) over the US for 1989-2013.Larger colored points with black outline indicate trends with at least 95 % significance, smaller colored points with black outline represent trends with 90-95 % significance, and colored points without outline indicate insignificant trends.Stations with cross and circle symbols are collocated between the ISD and IMPROVE networks over 1989-1996 for composite time series analysis in Fig. 6.

Figure 5 .
Figure 5. Spatial distribution of relative trends in 1/Vis and PM 2.5 over Europe for 2006-2013.Larger colored points with black outline indicate trends with at least 95 % significance, smaller colored points with black outline represent trends with 90-95 % significance, and colored points without outline indicate insignificant trends.Stations with cross and circle symbols are collocated between the ISD and EMEP networks for composite time series analysis in Fig. 6.

Figure 6 .
Figure 6.Composite time series and trends of (top) 1/Vis and b ext for collocated ISD and IMPROVE stations (Fig. 4) over 1989-1996 and (bottom) 1/Vis and PM 2.5 for collocated ISD and EMEP stations (Fig. 5) over 2006-2013.Only stations with significant trends of > 90 % confidence are collocated.The long ticks on the horizontal axis indicate the January of the year.Data gaps represent months with less than 75 % of the total grids.Error bars show the 25th and 75th percentile of all monthly values of collocated stations.

Figure 7 .
Figure 7. Spatial distribution of relative trends in 1/Vis over the US for 1945-1988.Larger colored points with black outline indicate trends with at least 95 % significance, smaller colored points with black outline represent trends with 90-95% significance, and colored points without outline indicate insignificant trends.The red rectangle defines the eastern US region for composite time series analysis in Fig. 8.

Figure 8 .
Figure 8. Composite time series of 1/Vis and SO 2 emission over the eastern US region.The long ticks on the horizontal axis indicateof the year, where all annual values are plotted.Light green dots represent the average monthly 1/Vis of all qualified stations (error bars showing the 25th and 75th percentile) in the defined region.Red dots show the number of grid cells for averaging, and data gaps indicate months with less than 75 % of the total grids for each period.Blue lines and text represent the 1/Vis trends calculated using the monthly anomalies for each period.Trends in parentheses are the 95 % confidence intervals.Black lines are the annual 1/Vis averaged from at least 8 monthly values.SO 2 emissions for the entire US from S. J.Smith et al. (2011a) are in orange.Purple indicates EDGAR SO 2 emissions for the entire US (dashed) and for the defined region (solid) in Fig.7.

Figure 9 .
Figure 9. Spatial distribution of relative trends in 1/Vis over Europe for 1973-2005.Larger colored points with black outline indicate trends with at least 95 % significance, smaller colored points with black outline represent trends with 90-95% significance, and colored points without outline indicate insignificant trends.Red rectangles define the eastern and western Europe regions for composite time series analysis in Fig. 10.

Figure 10 .
Figure 10.Regional time series analysis of 1/Vis and SO 2 emission over western and eastern Europe.The long ticks on the horizontal axis indicate January of the year, where all annual values are plotted.Light green dots represent the average monthly 1/Vis of all qualified stations (error bars showing the 25th and 75th percentile) in the defined region.Red dots show the number of grid cells for averaging, and data gaps indicate months with less than 75 % of the total grids for each period.Blue lines and text represent the 1/Vis trends calculated using the monthly anomalies for each period.Trends in parentheses are the 95 % confidence intervals.Black lines are the annual 1/Vis averaged from at least 8 monthly values.The Smith SO 2 emissions in orange are the total emission of all countries listed in Table2for each region.The EDGAR SO 2 emissions in purple are summed from all pixels inside the defined region (Fig.9).

Figure 11 .
Figure 11.Spatial distribution of relative trends in 1/Vis over eastern Asia for 1973-2013.Larger colored points with black outline indicate trends with at least 95 % significance, smaller colored points with black outline represent trends with 90-95% significance, and colored points without outline indicate insignificant trends.Red rectangles define the northern and southern China regions for composite time series analysis in Fig. 12.

Figure A2 .
Figure A2.Scatter plot of monthly b sp (measured by nephelometers) and b ext (estimated from aerosol speciation data) from all IMPROVE stations with b sp measurements for 56 IMPROVE sites over 1993-2013.The intercept of ∼ 12 Mm −1 corresponds to Reyleigh scattering.

Figure A3 .
Figure A3.SO 2 emission for several major countries.Data are from S. J. Smith et al. (2011a).The top left and top right panels include major countries of western and eastern Europe, respectively.Vertical lines represent division years of the study periods that roughly indicate transition points of emission trend.

Figure A4 .
Figure A4.Spatial distribution of relative trends in 1/Vis over the US for 1929-1944.Larger colored points with black outline indicate trends with at least 95 % significance, smaller colored points with black outline represent trends with 90-95 % significance, and colored points without outline indicate insignificant trends.
Wang et al. (2012)creases significantly with the mean b ext , indicating the tendency for better 1/Vis representativeness in more polluted regions.As previously discussed, at lower b ext more reported Vis are close to the threshold Vis, thus the true 1/Vis tends to be less well resolved.Also, because the Vis data are reported in discrete values, clean stations with a narrow dynamic range of b ext have few reportable Vis to capture the continuous b ext variation.Moreover, the increment of adjacent reportable Vis is relatively coarse in cleaner conditions(WMO, 2008), and atmospheric homogeneity might break down for longer distances.All these factors weaken the ability of Vis to capture b ext variation in clean environments.Wang et al. (2012)found low correlation of 1/Vis with PM 10 over the US and Canada, and similarly attributed this to low aerosol concentrations and higher Vis uncertainty over North America.Thus we apply the 40 Mm −1 threshold of annual 1/Vis to further filter the data as introduced inSect.3.1.4.Without this screening, seven of eight stations with mean 1/Vis < 40 Mm −1 were found to exhibit low correlations (r < 0.25) with collocated b ext .Different thresholds from 10 to 70 Mm −1 were tested, and thresholds above 40 Mm −1 ceased to improve the consistency with the few sites reporting b ext .

Table 1 .
Summary of Pearson correlation coefficients (r) between annual 1/Vis and SO 2 emissions for five regions. ).