key: cord-0737643-xz5klqe5
authors: Stieb, David M.; Evans, Greg J.; To, Teresa M.; Brook, Jeffrey R.; Burnett, Richard T.
title: An ecological analysis of long-term exposure to PM(2.5) and incidence of COVID-19 in Canadian Health Regions
date: 2020-08-26
journal: Environ Res
DOI: 10.1016/j.envres.2020.110052
sha: cfada15663b9493468f03a06e4db4558f86c58be
doc_id: 737643
cord_uid: xz5klqe5

BACKGROUND: Ambient fine particulate matter (PM(2.5)) is associated with a wide range of acute and chronic health effects, including increased risk of respiratory infection. However, evidence specifically related to novel coronavirus disease (COVID-19) is limited. METHODS: COVID-19 case counts for 111 Canadian health regions were obtained from the COVID-19 Canada Open Data portal. Annual PM(2.5) data for 2000-2016 were estimated from a national exposure surface based on remote sensing, chemical transport modelling and ground observations, and minimum and maximum temperature data for 2000-2015 were based on a national interpolated surface derived from thin-plate smoothing splines. Population counts and sociodemographic data by health region were obtained from the 2016 census, and health data (self-rated health and prevalence of smoking, obesity, and selected chronic diseases) by health region, were obtained from the Canadian Community Health Survey. Data on total number of COVID-19 tests and changes in mobility comparing post-vs. pre-introduction of social distancing measures were available by province. Data were analyzed using negative binomial regression models. RESULTS: After controlling for province, temperature, demographic and health characteristics and days from peak incidence by health region, long-term PM(2.5) exposure exhibited a positive association with COVID-19 incidence (incidence rate ratio 1.07, 95% confidence interval 0.97-1.18 per μg/m(3)). This association was larger in magnitude and statistically significant in analyses excluding provinces that reported cases only for aggregated health regions, excluding health regions with less than median population density, and restricted to the most highly affected provinces (Quebec and Ontario). CONCLUSIONS: We observed a positive association between COVID-19 incidence and long-term PM(2.5) exposure in Canadian health regions. The association was larger in magnitude and statistically significant in more highly affected health regions and those with potentially less exposure measurement error. While our results generate hypotheses for further testing, they should be interpreted with caution and require further examination using study designs less prone to bias.

There is considerable evidence from multiple lines of research including toxicology, human clinical studies and epidemiological studies that air pollution in general and fine particulate matter (PM 2.5 ) more specifically are associated with a wide range of acute and chronic health effects (Thurston et al., 2017) . Based on estimates from the Global Burden of Disease initiative, PM 2.5 accounts for the greatest burden of mortality of any environmental exposure, accounting for approximately 3 million worldwide deaths annually (GBD 2017 Risk Factor Collaborators, 2018 . It is well established that acute exposure increases the risk of emergency visits and hospital admissions for respiratory infections including pneumonia (Atkinson et al., 2014; Domingo and Rovira, 2020; Nhung et al., 2017) . There is also a growing body of evidence that long term exposure increases the risk of morbidity and mortality from respiratory infection (Mehta et al., 2013; Neupane et al., 2010) .

Evidence specifically related to PM 2.5 and novel coronaviruses such as severe acute respiratory syndrome coronavirus (SARS-CoV-1) and middle east respiratory syndrome (MERS) is limited. Studies based on the SARS-CoV-1 outbreak suggest that meteorology and exposure to air pollution increased transmission (Cai et al., 2007) and worsened patient prognosis (Kan et al., 2005) . Notably SARS-CoV-1 patients from more polluted regions were twice as likely to die as those in less polluted places (Cui et al., 2003) . There is also evidence that air pollution exposure more generally adversely affects respiratory immune defences (Domingo and Rovira, 2020; Yang et al., 2020) , and emerging evidence suggesting that novel corona virus disease incidence and mortality may be increased in relation to both acute (Zhu et al., 2020) and chronic exposure (Andree, 2020; Liang et al., 2020; Ogen, 2020; Wu et al., 2020) . This evidence implies that deterioration in air quality over short time periods (e.g. from wildfire smoke, other local burning, specific meteorological events such as temperature inversions) may lead to more cases of severe COVID-19 infections, adding further demand to the healthcare system. Conversely, improving air quality by reducing both the occurrence of acute events and long term average concentrations, may help to protect communities from COVID-19 and reduce the burden on hospitals.

In this study, we conduct an ecological analysis of COVID-19 cases and 17 year average PM 2.5 concentrations among Canadian health regions. While ecological analyses have many limitations which preclude attribution of cause and effect, they can be readily conducted once data are available and permit the generation of hypotheses to be more rigorously examined in subsequent studies.

COVID-19 case counts compiled from publicly available reports for 111 health regions were obtained from the COVID-19 Canada Open Data portal (Berry et al., 2020) . Health regions are defined by provincial ministries of health; in some jurisdictions, they correspond to areas served by local public health departments or authorities (Statistics Canada, 2018a). The provinces of British Columbia and Saskatchewan reported counts for groups of health regions (up to four per group), which we distributed to individual health regions in proportion to population. Similarly, health region was not reported for many cases in Nova Scotia. These were distributed among four health regions in proportion to population. In an attempt to account for stage of outbreak, response to distancing measures, and correlation between case counts and deaths, we also obtained data on days elapsed since the first case, days since peak daily incidence of new cases, and deaths at the health region level (Berry et al., 2020) , as well as date of declaration of public health emergency or state of emergency at the provincial level (Boire-Schwab et al., 2020) . Total number of COVID-19 tests was only available by province. Annual (Boys et al., 2014) . These data have been used extensively in air pollution epidemiology studies in Canada (Pappin et al., 2019) . Greenness data for 2000-2019 were based on growing season maximum normalized difference vegetation index (NDVI) data from the MODIS onboard the Terra satellites (Didan, 2015; Gorelick et al., 2017) . NDVI values vary from -1 to 1, negative values indicating features such as water and positive values indicating vegetation. We employed positive values only in calculating average greenness by health region. Temperature data for 2000-2015 (annual minimum of lowest monthly maximum temperature -henceforth referred to as "minimum temperature," and annual maximum of highest monthly minimum temperature -henceforth referred to as "maximum temperature") were based on a national interpolated surface of available observations derived using thin-plate smoothing splines, implemented in Australian National University Spline (ANUSPLIN) climate modeling software (Wang et al. 2018) . All exposure data were available by 6 character postal codes (CANMAP 2015) , which were mapped to 2018 health region boundaries in R (GISTools (Brunsdon and Chen, 2014) , rgdal (Bivand et al., 2019) and raster (Hijmans, 2020) packages). Postal codes are used by Canada Post corporation for mail delivery, and are analogous to American zip codes -there are currently approximately 875,000 postal codes in Canada; as they are point locations, they do not have a specified surface area (CANUE, 2018) . Maps were generated using the ggplot2 (Wickham, 2016) and broom (Robinson and Hayes, 2020) packages. Population counts and sociodemographic data by health region on percent of population 65 or older, with income less than the low income cutoff (LICO), and Black (Adams et al., 2020; Dyer, 2020; Yancy, 2020) , were obtained from the 2016 census (Statistics Canada, 2018b) . LICOs are defined as J o u r n a l P r e -p r o o f income levels below which families spend a disproportionate share of their income on necessities, and are family-size and community-size specific (Statistics Canada 2017). Health data by health region on factors thought to increase susceptibility to COVID-19 (Adams et al., 2020) (percent of population who rate health as fair or poor, are daily or occasional smokers, overweight, obese, have asthma, chronic obstructive pulmonary disease (COPD), hypertension or diabetes), were based on data from the 2017 and 2018 Canadian Community Health Survey (CCHS) (Statistics Canada, 2020), which is an annual national cross-sectional survey of individuals 12 years of age and over. Changes in mobility by province comparing post-vs. pre-introduction of social distancing measures were based on aggregated data from Google Account users who opted-in to location history for their account (Google LLC, 2020) .

Data were analyzed using negative binomial regression models, specifying PM 2.5 and covariates (including province) as fixed effects. Log population was included as an offset. PM 2.5 and covariates were first regressed individually vs. case counts, then those exhibiting statistically significant associations, were included in multivariate models. Prevalence of asthma, COPD, hypertension and diabetes were excluded from multivariate models with PM 2.5 , since they could be intermediate in a putative causal pathway with COVID-19 incidence. The most parsimonious model was selected based on the Akaike Information Criterion (AIC) (Akaike, 1974) . Presence of residual spatial autocorrelation was examined by mapping model residuals and computing Moran's I (Bivand and Wong, 2018) . Sensitivity analyses were conducted by excluding Montréal (which accounted for 27.5% of cases but only 5.5% of population); excluding British Columbia, Saskatchewan and Nova Scotia (which reported cases counts only for aggregated health regions); excluding health regions with less than median population density (with presumably greater exposure measurement error over larger, more sparsely populated health regions); restricting the analysis to Ontario and Quebec, the two provinces with the highest incidence, attributed in part to provincial level policies on testing (Weeks, 2020) , and timing of school vacation periods (Perreaux, 2020) respectively; and specifying province as a random rather than fixed effect . Analysis J o u r n a l P r e -p r o o f was conducted in R (R Core Team, 2019) using the lme4 package (Bates et al., 2015) . Research ethics board approval was not required because all data were publicly available and aggregated at health region level. The work described here has been carried out in accordance with the Uniform Requirements for manuscripts submitted to Biomedical journals.

There were 73,390 cases up to May 13, 2020 and overall incidence was 208.8 cases/100,000 (2016 population). Incidence was highest in Quebec (489.0/100,000) and Ontario (166.5/100,000) and lowest in Nunavut (0/100,000) and the Northwest Territories (12.0/100,000) (Figure1). The ten health regions with the highest incidence were distributed among multiple provinces, and included three large Table S1 ). Deaths were very strongly correlated with case counts by health region (R 2 =0.95). PM 2.5 concentrations averaged 6.1 (standard deviation 2.1) and were highest in urban areas as well as more generally in southern Ontario and Quebec (Figure2) . Variability in other exposures, sociodemographic and health characteristics by health region are summarized in Table 1 . There was little variability in days since declaration of emergency. A correlogram is provided in Figure S2 . PM 2.5 was strongly positively correlated with percent Black, minimum and maximum temperature. Population density exhibited strong positive correlations with percent Black and low income, PM 2.5 , minimum and maximum temperature, as well as strong negative correlations with prevalence of obesity and smoking.

Prevalence of obesity was strongly positively correlated with prevalence of diabetes, hypertension and smoking, while percent age 65 and older was strongly positively correlated with prevalence of hypertension and COPD. Changes in mobility during the period March 13 (the date after which social distancing recommendations were introduced in most provinces)-April 26 vs. January 3-February 6, J o u r n a l P r e -p r o o f were generally comparable among provinces (Supplementary Table S3 ). Notably, changes in visits/length of stay were somewhat more modest in Nunavut, which had no cases, and visits/length of stay in parks increased substantially in British Columbia, the most populous regions of which typically experience milder weather during the March/April period compared to elsewhere in Canada.

Results of regression models are summarized in Table 2 . In bivariate models, PM 2.5 , minimum and maximum temperature, percent Black population, percent of population <LICO, population density, and days since first case were significantly positively associated with COVID-19 counts, while days since peak incidence, percent of population 65 or older, and prevalence of asthma, COPD, hypertension, overweight and obesity exhibited significant negative associations. COVID-19 counts were not significantly associated with total number of tests or percent change in mobility to work locations by province (not shown). In the best fitting multivariate model, minimum temperature exhibited a significant positive association, and percent age 65 and older and days since peak incidence exhibited a significant negative association with COVID-19 counts. PM 2.5 and percent low income remained positively associated with COVID-19 counts, and prevalence of overweight and obesity remained negatively associated, but the associations were no longer significant. Additional multivariate models are summarized in Table S4 . The association of PM 2.5 with COVID-19 counts was not sensitive to the sequential addition to the best-fitting model of population density, percent Black population, or days since first case; the PM 2.5 coefficient changed by <10% and additional covariates were not significantly associated with COVID-19 counts. However, the association of PM 2.5 with COVID-19 counts was sensitive to the removal of minimum temperature from the best-fitting model (IRR increased to 1.10, 95% CI 1.00-1.21, p<0.05). Mapping of model residuals and Moran's I did not indicate presence of residual spatial autocorrelation.

In sensitivity analyses (Table 3) , the association of PM 2.5 with COVID-19 incidence was unchanged when Montréal was excluded or when province was included as a random effect rather than J o u r n a l P r e -p r o o f a fixed effect. In models excluding British Columbia, Saskatchewan and Nova Scotia, excluding health regions with less than median population density, and restricted to Ontario and Quebec, the association of PM 2.5 with COVID-19 incidence was larger in magnitude than the national analysis, and statistically significant. These analyses accounted for 85-95% of cases and 62-82% of population. Notably, in the subset of health regions with population density greater than the median, greenness exhibited a significant negative association with case counts on its own, but was no longer significant in a multivariate model which had a larger AIC than the final model. Results were also similar when based on data from earlier in the course of the pandemic when there were 35,986 cases (April 19).

We found that after controlling for province, temperature, health and demographic characteristics, and time since peak incidence by health region, long-term PM 2.5 exposure exhibited a non-statistically significant positive association with COVID-19 incidence among 111 Canadian health regions. This association was larger in magnitude and statistically significant in models excluding provinces which reported cases only for aggregated health regions, excluding health regions with less than median population density, and restricted to Ontario and Quebec, in each case accounting for a majority of COVID-19 cases. The association we observed parallels that reported in an analysis of PM 2.5 and COVID-19 mortality in approximately 3,000 US counties (Wu et al., 2020) . While this is striking, especially given the lower COVID-19 incidence and narrower range of PM 2.5 exposure in Canada compared to the US, findings should be interpreted with caution; the analyses examined different outcomes (incidence vs. mortality) and associations with PM 2.5 in both studies could result from similar biases related to ecologic analyses. PM 2.5 exposures and other risk factors for COVID-19 incidence can vary over a smaller scale than health region, particularly for large health regions. Greenness in particular varies at a much smaller scale than health region, even in predominantly urban health regions. Our J o u r n a l P r e -p r o o f analysis entailed coarser spatial resolution and a much smaller number of health regions than US counties. However, analyses excluding provinces that reported cases only for aggregated health regions, excluding health regions with less than median population density, and restricted to Quebec and Ontario resulted in a larger magnitude and statistically significant association, which could reflect reduced exposure measurement error for these generally smaller, more densely populated health regions. The US study also found that the mortality risk ratio was elevated in an analysis restricted to urban areas, but it was no longer significant (Wu et al., 2020) . Owing to the smaller sample, we were also constrained in the number of covariates that could be included in multi-variate models. Still, we did account for unmeasured confounders by province by specifying it as a fixed or random effect, and found that the magnitude of the association was larger in more highly urban and more highly affected areas.

Mobility data indicated that effects of social distancing restrictions were generally similar among provinces. We employed publicly available laboratory confirmed case counts, thus our results are not generalizable to milder or asymptomatic cases which did not undergo laboratory confirmation, and case criteria may differ by province and change over time. Data on COVID-19 testing were only available at the provincial level. Although criteria for testing are standardized, access may differ by health region, which could increase the observed incidence in more highly urban areas with potentially better access and generally higher PM 2.5 concentrations, potentially resulting in a spurious association between PM 2.5 exposure and COVID-19 incidence. While COVID-19-related hospitalizations and deaths may be less affected by artefacts introduced by differential eligibility for and/or access to testing by health region, the number of hospitalizations and deaths in Canada is relatively small. At the health region level, a large number of small or zero counts would be expected, which introduces other limitations in the analysis. In any case, hospitalization data were not available nationally at the health region (or smaller) level. Our observation that deaths were very strongly correlated with case counts by health region provides evidence that differences in case numbers by health region were not driven by differences in testing, J o u r n a l P r e -p r o o f which would be expected to introduce more scatter in the relationship between deaths and cases.

Finally, associations observed based solely on area level measures may not exist at the individual level (ecological fallacy).

The magnitude of the association both in our study and the US study was several fold larger per unit PM 2.5 than hazard ratios typically observed in cohort studies of mortality (Crouse et al., 2015) . It is premature to speculate on pathophysiological mechanistic explanations for this observation. While there is existing evidence that PM 2.5 increases the risk of respiratory infections, including the SARS virus, the most likely explanation for the large magnitude association is residual confounding by unmeasured factors.

In addition to the Wu et al. (2020) study, another American ecological study based on approximately 3,000 counties found that nitrogen dioxide (NO 2 )was significantly associated with COVID-19 mortality and case-fatality rates, while PM 2.5 exhibited a marginally significant association with mortality (Liang et al., 2020) . A significant positive association of PM 2.5 with COVID-19 incidence and hospital admissions was also reported in the Netherlands (Andree, 2020) . Ogen et al. (2020) reported co-location of high NO 2 concentrations estimated from remote sensing and high COVID-19 mortality counts in northern Italy and Madrid. However, the analysis did not account for the underlying population at risk, population density, timing of onset of cases or introduction of control measures, or sociodemographic or health characteristics. In a time-series study in 120 Chinese cities, significant positive associations were observed between short term (two week) exposure to PM 2.5 , NO 2 and ozone and COVID-19 incidence (Zhu et al., 2020) .

We also found that percent Black population was positively associated with COVID-19 incidence (although it was not included in the final multivariate model). This is consistent with reports from the US of disproportionately high COVID-19 incidence and mortality among Black Americans (Dyer, 2020; Yancy, 2020) . Percent of population less than LICO was also positively associated with COVID-19 incidence.

Those with lower income have more barriers to self isolating and social distancing. Inequality in access to care affecting access to COVID-19 testing and outcome following infection would be expected to be less in Canada with universal healthcare. Negative associations of percent of population age 65 and older, prevalence of asthma, COPD, hypertension, overweight and obesity with COVID-19 incidence may reflect a higher prevalence in predominantly rural health regions with few cases.

Subsequent studies should examine effects of both acute and chronic exposure to PM 2.5 , and air pollution more generally, on COVID-19 morbidity and mortality using more highly spatially resolved area level data, as well as individual level data to determine whether results from ecological analysis are borne out.

We report an analysis of COVID-19 incidence by Canadian health region. After controlling for temperature, demographic and health characteristics and days from peak incidence by health region, long-term PM 2.5 exposure exhibited a positive association with COVID-19 incidence, paralleling results of a recent American analysis of mortality. The association was larger in magnitude and statistically significant in more highly affected health regions and those with potentially less exposure measurement error. While parallels with US findings are noteworthy, and generate hypotheses for further testing, results should be interpreted with caution and require further examination using study designs less prone to bias.

J o u r n a l P r e -p r o o f Figure S2 . Correlogram. Size and colour intensity of markers reflect magnitude and statistical significance of Spearman correlation. Correlations with p>0.01 are blank. Figure generated using Hmisc (Harrell et al., 2020) and Corrplot (Wei and Simko, 2017) packages for R. 

Population-Based Estimates of Chronic Conditions Affecting Risk for Complications from Coronavirus Disease, United States

A new look at the statistical model identification

Incidence of COVID-19 and Connections with Air Pollution Exposure: Evidence from the Netherlands (preprint)

Epidemiological time series studies of PM2.5 and daily mortality and hospital admissions: a systematic review and metaanalysis

Fitting Linear Mixed-Effects Models Using lme4

Open access epidemiologic data and an interactive dashboard to monitor the COVID-19 outbreak in Canada

rgdal: Bindings for the "Geospatial

Comparing implementations of global and local indicators of spatial association 27

COVID-19 Emergency Measures Tracker

GLMM FAQ [WWW Document

Fifteen-year global time series of satellite-derived fine particulate matter

GISTools: Some further GIS capabilities for R

Influence of meteorological factors and air pollution on the outbreak of severe acute respiratory syndrome

CanMap Postal Code

Postal Code User Guide. Canadian Urban Environmental Health Research Consortium

Ambient PM2.5, O3, and NO2 exposures and associations with mortality over 16 years of follow-up in the canadian census health and environment cohort (CanCHEC)

Air pollution and case fatality of SARS in the People's Republic of China: an ecologic study

MOD13Q1 MODIS/Terra Vegetation Indices 16-Day L3 Global 250m SIN Grid V006 [Data set

Effects of air pollutants on the transmission and severity of respiratory viral infections

Covid-19: Black people and other minorities are hardest hit in US

Global, regional, and national comparative risk assessment of 84 behavioural, environmental and occupational, and metabolic risks or clusters of risks for 195 countries and territories, 1990-2017: a systematic analysis for the Global Burden of Disease Study

Google COVID-19 Community Mobility Reports

Google Earth Engine: Planetary-scale geospatial analysis for everyone. Big Remote. Sensed Data Tools Appl

Hmisc: Harrell Miscellaneous

raster: Geographic Data Analysis and Modeling. R package version 3

Relationship between ambient air pollution and daily mortality of SARS in Beijing

Urban Air Pollution May Enhance COVID-19 Case-Fatality and Mortality Rates in the United States

Ambient particulate air pollution and acute lower respiratory infections: A systematic review and implications for estimating the global burden of disease

Long-term exposure to ambient air pollution and risk of hospitalization with community-acquired pneumonia in older adults

Short-term association between ambient air pollution and pneumonia in children: A systematic review and meta-analysis of time-series and case-crossover studies

Assessing nitrogen dioxide (NO2) levels as a contributing factor to coronavirus (COVID-19) fatality

Examining the Shape of the Association between Low Levels of Fine Particulate Matter and Mortality across Three Cycles of the Canadian Census Health and Environment Cohort

Why Quebec's coronavirus cases have skyrocketed

R: A language and environment for statistical computing

Convert Statistical Analysis Objects into Tidy Tibbles

Table 13-10-0113-01 Health characteristics, two-year period estimates [WWW Document

Health Regions: Boundaries and Correspondence with Census Geography

Peer Groups

Housing, Aboriginal Peoples, Education, Labour, Journey to Work, Mobility and Migration, and Language of Work for Canada, Provinces and Territories, and Health Regions

A joint ERS/ATS policy statement: what constitutes an adverse health effect of air pollution? An analytical framework

Regional Estimates of Chemical Composition of Fine Particulate Matter Using a Combined Geoscience-Statistical Method with Information from Satellites, Models, and Monitors

Water balance model for Canadian Urban Environment (CANUE) and Health Research Consortium. Please contact CANUE (info@canue.ca) for additional information

Ottawa to intervene in lagging regions amid COVID-19 test backlog

Visualization of a Correlation Matrix

ggplot2: Elegant Graphics for Data Analysis

Exposure to air pollution and COVID-19 mortality in the United States

COVID-19 and African Americans

The Impact of PM2.5 on the Host Defense of Respiratory System

Association between short-term exposure to air pollution and COVID-19 infection: Evidence from China