A Bayesian Hierarchical Approach for Relating PM2.5 Exposure to Cardiovascular Mortality in North Carolina

Considerable attention has been given to the relationship between levels of fine particulate matter (particulate matter ≤ 2.5 μm in aerodynamic diameter; PM2.5) in the atmosphere and health effects in human populations. Since the U.S. Environmental Protection Agency began widespread monitoring of PM2.5 levels in 1999, the epidemiologic community has performed numerous observational studies modeling mortality and morbidity responses to PM2.5 levels using Poisson generalized additive models (GAMs). Although these models are useful for relating ambient PM2.5 levels to mortality, they cannot directly measure the strength of the effect of exposure to PM2.5 on mortality. In order to assess this effect, we propose a three-stage Bayesian hierarchical model as an alternative to the classical Poisson GAM. Fitting our model to data collected in seven North Carolina counties from 1999 through 2001, we found that an increase in PM2.5 exposure is linked to increased risk of cardiovascular mortality in the same day and next 2 days. Specifically, a 10-μg/m3 increase in average PM2.5 exposure is associated with a 2.5% increase in the relative risk of current-day cardiovascular mortality, a 4.0% increase in the relative risk of cardiovascular mortality the next day, and an 11.4% increase in the relative risk of cardiovascular mortality 2 days later. Because of the small sample size of our study, only the third effect was found to have > 95% posterior probability of being > 0. In addition, we compared the results obtained from our model to those obtained by applying frequentist (or classical, repeated sampling-based) and Bayesian versions of the classical Poisson GAM to our study population.

VOLUME 112 | NUMBER 13 | September 2004 • Environmental Health Perspectives Research | Article Researchers have found that acute episodes of increased particulate matter (PM) are associated with nonaccidental mortality (Goldberg et al. 2001), total mortality (Katsouyanni et al. 2001;Laden et al. 2000;Mar et al. 2000;Wichmann et al. 2000), cardiovascular deaths (Hoek et al. 2001;Ostro et al. 2000), respiratory deaths (Braga et al. 2001;Hoek et al. 2001), elderly deaths (Katsouyanni et al. 2001), asthma in children and the nonelderly (Lin et al. 2002;Norris et al. 1999;Sheppard et al. 1999), and morbidity (Schwartz 1999;Zanobetti et al. 2000). In all of these studies, the approach taken by the researchers to establish a connection between ambient PM levels and health end points consists of relating measured PM levels on a given day to mortality or morbidity rates on the same or closely following days while adjusting for possible confounding factors such as weather, day of the week, and long-term trends in mortality rates. By far, the most common model used to establish this relationship is the Poisson generalized additive model (GAM). Poisson GAMs are well suited for addressing the question of whether levels of ambient PM in the outdoor environment are associated with health end points, but they may not be the best approach for quantifying the relationship between PM exposure and health end points because direct exposure data cannot be collected for large populations over long periods of time. As a result, Poisson GAMs cannot give direct estimates of increases in the relative risk of morbidity and mortality as a result of exposure to PM.
In attempting to explore the relationship between PM exposure and morbidity or mortality, care should be taken not to assume that the relationship between ambient levels and mortality implies a similar connection between exposure and mortality. It is well documented that ambient levels poorly approximate true exposure (Dockery and Spengler 1981;Lioy et al. 1990;Spengler et al. 1985;Tamura and Ando, unpublished data), and ignoring the discrepancy between exposure and ambient levels in investigations of health effects can lead to biases and underestimation or overestimation of the uncertainty about effects even in simple models (Armstrong et al. 1992). One recent study from the Health Effects Institute (HEI; Cambridge, MA) shows that PM studies are no different: ignoring exposure information can result in biases and misrepresentation of uncertainty when linking PM to health effects (Samet et al. 2000).
In an effort to include exposure information in a model linking levels of PM ≤ 10 µm in aerodynamic diameter (PM 10 ) and mortality, an HEI study (Samet et al. 2000) proposed a multistage Bayesian Poisson regression model, a generalization of the GAM, that includes exposure information. The focus of the HEI study was on Baltimore, Maryland, where daily mortality, PM 10 , and weather variables were collected from 1987 through 1994. Within Baltimore, Samet et al. used the Poisson GAM form to relate PM 10 exposure (instead of ambient levels) to mortality. At the next stage of the hierarchy, the latent exposure is related to ambient PM levels using a linear regression form. To provide information about the coefficients of the regression relating the latent exposure to ambient levels, Samet et al. hypothesized that the same linear form is appropriate for each of five exposure studies and linked the coefficients in each study and the Baltimore population together through another level in the hierarchy.
Although the approach of Samet et al. (2000) takes an important step forward by including exposure information in an epidemiologic model, the method of relating ambient levels to exposure levels could be improved. The assumption that the linear relationship between PM 10 levels and true exposure is similar between the Baltimore population and the populations in the five exposure studies may be unwarranted. In contrast to this HEI approach, an alternative approach for relating ambient pollutant levels to true personal exposure that has gained acceptance more recently is the use of computer exposure simulators. Zidek et al. (2003) presented a general statistical framework for the construction of these simulators. Exposure simulators use activity data and microenvironment pollutant-level data to estimate pollutant exposure levels for individuals. One of the most sophisticated exposure simulators to date for PM is the Stochastic Human Exposure and Dose Simulation (SHEDS-PM) (Burke et al. 2001). For a single individual, SHEDS-PM stochastically simulates a PM level for each of the environments in which the individual spends time. Once SHEDS-PM has defined the microenvironmental levels, the total PM exposure for the individual is estimated by weighting the PM levels in the various environments by the amount of time the individual spends in each of those environments. By examining the estimated PM exposure levels of several individuals created in this manner, the distribution of exposure levels for a population can be characterized.
Building upon the Bayesian model used in the HEI study (Samet et al. 2000), we propose a Bayesian hierarchical model for modeling the relationships among levels of ambient fine PM (particulate matter ≤ 2.5 µm in aerodynamic diameter; PM 2.5 ), average exposure to PM 2.5 , and cardiovascular mortality that incorporates an exposure simulator similar to SHEDS-PM. Unlike most studies, our model allows us to directly quantify the effect of exposure to PM 2.5 on cardiovascular mortality. Bayesian hierarchical modeling is a framework that allows multiple data sources and statistical modeling techniques to be incorporated into a single coherent statistical model (Gelman et al. 1995). In contrast to the Poisson GAM, our model describes the hierarchical nature of the process that connects monitor readings of PM 2.5 to cardiovascular mortality by using a three-level hierarchy. The hierarchy is summarized in Table 1. At the first level, we describe the relationship between PM 2.5 monitors and a continuous surface of ambient PM 2.5 concentrations by allowing for monitor error and considering the spatial properties of PM 2.5 . At the next level, we link average ambient PM 2.5 concentrations at the county level to average population exposure at the county level using an exposure simulator similar to SHEDS-PM. Finally, the third level links average exposure levels to daily cardiovascular mortality counts using the Poisson GAM form. By incorporating all of these levels into a single Bayesian hierarchical model, we are able to estimate the effect of PM 2.5 exposure on cardiovascular mortality and to combine several disparate sources of data in a meaningful way. Although not clearly marked in Table 1, note that the modeled process from level 1 feeds into the modeling technique for level 2, and the modeled process from level 2 feeds into the modeling technique for level 3. By fitting our model using 3 years of data in seven counties in North Carolina (Alamance, Chatham, Durham, Guilford, Johnston, Randolph, and Wake), we found that increased PM 2.5 exposure is related to increased risk of cardiovascular mortality on the same day and the next 2 days. The size of the observed effect is greater than that observed between ambient PM 2.5 levels and cardiovascular mortality, although similar patterns in the effects appear. Each monitor in North Carolina takes readings on a daily, 1-in-3-day, or 1-in-6-day schedule. Daily meteorologic data across North Carolina were obtained from the National Oceanographic and Atmospheric Association's (NOAA) National Climatic Data Center (Asheville NC) via online subscription (NOAA 2003). For each county, the values of the three variables of interest (daily maximum temperature, average wind speed, and relative humidity) were assumed to be equal to the values of those variables reported by the weather station closest to the centroid of the county. We imputed missing meteorologic data (~2% missing overall) by calculating the average value for all other counties with complete data on the same day and substituting that average value for the missing value. Data on human activity patterns were obtained from the Consolidated Human Activities Database (CHAD; U.S. EPA 2003a). This database contains the results of 12 studies in which individual 24-hr details of activities and the environments in which those activities took place were recorded. We restricted our use of the database to records contained in the National Human Activity Pattern Survey (NHAPS) portion of the CHAD and to records of individuals > 20 years of age. Demographic data on the county level were obtained from the U.S. Census Bureau (2003). The population counts for the 2000 census were assumed to be representative of the population counts across the time period studied (1999)(2000)(2001). We used two level-3 summary files in our analysis, P1 and PCT35, which include total population counts by county and the number of individuals > 16 years of age in each county by sex, age, and employment status, respectively.

Materials and Methods
The model that we propose for relating PM 2.5 readings at monitors to daily cardiovascular mortality counts is a three-level hierarchical Bayesian model. The three levels in our model are as follows: a) linking monitor readings to ambient levels over the study region, b) linking ambient levels to exposure levels, and c) linking exposure levels to mortality (Table 1).
Level 1. Central to our model relating PM levels to mortality is that, for any given day, a continuous surface of ambient PM 2.5 levels exists over the study region; this is what would be measured if we obtained an infinite number of monitor readings (spatially dense) without error each day. The first level of our model specifies the spatial distribution of PM 2.5 and relates that distribution to readings taken at monitors on a single day.
We conducted a spatial analysis of PM 2.5 and determined that PM 2.5 exhibits strong spatial correlation over the region of interest [details reported by Calder et al. (2003)]. In order to incorporate this information into a statistical model, we assigned a joint multivariate normal distribution to any set of observations of the PM 2.5 surface. Although we acknowledge that PM 2.5 readings tend to be right-skewed rather than normally distributed, this simplification is not expected to have a strong impact on the overall model fit and simplifies model fitting considerably. On any day t and for any set of sites s(1), … , s(n Ψ ), the distribution of the PM 2.5 surface Ψ t at those is a design matrix of covariates, θ is a parameter vector, and Σ is an n Ψ × n Ψ spatial covariance matrix constructed using information from our exploratory spatial analysis of outdoor PM 2.5 levels. For each site, s(1), … , s(n Ψ ), M t includes a row with elements representing an overall mean, maximum temperature, average wind speed, and two sinusoidal terms that capture seasonal cycles. We considered the corresponding five regression coefficients, θ = (θ 0 , … , θ 4 ), to be unknown, and we minimized prior influence by placing vague N(0, 100) priors on these parameters.
The sites s(1), … , s(n Ψ ) for which the spatial distribution of PM 2.5 is estimated need not be locations with monitors. The matrices M t and Σ are defined for any location in our modeled domain. In fact, in our implementation we modeled the spatial process at several locations that do not have monitors to better characterize the average ambient level over the entire spatial area of each county.
In relating monitor readings to the ambient surface we have defined, we assumed that the PM 2.5 monitors measure the ambient PM 2.5 surface with some error (measurement Article | Relating PM 2.5 exposure to mortality Environmental Health Perspectives • VOLUME 112 | NUMBER 13 | September 2004 error and other random sources of error) at their locations: , where X t (s) is the monitor reading at monitoring site s at time t, Ψ t (s) is the value of the ambient surface at the location of monitoring site s at time t, and σ x 2 is the variance of the measurement error. This construction automatically incorporates the additional uncertainty about the ambient PM 2.5 surface on days when fewer monitors take readings. Days when more monitors take readings (every third or sixth day) will carry more information about the ambient surface than will days when only a subset of daily monitors takes readings, so our uncertainty about the ambient surface will be smaller on these days.
In order to construct a prior distribution for σ x 2 , the variance of the measurement error at the PM 2.5 monitors, precision and accuracy data were downloaded from the AIRS/AQS database (U.S. EPA 2003b). Using these data, we developed an inverse-gamma (649, 1433.405) prior distribution (mean = 2.2, variance = 7.5 × 10 -3 ) for σ x 2 . This prior was developed using a simple conjugate inversegamma/normal model [e.g., Gelman et al. (1995)] with an inverse-gamma (1, 1) prior on σ x 2 before observing data. By creating a continuous surface of ambient PM 2.5 levels, we gained several advantages over the more common "monitor averaging" approach. First, information on the ambient PM 2.5 level on any given day is shared across counties, allowing more accurate characterization of ambient levels in all locations. Second, the interpolation of a continuous ambient surface allows inference about the ambient level in counties that do not contain any PM 2.5 monitors, thereby giving better representation to rural counties. Third, the Bayesian specification of the prior distribution on the ambient level allows natural incorporation of seasonal cycles and meteorologic effects on PM 2.5 levels. Finally, we can characterize the average ambient level in any county on any day by averaging the spatial surface over the county.
Level 2. Level 2 of our model links average ambient PM 2.5 levels in a county to the average exposure level within that county. In this level of the model, we used a deterministic population-level exposure simulator to assist in relating ambient levels to true exposure. Our simulator uses human activity data, information about PM 2.5 levels in indoor environments, and the average ambient concentration on a given day to approximate the exposure level of several individuals in a county on that day. Then, the exposure levels for these individuals are averaged to estimate an average exposure level for all individuals in the county on that day. The population-level exposure simulator used in our model is an adaptation of the SHEDS-PM simulator proposed by Burke et al. (2001). Like SHEDS-PM, our simulator calculates exposure for an individual person using an activity diary and ambient PM 2.5 levels as inputs. This process is repeated for several individuals, and the resulting average exposure is estimated as the mean of the individual exposure levels.
Assuming that the outdoor PM 2.5 level is known and the activity pattern of an individual is known, our simulator calculates individual exposure as follows: where ζ ict is the exposure level for individual i in county c on day t, m ico is the number of minutes the individual spends outdoors, m ice is the number of minutes the individual spends in indoor microenvironment e (residential, office, school, store, vehicle, restaurant, and bar), m ic,smoke is the number of minutes the individual spends with smokers present, m ic,cook is the number of minutes the individual spends cooking, L oct is the ambient PM 2.5 level in county c on day t, L ect is the PM 2.5 level in indoor microenvironment e in county c on day t, L smoke is the addition to the PM 2.5 level in the current microenvironment when smokers are present, L cook is the addition to the PM 2.5 level in the current microenvironment when the individual is cooking, and 1,440 is the number of minutes in a day. When the simulator is implemented in our statistical model, L oct is set equal to the average ambient level in the county at time t, Ψ ct . Additional PM 2.5 measures from smoking and cooking are fixed at 10 µg/m 3 [based on values reported by Burke et al. (2001)] and 5 µg/m 3 [based on findings of Wallace et al. (2003)]. We kept these values constant to simplify computation; a more accurate approach would be to account for the brief shock these activities give to indoor PM 2.5 levels stochastically. Note that this equation makes no distinction between the toxicity of indoor and outdoor particles in our model. The values of L ect for indoor microenvironments are calculated as linear functions of the outdoor level: L ect = a e + b e L oct for e in the set {residential, office, school, store, vehicle, restaurant, bar}.
Values of a e and b e are shown in Table 2. These values were calculated using simplifications of values reported by Burke et al. (2001) for SHEDS-PM.
In each of the counties in which we hope to model the relationship between exposure and cardiovascular mortality, we applied the exposure simulator to several individuals to estimate an average exposure value. In order to apply the simulator, we used activity data that are representative of the true activity patterns in each county in which we modeled the mortality/exposure link. We simulated the activity data by randomly sampling 100 individuals from the county of interest using census demographic information (U.S. Census Bureau 2003) and matching each individual with an activity record from the CHAD (U.S. EPA 2003a). These activity records are drawn from diaries kept across the entire country. Despite possible geographic mismatches, this method of obtaining activity information is usually sufficient for obtaining representative activity information (Özkaynak H, personal communication). To simplify model implementation, a single activity pattern was associated with each individual, and no adjustments were made for different times of the year (i.e., winter vs. summer activity patterns).
To account for possible discrepancy between the simulator predicted value of exposure and true exposure levels, we specified that the average exposure level in a given county is normally distributed around the value predicted by the simulator: where Z ct is the average exposure level in county c at time t, Ψ ct is the average ambient level in county c at time t, ξ(Ψ ct ) is the average exposure level predicted by the simulator in county c at time t as a function of the average ambient level, and σ z 2 is the variance of the error in the simulator. We place a uniform (0, 25) prior on σ z 2 . Although there is not enough information in the data to estimate σ z 2 accurately, allowing it to be random incorporates our uncertainty in the simulator into the model resulting in more accurate uncertainty estimates at the third level.
Level 3. In the third level of the model, we linked exposure directly to mortality using the Poisson GAM form commonly used in studies of the link between PM 2.5 and mortality. Mortality was assumed to be Poisson distributed with a mean that depends on average PM 2.5 exposure in the current and 3 previous days as well as the values of several confounders: where Y ct is the mortality in county c on day t, E c is the expected daily mortality rate in county c (necessary for adjusting the mean level so that the β and η parameters have the same interpretation in all counties), λ ct may be interpreted as a relative risk of death in county c on day t, µ is an overall baseline relative risk of death in the study region over the time period studied, β 0 , … , β 3 are parameters describing the influence of county-level average exposure on mortality rate, f p (C pct ) are transformations of confounding variables, and η 1 , … , η P are parameters describing the influence of confounding variables on mortality. For our data set, confounding variables included a factor variable for the day of the week, a cubic spline transformation of time to account for long-term trends in cardiovascular mortality, a cubic spline transformation of maximum temperature, a cubic spline transformation of relative humidity, and cubic spline transformations of 1-to 3-day lagged values of maximum temperature and relative humidity. The cubic spline transformation of time included 21 evenly spaced knots, and the cubic spline transformations of maximum temperature and relative humidity each included five evenly spaced knots. The model was not assessed for sensitivity to the placement of these knot locations. We reparameterized the confounding variable term into a design matrix (C~) and coefficient vector (γ), and we placed vague N(0, 100) priors on the coefficients. We also placed vague N(0, 100) priors on all of the β-parameters describing the strength of the relationship between PM 2.5 exposure and cardiovascular mortality at different lags as well as on the overall mean relative risk parameter, µ. Summary. Although we have introduced a three-level model, we emphasize that the three levels of the model are all fitted simultaneously as a single coherent statistical model. There are three main advantages to creating a hierarchical Bayesian model for solving such a complex problem. The most important advantage is that uncertainty in parameters is propagated throughout the model. For example, our uncertainty about the true ambient surface (due to errors in the monitors and the necessity of spatial interpolation) carries through to result in a corresponding level of uncertainty about the effect of exposure on cardiovascular mortality. The second important advantage of hierarchical Bayesian modeling is that it is simple to specify large, complex models using simpler statements about conditionally independent parameters. It would be impossible to specify the joint distribution of the thousands of parameters involved in our model if we tried to model the spatial properties of PM 2.5 , the relationship between exposure and ambient levels, and the relationship between exposure and cardiovascular mortality simultaneously. In contrast, the hierarchical approach allows us to specify each level of the model conditionally independent of other levels and to combine the information at the end to obtain a joint distribution of all parameters. The third advantage is that elements of the hierarchy can be substituted without changing the overall form of the model. For instance, we could substitute a different exposure simulator in the second level of the model.

Results
Model fitting was performed using a Markov chain Monte Carlo algorithm (Gelfand and Smith 1990;Geman and Geman 1984;Hastings 1970). The algorithm was implemented with custom C++ software developed using Microsoft Visual Studio (Microsoft Corporation, Redmond, WA). Random number generation was performed using functions from the Numerical Algorithms Group library (NAG, Ltd, Oxford, UK). The algorithm was run for 200,000 iterations, 50,000 of which were discarded as "burn-in" iterations. To reduce the storage space for the samples, the remaining 150,000 samples were thinned by a factor of 50, resulting in a total of 3,000 draws from the joint posterior distribution.
The marginal posterior distributions of several important parameters are summarized in Table 3. For each of the parameters, we include an estimate of the posterior mean (calculated by averaging samples from the posterior distribution) and posterior median (calculated as the median of the sample), a Monte Carlo error for the mean, and a posterior 95% credible interval. The Monte Carlo error for the mean describes how far off our estimate of the true posterior mean is as a result of using a Monte Carlo method for exploring the posterior; it does not describe the uncertainty in the actual parameter. The 95% credible interval does describe the uncertainty in the parameter; it is an equal-tail interval such that the posterior probability that the parameter falls within the interval is 95%. Credible intervals are the Bayesian analogue of the confidence interval but are much easier to interpret because they give direct information about the probability of a parameter falling within certain bounds.
The posterior analysis indicates a positive effect of PM 2.5 exposure on the relative risk of cardiovascular mortality. The posterior marginal expectations of the parameters indicate that a 10-µg/m 3 increase in average PM 2.5 exposure is associated with a 2.5% increase (95% credible interval, -3.9 to 9.6) in the relative risk of current day cardiovascular mortality, a 4.0% increase (-3.3 to 12.2) in the relative risk of cardiovascular mortality the next day, an 11.4% increase (2.8 to 19.8) in the relative risk of cardiovascular mortality 2 days later, and a 1.1% decrease (-7.5 to 5.2) in the relative risk of cardiovascular mortality 3 days later. These rates were calculated by multiplying the β-value corresponding to the effect by 10 and exponentiating. Only the effect on the second day after exposure has a > 95% posterior probability of exceeding zero. Note that the estimates presented are marginal expectations and therefore cannot be added together (e.g., to get an overall risk of cardiovascular mortality from exposure to PM 2.5 ) in a meaningful way. The negative estimate on the third day might be considered an unexpected effect, but it does lend some support to the theory of harvesting (Schwartz 2000). This theory hypothesizes that individuals close to dying of cardiovascular-related causes may die soon after a spike in PM 2.5 exposure, leaving only healthier individuals and consequently decreasing the overall risk of cardiovascular mortality in the total population.
We are unaware of any other study that has attempted to directly estimate the effect of PM 2.5 exposure on mortality, but some related estimates for PM 10 are available from the HEI study (Samet et al. 2000). In that study, a 10-µg/m 3 increase in PM 10 exposure is associated with a 1.4% increase in same-day relative risk of mortality. Although the uncertainty about the HEI estimate is much smaller (probably as the result of a longer time period of study), the point estimate is similar to the one obtained in our analysis.
Although our main goal in this analysis was to demonstrate the effect of PM 2.5 exposure on Article | Relating PM 2.5 exposure to mortality Environmental Health Perspectives • VOLUME 112 | NUMBER 13 | September 2004 cardiovascular mortality, we can also address the effect of changes in the ambient level on the relative risk of cardiovascular mortality. To determine the relationship between ambient levels and relative risk induced by our model, we examined the joint posterior distribution of average ambient levels, Ψ ct , and log relative risk, λ ct , on the same and closely following days. Figure 1 shows smoothed images of the joint distributions combining information across counties. Lines have been added to the figures to illustrate the overall direction of the effect; the line is chosen to minimize the sum of squared distances between samples from the distribution (not shown) and the line. The slope of the line is a summary of the effect of an increase in average ambient level on the log relative risk of cardiovascular mortality, although it is not a parameter in the model. By exponentiating the slope of the line, we obtain an estimate of the proportional increase in relative risk associated with a unit change in ambient level. The lines imply that a 10-µg/m 3 increase in ambient level is associated with a 0.09% increase in the relative risk of cardiovascular mortality on the same day, a 0.2% increase the next day, a 1.0% increase 2 days later, and a 1.4% decrease 3 days later. As with the estimates of effect of exposure on cardiovascular mortality, these estimates are marginal effects and should be interpreted individually; they should not be combined to find an overall effect. These estimates tend to be lower than some comparable estimates reported in the epidemiologic literature. The effect of 2-day mean ambient levels on total mortality has been estimated at 3.3% for chronic obstructive pulmonary disease, 2.1% for ischemic heart disease [both estimates from Schwartz et al. (1996)], and 1.5% for total mortality from natural causes (Klemm et al. 2000), all higher than our largest estimate. This result is not surprising because the inclusion of an exposure link in our model should weaken the direct relationship between ambient levels and mortality. The trend of a weaker association between ambient levels and mortality than between exposure and mortality is similar to the trend reported in the HEI study (Samet et al. 2000).
Although the assessment of the relationship between PM 2.5 and cardiovascular mortality is the main focus of this analysis, estimates of other parameters provide insights into some components of the model. For instance, the estimate of θ 0 , the baseline average ambient PM 2.5 level over all days examined (temperature at 0°F, wind speed at 0 miles/hr), indicates that baseline ambient PM 2.5 levels averaged approximately 9.7 µg/m 3 over the study region from January 1999 through December 2001. The Bayesian model provides an uncertainty estimate for this parameter as well; the baseline ambient PM 2.5 level averaged between 6.1 µg/m 3 and 13.2 µg/m 3 with 95% posterior probability. Some other effects to note are a positive relationship between maximum daily temperature and ambient PM 2.5 levels (an increase of 1°F in maximum temperature is associated with an increase of 0.09 µg/m 3 in daily average ambient PM 2.5 level) and a negative relationship between daily average wind speed and ambient PM 2.5 level (an increase of 1 mile/hr in average daily wind speed is associated with a decrease of 0.08 µg/m 3 in daily average ambient PM 2.5 level). Finally, it is of interest to examine the relationship between average ambient levels and average exposure levels in the counties of interest. The estimates of these values are presented in Table 4 along with some demographic information that was used to choose individuals for the simulator. No correlation between the demographic data and posterior mean exposure levels was observed for the seven counties in our study.
Another interesting parameter estimated in our model is the relative risk of cardiovascular mortality in each county at each time step, λ ct . Examining the relative risk of cardiovascular mortality over the time period studied reveals some interesting patterns. All counties showed similar patterns, so we only present the results for Alamance County (Figure 2). The relative risk of cardiovascular mortality in each county follows a sinusoidal pattern that peaks when the seasonal cycle for PM 2.5 is at its lowest Article | Holloman et al.

1286
VOLUME 112 | NUMBER 13 | September 2004 • Environmental Health Perspectives  point (as implied by the estimates of θ 3 and θ 4 ). The relative risk includes the influence of all of the confounding variables (maximum temperature, relative humidity, long-term cardiovascular mortality trend, and day of the week) in addition to the effect of PM 2.5 exposure on cardiovascular mortality. Therefore, we conclude that overall cardiovascular mortality is significantly affected by numerous factors other than PM 2.5 ; however, our analysis shows that PM 2.5 exposure plays an important role in determining the relative risk of cardiovascular mortality.

Model validation and comparison.
In order to assess whether our model gives reasonable results, we fitted different forms of the model and compared the results obtained in each case. We first considered the effect of eliminating both the spatial interpolation of ambient levels (level 1) and removing the exposure link (level 2 of our model). We call this alternate model 1. We can only fit this model in three of the seven original counties (Durham, Guilford, and Wake) because only these three counties contain at least one daily PM 2.5 monitor. In each county, we first obtained a PM 2.5 reading on each day by averaging the PM 2.5 readings from all monitors that took readings on that day in the county. Prior distributions for all parameters that remain in the model (µ, β-parameters, and γ-parameters) are the same as in our full Bayesian model. We compared the results of this model with results obtained by fitting Poisson GAMs in each of the three counties individually.
The second alternate model that we fitted replaces level 2 of our Bayesian model with a simplified exposure link. Rather than including an exposure simulator, we constructed alternate model 2 by hypothesizing that exposure is equal to the ambient level plus some error [i.e., The remainder of the model is specified exactly as in our original Bayesian model. Summaries of the parameters of most interest, the β-parameters, appear in Table 5, which reports marginal posterior means and 95% credible intervals for the Bayesian models (alternate models 1 and 2) and maximum likelihood estimates with 95% confidence intervals for the classical Poisson GAMs. Note that the parameters for alternate model 2 are interpreted as the effect of a oneunit increase in PM 2.5 exposure on the log relative risk of cardiovascular mortality, whereas the parameters in the other models relate ambient PM 2.5 levels to the log relative risk of cardiovascular mortality.
The results from alternate model 1, the Bayesian model with no spatial interpolation or exposure link, are comparable with the results obtained by fitting the classical Poisson GAM in each of the three counties. This similarity gives evidence that the Bayesian approach produces results similar to those ordinarily obtained using the classical Poisson GAM approach. However, using a Bayesian model allows the incorporation of additional data sources and levels into the hierarchy, so the Bayesian model is more readily expanded.
As expected, the results from alternate model 2 are different from the results obtained from the classical models and alternate model 1; alternate model 2 summarizes the effect of PM 2.5 exposure, not ambient level, on mortality. The results from alternate model 2 are more comparable with those obtained from our full Bayesian model. This similarity indicates that our model is robust to our choice of exposure simulator. However, we do not conclude that the exposure simulator is unnecessary because increased accuracy of simulated exposures will lead to more accurate estimates of the effect of exposure on mortality.

Conclusions
By constructing a hierarchical Bayesian model that divides the process linking PM 2.5 monitor readings and mortality into three intuitive levels, we have shown that elevated PM 2.5 exposure is related to increased risk of cardiovascular mortality in the closely following days. We found that increases in the level of PM 2.5 exposure are most closely related to increased relative risk of cardiovascular mortality 2 days later. In addition, we have demonstrated that the effect of increased levels of exposure on cardiovascular mortality is not equivalent to the effect of increased levels of ambient PM 2.5 on cardiovascular mortality. Our results are similar to those reported in several studies lending additional support to our findings. In addition, we estimate that the association between ambient levels and relative risk of cardiovascular mortality on closely following days is lower than what has been previously reported in the literature.
Despite the sophistication of our model, the second level of the model leaves room for improvement. A deficiency of the second level is the absence of real exposure data. Another limitation of the second level is the simplicity of our exposure simulator; our exposure simulator ignores changes in people's activity patterns over different days of the week and different seasons, uses fixed values to relate indoor and outdoor PM 2.5 values, and may introduce biases in estimation by assuming that the outdoor level is the same for each individual, calculating individual exposures, and then averaging across individuals (Freedman 1999).
Future work on this type of model might focus on addressing the weaknesses in the second level of our model. For example, if real exposure data can be acquired, a data-driven version could be substituted without substantially changing the structure of the model. Similarly, a more complex exposure simulator that takes seasons and the day of the week into account could be substituted to improve the reliability of the results. Nonetheless, the results obtained by incorporating a simple exposure simulator into the model provide valuable insight into the relationship between PM 2.5 exposure and cardiovascular mortality. Table 5. Estimates of the β-parameters (credible intervals) in alternative models.