Introduction

Extreme temperatures are associated with increased daily mortality in many regions of the world (Patz and Khaliq 2002). Because human activity is likely to increase overall global average temperatures, research efforts have focused on the health effects of exposure to high temperatures and heat waves in summer. In the United States, increased mortality during high-temperature days has been extensively investigated. Semenza et al. (1996) reported that a heat wave in Chicago, Illinois, in 1995 was associated with an increase in the death rate among socially isolated people who had no air conditioning. In studies of multiple U.S. cities, similar results were reported (Curriero et al. 2002). In Europe, excess mortality during high-temperature days has also been noted. Le Tertre et al. (2006) also reported an association between the 2003 heat wave in France and increases in all causes of mortality in nine French cities. Stafoggia et al. (2006) explored vulnerability to heat-related mortality in four Italian cities: Bologna, Milan, Rome, and Turin. The populations particularly vulnerable to high summer temperatures were the elderly, women, widows and widowers, those with particular medical conditions, and those in nursing homes and health care facilities.

Air pollution is also associated with increased daily mortality (Pope 2000). A large number of daily mortality time-series analyses have provided sufficiently convincing evidence that nonaccidental mortality, including cardiopulmonary mortality, is associated with ambient particulate matter (PM) exposure in the United States (Ostro et al. 2007), Canada (Burnett et al. 2000), Rome (Forastiere et al. 2007), China (Kan et al. 2007), Korea (Lee et al. 2000), Greece (Katsouyanni et al. 1997), and Chile (Cakmak et al. 2007). The estimated effect is generally in the range of 1.0–8.0% excess deaths per 50-µg/m³ increments in 24-hr average concentrations of particulate matter ≤ 10 µm in aerodynamic diameter (PM₁₀) (Schwartz and Zanobetti 2000).

Although the independent impacts of high temperature and air pollution on daily mortality have been widely explored, few studies have examined the interaction between high temperature and air pollution (Samet et al. 1998). Investigating the effects of the synergy between air pollution and high temperature on mortality, although desirable, is difficult, because a suitable study site is not easily available. The Chinese city of Wuhan, however, provides an opportunity to examine these synergistic effects; it has been called an "oven" city because of its extremely hot summers. Previous studies in Wuhan (He et al. 1993; Qian et al. 2004) have shown high air pollution levels, with concentration ranges wider than those reported in the published literature for other locations. Therefore, we tested the hypothesis that temperature extremes modify the mortality effects of air pollution.

Methods

Study area and population. Wuhan is the capital of Hubei Province, which is located in the middle of the Yangzi River delta, at 29°58´–31°22´ north latitude and 113°41´–115°05´ east longitude. Its population is approximately 7.5 million people, of whom approximately 4.5 million reside in nine urban core districts within an area of 201 km². Wuhan has a subtropical, humid, monsoon climate with a distinct pattern of four seasons. Its average daily temperature in July is 37.2°C, and the maximum daily temperature often exceeds 40°C. The major industries in Wuhan include ferrous smelters, chemical plants, power plants, and machinery plants. The major sources of air pollution in the city are motor vehicles and the burning of coal for domestic cooking, heating, and industrial processes.

Data sources. Mortality data from 1 July 2000 to 30 June 2004 were obtained from the Wuhan Centres for Disease Prevention and Control (WCDC). The government requires that a decedent's family obtain a death certificate from a hospital or a local community clinic to remove the deceased person from the government-controlled household registration. The local WCDC issues two copies of the death certificate according to the certificate from the hospital or the clinic. One copy is submitted to the public safety department to stop the decedent's address registration, and the other copy is used for the cremation.

The WCDC electronically archives all death certificates. In 1992, the WCDC became the first center in China to standardize its system for mortality data collection. The system's requirements are as follows: a) mortality data must be validated four times per year; b) death events collected from the WCDC must conform with those collected from the Wuhan Police Department; c) no data may be missing from any death certificate; d) unclear causes and diagnosis may not constitute > 2% of deaths in urban districts; and e) a correct coding rate of > 98% must be achieved for cause-specific deaths. For deaths that occurred before 1 January 2003, the International Classification of Disease, Ninth Revision [ICD-9; World Health Organization (WHO) 1978] codes were applied; for deaths that occurred after 31 December 2002, ICD-10 (WHO 1993) codes were applied. Total mortality was divided into the following major causes: nonaccidental mortality (ICD-9 codes 1–799; ICD-10 codes A00–R99), cardiovascular diseases (ICD-9 codes 390–459; ICD-10 codes I00–I99), stroke (ICD-9 codes 430–438; ICD-10 codes I60–I69), cardiac diseases (ICD-9 codes 390–398 and 410–429; ICD-10 codes I00–I09 and I20–I52), respiratory diseases (ICD-9 codes 460–519; ICD-10 codes J00–J98), and cardiopulmonary diseases (ICD-9 codes 390–459 and 460–519; ICD-10 codes I00–I99 and ICD-10 J00–J98). The Human Subject Protection Office of the Penn State College of Medicine approved the current study protocol.

Pollution data were collected by the Wuhan Environmental Monitoring Center (WEMC) and certified by the U.S. Environmental Protection Agency. Daily concentrations of PM₁₀, sulfur dioxide, nitrogen dioxide, and ozone (8-hr mean concentrations, 1000–1800 hours) were collected for the study period. The monitoring system strictly followed the quality assurance/quality control procedure set by the State Environmental Protection Administration of China (1992). Briefly, the WEMC conducts regularly scheduled performance audits and precision checks on the air-monitoring equipment. Quarterly performance audits are also conducted to assess data accuracy. PM₁₀ measurements were collected using PM₁₀ beta attenuation mass monitors, (model 7001); SO₂ measurements were collected using an ultraviolet fluorescence SO₂ analyzer (model 4108); NO₂ measurements were collected using a chemiluminescent NO₂ analyzer (model 2108); and O₃ measurements were collected using an ultraviolet photometry O₃ analyzer (model 1008), all from Dasibi Environmental Corporation (Glendale, CA, USA). Meteorologic data were provided by the Wuhan Meteorological Administration.

Statistical methods. We used quasi-likelihood estimation within the context of the generalized additive models (GAMs) to model the natural logarithm of the expected daily death counts as a function of the predictor variables (Hastie and Tibshirani 1990). We examined the effect estimates for each pollutant at 0-, 1-, 2-, 3-, and 4-day lags, and at lag 0–1 day and lag 0–4 day average concentrations prior to the death events. In general, the largest pollutant effects were observed at the lag 0–1, where pollution concentrations were evaluated at the average of the day of death (lag 0) and 1 day before death (lag 1). Therefore, for purposes of this study we focused on the results of the lag 0–1 model. All model analyses were performed using R, version 2.5.0, using the mgcv package, 1.3-24 (The R Foundation for Statistical Computing 2007).

There were two steps in the model building and fit: development of the best base model (without a pollutant) and development of the main model (with a pollutant). The latter was achieved by adding the air pollution variable(s) to the final and best cause-specific base model, assuming a linear relationship between the logarithmic mortality count and the air pollutant concentration. To obtain the best base model, the GAM analyses were performed covering two major areas. First, we controlled for potential confounding of yearly, seasonal, and subseasonal variations and for other time-varying influences on mortality. To begin, we included indicators for days of the week to take into account the change in traffic volume between workdays and weekends. We then regressed the natural logarithm of the daily death counts on a day sequence to adjust for time trends using either natural splines (ns) or penalized splines (ps). Furthermore, visual inspection of the mortality time-series showed two peaks of death counts over the two periods 28 July–3 August 2003 (sum03) and 1 December–31 December 2003 (win03). We added a factor variable for the three periods (sum03, win03, and others) and performed local smoothing by specifying the "by" option for these three periods to control for the extreme peaks of death counts. Second, we controlled for potential confounding of relevant weather variables, which is important during unusually high and low temperatures in Wuhan. We controlled for weather variables using a) indicator variables for extremely hot days, cold days, and humid days; and b) ns or ps for the temperature and humidity, respectively. The extremely hot and cold days were defined as those days on which the highest or lowest daily average temperatures were > 95th percentile or < 5th percentile of the 4 years of data, respectively (Dockery et al. 1992). The 5th and 95th percentiles for temperature were 3.6 and 31.7°C.. Similarly, the extremely humid days were those days with daily average relative humidity > 95th percentile of the 4 years of data. The goal in the previous two steps was to obtain conservative estimates on the subsequent pollution mortality associations.

Taking into account the literature review and the common protocol of the Health Effects Institute's program of the Public Health and Air Pollution in Asia, we used four competing approaches to determine the appropriate degrees of freedom (df) for the time and weather in developing the best base model for each cause-specific mortality model (Curriero et al. 2002; Dominici et al. 2003). These include two ns methods that used the fixed df, the sequential ns method, and the ps method, where the former three ns methods were parametric-based regression splines and always used 2 df and 3 df for the local smoothers for sum03 and win03, respectively. For the two fixed-df models, we considered 6 and 8 df/year for time, 3 and 4 df for temperature, and 3 and 4 df for humidity over the entire 4-year study period. For the sequential method, we started with a reduced model (only days of week, extreme weather indicators, and local smoothing terms). We tried 3–8 df/year for the time and then chose the df that had the smallest sum of the absolute partial autocorrelation values over a 30-day lag period. Next we added temperature to the above model using 2–4 df. We repeated this process for relative humidity after including temperature, time trend, days of week, and extreme weather indicators. We ran the ps model to select the optimal df for overall time trend, local time intervals, temperature, and relative humidity. We initialized the df as 8 df/year for time, 3 df for sum03, 3 df for win03, and 3 df for both temperature and relative humidity. We observed that the local smoothing df remains the same or within 1–2 df differences from the dfs used in the sequential method for various cause-specific mortality. The criteria for selecting the best-fitting model are as follows: a) the absolute value of the partial autocorrelation < 0.1 for all 30-day lags; and b) the smallest sum of the absolute partial autocorrelation values over a 30-day lag period.

To address whether estimated effects are valid and whether they are strongly influenced by different model specifications during the modeling process, we conducted a series of sensitivity analyses in two areas. The first area concerns different smoothing approaches for time, temperature, and humidity. These included a) alternating smoothing order in the sequential method from time, temperature, and humidity to temperature, humidity, and time; b) using fixed df for time, temperature, and humidity (e.g., 6 df for time/year, 3 df for temperature, and 3 df for humidity; and 8 df for time/year, 4 df for temperature, and 4 df for humidity); and c) using the ps approach. The second area concerns model specifications, where the best main models were fitted alternatively by a) adding influenza epidemics; b) adding an indicator for the period of ICD-10 use; c) removing Wuhan, the most industrialized district; d) removing extreme temperature data; e) redefining extreme temperature; and f ) adding the lag climate variable

Last, we redefined the temperature groups using different percentile cutoffs of the temperature ranges (3rd, 7th, 10th, and 15th percentiles) to determine whether the effects observed using the 5th percentiles were significantly changed.

We used several approaches to investigate the validity of the linearity assumption for each air pollutant. First, we replaced the linear term of the pollutant concentrations with a smooth function with 3 df using ns. Both the likelihood ratio test with 2 df (which compares the original main model with the smoothed model) and the visual inspection approach were used to assess whether the smoothed exposure–response curve resembles a straight line. Next, we performed piecewise regression models by allowing different slopes of pollutant concentrations before and after a cutoff point. The cutoff points of PM₁₀ were tested from zero to 150 µg/m³ in 25-µg/m³ increments. The best piecewise regression model was the one in which the cutoff point minimized the generalized cross-validation value. In general, assuming the linearity of air pollution effects on the logarithm of mortality appears to be appropriate.

To investigate the synergetic effects between air pollution and temperature, our main models were built to include additional season indicators and two interaction terms between a linear term of air pollution and an indicator of either extreme high temperature or extreme low temperature (the normal temperature serves as the reference). The effect estimates were expressed using a percentage change in the mean number of daily deaths per 10-µg/m³ increments in 24-hr mean concentrations of a pollutant (8-hr mean concentrations for O₃). The associated upper and lower 95% confidence limits by weather condition were obtained by taking the exponential of the upper and lower 95% confidence limits of the estimated βs. The overall test of the interaction effects between extreme high and low temperatures and air pollution was performed using the likelihood ratio test with 2 df.

Results

The daily mean concentrations of PM₁₀, SO₂, and NO₂ were much lower during high-temperature days than during low-temperature and normal-temperature days (Table 1). The 8-hr mean concentrations of O₃, as expected, were highest during the high-temperature days. There was great variation in the daily average temperature (33.1°C vs. 2.2°C) but small variation in the daily average relative humidity among the three temperature groups.

Table 1.

Table 2.

Table 3.

Table 4.

Table 5.

Figure 1. Cause-specific mortality plots for PM10 stratified by varying percentiles of temperature cutoff points (3, 5, 7, 10, 15) at lag 0-1 day. Values shown are the mean percentage of change in daily mortality per 10-µ/m3 increase in PM10 concentration and 95% Cl.

Table 6.

Table 7.

Table 8.

There were considerable variations in mean daily levels of pollutants (Table 2). The mean daily concentrations of SO₂ and NO₂ generally increased during the study period across the three temperature groups. Despite spatial variations in the daily mean concentrations, which were mainly driven by the highest PM₁₀ and SO₂ concentrations measured at the Wugan station located near a smelter, we found that the distributions of PM₁₀ over distances were fairly homogeneous, as shown by the high Pearson correlation coefficients between measurements from the monitoring stations (0.50–0.97). SO₂ and NO₂ were similarly homogeneously distributed except during the high-temperature days.

We collected information on a total of 89,131 nonaccidental death cases. The daily mean number of nonaccidental deaths was 61, with a maximum of 213 and with a main contribution of cardiopulmonary mortality (daily mean of 35). The majority of individuals died when they were ≥ 65 years of age (71.9%). The mean age of nonaccidental deaths was 69 years, with a range of 0–106 years. Persons ≥ 65 years of age contributed to more than half of the daily deaths for each of the underlying causes of death. The percentage of deaths in the 0–4 year age group was 1.5%. There were only 11 no-death days, all with normal temperature (Table 3). Each variance was greater than the mean, indicating that the mortality data followed the overdispersed Poisson distributions across the three temperature groups, which warrant additional control for weather and temporal trends in the data.

We observed consistent associations between daily mortality and PM₁₀, NO₂, and SO₂ (Qian et al. 2007a, 2007b). In general, using different smoothing approaches did not change the effect estimates significantly, nor did using different model specifications. We also observed a consistent interaction of PM₁₀ with temperature (Table 4). The PM₁₀ effects were strongest on extremely high-temperature days (daily average temperature, 33.1°C), less strong on extremely low-temperature days (2.2°C), and weakest on normal-temperature days (18.0°C). The estimates of the mean percentage of change in daily mortality per 10-µg/m³ increase in PM₁₀ concentrations at the average of lags 0 and 1 day during high temperature were 2.20% [95% confidence interval (CI), 0.74–3.68] for nonaccidental; 3.28% (1.24–5.37) for cardiovascular; 2.35% (–0.03 to 4.78) for stroke; 3.31% (–0.22 to 6.97) for cardiac; 1.15% (–3.54 to 6.07) for respiratory; and 3.02% (1.03–5.04) for cardiopulmonary mortality. Interestingly, we did not observe consistent stronger temperature effects of modification for the majority of outcomes in the elderly (Table 5). One possible explanation might be that the elderly were more likely to stay inside the house on hot days, avoiding exposure to extreme temperature. For the gaseous pollutants, the only interaction observed was that of O₃ on nonaccidental mortality. We found that the estimated PM₁₀ effects using the 5th percentile cutoff were generally similar to the effects estimated using the 3rd percentile (Figure 1). Except for respiratory mortality, we observed that the estimated PM₁₀ effect decreased with increasing percentile on the high-temperature days. Figure 1 also shows that the relationship of daily mortality with temperature is U-shaped, which is consistent with other studies (Gouveia et al. 2003).

The estimated PM₁₀ effects were attenuated in the two pollutant models (Table 6). For example, inclusion of NO₂ in the model substantially reduced the PM₁₀ effect for nonaccidental mortality at normal temperature, whereas the inclusion of SO₂ had less influence. These relationships were also present at low temperatures. Conversely, at high temperatures, the inclusion of either NO₂ or SO₂ had little influence on the association of PM₁₀ with nonaccidental mortality. Although PM₁₀ was correlated with both NO₂ and SO₂ (Table 7), the attenuation of the estimated effects in two-pollutant models might not be due simply to confounding, but rather an indicator of the source-related component of PM responsible for the adverse health effect. The sources and composition of PM₁₀, and hence the toxicity, vary with temperature. Thus, temperature may be serving as an indicator of PM₁₀ composition. The interaction of O₃ on nonaccidental mortality was attenuated but remained significant after controlling for PM₁₀ and SO₂ in the copollutant models (Table 8). Because temperature was positively correlated with O₃ (r = 0.52), part of the interaction between PM₁₀ and high temperature might be due to O₃.

Discussion and Conclusion

We observed that high temperatures enhanced PM₁₀ mortality effects, even though PM₁₀ daily concentrations were lower on the extremely high-temperature days than on the normal-temperature and low-temperature days.

The small number of previous relevant studies reported conflicting results on this interaction. Samet et al. (1998) found no significant evidence that weather variables modified the pollution–mortality relationship. However, Katsouyanni et al. (1993) found a significant effect of the interaction between SO₂ and high temperature on total mortality but no significant interactions between high temperature and either smoke or O₃. We speculate that the following environmental features are related to the significant synergistic effects of PM₁₀ and high temperature in Wuhan. First, the maximum summer temperature often exceeded 40°C and lasted about 2 weeks. Wuhan's special topography causes narrow differences in daily high and low temperatures. Even around midnight in the summer, indoor air temperatures > 32°C are not uncommon. Thus, the city residents were exposed to high temperatures for longer periods than residents of many other cities. Second, few residences in Wuhan were built with energy conservation in mind; a vast amount of radiant energy can easily infiltrate buildings and be absorbed, even when all windows are closed. The temperature inside is commonly comparable to the temperature in the shade outside. In addition, air conditioners have seldom been used because of the high cost of electricity. Third, the most commonly used means for cooling are fans, which can be effective in protecting against heat stress in areas without extremely high temperatures. However, with the temperatures in Wuhan, the use of fans could contribute to heat stress by exacerbating dehydration (Centers for Disease Control and Prevention 1995). Finally, approximately 4.5 million permanent residents plus approximately 1 million transients live in the urban core districts with an area of 201 km². This high population density adds to the urban "heat island" effect, which would make the temperature somewhat higher in the urban core areas than in the suburban areas.

The mechanism underlying the synergistic effects of ambient particle pollution and extremely high temperatures on daily mortality is not yet clear. Some potential explanations have been proposed, especially for the elderly (Easterling et al. 2000). Brunekreef and Holgate (2002) hypothesized that air particles increase the risk of cardiopulmonary mortality through direct and indirect pathophysiologic mechanisms, including pulmonary and systemic inflammation, accelerated atherosclerosis, altered cardiac autonomic function, and increase of inflammatory cytokines in the heart. Many studies have addressed the mechanisms by which high temperature is associated with increased mortality. In animal studies, Keatinge et al. (1986) observed dehydration, increased intracranial and arterial hypertension, endothelial cell damage, and cerebral ischemia during the onset of heat stroke in animals exposed to high temperatures. In a clinical trial study, Gordon et al. (1988) found that exposure to high temperatures increased plasma viscosity and serum cholesterol level. Tsai et al. (2003) suggested that high temperature may help precipitate coronary artery disease and cerebral infarction. Flynn et al. (2005) observed that many of the elderly who died in the heat wave in France during the first 2 weeks of August 2003 were dehydrated, hypernatremic, and hyperkalemic, with evidence of renal failure (Vanhems et al. 2003). The investigators postulated that the most probable causes of death during the heat wave were thromboembolic disease and malignant cardiac arrhythmias as well as heat-induced sepsislike shock (Flynn et al. 2005).

Our study has several limitations. First, both ICD-9 and ICD-10 codes were used. The change in ICD coding might produce misclassification in cause-specific mortality. To address this uncertainty, we examined daily death counts between ICD-9 and ICD-10 mortality data in 2002. We found high concordance rates between the two-coded mortality data, and the maximum change in the estimated pollution mortality effect was 0.09%. These results support our contention that the change in the ICD coding system did not significantly affect the associations identified in this study. Second, there might be other important unknown and unmeasured factors. For example, socioeconomic status can play an important role as an effect modifier. Unfortunately, we do not currently have data on hand to explore the effects of these factors. Third, interpretation of the effects of interaction between O₃ and temperature requires caution, because O₃ data were obtained from only one monitoring station. The limited O₃ data may also restrict our ability to reach any reliable conclusion. Last, measurement errors in exposure are clearly applicable to this study. However, this measurement error generally belongs to the Berkson type and thus is nondifferential in nature, which is likely to cause a bias toward the null and lead to underestimated associations (Armstrong 1998).

In conclusion, we found synergistic effects between PM₁₀ and extremely high temperature on daily mortality in this highly polluted city. Further studies are needed to confirm these findings.