Associations between long-term exposure to chemical constituents of fine particulate matter (PM2.5) and mortality in Medicare enrollees in the eastern United States.

BACKGROUND: Several epidemiological studies have reported that long-term exposure to fine particulate matter (PM2.5) is associated with higher mortality. Evidence regarding contributions of PM2.5 constituents is inconclusive.
OBJECTIVES: We assembled a data set of 12.5 million Medicare enrollees (≥ 65 years of age) to determine which PM2.5 constituents are a) associated with mortality controlling for previous-year PM2.5 total mass (main effect); and b) elevated in locations exhibiting stronger associations between previous-year PM2.5 and mortality (effect modification).
METHODS: For 518 PM2.5 monitoring locations (eastern United States, 2000-2006), we calculated monthly mortality rates, monthly long-term (previous 1-year average) PM2.5, and 7-year averages (2000-2006) of major PM2.5 constituents [elemental carbon (EC), organic carbon matter (OCM), sulfate (SO42-), silicon (Si), nitrate (NO3-), and sodium (Na)] and community-level variables. We applied a Bayesian hierarchical model to estimate location-specific mortality rates associated with previous-year PM2.5 (model level 1) and identify constituents that contributed to the spatial variability of mortality, and constituents that modified associations between previous-year PM2.5 and mortality (model level 2), controlling for community-level confounders.
RESULTS: One-standard deviation (SD) increases in 7-year average EC, Si, and NO3- concentrations were associated with 1.3% [95% posterior interval (PI): 0.3, 2.2], 1.4% (95% PI: 0.6, 2.4), and 1.2% (95% PI: 0.4, 2.1) increases in monthly mortality, controlling for previous-year PM2.5. Associations between previous-year PM2.5 and mortality were stronger in combination with 1-SD increases in SO42- and Na.
CONCLUSIONS: Long-term exposures to PM2.5 and several constituents were associated with mortality in the elderly population of the eastern United States. Moreover, some constituents increased the association between long-term exposure to PM2.5 and mortality. These results provide new evidence that chemical composition can partly explain the differential toxicity of PM2.5.

Introduction

Regulatory control of particulate matter (PM) could be dramatically improved with robust quantification of the evidence regarding the toxicity of various constituents of the PM mixture and of their sources [U.S. Environmental Protection Agency (EPA) 2004]. Currently, PM is regulated based on the total mass concentration without regard to its chemical composition, but scientific evidence on which types of constituents are most harmful could result in more effective regulations. However, the knowledge regarding differential toxicities of PM constituents has been identified as a crucial research gap (National Research Council 2004).

Responding to the need for such evidence, for the last decade, data have been accumulated for the constituents of fine particulate matter (≤ 2.5 μm in aerodynamic diameter; PM2.5) nationwide in the United States and have provided opportunities for studying the association between morbidity/mortality risk and PM2.5 constituents. Using these data, numerous epidemiological studies have reported evidence of health effects associated with PM2.5 constituents, both in short-term (a few days previous) and long-term (a few years previous) exposure time frames. Studies focusing on the short-term health effects include those by Bell et al. (2014), Cao et al. (2012), Ito et al. (2011), Kim et al. (2012), Levy et al. (2012), Ostro et al. (2009), Peng et al. (2009), and Zhou et al. (2011). Fewer studies have investigated the long-term health effects of different PM2.5 constituents, including those by Dockery et al. (1993), Ostro et al. (2010), and Pope et al. (1995, 2002). However, studies have not reported consistent findings regarding associations with specific constituents. Such discrepancies may derive from different aspects of the study design (e.g., population, confounding control, time frame, and statistical analysis), and the U.S. EPA (2009) has called for further research.

To fill this research gap, we investigated the differential toxicity of long-term PM2.5 exposure according to its chemical composition, based on a large-scale national database including approximately 12.5 million Medicare enrollees (≥ 65 years of age). Combining several sources of data, we constructed a monthly multi-site time-series data set for 518 PM2.5 monitoring locations in the eastern region of the United States during 2000–2006. The data include, for each monitoring location, monthly mortality rates, monthly values of the average PM2.5 concentration over the previous 12 months, 7-year average concentrations of PM2.5 constituents, and community-level confounding variables on socioeconomic status (SES) and racial composition. Using a Bayesian hierarchical (BH) regression model, we estimated spatially varying (SV) mortality rates associated with previous-year PM2.5 and identified chemical constituents that explained the spatial variability of the mortality rates, controlling for PM2.5 and community-level characteristics.

Methods

Data description. PM2.5 total mass. We obtained daily (24-hr average) concentrations of PM2.5 at 518 monitors in the eastern United States (Figure 1) for 2000–2006 from the U.S. EPA Air Quality System (AQS) database (U.S. EPA 2014). Using the daily PM2.5 data, we calculated monthly long-term exposure to PM2.5 as described in detail by Greven et al. (2011). In brief, for the first day of every month and at each of the 518 monitor locations, we calculated the previous 1-year average of daily PM2.5 concentrations (x) for ith monitor at jth month, with i = 1,…, n and j = 1,…, n. Because not all monitors had valid measurements for the entire study period, the number of monthly PM2.5 values at a given monitoring location (n) varied from 33 to 70.

Figure 1

Map of 518 PM2.5 monitor locations and 174 PM2.5 chemical constituent monitor locations in the eastern region of the United States.

PM2.5 chemical constituents. Although the U.S. EPA measures > 50 PM2.5 chemical constituents, we focused on 6 identified in previous research as contributing substantially to PM2.5 total mass (Bell et al. 2007): elemental carbon (EC), organic carbon matter (OCM), sulfate (SO42–), silicon (Si), nitrate (NO3–), and sodium (Na). We obtained 7-year averages of the 6 chemical constituents at 174 monitors in the eastern region for 2000–2006 from the U.S. EPA AQS database, as described in detail elsewhere (Bell et al. 2007).

One challenge in our study is that the PM2.5 constituents are measured at monitors (174 monitors) that are not collocated with the PM2.5 monitors (518 monitors) (Figure 1). We assumed that the levels of constituents are spatially homogenous within a 6-mile radius, and therefore linked PM2.5 monitors to PM2.5 constituent monitors within 6 miles. Of the 518 PM2.5 monitors, 241 had PM2.5 constituent monitors within 6 miles, and we assigned 7-year averages of each PM2.5 constituent of the closest constituent monitor to each of 241 PM2.5 monitors. For the remaining 277 monitoring locations, we treated the levels of PM2.5 chemical constituents as missing and applied a statistical approach to impute the missing data, as described in “Statistical methods.” For either measured or imputed values, we let z = (z1,…z6)´ denote the 7-year average concentrations of the six chemical constituents for ith PM2.5 monitor location.

Mortality count and total number of people at risk. Mortality counts and the total number of people at risk were obtained at ZIP-code level from billing claims of Medicare enrollees who are fee-for-service Medicare beneficiaries (≥ 65 years of age) (Greven et al. 2011). For each of the 518 PM2.5 monitor locations, we calculated monthly numbers of deaths and people at risk among the Medicare enrollees residing in each ZIP code with a centroid < 6 miles from a PM2.5 monitor location. Depending on the location, 6-mile buffers around the monitors included the centroids of at least 3 and up to 20 different ZIP codes, and the data were aggregated over 3–20 ZIP codes. We let Y and N denote the number of deaths and the people at risk for ith monitor location at jth month. For the whole study period across all 518 locations, the total size of the study population was 12.5 million enrollees, with the total number of deaths equal to 2.2 million approximately residing in 4,974 ZIP codes. For the 241 locations with PM2.5 constituent data available, 1.2 million deaths occurred among 7.5 million enrollees approximately living in 3,425 ZIP codes.

Community-level confounders. We obtained ZIP code–level data on community-level confounding variables including SES and racial composition from the U.S. Census 2000 (U.S. Census Bureau 2000). We averaged values over all ZIP codes with centroids within 6 miles of each PM2.5 monitor and assigned the averaged value to each monitor. We let w = (w1,…w5)´ denote the five community-level confounders: median family income, proportion of people with high-school diploma or equivalent, proportion of residents in urban environment, proportion of white residents, and proportion of black residents.

Statistical methods. We analyzed the linked data using a BH Poisson regression model. The first level, a Poisson regression model with SV random effects, was used to estimate the association between month-to-month variation in mortality rate and month-to-month variation in long-term (previous 1-year average) PM2.5:

Y ~ Poisson(λ), i = 1,…, n and j = 1,…, n,Log(λ) = log(N) + α0 + α1x, [1]

where Y and N are the number of deaths and the size of the population at risk for the ith monitoring location and jth month, x is the previous-year average PM2.5 centered around the location-specific average (i.e., x), and α0 and α1 are the location-specific (SV) random intercepts and slopes. The parameter α0 represents the SV baseline mortality rate when the previous-year average PM2.5 is equal to its location-specific average (i.e., x 0). The parameter αi1 represents the SV association between month-to-month variation in mortality rate and month-to-month variation in previous-year PM2.5.

The second level of the BH model regresses the location-specific 7-year averages of PM2.5 constituents and community-level confounders on the SV intercept and slope, α0 and α1:

α0 = β0 + Σ β* + Σ β6* + β12 + ε0 [2]

α1 = γ0 + Σ γ* + Σ γ6* + ε1, [3]

where z* is the level of the kth chemical constituent and w* is lth community-level confounder at ith location, and ε0 and ε1 are random errors. We centered and scaled all explanatory variables to simplify interpretation and reduce multicollinearity. Note that is included in the SV intercept in Equation 2 to control for total PM2.5 concentration when estimating the effects of constituents on the SV mortality rate (Mostofsky et al. 2012).

To account for potential residual spatial correlation in the second level, we assumed the error terms could be spatially correlated using a standard approach (Gelfand et al. 2003) (see Supplemental Material, “Accounting for residual spatial correlation”). We fit our BH model using a Monte Carlo Markov chain (MCMC) method [see Supplemental Material, “Two-stage estimation and the Markov Chain Monte Carlo (MCMC) algorithm”]. All computations were conducted using R statistical software (R Core Team 2013).

There are four sets of parameters of interest. From Equation 1, we obtained a) the SV (i.e., monitor-specific) baseline mortality rates when the previous-year PM2.5 was equal to its monitor-specific overall average (SV intercepts α0 for each location i, expressed as deaths/month/1,000 persons); and b) the SV association between month-to-month variation in mortality rate and month-to-month variation in previous-year PM2.5 (SV slopes α1 for each location i, expressed as the percentage increase in the mortality rate associated with a 1-μg/m3 increase in previous-year PM2.5). From Equations 2 and 3, we obtained c) the association between the SV intercepts and the monitor-specific 7-year averages of PM2.5 constituents, adjusted by community-level confounders and previous-year PM2.5 (the β coefficients from Equation 2, expressed as the percentage increase in the mortality rate associated with a 1-SD increase in the 7-year average concentration of each constituent), and d) the association between the SV slopes and the monitor-specific 7-year averages of PM2.5 constituents, adjusted by community-level confounders (the γ coefficients from Equation 3, expressed as the percentage increase in the mortality rate ratio for previous year PM2.5 associated with a 1-SD increase in the 7-year average concentration of each constituent).

To find the best fit for Equation 2, we conducted an extensive sensitivity analysis. We considered the following eight models: no explanatory variable, constituents only, community-level confounders only, or both constituents and community-level confounders as explanatory variables, all with and without spatially correlated errors. Among the eight options, we chose the best fit based on the Deviance Information Criteria (DIC) (Spiegelhalter et al. 2002).

There were 277 PM2.5 monitoring locations with missing values for the constituents. Separately for each constituent, we fit a Bayesian spatial Gaussian process (GP) model based on the observed data (i.e., 241 locations) and estimated a spatial correlation using the spBayes R package (Finley et al. 2007) and imputed the missing values based on the posterior predictive sample means for the 277 PM2.5 monitors (see Supplemental Material, “Bayesian spatial Gaussian process (GP) for missing imputation”). Before using the imputed constituent levels in the analysis, we confirmed that the Bayesian spatial GP modeling was appropriate for imputation via a cross-validation (CV) study (see Supplemental Material, “Cross validation study”).

We conducted the analysis for the complete-case data (n = 241 monitoring locations with the data available for both PM2.5 total mass and the chemical constituents) and for the all-sites data (n = 518 monitoring locations using imputed values for the 277 locations without measurements for PM2.5 chemical constituents). Also, we analyzed the data for the entire elderly population (≥ 65 years) and stratified by two age groups (65–74 vs. ≥ 75 years).

Results

Table 1 reports summary statistics for each variable for the complete-case data (n = 241) and for the all-sites data (n = 518). Figure 2A displays maps of 7-year averages of PM2.5 exposure levels (micrograms per cubic meter) and Figure 2B presents maps of 7-year averages of mortality rates (deaths/month/1,000 persons) for 518 monitoring locations. Figure 3 shows maps of 7-year averages of each chemical constituent (micrograms per cubic meter) for the 241 locations with available data. SO42– and NO3– levels seem to exhibit strong spatial correlations; OCM, Si, and Na levels moderate spatial correlations; and EC levels weak spatial correlations, with high values only in a few locations. Estimated spatial correlations obtained from the spatial GP model between pairs of monitors with a distance of about 40 miles are 0.05, 0.20, 0.21, 0.19, 0.21, and 0.20 for EC, OCM, SO42–, Si, NO3–, and Na, respectively. All five community-level confounders are also spatially mapped (see Supplemental Material, Figure S1) over the 518 locations.

Table 1

Summary statistics for each variable (mean ± SD).

Variable	Complete-case data (n = 241)	All-sites data (n = 518)
Population size at risk (n/month)	16901.00 ± 17307.31	14538.92 ± 15453.71
Mortality count (n/month)	76.01 ± 74.04	64.93 ± 66.7
Long-term PM2.5 exposure level (μg/m3)	14.56 ± 1.88	13.7 ± 2.13
PM2.5 chemical constituents (μg/m3)
Elemental carbon (EC)	0.71 ± 0.33	0.68 ± 0.24a
Organic carbon matter (OCM)	4.1 ± 1.06	4.05 ± 0.90a
Sulfate (SO42–)	4.22 ± 0.81	4.14 ± 0.80a
Silicon (Si)	0.09 ± 0.03	0.09 ± 0.03a
Nitrate (NO3–)	1.86 ± 0.86	1.68 ± 0.85a
Sodium (Na)	0.16 ± 0.08	0.17 ± 0.07a
Community-level confounders
Family income ($)	38247.39 ± 11149.87	40305.00 ± 12461.21
Percent high school graduate	0.78 ± 0.06	0.79 ± 0.08
Percent urban	0.93 ± 0.16	0.86 ± 0.23
Percent white	0.64 ± 0.20	0.70 ± 0.20
Percent black	0.25 ± 0.18	0.19 ± 0.17
For population size at risk, mortality count and long-term (previous 1-year average) PM2.5, location-specific monthly values are averaged across locations for the whole study period (2000–2006). For chemical constituents, location-specific 7-year averages are averaged across locations. For community-level confounders, location-specific values are averaged across locations. aNumbers were calculated including the imputed PM2.5 constituent levels for the 277 PM2.5 monitoring locations with missing constituent levels.

Figure 2

(A) Map of 7-year (2000–2006) averages of monthly long-term (previous 1-year average) PM2.5 exposure (μg/m3) for all PM2.5 monitor locations (n = 518). (B) Map of 7-year (2000–2006) averages of monthly mortality rate (deaths/month/1,000 persons) for all PM2.5 monitor locations (n = 518).

Figure 3

Maps of 7-year (2000–2006) averages of each of the six PM2.5 chemical constituents (μg/m3) for 241 monitor locations with available data.

The correlations among the 7-year averages of PM2.5 total mass, the PM2.5 chemical constituents, and community-level confounders are summarized in Supplemental Material, Table S1, for the complete-case data (n = 241). We observed that PM2.5 is correlated positively with OCM (0.43), SO42– (0.61), and the proportion of black residents (0.32) and inversely with Na (–0.41). The highest positive correlations among the constituents were observed between EC and OCM (0.44), SO42– and OCM (0.41), and Si and OCM (0.43). Among the community-level confounders, strong positive correlations were observed as 0.62 between median family income and the proportion of people with high school diploma or equivalent, and 0.50 between the proportions of white residents and high school graduates, whereas the strongest negative correlation was –0.84 between the proportions of white and black residents. Between the constituents and the community-level confounders, the highest correlations were observed for OCM at 0.39 and for Si at 0.34 with the proportion of black residents.

Before the BH regression modeling, we conducted CV studies for our imputation method for the missing constituent levels. The sample correlation coefficients between the observed and predicted values for the test data are 0.64–0.94 for all constituents averaged over 5 CV data sets (see Supplemental Material, Table S2). The root mean square error (RMSE) for prediction for each constituent and the average RMSE over five CV data sets is about half of the sample standard deviation for all constituents (see Supplemental Material, Table S3). Scatter plots for the observed versus predicted data show that the points generally follow the reference line (meaning observed values = predicted values) (see Supplemental Material, Figure S2). Based on the CV study results, we concluded that the Bayesian spatial GP method was appropriate for imputing the missing constituents in our study.

We analyzed the complete-case data (n = 241) and the all-sites data (n = 518), separately. We fit Equations 1, 2, and 3 and Supplemental Material, Equation S1 with eight different options and chose the best fit based on the smallest DIC (see Supplemental Material, Table S4). We obtained the smallest DIC for the model including both chemical constituents and community-level confounders as explanatory variables and with spatially independent errors both for the SV intercept and slope model in complete-case data as well as in the all-sites data.

Figure 4 displays results from the first level of the BH model for complete-case data (left panels) and all-sites data (right panels), respectively. Both sets of data showed similar results. Figure 4A shows the estimated monthly mortality rate when the previous-year PM2.5 is equal to its monitor-specific overall average ranges from 3.37 to 6.15 (deaths/month/1,000 persons) over the study region from all-sites data analysis. Also, Figure 4B shows that the estimated association of mortality rate with a 1-μg/m3 increase in the previous-year PM2.5 is from –1.0 to 4.6 (percent increase in mortality rate).

Figure 4

(A) Maps of the estimates (posterior means) of the SV intercept from the complete-case data analysis (n = 241, left) and the all-sites data analysis (n = 518, right). The values represent the monthly mortality rate (deaths/month/1,000 persons) when previous-year PM2.5 is at location-specific average. (B) Maps of the estimates (posterior means) of the SV slope from the complete-case data analysis (n = 241, left) and the all-sites data analysis (n = 518, right). The values represent the percent increase in the monthly mortality rate associated with a 1-μg/m3 increase in previous-year PM2.5.

Figure 5 reports the results from the second level of the BH model for complete-case data (left-solid bars) and all-sites data (right-dashed bars), respectively. Results were similar between complete-case data and all-sites data, but the all-sites estimates were somewhat smaller and their confidence intervals are narrower. In both analyses, we observed that adjusting for the community-level confounders and PM2.5 total mass, EC, Si, and NO3– were positively associated with mortality rate (the SV intercept, β), whereas SO42– was inversely related to mortality (Figure 5A). Meanwhile, SV slope estimates (γ) indicated that the percentage increase in mortality rate with a 1-unit increase in average previous-year PM2.5 was greater than expected when combined with a 1-SD increase in SO42– and Na (Figure 5B).

Figure 5

(A) Posterior estimates with 95% posterior intervals for the β regression coefficients in the second-level SV intercept model. Solid error bars are for the complete-case data (n = 241) and dashed error bars are for the all-sites data (n = 518). Values correspond to the estimated percentage increase in monthly mortality rate associated with a 1-SD increase in each constituent, adjusted for previous-year average of PM2.5 total mass and for community-level covariates. (B) Posterior estimates with 95% posterior intervals for the γ regression coefficients in the second-level SV slope model. Solid error bars are for the complete-case data (n = 241) and dashed error bars are for the all-sites data (n = 518). Values correspond to the estimated percentage increase in the association between previous-year average of PM2.5 and mortality when combined with a 1-SD increase in each constituent, adjusted for community level covariates.

Finally, we analyzed the data stratified by two age groups (65–74 vs. ≥ 75 years). For the SV intercept model (see Supplemental Material, Figure S3), results for the younger age group were similar to those for the all-age analysis (i.e., main effect estimates for EC, Si, and NO3– were significant), whereas for the SV slope model (see Supplemental Material, Figure S4), both age groups (65–74 and ≥ 75 years) had results similar to those for the all-age analysis (≥ 65 years).

Discussion

We investigated a) whether month-to-month changes in mortality rates were associated with month-to-month changes in the previous-year average exposure to PM2.5, and b) whether 7-year average levels of PM2.5 chemical constituents modified this association.

For the SV intercept, EC, Si, and NO3– were positively associated with mortality rates after adjusting for PM2.5 total mass and the community-level confounders. For EC, our results are consistent with those of a previous cohort study of female public school professionals in California (Ostro et al. 2010). Evidence for the toxicity of Si was also found in other studies. Ostro et al. (2010) reported that long-term exposure to Si was positively associated with pulmonary mortality. Si may serve as a surrogate for toxic constituents found in mineral dust associated with traffic (Reff et al. 2009). A previous study reported that long-term exposure to traffic or traffic-related pollutants is associated with cardiopulmonary mortality (Jerrett et al. 2005). Few studies examined associations for NO3– with mortality in a long-term framework. Ostro et al. (2010) reported that long-term exposure to NO3– was significant for mortality in a single-pollutant model but not in a multipollutant model. Additionally, SO42– was found to be inversely associated with mortality rates, which is inconsistent with previous studies where positive associations were found (Dockery et al. 1993; Ostro et al. 2010; Pope et al. 1995, 2002). Although the observed positive associations were from single-pollutant approaches, our study used a multipollutant analysis that also included adjustment for the PM2.5 total mass and community-level confounders. Therefore, the inverse associations that we observed should be interpreted with caution, because they may be an artifact of multicollinearity resulting from correlations between SO42– and other constituents, PM2.5 total mass, and community-level confounders.

SO42– and Na were significant modifiers of monitor-specific associations between previous-year average PM2.5 and mortality rates. Previous long-term exposure studies for PM2.5 constituents have reported significant positive associations of SO42– with all-cause mortality (Dockery et al. 1993; Pope et al. 2002) or with cardiopulmonary mortality (Ostro et al. 2010; Pope et al. 2002). Na was also a significant modifier that strengthened the association between long-term PM2.5 and mortality. Few studies have estimated associations between Na and health outcomes, with some showing evidence of associations with mortality (Krall et al. 2013) or hospital admission (Zanobetti et al. 2009).

Several possible mechanisms have been proposed in human subject studies linking constituents to biomarkers: systemic inflammation and oxidative stress associated with EC (Neophytou et al. 2013), altered DNA methylation related to Si (Hou et al. 2014), and inflammation related to NO3– and SO42– (Wu et al. 2012). However, because of limited evidences from experimental/toxicological studies, the biological pathway through which short-term exposure to PM2.5 and its components affects health is still an area of active investigation, and the mechanisms for long-term exposure are less understood.

The U.S. EPA measures > 50 different chemical constituents. Analyzing all available constituents would present problems of multiple comparisons. We selected the six constituents that were previously shown to be the largest contributors to PM2.5 total mass and/or co-vary with PM2.5 total mass (Bell et al. 2007). Also, the reliability of a community-level average of PM2.5 constituent exposure varies by constituent. For the six constituents investigated in the present study, the average correlation of monitors in close proximity (< 5 km) ranges from 0.60 to 0.93 and for larger distances (20–50 km) ranges from 0.46 to 0.88 (Bell et al. 2011). The spatial heterogeneity of many other constituents may be larger, limiting the interpretation of community-level exposures. However, we recognize that other constituents have also been found to be associated with human health. In particular, associations of health outcomes with PM2.5 metal constituents that were not included in our analysis, such as aluminum, calcium, chromium, lead, manganese, nickel, titanium, vanadium, and zinc, have been reported in previous studies (Bell et al. 2014; Cavallari et al. 2008; Hsu et al. 2011; Lippmann et al. 2006; Wu et al. 2012).

One limitation of the available air pollution data is that monitors that measure PM2.5 total mass and monitors that measure the PM2.5 chemical constituents are misaligned (Figure 1). We addressed this limitation by assuming that ambient levels of PM2.5 constituents were homogeneous within a 6-mile radius. This spatial homogeneity assumption for air pollutants builds on previous research (Bell et al. 2011) where 6 miles (about 10 km) in radius is a reasonable buffer size for the homogeneity assumption. In Bell et al. (2011), the estimated spatial correlations between pairs of monitors with distances of 5–10 km are 0.67, 0.85, 0.95, 0.62, 0.95, and 0.59 for EC, OCM, SO42–, Si, NO3–, and Na, respectively. However, spatial variability varies by constituent—for example, with more heterogeneity for Si or Na than for SO42– or NO3–—and different buffer sizes may be applied for different constituents when aligning various sources of data.

Another limitation for air pollution data is that monitors that measure PM2.5 total mass are much denser than monitors that measure PM2.5 chemical constituents (Figure 1). When we aligned the two kinds of monitors, missing data occurred for almost half of the PM2.5 monitors. To avoid simply removing the observations with missing values and reducing the sample size to half, we adopted a Bayesian spatial GP modeling and conducted a single value imputation for the missing data separately for each constituent. We compared the results between the complete-case data analysis and the all-sites data analysis with the imputed values. Although imputation did not change our primary conclusions, results based on the imputed data should be interpreted with caution. Specifically, using a single-value imputation does not incorporate uncertainty for prediction, and measurement error can occur for the explanatory variables in regression modeling (Gryparis et al. 2009).

In our study, the PM data are the ambient levels, which we use to approximate the actual human exposure. The ambient level of a given pollutant is not a perfect surrogate of personal exposure to that pollutant, which can induce exposure measurement error into the analysis with variations in error by constituent. In a multipollutant analysis such as the present study, this type of error may induce upward bias in regression coefficient estimates, resulting in anticonservative inference on health effects. However, several authors have shown that this type of bias barely occurs in situations in which the amount of error or the correlations among pollutants in analysis are extremely large (Schwartz and Coull 2003). Therefore, it is unlikely that differences between ambient levels and personal exposures explain the observed associations in our study.

Our analysis is based on multi-site time-series data where long-term exposure was estimated by calculating previous 1-year average of daily exposure values at each temporal point (i.e., first day of each month). However, results may be sensitive to different choices of time frames. Kim et al. (2012) reported that different lag values should be selected for the short-term effects of PM2.5 constituents depending on health outcomes. Shorter or longer time frames than a year could be considered for examining long-term health effects of PM2.5.

Our study focused on the eastern region of the United States, and our findings may not be generalizable to other areas because the characteristics of PM mixtures and populations are quite different across the United States (Bell et al. 2007), and effect modification by the chemical composition of PM2.5 may vary among regions. Also, we focused on the elderly population, which may be more susceptible to effects of exposure than other age groups. In our study, slight differences in results were found between two age groups (65–74 and ≥ 75 years).

To our knowledge, this is the first large-scale study (covering the eastern United States) to investigate the association between long-term exposure to PM2.5 and mortality rate and effect modification by the chemical constituents of PM2.5. Unlike previous studies of PM2.5 constituents, we used a BH regression approach, where PM2.5 constituents were modeled as potential modifiers of the main effect of PM2.5 on health outcomes. Despite limitations, our findings add new evidence regarding the differential toxicity of PM2.5 constituents and their potential influence on the long-term health effects of PM2.5.

Click here for additional data file.