Low concentrations of fine particle air pollution and mortality in the Canadian Community Health Survey cohort

Background Approximately 2.9 million deaths are attributed to ambient fine particle air pollution around the world each year (PM2.5). In general, cohort studies of mortality and outdoor PM2.5 concentrations have limited information on individuals exposed to low levels of PM2.5 as well as covariates such as smoking behaviours, alcohol consumption, and diet which may confound relationships with mortality. This study provides an updated and extended analysis of the Canadian Community Health Survey-Mortality cohort: a population-based cohort with detailed PM2.5 exposure data and information on a number of important individual-level behavioural risk factors. We also used this rich dataset to provide insight into the shape of the concentration-response curve for mortality at low levels of PM2.5. Methods Respondents to the Canadian Community Health Survey from 2000 to 2012 were linked by postal code history from 1981 to 2016 to high resolution PM2.5 exposure estimates, and mortality incidence to 2016. Cox proportional hazard models were used to estimate the relationship between non-accidental mortality and ambient PM2.5 concentrations (measured as a three-year average with a one-year lag) adjusted for socio-economic, behavioural, and time-varying contextual covariates. Results In total, 50,700 deaths from non-accidental causes occurred in the cohort over the follow-up period. Annual average ambient PM2.5 concentrations were low (i.e. 5.9 μg/m3, s.d. 2.0) and each 10 μg/m3 increase in exposure was associated with an increase in non-accidental mortality (HR = 1.11; 95% CI 1.04–1.18). Adjustment for behavioural covariates did not materially change this relationship. We estimated a supra-linear concentration-response curve extending to concentrations below 2 μg/m3 using a shape constrained health impact function. Mortality risks associated with exposure to PM2.5 were increased for males, those under age 65, and non-immigrants. Hazard ratios for PM2.5 and mortality were attenuated when gaseous pollutants were included in models. Conclusions Outdoor PM2.5 concentrations were associated with non-accidental mortality and adjusting for individual-level behavioural covariates did not materially change this relationship. The concentration-response curve was supra-linear with increased mortality risks extending to low outdoor PM2.5 concentrations.

Background Exposure to ambient fine particle air pollution (PM 2.5 ) is responsible for an estimated 2.9 million deaths annually and 83 million disability-adjusted life years lost [1], with several large epidemiological cohort studies linking long-term exposure to PM 2.5 to all-cause and cause-specific mortality [2][3][4]. Even in settings with relatively low concentrations of air pollution, such as Canada, the relationships persist [5,6]. Despite these findings, there remain two key areas of potential bias and uncertainty that past work has been unable to address simultaneously. The first is the inability to directly adjust for individual-level behavioural risk factors associated with chronic disease mortality, such as smoking, diet, and exercise, or health measures such as body mass index; various indirect methods for adjustment have been applied elsewhere [7,8]. The second regards the shape of the concentration-response curve for PM 2.5 and mortality. This issue has become increasingly pertinent as clean air regulations have succeeded in reducing PM 2.5 concentrations across North America and elsewhere, and thus a better understanding of the shape of the PM 2.5 -mortality associations at low concentrations are required for cost-benefit assessments of future reduction efforts.
The purpose of this study was to provide an updated and extended analysis of the Canadian Community Health Survey-Mortality cohort [9] including [1]: additional years of follow-up to 2016 [2]; improvements in the resolution of PM 2.5 exposure (approximately 1km 2 grid) [3]; annual residential history from 1981 to 2016 for all cohort members from a linkage to postal code records [4]; time-varying contextual covariates [5]; inclusion of immigrants to Canada, and [6] an improved linkage between survey respondents and death records. We examine the shape of the concentration-response curves using a Shape Constrained Health Impact Function (SCHIF) [10] and perform sensitivity analyses.

CCHS-mortality cohort
The Canadian Community Health Survey (CCHS) is a national cross-sectional survey of the Canadian population that collects information related to health status, health care utilization, and health determinants. From 2000 to 2007 the survey was administered every 2 years to approximately 130,000 respondents; from 2007 onwards, data has been collected on an ongoing basis from 65,000 respondents per year and released annually with response rates declining over time ( Fig. 1) [11][12][13][14][15][16]. The .8%. f see methods for list of exclusion criteria, totals will exceed number of deleted person-years given that more than one exclusion criteria may apply to a single person-year; immigrated to Canada less than 10 years before survey date n = 541,600, age during follow-up period exceeds 89 years n = 161,000, no postal code n = 5,009,900, could not be linked to air pollution values n = 5,711,600, could not be linked to Can-MARG values n = 7,668,000, could not be linked to Census Metropolitan Area/Census Agglomeration size n = 4,800,600, could not be linked to airshed n = 3500, 3-year moving average being informed by only 1 year of exposure n = 4,321,500, year after subject death n = 343,600, year before survey interview date n = 13,570,300. g from 452,700 unique individuals CCHS data are sampled from approximately 98% of the Canadian population aged 12 and older living in private dwellings within the 115 Heath Regions covering all provinces and territories. Individuals living on Indian Reserves and on Crown Land, institutional residents, full-time members of the Canadian Forces, and residents of certain remote regions are excluded.
Consent to record linkage and data sharing was obtained at time of survey ( Fig. 1) and only those CCHS respondents who agreed were linked to death records and residential history through Statistics Canada's Social Data Linkage Environment (SDLE) [17]. The linkage was approved by Statistics Canada [18] and is governed by the Directive on Microdata Linkage [19]. The linkage occurred within the Derived Record Depository, a highly secure linkage environment comprised of a national dynamic relational database of basic personal identifiers. Survey and administrative data are linked to the Derived Record Depository using a SAS-based generalized record linkage software (G-link) that supports deterministic and probabilistic linkage based on the mathematical theory of record linkage developed at Statistics Canada [20]. Mortality linkage to the Derived Record Depository between 2000 and 2016 was 99.8% [21]. A list of linked unique individuals is created through linkages that are deterministic (matching records based on unique identifiers) and probabilistic (matching records based on nonunique identifiers such as names, sex, date of birth, and postal code and estimating the likelihood that records are referring to the same entity). For the CCHS cycles considered, there was a linkage rate to the Derived Record Depository of 95.2% and a false error rate for the CCHS to SDLE linkage of 0.4% [22].
There were 666,000 CCHS respondents who agreed to data linkage ( Fig. 1), reduced to 540,900 after excluding subjects with death dates prior to survey response dates (i.e. either the death record or linkage must be incorrect) or who were below the age of 25 or above 89 at time of survey as they are less likely to reside at the same postal code as their income tax mailing address [23]. The CCHS to SDLE linkage rates across key indicators were consistently high, ranging from 94.4% for the 20-29 age group to 96.2% for the 80-89 age group, 95.5 and 95.3% for males and females respectively, and by province/territory from 91.8% for the Yukon to 96.7% for Newfoundland and Labrador [22].
After linkage, we stacked the CCHS cycles into one data file. We standardized the variable categorizations when discrepancies between cycles existed. The covariates (listed, along with categorizations in Table 1) included socio-economic, behavioural, and contextual measures. More information about the definitions and classifications of these variables can be found elsewhere [9]. Provincially standardized deciles were calculated according to the distribution of residents in each province and the ratio of their total household income to the low-income cut-off for their corresponding household and community size. As this measure excluded subjects living in territories, we took the mean income within each decile and used these as cut-offs to categorize subjects living in territories by income into deciles. Once all subjects were placed in deciles, we merged groups to create quintiles.
Postal code history was complete from 1981 to 2016 for 35.0% of respondents and 12.6% of respondents had no postal code history. There were gaps in postal code histories for 52.4% of respondents, which is to be expected, as taxes may not have been filed for various reasons (e.g immigration, death, or age). We imputed complete or partial postal codes only when bookended by postal codes with sufficient similarity before and after the gap [24]. For example, if a postal code in 2008 was K1A 0T6 and then K1A 0K9 in 2012, a partial postal code of K1A 0** would be imputed for the four missing years from 2009 to 2011. We did not impute postal codes if a gap existed at the beginning or end of the follow-up period or after a person's death; full or partial postal codes (two to five digits) were imputed for 1.5% of person-years.
We organised the cohort into a person-year file with each year of exposure (1981-2016) per person representing a row of data. Subsequently we excluded specific person-years [1] once they turned age 90 during followup [2], if the person had immigrated to Canada less than 10 years prior to survey interview [3], if there was no postal code [4], if the postal code could not be linked to air pollution or contextual covariates [5], if the PM 2.5 three-year moving average with a one-year lag was calculated by fewer than 2 years of exposure data, or [6] if the person-year was before survey interview date or after a person's death (Fig. 1). We excluded recent immigrants to Canada (10 years or less) since they have spent the majority of their lives outside of Canada with unknown exposure, and this time exceeds the number of years in Canada where exposure can be estimated.

Exposure file and analytical file
The task of linking contextual covariates and air pollution values to the cohort required the creation of a master list of postal codes with their respective points of latitude and longitude and census geography. We pro-  We developed and used annual exposure estimates of PM 2.5 from 1998 to 2012 by relating satellite retrievals of aerosol optical depth (AOD) to near-surface PM 2.5 concentrations using the geophysically-based relationship simulated by a chemical transport model [28]. Ground monitoring data from the National Air Pollution Surveillance (NAPS) network were then incorporated, along with other North American-based measurements, to constrain these PM 2.5 estimates with geographically weighted regression. The resulting ambient PM 2.5 surface provided estimates for North America at about a 1km 2 resolution [28]. Spatial variation from this surface was used with simulated PM 2.5 and consistently constrained with local ground-based monitors to extend our PM 2.5 coverage to 2015 [29].
The ambient warm season daily-maximum eight-hour average O 3 surfaces were developed by Environment and Climate Change Canada for 2002-2015 using chemical transport modelling informed by surface observations as hourly estimates from 2002 to 2015 [30][31][32]. Estimates of NO 2 were created using a national land use regression model (LUR) informed by on satellite-derived NO 2 (10 km resolution), distances to highways and major roads, and roadway kernel density gradients [33]. Ozone and NO 2 values were back-casted to obtain exposures for 1981-2015 using ground monitoring data from the Canadian National Air Pollution Surveillance program. Annual adjustment factors were calculated at a census division level from the ratio of observed concentration to the values in the surface for the reference year (see Pinault et al. for more detail [9]).
We linked postal codes to PM 2.5 in ArcGIS Desktop 10.5.1 using the points of latitude and longitude from the master postal code list and the air pollutant surfaces. In cases where there were multiple points of latitude and longitude for a single urban postal code, we used equal weighting of the multiple air pollutant values to provide a singular value. In rural communities, we took the population-weighted average of the values associated with duplicate postal codes. We used population-weighing to average multiple values to create inputs for partial postal codes (2 to 5 digits).

Covariates
Contextual covariates were available at various census geographies and we merged these to individual personyears via postal codes (as described below). We created historic measures when possible to reflect neighbourhoodlevel changes over the time.
Regions of Canada that share air quality characteristics and movement patterns have been defined by the Air Quality Management System (AQMS) as six distinct airsheds [34]. By subdividing the country into large geographic areas, adjustment for the broad spatial variation in mortality rates can be performed [9,34]. We assigned airshed to the cohort by postal code. We used a population weighted mode in cases where there were multiple points of latitude and longitude for a single postal code.
We developed a historic community size variable to account for different sizes of metropolitan regions and changes in population over time, classifying Census Metropolitan Areas (CMAs: major urban core, 100,000+ residents) and Census Agglomerations (CAs: smaller urban cores, 10,000+ residents) by population counts [35]. Since CMA/CAs cover large areas that can include farmland near the urban-rural fringe and residential enclaves of commuters to the city, we created a measure that accounts for differences in urban form within these CMA/CAs. We used population density measures (1991-2016) and frequencies for different modes of transportation at the neighbourhood level (1996-2016) to categorize census tracts as active urban core, transit-reliant suburb, car-reliant suburb, exurban, and non-CMA/CA [36]. Both CMA/CA size and urban form were attached to the postal code list via census geography before merging with the cohort. In cases where there were multiple points of latitude and longitude representing a postal code, we used a population-weighted mode to assign categories.
The Canadian Marginalization Index (Can-MARG) is a measure of community-level marginalization comprised of four factors: material deprivation (e.g. proportion of people living in dwellings in need of repair), residential instability (e.g. proportion of people who live in a dwelling that they do not own) dependency (e.g. proportion of seniors and youth compared to those who are not), and ethnic concentration (e.g. proportion of recent immigrants and self-reported visible minorities) [37]. We used historic census tract-level Can-MARG values in CMA/CAs, and a population-weighted aggregation of the dissemination area-level Can-MARG values at the census subdivision level in rural areas outside of CMA/CAs that are not covered by census tracts. We assigned Can-MARG values to points of latitude and longitude before quintiles were assigned.

Statistical analysis
We calculated for each individual and year of follow-up a three-year moving average for PM 2.5 with a one-year lag, (e.g. the exposure in 2002 is the average of exposures in 1999, 2000, and 2001).
We ran standard Cox proportional hazard models to assess the relationship between PM 2.5 exposure and non-accidental death (ICD-10 codes A to R) from survey interview year to the end of follow-up period or year of death. We started model building with a baseline hazard function for PM 2.5 stratified by five-year age groups, sex, and survey cycle to ensure that respondents within these strata would be broadly comparable. We calculated new hazard ratios for models that included each socio-economic and behavioural covariate individually. We included covariates in the partially-adjusted model if the log difference between the new hazard ratio and the baseline was more than 10%. Subsequently, we added contextual covariates individually to the partially-adjusted model and included them in the final model using the same criteria (comparing to the partially-adjusted model that included socio-economic and behavioural covariates). All covariates considered for inclusion in the final model and the associated hazard ratios are found in Table 1.
We examined the shape of the association between PM 2.5 and mortality with a SCHIF [10]. This method is based on a construction of several transformations of concentration and fitting the transformed variable in a Cox model, estimating the log-hazard ratio for a unit change in the transformed variable and its standard error. An ensemble of all models examined was then constructed using a weighted average of the predicted log-hazard ratio and any concentration, with weights defined by the AIC of each model. The transformations are variations on a sigmoidal function which yields supralinear, near linear, and sub-linear shapes.

Sensitivity analyses
We examined effect modification by select socio-economic and behavioural covariates, and by high-and lowexposure groups to the combined oxidant capacity of NO 2 and O 3 (henceforth: O X ) which is calculated as the redox-weighted oxidant capacity [38] i.e. a weighted average using redox potentials as the weights (O x wt = [(1.07 V × NO 2 ) + (2.075 V × O 3 )]/3.145 V) ( Table 4). We examined multiple pollutant models to investigate whether the inclusion of other common pollutants (NO 2 , O 3 , and O X ) in the model may modify the PM 2.5mortality relationship [5,39].

Results
There were 4,452,700 person-years in the cohort after exclusion criteria were applied (Fig. 1) from 452,700 unique individuals. Entry into the cohort and length of the follow-up period varied by survey cycle, with the first cohort having up to 15 years of follow-up. For those who died, the average follow-up period was 5.1 years (s.d. 3.4); it was 6.5 years (s.d. 4.1) for those who survived the follow-up period. There were 50,700 non-accidental deaths. Of these, there were 7900 deaths from ischemic heart disease, 2800 from cerebrovascular disease, and 4300 from other cardiovascular diseases; 900 from pneumonia, 2800 from COPD, and 1100 from other respiratory diseases; 5500 from lung cancer, 1300 from colon cancer, 1300 from breast cancer, 1100 from pancreatic cancer, and 9900 from all other cancers. Further, there were 1700 deaths from diabetes, 3900 deaths from neuropsychiatric conditions, 2200 from digestive diseases, 1100 from genitourinary diseases and 3000 from all other non-accidental causes.
Exposure to PM 2.5 was higher in women, more recent immigrants, and non-Indigenous people. Being single, university educated, and in the poorest income quintile were also associated with higher exposures (Table 1). We observed higher exposure to PM 2.5 in people living in the largest CMAs and in the East Central airshed (which includes Toronto and Montreal). The distribution of exposure estimates for PM 2.5 , NO 2 , O 3 , and O X is found in Table 2.
The cohort was generally representative of the Canadian population, as seen through their mortality rates by subgroup (Table 1). Immigrants and non-Indigenous people had lower mortality rates compared to their counterparts. Being married, holding a university degree, and being employed were associated with a lower risk of mortality. As expected, there were clear trends in mortality risk with income, education, and immigrant status.
The unadjusted model had a hazard ratio of 0.96 (95% CI 0.92-1.00) which increased to 1.11 (95% CI 1.04-1.18) when adjusted by the socio-economic, behavioural, and contextual covariates that met the inclusion threshold ( Table 3). All covariates except for body mass index (BMI), employment status, and urban form met the criteria and were included in the final model. When we added the behavioural covariates to a model that included only socio-economic covariates the hazard ratio increased from 1.05 (95% CI 1.00-1.09) to 1.09 (95% CI 1.05-1.15). Conversely, when we added the behavioural covariates to a model that included both the socio-economic and contextual covariates, they lowered the PM 2.5 hazard ratio from 1.13 (95% CI 1.06-1.21) to 1.11 (95% CI 1.04-1.18).
The SCHIF characterisation of the PM 2.5 -mortality association (for all cohort members) displayed a supra-linear shape that rises in a steeper fashion compared to the standard log-linear model prediction for lower concentrations and changes in a more moderate manner for higher levels (Fig. 2). Note that the SCHIF displays wider uncertainty intervals compared to the log-linear model at low concentrations, in part due to the additional variation associated with model shape, a feature captured by the SCHIF but not the log-linear model. We observed a positive and statistically significant (p < 0.05) association between PM 2.5 and non-accidental mortality for all concentrations examined as indicated by the SCHIF hazard ratio predictions. We assessed effect modification within the PM 2.5 -mortality relationship by separating the cohort by age, sex, immigrant status (i.e. immigrants who had been in Canada for 10 or more years vs. non-immigrants), and educational attainment, and comparing resulting hazard ratios with Cochrane's Q (

Discussion
Using a cohort comprised of several cycles of a health survey with up to a 15-year follow-up period and high resolution exposure estimates, we found that exposure to PM 2.5 was associated with an 11% increase in non-accidental mortality per 10 μg/m 3 after extensive adjustment for socio-economic, behavioural, and contextual covariates. The hazard ratio for the full cohort was similar to that of the Nurse's Health Study in the United States (1.13 95% CI 1.05-1.22) that adjusted for individual-level socio-economic and behavioural covariates [40] and a cohort from England (1.13 95% CI 1.00-1.25) that controlled for smoking, BMI, income, age, and sex [41]. Burnett and colleagues [42] report hazard ratio estimates for a 10 μg/m 3 change in long-term exposure to PM 2.5 and non-accidental mortality in 41 cohorts conducted globally, 36 of which included adjustment for behavioural risk factors. The pooled hazard ratio among these 36 cohorts was 1.09 (95% CI 1.05-1.12), a value similar to that observed in our current study (1.11 95% CI 1.04-1.18). A version of the 2001 CanCHEC census-based cohort produced a hazard ratio that is similar to this work (1.09 95% CI 1.07-1.11) [6]. The impact of individual-level behavioural risk factors on the PM 2.5 -mortality association was assessed to address a common critique of many large administrative cohort studies examining the air pollution-mortality relationship. The inclusion of behavioural covariates to a model including socio-economic and ecological covariates lowered the PM 2.5 hazard ratio 2% (from 1.13 to 1.11). This modest change in the hazard ratio can be interpreted to indicate that the behavioural covariates were being adequately controlled for by the socio-economic and ecological covariates in the established relationship between PM 2.5 exposure and non-accidental mortality. This finding is similar to the previous CCHS cohort analysis and analysis of a Medicare-based cohort; both reported that adjustment for behavioural covariates had a minimal effect on hazard ratios [3,7]. There is evidence (Tables 3 and 4) for a small increase in risk of PM 2.5 -related mortality in occasional or regular drinkers but this may be masked by null effects from the inclusion of other behavioural covariates (fruit and vegetable consumption, smoking behaviours) and this confounding is likely the result of the spatial distribution of drinking behaviours, with binge drinkers having the largest mortality risk but lower PM 2.5 exposures. This study, through its inclusion of multiple covariates and an explicit a priori analysis approach for model building therefore provides the most extensive evidence to date that, in the Canadian context, missing data on behavioural risk factors for mortality have a minimal confounding bias on the PM 2.5 -mortality association.
The increase in the PM 2.5 hazard ratio with the addition of the ecological covariates was largely driven by the addition of airsheds. Not only do these airsheds characterize broad air movement patterns, they also capture areas with similar composition of PM 2.5 (e.g., proportion of PM 2.5 composed of nitrate is highest in the Prairie airshed, whereas the Southern Atlantic airshed is composed of a notably higher proportion of black carbon) [34]. They also delineate general socio-cultural groups with distinct mortality risk factors beyond those captured by the typical socioeconomic census variables included in our survival models. The three airsheds with the largest hazard ratios, along with high material deprivation, all have the lowest levels of air pollutants which would account for the negative confounding effect observed in Table 3. Further, the largest airshed (East Central) contains both Toronto and Montreal, the two largest CMAs in Canada and significant population hubs. High PM 2.5 exposure and related mortality are largely driven by the population of Toronto (21% of the national population in 2006) where the mean PM 2.5 exposure is 9.33 μg/m 3 whereas the mean in the rest of the country is 7.68 μg/m 3 [43]. These results are consistent with a descriptive analysis of PM 2.5 exposure in 2006 long-form census respondents [9]. Although urban areas are the most common residence for both high income and highly educated Canadians, rural residences are more common among the high income earners than university graduates (i.e within the highest income quintile, 73.7% urban vs 26.3% urban fringe or rural; among those who are university educated, 82.6% urban vs. 17.3% urban fringe or rural). The greater tendency for high-income Canadians to live in rural areas is consistent with the findings in this paper. As a result, PM 2.5 exposure by income categories is a slightly more linear pattern than education in both of these studies. We estimated the shape of the concentration-response (CR) function for the PM 2.5 -mortality association. A slight supra-linear association (Fig. 2) was found, with a steep CR function at the low to median PM 2.5 range which levelled off slightly after approximately 10 μg/m 3 . The SCHIF hazard ratio predictions indicated a positive and significant association between PM 2.5 and non-accidental mortality for all concentrations, suggesting risks to concentrations below 2 μg/m 3 . Previous work using a CCHS-based cohort used a spline-based procedure and found that the shape of the relationship between non-accidental mortality and PM 2.5 was supra-linear in shape with a threshold of 4.5 μg/m 3 , but was limited due to wide confidence intervals [9]. A study in China using a SCHIF function found non-linear relationships for multiple causes of death [44]. Such a relationship, when applied in a health impact framework, as in the Global Burden of Disease [45,46] and in the recent Global Exposure Mortality Model [42] suggest benefits both from reducing PM 2.5 concentrations areas with the highest concentrations and from continuing to reduce them in relatively cleaner areas, including Canada, where it is estimated that the entire population now lives in areas with ambient PM 2.5 concentrations below the current WHO Guideline [47]. Worldwide it is estimated that small absolute reductions under 3 μg/m 3 could prevent hundreds of thousands of deaths in areas that comparatively have low levels of PM 2.5 [48].
The risk of non-accidental mortality from exposure to PM 2.5 was 4% higher in males over females (males 1.13, females 1.09), a pattern that has emerged in similar work. The hazard ratios from the current study are more aligned with the ESCAPE European pooled cohort (males 1.14 95% CI 1.04-1.24; females 0.99 95% CI 0.92-1.07) [2] albeit with a higher hazard ratio for women when compared to the previous version of CCHS-based cohort (males 1.34 95% CI 1.24-1.46; females 1.18 95% CI 1.09-1.28) [9]. Hazard ratios were lowest for members of the cohort aged 75 and older (1.04) and were similar for those aged 65 and under (1.14) and 65 to 75 (1.13); this is similar to the European study which found that risk decreases with age (< 60 years 1.16 95% CI 1.00-1.34; 60-75 years 1.10 95% CI 1.00-1.20; ≥75 years 1.03 95% CI 0.95-1.11). When we divided the cohort into immigrants (in Canada for 10 years or more) and non-immigrants, the PM 2.5 -mortality association increased for non-immigrants and was null among the immigrant population. This result is consistent with prior Canadian census-based cohort studies [5] and is possibly the result of what is termed the "healthy immigrant effect" [49][50][51][52][53], likely intensified by the preferential settlement of immigrants into the largest cities which have higher PM 2.5 exposure. The hazard ratio for high school graduates (1.14) was higher than for those without a diploma (1.08) which is to be expected given that the latter is more likely to live in rural areas [43], and have a mean PM 2.5 exposure that is lower than other educational groups (Table 1) [43]. We examined effect modification by behavioural covariates (i.e., fruit and vegetable consumption, smoking behaviour, and alcohol consumption) and found significant difference in the resulting hazard ratios only in the case of alcohol consumption. Effect modification analyses on the ESCAPE cohort also found no effect modification by fruit and vegetable consumption or smoking behaviour, but did not consider alcohol consumption [2]. The multiple pollutant models indicated that the relationship between non-accidental mortality and PM 2.5 exposure are attenuated when we included other pollutants (NO 2 , O 3 , and O X ) in the models. These findings indicate both that PM 2.5 is associated with mortality and that the inclusion of gaseous co-pollutants, O x in particular, may better characterize the biologically active aspects of PM 2.5 and the overall air pollution mixture compared to the PM 2.5 mass concentration [5]. Weichenthal et al. looked at the effect modification of oxidant gases on PM 2.5 more specifically and found that spatial variations in O x concentrations may act as surrogates for the presence or absence of harmful air pollutant mixtures that enhance PM 2.5 toxicity [42]. We examined the PM 2.5 -mortality association in both lowand high-O x person-years and found a 24% difference in risk. Our findings support these previous studies using different longitudinal Canadian cohorts and that knowledge of interactions between PM 2.5 and oxidant gases leading to adverse health will improve risk management activities and public health.
We performed this analysis on an extended and updated version of a cohort described in a previous study by Pinault et al. [9] with improvements to the exposure assessment and linkage to death, postal code history, environmental exposures, and contextual covariates. While some of the results are comparable to the previous cohort (e.g. socio-economic + behavioural covariate models are within a 1% margin), there are differences in the covariates included in the final models and the resulting hazard ratios. This is not unexpected since the contextual covariates addressing area-level marginalization in the two studies were created differently (area-level proportions of specific variables vs. a principle component analysis which resulted in four factors), and measured at different geographical units (census divisions vs. census tracts and census subdivisions). Another difference is that the updated cohort and current work includes immigrants who have lived in Canada for ten or more years whereas the previous work only included those who had been in Canada for 20 or more years. This newly included group of semi-recent immigrants (10-19 years in Canada) have substantially lower hazard ratios of mortality compared to the non-immigrant population (Table 1). Their inclusion in the current study acts to reduce the overall PM 2.5 hazard ratio (Table 4).
This large, national cohort is an extension and improvement to the previous CCHS-Mortality cohort, with an updated linkage and extended follow-up period for mortality and postal code history which now spans 36 years (1981 to 2016). More broadly the cohort has many strengths, including the fine resolution of the PM 2.5 estimates (1km 2 ), the ability to incorporate mobility across the follow-up years, an explicit a priori model building strategy, the inclusion of multiple time-varying contextual covariates to address spatial, neighbourhood-and city-level characteristics, and most uniquely the behavioural covariates such as smoking behaviours, alcohol consumption, diet, and exercise to control for health behaviours related to mortality that are not typically found on cohorts of this size.
This cohort and the analysis are limited by the data available. First, postal code history was derived from tax and administrative data. Historical postal codes reflect the mailing address as reported on a tax return and not necessarily a person's residence; in 92.9% of cases the postal code reflects the person's residence at time of survey [23]. Similarly, outdoor ambient levels of PM 2.5 at a person's residence may not reflect their actual exposure. Sensitivity analysis performed with the 2001 CanCHEC found that finer scale resolution (1km 2 ) estimates of PM 2.5 resulted in lower AIC values and higher hazard ratios in the PM 2.5 -mortality model for non-accidental death compared to a 10km 2 or 5km 2 grid indicating that exposure estimates that are more specific to a person's residence are appropriate [54]. Gaps in postal code history are imputed under the assumption that the person did not leave the country or community during that time. In assigning contextual covariates by postal code, misclassification may occur from taking the mode or mean when estimating a single value to represent multiple points of latitude and longitude for a single postal code. Second, in contrast to the CanCHEC cohorts (Pappin AJ, Crouse DL, Christidis T, Pinault LL, Tjepkema M, Erickson A, Brauer M, Weichenthal S, van Donkelaar A, Martin RV, Brook J, Hystad P. Burnett RT. Associations between low levels of fine particulate matter andmortality within Canadian cohorts. Environ Health Persp., under review), this cohort does not completely represent the full Canadian population; the Canadian Community Health Survey is not a census of the population and survey weights were not used in this analysis. Further, in creating this cohort persons were removed if they did not consent to data linkage or if they could not be linked to the SDLE. The CCHS over-samples rural communities [55] which results in a disproportionate sample in areas with low levels of PM 2.5 and higher rates of mortality. The sampling framework and un-weighted analysis likely caused the null unadjusted hazard ratio which became positive as covariates were added to the model to address confounding. These results are consistent with the Agricultural Health study which examined non-accidental death related to PM 2.5 in rural communities in two American states (Iowa and North Carolina) and found a protective hazard ratio in minimally and fully adjusted models [56]. Regardless, the protective unadjusted hazard ratio should not come as a surprise as contextual and socio-economic covariates are included in models because we know that they are related to both PM 2.5 and mortality and can act as confounders (see Table 1 for the mortality Hazard Ratios by individual covariates). Given that these factors covary with both mortality and PM 2.5 their inclusion in the models is crucial. We suggest that the unadjusted model is not reflective of the PM 2.5 -mortality relationship and that the direction or magnitude should not be over-interpreted. Third, although this cohort includes behavioural covariates these are selfreported and in some cases there are missing responses. To avoid introducing bias into the cohort, we used dummy variables to code missing information rather than excluding non-respondents outright. Finally, the cohort itself is limited by follow-up and some persons have as few as 4 years of follow-up (with a maximum follow-up of 15 years).

Conclusions
We provided an update to the Canadian Community Health Survey-Mortality cohort, with a new linkage of the survey respondents to death records, inclusion of additional survey cycles, an extension of the annual residential history and mortality follow-up period, a finer scale of air pollution exposure, time-varying contextual covariates, and the inclusion of immigrants who have lived in Canada for 10-20 years (rather than only those who have been in Canada for 20+ years). The risk of non-accidental mortality from ambient PM 2.5 was found even at low levels although the hazard ratio was attenuated in models that included other pollutants (NO 2 , O 3 , and O X ). The PM 2.5 -mortality association displayed a supra-linear concentration-response curve. The inclusion of behavioural covariates that could confound the PM 2.5 -mortality association (fruit and vegetable consumption, leisure exercise frequency, alcohol consumption, and smoking behaviours) did not appear to impact hazard ratios. Hazard ratios were higher for males, those aged 65 or less, and non-immigrants.