Long-term Fine Particulate Matter Exposure and Nonaccidental and Cause-specific Mortality in a Large National Cohort of Chinese Men.

BACKGROUND: Cohort studies in North America and western Europe have reported increased risk of mortality associated with long-term exposure to fine particles (PM2.5), but to date, no such studies have been reported in China, where higher levels of exposure are experienced.
OBJECTIVES: We estimated the association between long-term exposure to PM2.5 with nonaccidental and cause-specific mortality in a cohort of Chinese men.
METHODS: We conducted a prospective cohort study of 189,793 men 40 y old or older during 1990-91 from 45 areas in China. Annual average PM2.5 levels for the years 1990, 1995, 2000, and 2005 were estimated for each cohort location using a combination of satellite-based estimates, chemical transport model simulations, and ground-level measurements developed for the Global Burden of Disease (GBD) 2013 study. A Cox proportional hazards regression model was used to estimate hazard ratios (HR) for nonaccidental cardiovascular disease (CVD), chronic obstructive pulmonary disease (COPD), and lung-cancer mortality. We also assessed the shape of the concentration-response relationship and compared the risk estimates with those predicted by Integrated Exposure-Response (IER) function, which incorporated estimates of mortality risk from previous cohort studies in western Europe and North America.
RESULTS: The mean level of PM2.5 exposure during 2000-2005 was 43.7 μg/m3 (ranging from 4.2 to 83.8 μg/m3). Mortality HRs (95% CI) per 10-μg/m3 increase in PM2.5 were 1.09 (1.08, 1.09) for nonaccidental causes; 1.09 (1.08, 1.10) for CVD, 1.12 (1.10, 1.13) for COPD; and 1.12 (1.07, 1.14) for lung cancer. The HR estimate from our cohort was consistently higher than IER predictions.
CONCLUSIONS: Long-term exposure to PM2.5 was associated with nonaccidental, CVD, lung cancer, and COPD mortality in China. The IER estimator may underestimate the excess relative risk of cause-specific mortality due to long-term exposure to PM2.5 over the exposure range experienced in China and other low- and middle-income countries. https://doi.org/10.1289/EHP1673.

Introduction

Air pollution and especially particulate matter (PM) with aerodynamic diameter less than () became a focus for both the Chinese government and the Chinese public since the occurrence of heavy smog episodes across most of the areas in China in early 2013 (Chen et al. 2013). In addition, the Global Burden of Disease (GBD) 2010 identified outdoor air pollution as the fourth highest modifiable risk factor for disease burden in China, responsible for more than 1.2 million deaths in 2010. The long-term health effects of ambient air pollution have been studied in many high-income countries with evidence of relatively consistent associations between long-term exposure to ambient air pollution and nonaccidental and cardiovascular (CVD) mortality (Crouse et al. 2012; Gold and Mittleman 2013; Hoek et al. 2013; Beelen et al. 2014a, 2014b).

To date, most studies in China have focused on the effects of short-term exposure and have reported increased rates of daily mortality from natural causes, CVD, and respiratory disease that are consistent with observations in other global regions (Chen et al. 2012, 2013; HEI 2010; Zhou et al. 2015a, 2015b). The effects of long-term exposure have not been adequately investigated. Associations between long-term exposure to ambient air pollution and mortality in China were first reported in the late 1990s. Two aggregate-level (ecologic) studies estimated increased mortality from natural causes, CVD, respiratory disease, and lung cancer associated with long-term average exposure to total suspended particles (TSP) (Jin et al. 1999; Xu et al. 1996). A newer aggregate-level study also reported decreased life-expectancy at birth in northern China associated with coal use and with exposure to TSP (Chen et al. 2013). Recently, several studies have reported increased mortality using individual-level data (Cao et al. 2011; Dong et al. 2012; Wong et al. 2015; Zhang et al. 2011; Zhou et al. 2014). Among them, only a single study assessed the effects of , in this case based on high-resolution satellite data and applied to an elderly cohort in Hong Kong (Wong et al. 2015). Several cohort studies, including some with small sample sizes focused on PM with aerodynamic diameter less than () or TSP (Zhou et al. 2014). The other studies were conducted in a single city (Dong et al. 2012; Zhang et al. 2011). More data from large prospective studies are needed to address the association between and mortality among the Chinese population. Furthermore, almost all the studies so far have focused on urban areas, which may have higher pollution levels. Much of the populations in rural areas in China, however, are also exposed to higher levels of compared to urban and rural areas in most high-income countries (Brauer et al. 2012).

Previous estimates of disease burden due to exposure to in China and in other high-pollution locations have used the Integrated Exposure-Response (IER) function (Cohen et al. 2017). The IER incorporated estimates of mortality risk from cohort studies in western Europe and North America, where exposures were generally lower than exposures observed in China. Information from other types of exposure, such as active and passive smoking, were used to extrapolate risk to these higher levels. Direct evidence of the magnitude of mortality risk over the complete global exposure range is therefore lacking. This study was conducted to both obtain direct evidence of the association between exposure to in outdoor air and mortality in China, and to fill the exposure gap between concentrations observed in Western countries and globally.

In the present study, we estimated exposure levels used for the Global Burden of Disease (GBD) 2013 study and assessed the long-term effects on nonaccidental, CVD, COPD, and lung-cancer mortality in a cohort of 222,000 men. We also assessed the shape of the concentration–response relationship in this setting across the range of concentrations in China, including high levels.

Methods

Study Population

The design of the cohort has been described elsewhere (Chen et al. 2012; Yang et al. 2010; Zhou et al. 2008). Briefly, the participants for this study were recruited from 45 districts/counties across China in 1990–1991, randomly selected from China’s 145 Disease Surveillance Points (DSPs). DSPs are nationally representative and cover about 1% of China’s total population. A DSP site was either a district in an urban area or a county in a rural area. In each district/county, two or three residential units were randomly selected, and all resident men older than 40 y were invited to participate in the study, of whom about 80% accepted. As a result, the original study population consisted of 224,064 men. The baseline survey included a standardized questionnaire and physical measurements administered by trained health workers to obtain information on demographics, lifestyle factors (smoking, alcohol drinking, and household solid-fuel use), personal medical history, height, and weight. Vital status was monitored by local DSP staff through death registries, with regular crosschecks from local residential records kept at the Public Security Department and the Social Welfare Department, supplemented by annual active confirmation by local residential committees. The underlying causes of death were ascertained from official death certificates, supplemented with information from medical records and coded according to International Classification of Diseases version 9 (ICD-9). Cohort follow-up was from 1990/1991 to 2006. Deaths from CVD (390–414, 420–459), including ischemic heart disease (IHD) (410–414) and stroke (431–439), COPD (490–496), and lung cancer (162) were specifically examined in the present study. This study received approval from the ethics committee of the Chinese Center for Disease Control and Prevention. All participants gave oral informed consent.

Air Pollution Exposure

Each participant was assigned an area of residence based on DSP site as determined by the baseline survey. Based on the estimates developed for the GBD 2013 study, the exposure level for all 45 district/counties was estimated for 1990, 1995, 2000, and 2005, as described in detail elsewhere (Brauer et al. 2016).

Briefly, exposures were estimated using a combination of satellite-derived and chemical transport model estimates, calibrated to surface measurements. The satellite-based estimates use aerosol optical depth (AOD) retrievals from satellites to estimate near-surface at ( at the equator) resolution by applying the relationship of to AOD simulated by the GEOS-Chem chemical transport model (van Donkelaar et al. 2015). As satellite retrievals were not available prior to 1998, for 1990 and 1995 estimates were based on the ratio of GEOS-Chem simulations using anthropogenic emissions between 2005 and the respective year of interest, but constant meteorology. These estimates were then combined with simulations from the TM5-FASST (FAst Scenario Screening Tool, a reduced-form version of the TM5 chemical transport model) (van Dingenen and Leitao 2014) based on internally consistent emissions inventories for 1990, 2000, and 2010, and “typical” (year 2001) meteorology, as well as dust and sea salt contributions. TM5-FASST estimates for 1995 and 2005 were generated by fitting a natural cubic spline. Incorporating these separate native estimates (satellite-based and TM5-FASST) helps characterize uncertainty (and how it varies globally) in the final estimate. Finally, we included ground measurements of PM collected from a variety of sources, which were used in a regression calibration approach to combine the mean of the satellite-based estimates and the TM5-FASST simulations with the measurements to produce final estimates at resolution. A total of 4,073 data points from 3,387 unique locations (46% from direct measurements of ) in 79 countries were included in the global calibration. Among these were 394 measurements from China (90 direct measurements of and 304 where was estimated from measurements). The measured concentrations in China (during the 2010–2013 period) ranged from annual average .

Statistical Analysis

We categorized the participants into four quartiles (Q1–Q4) of exposure and presented the baseline characteristics according to quartiles. We calculated age-adjusted nonaccidental mortality rate for participants of these four categories using the total population in the cohort as reference population. We used Cox proportional hazard models for survival, with the time-scale setting as duration from year of recruitment to the year of death or censored at the year of the follow-up in 2006. We assigned to each subject the average ambient concentration at their DSP site, based on the years 2000 and 2005, and the mean concentrations between 2000 and 2005 were used as the exposure level for the main analysis. This period corresponds to the availability of satellite retrievals. The basic model included age at entry as strata variable by 1 y. In model 1, individual risk factors, such as marital status; educational level ( or of education): smoking-related variables; passive smoking; occupational exposure to dust, asbestos dust, coal dust, or stone dust; occupational exposure to coal tar or diesel engine exhaust; alcohol drinking (regular drinkers, nonregular drinkers) and units of alcohol per week (liang/week); body mass index (BMI) (); and diet (consumption of fresh fruit and vegetables). In addition, indoor air pollution defined as household solid-fuel use for heating or cooking (no coal or wood) was included. As more than 70% of the participants were smokers, we made comprehensive adjustments for smoking, and the smoking-related variables included smoking status (never, former, and current), current smoker's years of smoking, current smoker's years of smoking squared, current smoker's cigarettes per day, current smoker's cigarettes per day squared, former smoker's years of smoking, former smoker's years of smoking squared, former smoker's cigarettes per day, and former smoker's cigarettes per day squared. BMI and BMI squared were included as continuous variables. As hypertension might be mediators of the relationship between particle exposure and mortality, we also examined a separate model with adjustment for hypertension. Hypertension was defined if the participant was diagnosed by a health professional or took any antihypertensive drugs in the last two days before the interview or the measured or for those who did not take any antihypertensive drugs. A total of 189,793 men were included in the analysis after exclusion of participants with missing data for the covariates in model 1.

In model 2, we added urban/rural, region, and area-level mean years of education as proxy for contextual socioeconomic status (SES). The term urban referred to a district in a city, and rural meant a county in the baseline recruitment, according to criteria by China National Statistics Bureau. Six geographical regions were included, i.e.,Northeast, North, East, South-Central, Southwest, and Northwest (Figure 1), to adjust for broad-scale spatial variation in mortality not accounted for by the risk factors included in the survival model. The DSP-level mean years of education was obtained from the China Census 2000. For each model, the mortality hazard ratio (HR) associated with a increase in mean concentration during 2000–2005 at each participant’s baseline DSP site was estimated. To assess the effect of individual area-level confounders, we performed incremental adjustment by adding the three area-level covariates, individually and in pairs, to the model adjusted for individual-level covariates only. We also calculated the percent changes in values in comparison with the age-adjusted basic model to indicate the effect of adjustment on the model’s degree of fit.

Figure 1.

Mean concentrations (2000–2005 mean) of in the 45 urban and rural DSP sites in the cohort. A: North; B: Northeast; C: East; D: Southwest; E: South Central; F: Northwest; G: Hong Kong Special Administrative Region; H: Macao Special Administrative Region. The software used to create the figure is QGIS (version 2.18; QGIS Community, open source).

We performed separate analyses stratified by urban/rural, household solid-fuel use, and smoking status to explore potential effect modification. We also combined the three northern regions and three southern regions and stratified the two large geographical regions. We performed sensitivity analyses with exposures to average for the year 1990 (baseline) and average exposure from 1990 to 2005, and exclusion of participants who died in the first 1–3 y.

To examine the shape of the concentration–response relationship between and mortality, we used an integrated modeling framework in which a class of flexible algebraic concentration–response functions was fit to survival models using standard computer software [Shape Constrained Health Impact Function, (SCHIF)]. We constructed the class by defining transformations of concentration as the product of either a linear or log-linear function of concentration multiplied by a logistic weighting function. This approach was recently illustrated with two large cohort studies: the American Cancer Society Cancer Prevention Study II (ACS) cohort and the Canadian Census Health and Environment Cohort (CanCHEC) (Nasari et al. 2016). The model is monotonically nondecreasing and suitable for health-effect assessment. The model incorporates both sampling and model shape uncertainty. In addition, we estimated the HR for four groupings of exposure based on quintiles relative to the first quintile.

We compared our risk estimates for IHD, stroke, COPD, and lung-cancer mortality to those used by GBD2015 (Cohen et al. 2017) in two ways: First, we used the R routine RMA from the package Metafor to conduct a meta-analysis (R Core Team, 2016). We included in the meta-analysis the logarithm of the HR estimates and their corresponding standard errors per μg/m3 of for each of the four causes of death (IHD, stroke, COPD, and lung cancer) reported in Cohen et al (2017); see Table 1. We specified the random effects model with the REML option to estimate the heterogeneity variance. The meta-analysis estimated a summary logarithm of the HR among the cohorts, examining the association between outdoor concentrations of and cause-specific mortality, in addition to a standard error of this estimate. Second, we calculated the relative risk estimate and confidence interval (CI) for the four causes of death based on IER functions as reported in Cohen et al. (2017) between the 5th () and 95th () percentile for the age range of 60–64 y, as this was the average age of follow-up in this cohort. We selected this exposure contrast because it represents the level of exposure experienced by the cohort subjects. We also selected a specific age range for comparison because the IER functions are age specific.A specific exposure contrast is required for this calculation because the IER is nonlinear in concentration, resulting in different relative risk estimates for different values of the exposure contrast. We then calculated the relative risk and corresponding 95% CI for each cause of death separately, using the meta-analysis summary estimate and standard error based on the same exposure contrast as was used for the IER. Here, the relative risk model assumed a linear association between and the logarithm of the HR. All statistical analyses were conducted using the SAS (version 9.3; SAS Institute Inc.) and R software (version 3.2.2; R Core Team).

Table 1

Characteristics of cohort participants at baseline by four categories by quartiles of concentration.

Variables	Total	PM_2.5 categorya
		Q1	Q2	Q3	Q4
Participants (n)	189,793	58,034	40,243	42,512	49,004
Individual-level
Age (years, mean±SD)	54.8±10.7	54.6±10.5	55.2±10.8	54.9±10.9	54.9±10.7
Marital status (% of married)	91.0	94.7	91.1	91.8	85.9
Smoking (%)
Never	25.9	22.9	30.0	21.7	29.8
Former	5.9	6.1	4.0	5.7	7.2
Current	68.2	71.0	66.0	72.5	63.0
Years of smoking for current smokers (mean±SD)	21.8±17.8	22.9±17.6	21.5±18.3	23.1±17.4	19.6±17.6
Years of smoking for former smokers (mean±SD)	1.7±7.6	1.7±7.7	1.2±6.3	1.6±7.5	2.1±8.5
Cigarettes per day in current smokers	13.5±13.2	13.5±12.7	16.5±15.7	13.4±12.8	11.2±11.4
Cigarettes per day in former smokers	1.1±5.1	1.1±5.4	0.8±4.6	1.0±5.1	1.2±5.2
Passive smoking (%)	79.7	83.2	77.1	83.1	74.9
Education, ≤6 years (%)	72.3	70.8	69.6	74.9	70.7
Indoor air pollution (%)
Coal	27.6	29.8	17.8	33.6	27.7
Biomass	35.9	25.0	43.5	45.6	34.1
No exposure	36.5	45.2	38.7	20.8	38.2
Alcohol drinking (%)	32.6	37.4	37.5	26.4	28.4
Amount of alcohol consumed (drinks/week)	10.7±21.9	12.1±22.5	16.2±28.4	8.0±18.9	7.0±15.2
Occupational exposure (%)	4.5	4.4	3.4	4.0	6.1
Fruit and vegetable consumption (%)	31.8	49.4	20.6	22.6	28.1
BMI (kg/m^2, mean±SD)	21.6±2.6	21.3±2.7	21.0±2.4	21.8±2.5	22.1±2.7
Hypertension (%)	26.9	24.9	22.5	31.3	29.2
Area-level factors
Urban (%)	20.2	13.7	19.3	17.2	31.2
Region
Northeast (%)	10.2	31.3	3.2	0.0	0.0
North (%)	8.8	2.1	6.6	24.5	4.6
Northwest (%)	11.3	5.7	0.0	38.4	3.6
Southwest (%)	20.7	26.6	22.9	4.2	26.2
South-central (%)	41.1	20.5	67.4	27.6	55.9
East (%)	7.9	13.8	0.0	5.4	9.7
Mean years of education (mean±SD)b	7.8±1.2	7.4±1.2	7.6±1.0	8.0±1.0	8.3±1.4
Number of deaths
Nonaccidental	50,022	15,466	10,348	13,263	13,779
CVD	18,859	4,797	3,567	5,433	5,062
COPD	11,989	4,231	2,043	2,754	2,961
Lung cancer	2,523	797	452	579	695
Age-adjusted mortality rate (per 100,000 person-years)		2,538	2,355	3,077	2,695

Quartile cutpoints of concentrations were based on mean exposures during 2000–2005 at each baseline cohort site. Minimum, ; 25th percentile, ; 50th percentile, ; 75th percentile, ; maximum, .

Area-level mean years of education.

Results

The geographic location and names of the 45 cohort sites are shown in Figure 1. Descriptive analysis of the participants at baseline is presented in Table 1. The mean age of the participants was (). In addition, 20.2% of the participants lived in urban areas, and 72.3% had an educational level of less than 6 y. Most participants (74.1%) were current or former smokers. The mean BMI was , and 63.5% reported household solid-fuel use (27.6% with coal and 35.9% with biomass). Participants in the highest exposure category tended to be in urban areas and from areas with higher educational levels. In total, 50,022 nonaccidental deaths occurred during the 15 y of follow-up. These deaths included 18,859 deaths from CVD, 11,989 deaths from COPD, and 2,523 deaths from lung cancer. The age-adjusted mortality rate was 2,538, 2,355, 3,077, and 2,695 per 100,000 person-years, respectively, for participants in the four exposure categories.

The average concentration in the 45 cohort sites increased from in 1990 to in 2005, with concentrations in urban areas consistently higher than concentrations in rural areas (Table 2). As shown in Table 3, a increase in based on average exposures during 2000–2005 at each participant’s cohort site was associated with increased risk of mortality from nonaccidental causes, CVD, IHD, stroke, and lung cancer in the model adjusted for individual-level covariates. In the fully adjusted model that includes both individual- and area-level covariates, we observed that a increase in was associated with a higher mortality risk for nonaccidental causes (; 95% CI: 1.08, 1.09), CVD (; 95% CI: 1.08, 1.10), IHD (; 95% CI: 1.06, 1.12), stroke (; 95% CI: 1.13, 1.16), COPD (; 95% CI: 1.10, 1.13), and lung cancer (; 95% CI: 1.07, 1.14). Further adjustment for hypertension only slightly reduced the effect size without major changes (Table S1). In the incremental adjustment, the effect estimates increased with inclusion of the three area-level covariates (Table S2). We observed generally larger reductions in the values when the contextual variables were added to the model in comparison with the individual risk factors, suggesting that there remains considerable spatial variation in mortality that is not accounted for by our individual risk factors (Table S2). We therefore included area/contextual variables to account for some of this residual variation. Furthermore, we included both region and urban/rural specification as stratification variables to adjust for potential confounding effects on the -mortality association of these two area factors. In this model specification, exposures are compared among subjects in the same region-urban/rural strata, and not between strata. We included these area-level factors as stratification variables and not as predictors in the survival model, as their baseline mortality rates vary substantially between strata (Table S3). We note that a single –mortality effect estimate is obtained from this model, summarizing associations among all strata.

Table 2

Estimated exposure level of () in 45 cohort sites from 1990 to 2005.

PM_2.5 category	1990	1995	2000	2005
Total (n=45)
Mean±SD	36.4±16.3	38.7±17.6	40.7±18.6	46.4±20.2
Min	3.5	3.8	4.0	4.4
Q1	23.8	25.0	26.0	31.4
Q2	31.3	33.6	35.2	42.5
Q3	49.5	52.9	55.7	63.4
Max	81.9	81.5	81.5	89.8
Urban (n=23)
Mean±SD	41.4±17.4	44.3±19.0	46.7±20.2	53.7±22.2
Min	3.5	3.8	4.0	4.4
Q1	27.5	27.4	28.5	36.5
Q2	36.4	38.6	40.7	51.6
Q3	59.0	65.4	70.6	79.6
Max	69.5	74.0	77.7	89.8
Rural (n=22)
Mean±SD	35.1±15.8	37.3±16.9	39.2±17.8	44.5±19.2
Min	13.2	13.7	14.2	16.4
Q1	22.6	24.0	25.2	29.1
Q2	30.2	32.5	34.2	41.1
Q3	49.5	52.9	55.7	60.7
Max	81.9	81.5	81.5	79.3

Table 3

Hazard ratios (HRs) (and 95% CI) of mortality associated with increase in levels.

Cause of death	No. of deaths	Age-adjusted basic model	Adjusted model (1)a	Adjusted model (2)b
Nonaccidental	50,022	1.03 (1.02–1.03)	1.04 (1.03–1.04)	1.09 (1.08–1.09)
Cardiovascular	18,859	1.07 (1.06–1.07)	1.06 (1.05–1.07)	1.09 (1.08–1.10)
IHD	3,752	1.14 (1.12–1.15)	1.11 (1.09–1.13)	1.09 (1.06–1.12)
Stroke	11,301	1.08 (1.07–1.09)	1.07 (1.06–1.08)	1.14 (1.13–1.16)
COPD	11,989	0.98 (0.97–0.98)	1.00 (0.99–1.01)	1.12 (1.10–1.13)
Lung cancer	2,523	1.03 (1.01–1.05)	1.04 (1.02–1.07)	1.12 (1.07–1.14)

Note: levels were based on mean exposures during 2000–2005 at each participant’s baseline cohort site.

Model 1: Age-adjusted basic covariates including marital status, educational level, smoking status, years of smoking, cigarettes per day, passive smoking, occupational exposure, alcohol drinking, units of alcohol per week, body mass index (BMI), consumption of fresh fruits and vegetables, and household solid-fuel use.

Model 2: covariates including urban/rural, region, and mean years of education.

In stratified analysis (Figure 2), the association was consistently significant for nonaccidental, CVD, IHD, stroke, and COPD mortality in all strata according to urban/rural, north/south, smoking status, or exposure to indoor air pollution. The effect size of nonaccidental causes, CVD, IHD, and COPD in the north was generally higher than the south. The association between and lung cancer was significantly positive for southern regions and negative for the northern regions. Table 4 presents the results of sensitivity analysis. HRs were slightly increased when we use the baseline level of exposure (1990) and mean levels of between 1990 and 2005. When we excluded participants who died 1–3 y after enrollment, the HRs slightly decreased but still remained significant.

Figure 2.

Hazard ratios (HRs) (and 95% CI) of mortality [(A): nonaccidental deaths; (B): CVD; (C): IHD; (D): Stroke; (E): COPD; (F): Lung cancer] associated with a increase in levels, stratified by urban/rural, region, smoking status, and indoor air pollution. distribution was based on mean concentrations during 2000–2005 for each participant’s cohort site at baseline. Note: Stratified models were adjusted for age, individual-level covariates including marital status, educational level, smoking status, years of smoking, cigarettes per day, passive smoking, occupational exposure, alcohol drinking, units of alcohol per week, body mass index (BMI), consumption of fresh fruits and vegetables, and indoor air pollution and area-level covariates, including urban/rural, region, and mean years of education. y-Axis scale was 1.5 for CVD and IHD and 1.3 for other causes of death.

Table 4

Hazard ratio (HR) (95% CI) per increase of for sensitivity analyses for exposure at baseline and for different inclusion criteria.

Cause of death	Exposure in 1990(n=189,793)	Average exposure between1990–2005(n=189,793)	Excluding deaths within 3 years, exposure 2000–2005(n=179,089)a
Nonaccidental	1.12 (1.11, 1.12)	1.10 (1.10, 1.11)	1.07 (1.06, 1.08)
Cardiovascular	1.10 (1.08, 1.11)	1.11 (1.10, 1.12)	1.08 (1.07, 1.09)
IHD	1.10 (1.07, 1.14)	1.10 (1.07, 1.13)	1.07 (1.04, 1.10)
Stroke	1.15 (1.14, 1.18)	1.18 (1.16, 1.19)	1.14 (1.12, 1.16)
COPD	1.19 (1.17, 1.21)	1.13 (1.12, 1.16)	1.12 (1.10, 1.13)
Lung cancer	1.10 (1.06, 1.14)	1.10 (1.06, 1.13)	1.12 (1.07, 1.14)

Excluded numbers of deaths: 10,083 for nonaccidental causes; 3,489 for cardiovascular diseases; 689 for IHD; 2,158 for stroke; 2,708 for COPD; and 454 for lung cancer.

The HR estimates by quintile of concentrations suggested nonlinear relationships between and mortality from all disease outcomes (Figure S1). Compared with the first quintile, the second quintile had a lower HR during the low concentrations for nonaccidental causes, CVD, stroke, and COPD. The HR estimate associated with the fourth quintile tended to have the highest value, and a lower estimate was observed for the fifth quintile. These plots suggest a sigmodal pattern with little change in the HR for lower concentrations and a diminishing change for the highest exposures.

The SCHIF concentration–response plots for the six causes of deaths display a pattern similar to that of the quintile-based analyses (Figure 3). The estimated curves are monotonic in concentration and display less curvature than the quintile-based analysis, a design feature of the SCHIF. In general, there is a greater change in relative risk for higher concentration in comparison with lower values. We note there is little change in the HR for concentrations below and little uncertainty in these estimates. This finding is due to the slightly negative association in this range, as evidenced by the quintile analyses and the requirement of the SCHIF to yield HR estimates greater than unity and increase with increasing concentration. We also note that the HR estimates increase more rapidly for stroke than for the other causes of death. The much wider CIs for lung cancer suggested larger uncertainty, likely due to the relatively smaller number of lung-cancer deaths in this cohort.

Figure 3.

Concentration–response curves and 95% confidence intervals (CIs) for the relationship between and mortality based on the model adjusted for individual-level and area-level covariates: (A) Nonaccidental causes; (B) CVD; (C) IHD; (D) stroke; (E) COPD; (F) lung cancer. concentration was based on mean exposures during 2000–2005 for each cohort site at baseline.

Figure 4 shows the comparison of HR estimates for the four causes of deaths among the present study, IER predictions, and meta-analysis, based on the same exposure contrast between the fifth () and 95th () percentile as was used for IER (Cohen et al, 2017) for the age range of 60–64 y. We found the HR estimate from our cohort based on exposure contrast was consistently higher than the corresponding IER estimates for all the four disease outcomes and very similar to the meta-analysis estimates for COPD and lung cancer. In comparison with the meta-analysis estimates, our HR was somewhat lower for IHD and slightly higher for stroke. Both this cohort and the meta-analysis estimates were much larger than the IER predictions for this particular exposure contrast experienced by our cohort subjects.

Figure 4.

Comparison of hazard ratio (HR) estimates and confidence interval (CI) for the four causes of deaths between the present study, IER predictions and meta-analysis for the outdoor air pollution (OAP) study, based on the same exposure contrast between the fifth () and 95th () percentile (as was used for IER) for the age range of 60–64 y. IER predictions: estimates based on IER functions as reported in Cohen et al, (2017) between the 5th () and 95th () percentile for the age range of 60–64 years. Meta-analysis: HRs (95% CI) presented was calculated from the meta-analysis summary estimate and standard error based on the same exposure contrast as was used for the IER (), using the cohorts [a list of the cohort was in Table 1 in Cohen et al. (2017)] examining the association between outdoor concentrations of and cause-specific mortality.

Discussion

In this Chinese cohort, we observed an elevated risk for mortality in association with long-term exposure to . An increase of of was associated with 9%, 9%, 9%, 14%, 12%, and 12% increases in the risk of mortality from nonaccidental causes, CVD, IHD, stroke, COPD, and lung cancer, respectively, in the fully adjusted model. The concentration response curve suggested a nonlinear relationship between and mortality in China, where the exposure is higher than exposure in developed countries.

To estimate exposure levels, we used a method combining satellite-based estimates, chemical transport models, and ground-level measurements developed for the GBD 2013 study. The use of a continuous spatial surface of estimates allows for assignment of exposure to both urban and rural participants using a common approach, reducing the potential for sampling bias in which only populations in the vicinity of ground measurements are included. The exposure estimates are also internally consistent throughout the period of follow-up. To provide full spatial and temporal coverage, these estimates include chemical transport model simulations and satellite retrievals, but the estimates also incorporate the substantial number of ground measurements from China that have recently become available.

We focus on the 2000–2005 average concentrations because this period included satellite retrievals, whereas the 1990 and 1995 estimates were based on scaling from the year 2000 spatial patterns with chemical transport model simulations, as described previously (Brauer et al. 2016). Although this back-casting does incorporate time-varying spatial variation from the chemical transport model simulations, these are at a lower spatial resolution than the satellite retrievals, and thus there is a likelihood of greater exposure misclassification for the 1990 and 1995 estimates.

We compared our effect estimates of the associations between and mortality with previous large cohort studies, including ACS (Krewski et al. 2009), Six City Study (SCS) (Lepeule et al. 2002), California Teachers Study (CTS) (Lipsett et al. 2011), Dutch Study of Diet and Cancer (DSDC) (Beelen et al. 2014a; Beelen et al. 2014b), Canadian Census Health and Environment Cohort (CanCHEC) (Crouse et al. 2012), and a Hong Kong cohort (Wong et al. 2015) (Figure S2). Previous studies have generally reported HRs in excess of one, although many lacked statistical significance. The HR estimate of nonaccidental and CVD mortality in our study is generally lower than the estimates reported in most of the other major cohort studies, including SCS, CanCHEC, and the Hong Kong Cohort, but higher than the ACS estimate. Our effect estimate of IHD is consistently lower than the estimates in listed cohort studies conducted in high-income countries (Crouse et al. 2012; Krewski et al. 2009; Lepeule et al. 2012; Lipsett et al. 2011; Villeneuve et al. 2015). In line with previous reviews (Hamra et al. 2014; IARC 2016), our study showed a significant effect of exposure on lung-cancer mortality. The observed increased risk of COPD mortality strengthens the current evidence and is consistent with the most recent results of the ACS and other cohort studies (Cohen et al. 2017; Schikowski et al. 2014). The differences of effect size may be explained by various characteristics of study contexts, such as the feature of cohort population (e.g., we included only men in the current study, and CTS included only women), exposure estimate methods, pollution levels, sources and components of air pollutants, and population susceptibility to air pollution (Hedley et al. 2002).

We found a generally higher effect size of nonaccidental causes, CVD, IHD, and COPD in the north than in the south. Different sources and mixtures of from different areas may have differential health impacts. in some areas may be dominated by industrial emissions and power generation, whereas diverse sources’ mixtures contribute to in the other areas across China. Chen et al. (2013), noted above, reported regional differences in life-expectancy associated with long-term exposure to air pollution, and lower life expectancy north of the Huai River in 2000 associated with a national program instituted in the 1950s that provided free coal to northern areas and with elevated TSP levels. The health differences between regions might also be caused by factors other than air pollution, such as health care and access to medical services, which could not be fully adjusted for due to data availability. Further studies are needed to address such differences, especially in the northern provinces, and although we included personal educational level as a proxy for individual-level SES and used mean years of education as area-level SES in the models, accounting for the potential effects of individual and contextual variables related to on mortality remains an important issue for future studies.

The association between air pollution and lung-cancer mortality was positive in the south but negative in the north, which may be due to residual confounding by factors such as occupational exposure, for which we had limited data and therefore was crudely characterized in our study. We also noted that the HR in smokers was higher than among participants who had never smoked, which is not in accordance with the European Study of Cohorts for Air Pollution Effects (ESCAPE), which reported elevated HRs for lung-cancer incidence in association with in individuals who had never smoked (Raaschou-Nielsen et al. 2013).

Trends in air pollution exposure are important in modeling the long-term health effects. Although air pollution levels may vary over time, the changes are usually gradual, and spatial patterns are usually preserved, affecting the regional pollution in a similar way in different time periods (Lipsett et al. 2011), as suggested by our analysis. However, we had limited ability to assess trends and spatial patterns prior to 2000, as estimated exposures for 1990 and 1995 were based on spatial patterns from satellite retrievals in the year 2000, scaled by chemical transport model simulations at a coarser spatial resolution. Accordingly, we used average exposure between 2000 and 2005 as the exposure level for the main analysis, with baseline exposure and average exposure between the years 1990 and 2005 for sensitivity analyses. The effect estimates were similar in all these analyses, which may reflect the similar spatial patterns in pollution over time during our study period.

The accuracy of the specific mortality outcomes is also essential in this study. The mortality outcome in our cohort was ascertained by local DSP staff members through the death surveillance system and confirmed and coded by staff from central DSP office in Beijing, based on official death certificates, medical records, or verbal autopsies when no medical records are available. The good quality and reliability of the DSP have been proved (Liu et al. 2016).

Estimating the exposure–response relationships is critical to assessing the effect of specific regulatory actions to improve air quality on population mortality rates (Cohen et al. 2004; GBD 2013 Risk Factor Collaborators et al. 2015). Previous cohort studies often fit natural, restricted, or smoothing splines with a few degrees of freedom or use simple algebraic nonlinear functions or threshold function to assess the shape of the association between ambient air-pollutant concentrations and mortality (Crouse et al. 2012; Jerrett et al. 2009; Krewski et al. 2009). We used a more complex method to identify the shape of the association between air pollution and mortality in cohort study designs using the Cox proportional hazards model. However, we restricted the shape of the concentration–mortality association such that the HR estimates are always greater than unity and are monotonically nondecreasing, with the amount of curvature that makes them suitable for health-effects assessment. Our SCHIF-estimated functions display a similar pattern to that of using the more unstructured quintiles of exposure. Our study provides much needed direct evidence of the shape of the exposure–response relationship for and mortality over the entire global range. For this reason, it also provides a test of predictions of relative effects based on the IER used in GBD (Burnett et al. 2014; Cohen et al. 2017). Our results suggest that the IER may currently underestimate the excess relative risk over the exposure range experienced in China and other countries with high exposure to . Although both the IER and the risk functions fit to the data in this study are constrained to be monotonically increasing, our results corroborate both the nonlinearity in CVD mortality relative risks predicted by the IER and the differences in cause-specific risk function shapes, which was also noted by Pope et al. (2011). They provide important new data points for future updates of the GBD risk functions, but, given the variation in cause-specific effect sizes observed in North American and European studies, it will be important to replicate this study in China and other high-exposure settings.

To our knowledge, this study is the first prospective cohort study to investigate the effect of long-term exposure to on CVD mortality and assess the concentration–response relationship in China. The large sample size and multicity coverage add strength and generalizability to our study. The estimates for exposure across all study areas from 1990 filled an important gap of previously unavailable historical air-pollutant exposure data. Furthermore, no previous cohort studies have examined the long-term health impacts of in rural areas in China. The findings from this study provide important estimates of long-term effects of air pollution for high-exposure settings typical in developing countries, and especially in China, where there is an urgent need to address the challenge of air-pollution exposure reduction for both urban and rural areas. Our estimates can lend support for decision makers to implement cost-effective air-quality improvement strategies and management plans.

Despite its novelty and many strengths, there are also some limitations to our study. First, we included only male participants, because this cohort was originally designed to investigate the effects of smoking (Niu et al. 1998). We therefore cannot assess the effects of long-term exposure to on nonaccidental and CVD mortality in women, although previous studies showed no systematic sex differences in the association between and mortality (Wong et al. 2015). Second, exposure was assessed at the center level and derived from gridded exposure estimates at resolution. As such, these estimates do not incor-porate finer-scale spatial variation such as that related to topography or proximity to roads and point sources (Beelen et al. 2008). The differences between the proxy estimates and true exposure values are a source of unavoidable measurement error, which may bias the interpretation of air-pollution health effects. In a cross-validation analysis in which 10% of the global measurement sites used for calibration were randomly selected for model evaluation and repeated for a total of three separate sets of 10% testing sites, we estimated a root mean squared error of (Brauer et al. 2016). Third, we were not able to adjust for multiple pollutants, such as sulfur dioxide, nitrogen dioxide, and ozone, due to data availability. Fourth, we cannot adjust for some important area-level SES variables, such as health care indicators, as the information was not available in the study period. It will be possible to consider this factor in the future follow-up of the cohort. However, more cohort studies of long-term health effects of exposure should be conducted in China to confirm the magnitude of effects found in the current study and address the limitations.

In conclusion, our study found increased nonaccidental, CVD, COPD, and lung-cancer mortality risks associated with high levels of exposure to ambient in China. The study provides information on exposure–response curves over a much broader range of exposure than previously studied specifically at high concentrations, contributing important new information for burden-of-disease assessments at the national and global levels. Our results suggest that the IER estimator used by the GBD project may underestimate the excess relative risk of cause-specific mortality due to long-term exposure to over the exposure range experienced in China and other low- and middle-income countries, but replication of these results is needed and should be a high priority for future research.

Click here for additional data file.