Estimated Acute Effects of Ozone on Mortality in a Rural District of Beijing, China, 2005–2013: A Time-Stratified Case-Crossover Study

Studies have shown that ozone (O3) has adverse impacts on human health. In China, O3 levels have continued to increase since 2010. When compared to the large number of studies concerning the health effects of PM2.5 in China, there have been limited explorations of the effects of O3. The Beijing region has one of the highest O3 concentrations in the country, but there appear to be no published studies regarding the health effects of O3 in Beijing. In this study, we applied a time-stratified case-crossover design to explore the effects of O3 on cause-specific mortality for a rural location near Beijing over the period 2005–2013. For year-round effects, we found that for all-causes mortality, with a 10-unit increase in O3 concentration, the odds ratios (ORs) were in the range of 1.009–1.020 for different lag days. The ORs for cardiovascular mortality with a 10-unit increase in O3 concentration were in the range of 1.011–1.017 for different lag days. For warm season effects, the ORs with a 10-unit increase in O3 concentration for all-cause mortality were in the range of 1.025–1.031 for different lag days. The ORs for cardiovascular mortality with a 10-unit increase of O3 concentration were in the range of 1.020–1.024 for different lag days. Our findings fill a knowledge gap that has hitherto existed in studies regarding O3 health impacts, and our results will strengthen the rationale for O3 control in China.


Introduction
In recent years, China has experienced increasing numbers of severe air pollution events. Studies have demonstrated adverse impacts of ozone on human health; among all of the air pollutants, ozone (O 3 ) and PM 2.5 (particles with aerodynamic diameters < 2.5 µm) are believed to have the most significant associations between cause-specific mortality and morbidity, especially cardiorespiratory morbidity [1,2]. Most studies have focused on the adverse effects of PM 2.5 ; based on the results of the research, policies regarding PM 2.5 control have been implemented. In China, strict policies for PM 10 (particles with aerodynamic diameters < 10.0 µm) and PM 2.5 control have been implemented. Due to these strategies, concentrations of atmospheric TSP (total suspended particles) and PM 10 have been decreasing since 1998. Moreover, levels of PM 2.5 continued to decrease during the period 2010-2016 [3][4][5][6][7][8][9]. However, O 3 continued to increase after 2010, especially in the three most developed areas in China: the Pearl River delta, Yangtze River delta and Beijing-Tianjin-Hebei [3][4][5][6][7][8][9]. There were no data published concerning O 3 before 2008. Since PM 2.5 control policies included reducing both reducing both NOx and VOCs (volatile organic compounds), O3 concentration increased with a more favorable ratio of NOx to VOCs [10].
Evidence concerning the adverse effects of O3 on human health, while being relatively limited, is nonetheless very convincing [11][12][13]. These studies have demonstrated links between short-term O3 exposure and adverse health effects, including respiratory illnesses, acute respiratory symptoms, emergency department visits, hospital admissions, and premature mortality. Based on this, the World Health Organization has set O3 standards and provided suggestions [14]. In 2012, the Ministry of Environmental Protection of the People's Republic of China set the following O3 air quality standards (GB 3095-2012): class 1 (remote) areas mandate daily 8-h and 1-h maxima of 100 and 160 μg/m 3 , respectively; in class 2 (urban/industrial and surrounding rural) areas, the corresponding values are 160 and 200 μg/m 3 . It is difficult to evaluate the O3 standards due to the lack of evidence in regard to the health effects of O3 in China. When compared to the large number of studies of PM2.5 health effects in China, there have been limited explorations of O3. There are only five such studies being reported for Mainland China. Moreover, because of the geographic diversity of the country, these studies are inconsistent in seasonal patterns and health outcomes [15][16][17][18]. Moreover, those studies were all performed in central and southern China, even though the Beijing region has one of the highest O3 concentrations in the entire country.
In the present study, we examined the acute effects of O3 on mortality in Miyun County, a suburban district of Beijing. Our purpose was to estimate the impacts of O3 on human mortality in northern China. To achieve this aim, we applied a time-stratified case-crossover design to explore the lag effects of O3 on cause-specific mortality while using a dataset from 2005-2013. We discuss the implications of seasonal modification on the effects of O3.

Study Area
Miyun County is located in the northeast of the Beijing urban area, about 30 km from the center of the city, and is regarded as a background area of Beijing. Air pollution and meteorological data were obtained from the Shangdianzi regional Global Atmosphere Watch (GAW) station (40.39° N, 117.07° E, 293.9 m a.s.l.). This station is located approximately 55 km northeast of the Beijing city area ( Figure 1). More details regarding this station can be found elsewhere [19,20].

Data Collection
The datasets consisted of daily mortality records from 1 January 2005 to 31 December 2013. Mortality data were obtained from the Chinese Center for Disease Control and Prevention. We used three death counts, i.e., cardiovascular diseases (ICD-10 code: I00-I99), respiratory diseases (J00-J99), and all-cause mortality (A00-R99).
Daily PM 2.5 and gaseous pollutants (SO 2 , NOx, and O 3 ) were obtained from Shangdianzi station. Shangdianzi station (SDZ, 40 • 39 N, 117 • 7 E, 293.9 m a.s.l.), is one of the regional Global Atmosphere Watch (GAW) stations in China. The station is located in the northern part of the North China Plain and in the Miyun County of Beijing, approximately 100 km and 55 km northeast of the urban area and the Miyun Township of Beijing, respectively (Figure 1b). Only sparsely populated small villages, and thus insignificant anthropogenic emission sources, lie within 30 km of the site. The station's instrument building is situated on the south slope of a hill surrounded by mountains in every direction, except the southwest. Due to the valley topography, the prevailing winds at SDZ are from the east-northeast and the west-southwest. Polluted air masses from urban areas and satellite towns of Beijing can therefore be easily transported to SDZ by southwesterly winds, while relatively clean air masses arrive from other wind directions.
Daily meteorological variables (mean and maximum temperature, relative humidity, pressure, wind speed, and direction) were recorded by China Meteorological Administration. From 1 January 2005 to 31 December 2013, 3287 days had data recorded. We used both 8-h maximum ground-level O 3 and 8-h maximum moving average O 3 .

Statistical Methods
A time-stratified case-crossover design was used to investigate the associations between O 3 and cause-specific mortality. In this design, "case" days when deaths occurred were compared with control days to assess the effects resulting from differences in exposure to O 3 . Control days were selected to be nearby to case days; in this way, only recent changes in exposure would be compared, and long-term or seasonal variation in exposure could be efficiently eliminated. Conditional logistic regression was used to calculate the odds ratio for cases as compared with controls for a unit increase in O 3 exposure.
We split our time series data into equally-sized, non-overlapping strata and then used a 35-day stratum length with an exclusion period of three days. The exposure of the case day (index day) was compared with the exposure of the control days, which were matched on the same day of the week within the same stratum. Both single pollutant models and multivariate models (containing all pollutants and meteorological factors) were calculated; separate models were used for all natural cause, and cardiovascular mortality. We also controlled for day of the week (DOW), with Sunday as the reference day. The estimates for O 3 were scaled to correspond to a 10 ppb increase.
As temperature may have larger effects on mortality than O 3 and is highly correlated with O 3 , we controlled for temperature by selecting control days within a similar temperature range as the case day.
Odds ratios and 95% confidence intervals (95% CIs) were estimated. The lag structure was an unconstrained distributed lag of the same-day 8-h maximum average ground-level O 3 concentration (lag 0) and ground-level ozone lag 0-3 days before the case-or control-day. To explore seasonal O 3 effects on mortality, we divided the data into separate datasets for the warm season (May-October) and the cold season (November-April).
We considered p < 0.05 as significant in our statistical tests (all were two-sided). We used R software (version 3.4.3) [21] and the "season" package [22,23] to perform the analysis.

Sensitivity Analysis
We also replaced the 8-h maximum O 3 Table 1 shows summary statistics for mortality data, air pollutants, and meteorological factors. The results show considerable variation in O 3 , temperature, relative humidity, and PM 2.5 , i.e., 2.10 to 200.60 µg/m 3 for O 3 , −15.9 to 32.8 • C for daily mean temperature, 8.0% to 98.0% for relative humidity, and 3.63 to 250.13 µg/m 3 for PM 2.5 . There were also ranges of 0.07-54.45 µg/m 3 for SO 2 and 0.70-90.51 µg/m 3 for NOx. There were a total of 21,941 all-cause deaths, 1858 respiratory deaths, and 12,275 cardiovascular deaths during the study period.  Figure 2 shows the correlation matrix of air pollutants and temperature. The correlation coefficient between PM 2.5 and O 3 was 0.26, much smaller than 0.4, indicating that PM 2.5 and O 3 had a weak linear correlation; the two pollutants could be incorporated into a regression model without causing model instability. In addition, there was a strong correlation between temperature and O 3 . The correlation coefficient between NOx and O 3 was −0.3, and that between SO 2 and O 3 was −0.19. Figure 3 presents boxplots for air pollutants during our study period. It could also be seen that the variation of temperature and O 3 concentrations were relatively larger than those of PM 2.5 , NOx, and SO 2 concentrations.    Figure 4 shows that PM 2.5 concentration did not fluctuate much between the four seasons, but there were significant seasonal variations (p < 0.01) in O 3 , NOx, and SO 2 concentrations. The highest concentration of O 3 was in summer, followed by spring and fall. The NOx and SO 2 patterns were opposite, highest in winter, followed by fall and spring.

Seasonal Characteristics of Health Outcomes, Temperature and Air Pollutants
In our study, the health outcomes also had seasonal patterns due to a complex array of causes, among which temperature, PM 2.5 , and O 3 were the most important.  Figure 4 shows that PM2.5 concentration did not fluctuate much between the four seasons, but there were significant seasonal variations (p < 0.01) in O3, NOx, and SO2 concentrations. The highest concentration of O3 was in summer, followed by spring and fall. The NOx and SO2 patterns were opposite, highest in winter, followed by fall and spring.

Seasonal Characteristics of Health Outcomes, Temperature and Air Pollutants
In our study, the health outcomes also had seasonal patterns due to a complex array of causes, among which temperature, PM2.5, and O3 were the most important.  Table 2 shows the lag structure of O3 effects on mortality throughout one year and during the warm season. A total of 3109 case-days and 53,285 control-days were included in the analysis. Significant (p < 0.05) associations between O3 and cause-specific mortalities on different lag days were observed, and odds ratios (ORs) increased with O3 concentration. Most of the estimates were statistically significant in both single-pollutant and multi-pollutants models, although there were differences between lag days. Estimated effects of O3 were moderately reduced but still significant after adjustment for PM2.5 and SO2 in two-pollutant models. This is probably because, in the atmosphere, O3 has a different formation path from those of PM2.5 and SO2 and thus does not typically covary with these pollutants. The weak correlations observed between O3 and these two pollutants indicate that the mortality effect of O3 exposure was at least partially independent. In addition, during the warm season, NO2 was an important confounder in the association between O3 and mortality.   Table 2 shows the lag structure of O 3 effects on mortality throughout one year and during the warm season. A total of 3109 case-days and 53,285 control-days were included in the analysis. Significant (p < 0.05) associations between O 3 and cause-specific mortalities on different lag days were observed, and odds ratios (ORs) increased with O 3 concentration. Most of the estimates were statistically significant in both single-pollutant and multi-pollutants models, although there were differences between lag days. Estimated effects of O 3 were moderately reduced but still significant after adjustment for PM 2.5 and SO 2 in two-pollutant models. This is probably because, in the atmosphere, O 3 has a different formation path from those of PM 2.5 and SO 2 and thus does not typically covary with these pollutants. The weak correlations observed between O 3 and these two pollutants indicate that the mortality effect of O 3 exposure was at least partially independent. In addition, during the warm season, NO 2 was an important confounder in the association between O 3 and mortality. In the single pollutant model (without adjusting for other pollutants) for all-cause mortality for a whole year, O 3 had significant effects on the current day, lag 1 day and lag 2 day; with a 10-unit increase in ambient O 3 concentration, the ORs were 1.021 (95% CI: 1.013-1.029), 1.010 (95% CI: 1.002-1.019), and 1.010 (95% CI: 1.001-1.018), respectively. As for cardiovascular mortality, in the single pollutant model, O 3 had significant effects on the current day, lag 1 day, lag 2 day, and lag 3 day; with a 10-unit increase in ambient O 3 concentration, the ORs were 1.017 (95% CI: 1.007-1.029) 1.013 (95% CI: 1.002-1.024), 1.011 (95% CI: 1.000-1.022), and 1.012 (95% CI: 1.001-1.023), respectively. There was no significant association observed between O 3 and respiratory mortality in our study.

Time-Stratified Case-Crossover
We included other pollutants in the two-pollutant models to estimate O 3 effects. Pearson correlation coefficients between any two pollutants were all < 0.4. Table 2 shows that estimated effects of O 3 were still significant with slight change after adjustment for PM 2.5 , NO 2 , and SO 2 for all-cause mortality. However, for cardiovascular mortality, after adjusting for PM 2.5 , the effects for lag 1 day and lag 2 day became insignificant, while those for the current day and lag 2 day were still significant though slightly decreased. After adjusting for SO 2 , O 3 effects for lag 1 day and lag 3 day became insignificant, while those effects for the current day and lag 2 day were still significant with only slight changes. After adjusting for NO 2 , the associations for different lags became insignificant. Table 2 also shows all of the significant (p < 0.05) effects of O 3 on health outcomes during the warm season for both single O 3 models and in the multiple-pollutants model matched for temperature. After matching for temperature to within one degree, we observed significant associations between O 3 and mortality, and the magnitude of the ORs during the warm season were larger than those for year-round estimates.

Sensitivity Analysis
We replaced 8-h maximum O 3 concentration with the 8-h moving average maximum O 3 concentration in all of the models and replaced daily mean temperature with daily maximum temperature and repeated the analysis. The results did not change significantly. We found that the model was robust to changes in stratum length.

Discussion
We found significant associations between cause-specific mortalities and ambient O 3 concentration increases. When O 3 concentration increased, the ORs of all-cause and cardiovascular mortality increased. For both mortalities, the estimated effects of O 3 were robust with adjustment for other pollutants (PM 2.5 , NO 2 , and SO 2 ), while that for respiratory mortality was not, and this is consistent with previous reports [24]. Larger estimates of O 3 appeared during the warm season for both all-cause and cardiovascular mortality. For year-round effects, the ORs with 10-unit increases of O 3 concentration for all-cause mortality were in the range of 1.009-1.020 for different lag days before controlling for other pollutants; the range changed to 1.009-1.025 after those controls. The ORs with 10-unit increases in O 3 concentration for cardiovascular mortality were in the range of 1.011-1.017 for different lag days before controlling for other pollutants; the range changed to 1.010-1.017 after those controls.
During the warm season, the ORs with 10-unit increases in O 3 concentration for all-cause mortality were in the range of 1.025-1.031 for different lag days before controlling for other pollutants; the range changed to 1.016-1.090 after those controls. The ORs with 10-unit increases of O 3 concentration in cardiovascular mortality were in the range of 1.020-1.024 for different lag days before controlling for other pollutants; the range changed to 1.021-1.078 after those controls.
Positive ORs estimates for O 3 in all-cause and cardiovascular mortalities became slightly larger when NO 2 or SO 2 were included in the model ( Table 2). The reason may be that O 3 and NO 2 or SO 2 are negatively correlated in the atmosphere. This negative correlation had an enhancement effect in the two-pollutant model.
We obtained larger estimates for lag 3 day exposure as compared with those of the current day for both all-cause mortality and cardiovascular mortality during the warm season (Table 2) while controlling for PM 2.5 . These observations are consistent with those from other cities (Shanghai, another metropolitan city in east China) [15] in China. The larger ORs estimated for lag 3 day suggest the accumulation of both acute and less acute health effects over longer periods. The reason why O 3 could affect the cardiovascular system might be that exposure to O 3 can induce inflammatory responses. As in vivo and in vitro experiments have demonstrated, O 3 may mediate a pulmonary inflammatory response; inflammation may subsequently activate hemostatic pathways, impairing vascular function and accelerating atherosclerosis. As cardiovascular mortality accounted for more than 60% of all-cause mortality, a similar result was observed in the estimates of O 3 on all-cause mortality.
The magnitude of O 3 estimated effect was much higher than those reported in the USA or Europe. According to a multisite time-series study in the USA, the pooled estimate for 95 urban communities was a 20 µg/m 3 increase of O 3 associated with approximately a 0.45-0.60% increase in mortality [12,25]. The concentrations of O 3 in Miyun County (annual mean 52-65 µg/m 3 and seasonal mean 35-85 µg/m 3 ) were much higher than those in North American cities (14-38 µg/m 3 ) [26]. In our smoothed plots of O 3 concentration against mortality risk, we observed a steeper slope in the high O 3 concentration range (Supplementary Material, Figure S1). It is worth investigating whether there is any association between long-term O 3 exposure and the acute effects of O 3 .
Although much higher than those reported in the USA and Europe, the magnitude of estimated O 3 effect in our study was similar to those of other cities in China [18,27,28]. In four cities in the Pearl River Delta (Guangdong Province, southern China), the strongest effect was on respiratory mortality, and the RR was 1.46~2.61% with a 10 µg/m 3 increase in O 3 concentration. In Hong Kong [29], Shanghai [15], and Jiangsu Province [30], the estimated effects (RR) for cardiovascular mortality were 1.31~1.75%. In Wuhan [31], the RR range was 1.03~1.64%. In Zhengzhou [32], the RR range was 1.28~1.79%.
Although the magnitudes of the estimated effects were similar, seasonal patterns were very different between our study and others in China. In our study, the most significant association was observed during the warm season, whereas in the other investigations [30,32], in cities of central-eastern (Zhengzhou, Wuhan, and Yangtze River Delta) and southern (Pearl River Delta and Hong Kong) China, there were significant associations between ambient O 3 and mortality in the cold season. In the latter studies, after adjusting for PM 10 , the estimated effects of O 3 on total and cardiovascular mortality increased from September through November, while those for respiratory mortality were only significant from January to August and in December.
The above differences might originate from geographic disparities. When compared to central or southern China, Beijing is a northern city with four distinct seasons as illustrated in Figure 4. Due to those seasons, air pollution in Beijing has a very distinct seasonal pattern. A major impact of seasons is peoples' exposure patterns to O 3 . In southern China, it is cooler and drier in the "cold" season compared to the "warm" season (higher temperatures and humidity), so people are more likely to open windows or stay outside, increasing their frequency of exposure to O 3 . In the warm season, people are more likely to close windows and use air conditioners. In Beijing, because of lower O 3 concentrations and lower temperatures in late fall and winter (with heating from 15 November to 15 March), people are exposed to very low levels of O 3 . Moreover, although some significant associations appeared in the single-pollutant O 3 model, they became insignificant after adjusting for PM 2.5 . Due to the heating supply, a much higher concentration of PM 2.5 appeared in the cold season; this may "cover" the effects of O 3 . The seasonal pattern in our study is similar to those of northern cities in the USA and Europe [33].
As noted above, we did not observe positive associations between O 3 and respiratory mortality. The first potential reason is the climate. Jerrett et al. [13] stated in their large cohort study across the USA that it was quite possible that no positive association between O 3 and respiratory mortality would be found in cool areas (cool in this case meaning a long period of average daily maximum temperature < 25.4 • C). Given Miyun's cooler climate (yearly mean maximum temperature 17-19 • C) relative to most of the USA, this might be one reason why we found no association between O 3 and respiratory mortality.
One strength of the present study is that we chose a rural district of Beijing as the study area. Beijing is one of the three highest O 3 pollution areas (the other two being the Yangtze and Pearl River deltas) in China, but we are unaware of any study regarding O 3 and health effects in Beijing. It is important to compare the results from various geographic regions in China for policymaking. Moreover, most health effect studies of O 3 have been conducted in urban areas, and very limited work has been done in suburban and rural areas. As O 3 and PM 2.5 are the two air pollutants that have the most significant associations with health effects, and people are often exposed to these two pollutants simultaneously, it is important to evaluate a relatively independent effect of O 3 on mortality. A rural district as a study object would be a good choice. The O 3 concentration in the Beijing city area was much lower than that in Miyun County, opposite to the PM 2.5 concentration. From 2005-2013, the mean concentration of O 3 in Miyun was 36.0 µg/m 3 , about 1.59 times higher than that in the city area (22.6 µg/m 3 ), whereas the mean concentration of PM 2.5 in Miyun was 44.0 µg/m 3 , only 60% of that in the city area (73.4 µg/m 3 ). Since people that were exposed to high concentrations of PM 2.5 may be more vulnerable to O 3 pollution, an area with relatively high O 3 and low PM 2.5 would be a good choice for exploring the health effects of ozone.
The second advantage is that Miyun County has a cooler summer compared to the city area. The annual daily mean temperature in Miyun was~1 • C lower than that in the Beijing city area. Moreover, there were only 57 high-temperature (>35 • C) days during 2005-2013 in Miyun, but 102 days in the Beijing city area. Also, there were 626 days with maximum temperature >30 • C in Miyun summers, and 702 d in the Beijing city area. A cooler summer and less developed economy suggest less air conditioner use and more open windows, increasing exposure to ambient O 3 . In our study, we found that the most significant association between O 3 and health outcomes was in summer. Almost all of the RRs in that season were significantly higher than those in the entire year or other seasons for both all-cause and cardiovascular mortalities. Since summer in Miyun has the least toxic PM 2.5 and maximum O 3 concentration, the effects on mortality may be the least confounded.
The third strength of our study is our design for controlling for temperature in the model. As O 3 correlates to sunlight, the most important confounder is temperature. It is well known that temperature plays an important role in the association between O 3 and mortality. In our study, the seasonal pattern in mortality was the reverse of that in O 3 : it was the lowest in summer when the highest concentration of O 3 appeared. This suggested that higher concentrations of O 3 were associated with lower mortality risks; this illustrated how temperature confounded the association between O 3 and mortality: higher temperature is associated with lower death risk. It is obvious that temperature had different effects among seasons; the most significant effects appeared in summer and winter with different lag days.
In our study, although stratification was designed to control for the collinearity between O 3 and temperature, within the strata there still might be a significant correlation between temperature and O 3 in a city with distinct seasons such as Beijing. In order to control more rigidly for the effects of temperature, we selected control days within a similar temperature range as the case days. An advantage of matching using a confounder is that the shape of the association between the confounder and the dependent variable is not important. This means that the association can take any shape (including non-linear forms), and the estimates would be robust.
One limitation of our study is that compared to the Beijing city area, there is a smaller resident population in Miyun County (470 thousand compared to nine million), and the number of respiratory deaths was relatively small. This small number limited our ability to detect a weak pollution association. Another limitation should be noted in interpreting the results of our study. The exposure data were obtained from only one air pollution monitoring station, and the pollutant measurements may differ from individual exposure levels. Therefore, further investigation is needed to explore this issue. In this case, several factors might be taken into account, such as the correlation between individual exposure or average population exposure and monitoring data. Therefore, the used O 3 exposure concentration should be most closely related to the individual exposure level when analyzing the health effects of ambient O 3 exposure.
A third limitation was that we did not include other confounders, such as socio-economic position or individual behavior, which could also influence the health effects of ozone. This was due to data limitation.

Conclusions
We used a time-stratified case-crossover model to account for the effects of O 3 on human health, and the analysis provided lag-specific estimates. We conclude that O 3 had major impacts on cause-specific mortalities in Beijing during the warm season. Our results suggest the need for further investigation of the pathophysiological mechanism of O 3 -associated cardiovascular impact in the northern city. Our work strengthens the evidence for the adverse impact of O 3 on human health, and our data should be helpful in disease prevention and policy development. Also, our findings fill a knowledge gap that has hitherto existed in studies regarding the health impacts of O 3 . The results will strengthen the rationale for O 3 control in China.
Author Contributions: Y.L. analyzing data and writing this manuscript; Y.S. helping to construct the case-crossover model; C.Z. interpretation the health outcome data; Z.M. interpretation the observing data and substantively revised this manuscript.

Conflicts of Interest:
The authors declare no conflict of interest.