Model framework
This study combines multiple datasets and models, as presented in Extended Data Fig. 1, to estimate the contribution of the 2023 Canadian wildfires to global PM2.5 exposure and health impacts under a near-real-time framework (http://tapdata.org.cn). We also analyse two additional years, 2021 and 2017, which were reported as the years with the second- and third-largest wildfire emissions in Canada since 2000 (Extended Data Fig. 2b), for comparison. We first used the GEOS-Chem chemical transport model36 at a spatial resolution of 2° × 2.5° and 3 near-real-time global fire emission inventories, that is, the Global Fire Emissions Database version 4 with small fires (GFEDv4.1s)37,65, the Quick Fire Emissions Dataset version 2.5 (QFEDv2.5r1)38 and the Global Fire Assimilation System version 1.2 (GFASv1.2)39, to simulate the global daily PM2.5 concentrations and the fractional shares in total PM2.5 concentrations contributed by wildfire emissions using a zero-out approach. Second, to improve the spatial resolution and accuracy of the global daily PM2.5 estimation, we developed a machine-learning-based PM2.5 retrieval model that combines data from multiple sources, including ground-monitoring measurements, satellite retrievals, reanalysis data and GEOS-Chem simulations, to estimate the global daily PM2.5 concentrations at a spatial resolution of 0.1° × 0.1°. The PM2.5 retrieval model was trained using GEOS-Chem simulations with the GFED, the QFED and the GFAS as a priori fire emissions, respectively, and three sets of global daily PM2.5 estimates were derived. Then the retrieved total PM2.5 concentrations based on the GFED, the QFED and the GFAS were multiplied by previously simulated fractional contributions with corresponding fire emissions, to obtain the PM2.5 exposure attributable to wildfires. The performance of the three fire emission inventories in estimating fire-related PM2.5 exposures was evaluated, and the GFED-based results were selected for presentation and further analysis to facilitate comparisons with other recent studies that use similar approaches9,10,11,44. To investigate the transboundary impact of Canadian wildfires globally and in North America, we further quantified the contributions of five regions’ wildfires (that is, Eastern Canada, Western Canada, Eastern USA, Western USA and other global regions; Extended Data Fig. 3b) to PM2.5 concentrations by conducting additional zero-out scenarios using the GFED emission inventory. Finally, we assessed the acute and chronic deaths attributable to PM2.5 pollution from Canadian wildfires using previously established exposure–response functions15,40. Further details of each analytical step are provided below.
Global fire emissions
A fire emission inventory is an essential input dataset for our analyses. We separately examined the impacts of fire emissions using three near-real-time fire emission inventories available for the year 2023, the GFED37,65, the QFED38 and the GFAS39. The model results and fire-related PM2.5 estimates using each of the three fire emissions datasets are discussed in ‘Model evaluation’.
GFEDv4.1s
The GFED inventory was developed for use in large-scale modelling studies. The latest GFEDv4.1s used in this study is archived at https://surfdrive.surf.nl/files/index.php/s/5y7TdE6ufwpkAW1. It is based on 500-m Moderate Resolution Imaging Spectroradiometer (MODIS) burned area maps66 supplemented with MODIS active fires converted to burned area67. After 2016, emissions are derived from MODIS active fires scaled to emissions based on the 2001–2016 period when both datasets overlapped. Emission factors, mostly from ref. 68, are used to convert fire carbon emissions to trace gases and aerosols. The emissions, including carbon, dry matter, CO2, CO, NOx, organic carbon, black carbon, PM2.5, total particulate matter and SO2 among others, are available from 1997 to 2023 at 0.25° × 0.25° globally37,65.
QFEDv2.5r1
The QFED inventory was developed by the National Aeronautics and Space Administration (NASA) and serves as the standard fire emissions in the GEOS data assimilation system and the Modern-Era Retrospective analysis for Research and Applications, version 2 (MERRA-2) reanalysis data products69. In this study, QFEDv2.5r1 is used, which is available at https://portal.nccs.nasa.gov/datashare/gmao/qfed/. On the basis of a top-down approach, QFED obtains the fire radiative power (FRP) and location from satellite observations from MODIS Level 2 fire products and MODIS Geolocation products and calculates the open combustion of non-fossilized vegetative or organic fuel38. It provides high spatiotemporal resolution and near-real-time global biomass burning emissions, including the 0.1° × 0.1° daily emissions of black carbon, organic carbon, SO2, CO, CO2, PM2.5, NH3, NOx and so on, from 2001 to present.
GFASv1.2
The GFASv1.2 data are used for the Copernicus Atmosphere Monitoring Service (CAMS) global atmospheric composition and regional air-quality forecasts, which can be found at https://ads.atmosphere.copernicus.eu/cdsapp#!/dataset/cams-global-fire-emissions-gfas?tab=overview/. Fire emissions are calculated based on FRP measurements from two MODIS instruments onboard NASA’s Terra and Aqua satellite that are first converted to estimates of the dry matter consumed by fire and then to emissions using biome-specific emission factors. The GFAS provides daily averaged biomass burning and vegetation fire emissions for 40 pyrogenic species (aerosols, reactive gases and greenhouse gases) from 2003 to the present, with a spatial resolution of 0.1° × 0.1° (ref. 39).
GEOS-Chem simulation
Using each of the three fire emission inventories (that is, the GFED, the QFED and the GFAS), we estimated the global PM2.5 concentrations using the three-dimensional global chemical transport model GEOS-Chem36. These concentrations are an important input to our PM2.5 retrieval model (‘Retrieval of global PM2.5 based on multi-source data fusion’) and are used to quantify the fractional contributions of wildfire emissions to total PM2.5 concentrations (‘Estimation of fire-specific PM2.5 exposure’ and ‘Fire source attribution’). GEOS-Chem has been used by numerous previous studies to simulate smoke pollution from wildfires8,43,70,71,72,73,74,75,76.
GEOS-Chem v.14.0.1 (https://zenodo.org/records/7271974/) is used in this study. The near-real-time meteorological data from the Goddard Earth Observation System-Forward Processing (GEOS-FP)77 of the NASA Global Modeling and Assimilation Office (GMAO) was used to drive the GEOS-Chem model. The GEOS-FP data span the time period from 2011 to the present, with a native resolution of 0.25° × 0.3125° and 72 vertical levels. We reduce the vertical levels to 47 and the spatial resolution to 2.0° × 2.5° to support the global chemical transport simulations. GEOS-Chem uses standard full chemistry with detailed oxidant–aerosol chemistry. Sulfate–nitrate–ammonium aerosol thermodynamics are computed with ISORROPIAv2.278.
All emissions in GEOS-Chem are configured by HEMCO (Harmonized Emissions Component) 3.079, to combine and regrid the different emissions. The global anthropogenic emissions (including shipping) of NOx, SO2, CO, NH3, black carbon, organic carbon and volatile organic compounds are provided by the Community Emissions Data System (CEDS) v2 inventory (https://data.pnnl.gov/dataset/CEDS-4-21-21/)80. Aircraft emissions are from the Aviation Emissions Inventory Code (AEIC)81 inventory. For global fire emissions, GFEDv4.1s, QFEDv2.5r1 and GFASv1.2 are used respectively. Dust82, sea salt83, lighting NOx (ref. 84), soil NOx (ref. 85) and biogenic volatile organic compounds (MEGAN v2.186) are calculated online in HEMCO. We use the non-local scheme implemented in ref. 87 for the boundary-layer mixing in GEOS-Chem, which emits all emissions into the atmospheric boundary layer, including wildfire emissions.
Retrieval of global PM2.5 based on multi-source data fusion
We estimate the global daily PM2.5 exposures from all sources at a 0.1° × 0.1° horizon resolution using a multilayer machine-learning retrieval model that fuses data from ground-monitoring measurements, satellite retrievals, GEOS-Chem model simulations, meteorological fields, reanalysis data and population distribution. The structure of the PM2.5 retrieval model is illustrated in Supplementary Fig. 1. We separately train our retrieval model with GEOS-Chem simulations based on the GFED, the QFED and the GFAS emissions (but with the same model structure) to derive three sets of global PM2.5 concentration estimates. Further details about the data sources and the model structure are provided below.
Data sources
Ground measurements
We collected PM2.5 surface monitoring data from different global regions as model input. We obtained surface PM2.5 measurements in Canada from Environmental Canada (https://data-donnees.az.ec.gc.ca/data/air/monitor/national-air-pollution-surveillance-naps-program/Data-Donnees/). We obtained PM2.5 measurements in the USA from the US Environmental Protection Agency (US EPA) AirNow (https://www.epa.gov/outdoor-air-quality-data/download-daily-data) and from the Interagency Monitoring of Protected Visual Environments (IMPROVE) (https://views.cira.colostate.edu/fed/QueryWizard/Default.aspx). To improve the representation of air pollution, we also collected PM2.5 monitoring data from: the AirFire programme of the US Forest Service (https://info.airfire.org/airmonitor-package) for the USA; the European Air Quality Portal (https://eeadmz1-cws-wp-air02.azurewebsites.net/) for Europe; the China National Environmental Monitoring Center (CNEMC; http://www.cnemc.cn/) for China; and the OpenAQ (https://openaq.org/) for other regions around the world. Hourly measurements were averaged as daily records and only daily records generated from at least 16 hourly data points were included. In summary, we collected approximately 453,000 valid daily records from about 1,610 monitors in North America, approximately 567,000 valid daily records from about 2,010 stations in Europe, approximately 612,000 valid daily records from about 1,720 stations in China, as well as approximately 68,000 valid daily records from about 850 stations in other regions for the year 2023.
Aerosol optical depth
Satellite aerosol optical depth (AOD) retrievals were extracted from the MODIS Level 2 aerosol products (MOD04 and MYD04) at a 0.1° spatial resolution88. To improve the data coverage and better reflect the aerosol loading during the day, we first fused the AOD retrievals from the Dark Target algorithm and the Deep Blue algorithm with daily linear regressions, and then fused the AOD from the Aqua and Terra satellites with daily linear regressions. As considerable gaps in the AOD data still existed after this data fusion, we used the CAMS modelling and reanalysis data (https://ads.atmosphere.copernicus.eu)89 with complete coverage to provide information on the spatial distribution of aerosols. CAMS parameters, including the AOD at 550 nm, black carbon AOD, organic carbon AOD, wildfire combustion rate, FRP and total column carbon monoxide, were adopted in the model. The MODIS AOD and various CAMS AODs were treated as separate predictors in the retrieval model.
GEOS-Chem simulations
The PM2.5 simulations from GEOS-Chem were used in this model. As three fire emission inventories were used in GEOS-Chem to simulate surface PM2.5 concentrations, we constructed three model training datasets with different GEOS-Chem simulations using the GFED, the QFED and the GFAS, respectively.
Other ancillary data
Smoke plume information in North America was collected from the Hazard Mapping System (HMS; https://www.ospo.noaa.gov/Products/land/hms.html#about)90, provided by the National Oceanic and Atmospheric Administration/National Environmental Satellite, Data, and Information Service (NOAA/NESDIS). The density-assigned and time-marked plumes polygons were manually generated from GOES-16 and GOES-17 ABI true-colour imagery. We assigned the time-specific plume data to daily plume density and included it in the model to provide valuable information on smoke plume distributions in North America. Meteorological fields, including daily average air temperature at 2 m, specific humidity at 2 m, relative humidity, surface pressure, boundary-layer height, total latent energy flux, evaporation from turbulence, U and V wind components at 10 m, and total precipitation, were extracted from GEOS-FP reanalysis data at a spatial resolution of 0.25° × 0.3125° and downscaled to the 0.1° modelling grid by the inverse distance weighting algorithm. The gridded population distribution data for 2020 were obtained from WorldPop (https://www.worldpop.org/) at a resolution of 30 arcseconds91,92 and we assumed a constant population distribution in 2017, 2021 and 2023. We constrained the gridded WorldPop population data with national total population for each year from the United Nations (https://population.un.org/wpp/)93, US Census Bureau, Population Division (https://www.census.gov/data/tables/time-series/demo/popest/2020s-state-total.html)94 and Statistics Canada (https://www150.statcan.gc.ca/)16.
Model structure
The total PM2.5 retrieval model was designed with a three-layer random forest structure following our previous work95 (Supplementary Fig. 1). The first-layer model predicts the high-pollution index, which is a binary variable indicating whether the station-day concentration is higher than the mean plus two standard deviations of PM2.5 concentrations of the corresponding station and month. In our previous work96, we found that a two-layer model including the high-pollution indicator and the Synthetic Minority Over-sampling Technique (SMOTE) resampling algorithm can correct the low- bias from the unbalanced training sample of high-pollution events. Following the approach, we applied the SMOTE resampling algorithm to increase the representation of high-pollution events in the model training samples. The second-layer model uses the prediction of high-pollution index as a predictor to predict the total PM2.5 concentrations. Previous studies have shown that the model trained with the residual can correct the systematic bias and improve the model prediction accuracy97,98. The third-layer model then predicts the residual between PM2.5 predictions and measurements and was trained by with-fire samples and no-fire samples separately, to highlight the potential differences in PM2.5 characteristics during fire events. Here the with-fire samples were defined as samples with CAMS combustion rate >0 or the HMS plume density >0 (refs. 8,9). The final PM2.5 concentration estimation is the prediction from the second-layer model plus the prediction from the third-layer model. The model was trained separately with data for 2023, 2021 and 2017 as well as with GEOS-Chem simulations with the GFED, the QFED and the GFAS emissions.
Estimation of fire-specific PM2.5 exposure
PM2.5 exposure attributable to wildfire emissions was then estimated using three fire emissions, respectively. In detail, the wildfire-related PM2.5 was quantified by multiplying the GEOS-Chem simulated fire contributions to total PM2.5 by the retrieved all-source PM2.5 concentrations. Supplementary Table 1 summarizes the GEOS-Chem simulations used in this study. We conducted the GEOS-Chem simulations with the GFED, the QFED and the GFAS inventories separately (that is, ‘base’ in Supplementary Table 1) as well as the no-fire GEOS-Chem simulation that turned off global fire emissions (that is, ‘nofire’ in Supplementary Table 1). All other emissions mentioned in ‘GEOS-Chem simulation’ are the same in the base and nofire simulations. Then the fraction of PM2.5 concentrations attributable to wildfires was calculated by equation (1) on a 2° × 2.5° grid, as determined by the GEOS-Chem model. We constructed three sets of wildfire-fraction data from the GEOS-Chem simulations driven by the three fire emissions. Likewise, the anthropogenic-related PM2.5 was calculated using equation (2) with similar model runs (turning off anthropogenic emissions, that is, ‘noanthro’ in Supplementary Table 1). The contributions from sources other than wildfires and anthropogenic activities were then obtained by subtracting their contributions from the total as in equation (3):
$${F}_{{\rm{f}}{\rm{i}}{\rm{r}}{\rm{e}},k}={{\rm{G}}{\rm{C}}}_{{\rm{f}}{\rm{i}}{\rm{r}}{\rm{e}},k}/{{\rm{G}}{\rm{C}}}_{{\rm{b}}{\rm{a}}{\rm{s}}{\rm{e}},k}=({{\rm{G}}{\rm{C}}}_{{\rm{b}}{\rm{a}}{\rm{s}}{\rm{e}},k}-{{\rm{G}}{\rm{C}}}_{{\rm{n}}{\rm{o}}{\rm{f}}{\rm{i}}{\rm{r}}{\rm{e}}})/{{\rm{G}}{\rm{C}}}_{{\rm{b}}{\rm{a}}{\rm{s}}{\rm{e}},k}$$
(1)
$${F}_{{\rm{a}}{\rm{n}}{\rm{t}}{\rm{h}}{\rm{r}}{\rm{o}},k}={{\rm{G}}{\rm{C}}}_{{\rm{a}}{\rm{n}}{\rm{t}}{\rm{h}}{\rm{r}}{\rm{o}},k}/{{\rm{G}}{\rm{C}}}_{{\rm{b}}{\rm{a}}{\rm{s}}{\rm{e}},k}=({{\rm{G}}{\rm{C}}}_{{\rm{b}}{\rm{a}}{\rm{s}}{\rm{e}},k}-{{\rm{G}}{\rm{C}}}_{{\rm{n}}{\rm{o}}{\rm{a}}{\rm{n}}{\rm{t}}{\rm{h}}{\rm{r}}{\rm{o}},k})/{{\rm{G}}{\rm{C}}}_{{\rm{b}}{\rm{a}}{\rm{s}}{\rm{e}},k}$$
(2)
$${F}_{{\rm{other}},k}=1-{F}_{{\rm{fire}},k}-{F}_{{\rm{anthro}},k}$$
(3)
where the subscript k represents the three fire emissions, that is, GFED, QFED and GFAS. GCbase,k and GCnoanthro,k represent the GEOS-Chem-simulated PM2.5 concentrations from the base and noanthro scenarios using fire emissions k, respectively. GCnofire represents the simulation with fire emissions turned off. Ffire,k, Fanthro,k and Fother,k represent the fractional contribution of wildfires, anthropogenic and other emissions to PM2.5 estimated from fire emissions k, respectively.
Then we spatially match the PM2.5 fractions from GEOS-Chem simulations (2° × 2.5°) with the 0.1° PM2.5 retrievals through bilinear interpolation. The PM2.5 retrievals were multiplied by the corresponding fractions to get the fire-related, anthropogenic-related and other-source-related PM2.5, as shown in equations (4)–(6):
$${C}_{{\rm{fire}},k}={C}_{{\rm{PM}},k}\times {F}_{{\rm{fire}},k}={C}_{{\rm{PM}},k}\times ({{\rm{GC}}}_{{\rm{fire}},k}/{{\rm{GC}}}_{{\rm{base}},k})$$
(4)
$${C}_{{\rm{anthro}},k}={C}_{{\rm{PM}},k}\times {F}_{{\rm{anthro}},k}={C}_{{\rm{PM}},k}\times ({{\rm{GC}}}_{{\rm{anthro}},k}/{{\rm{GC}}}_{{\rm{base}},k})$$
(5)
$${C}_{{\rm{other}},k}={C}_{{\rm{PM}},k}\times {F}_{{\rm{other}},k}$$
(6)
where CPM,k represents the total PM2.5 concentrations estimated from the machine-learning-based model using fire emissions k (k = GFED, QFED, GFAS) in ‘Retrieval of global PM2.5 based on multi-source data fusion’, Cfire,k, Canthro,k and Cother,k represent the fire-, anthropogenic- and other-source-related PM2.5 concentrations based on fire emissions k.
Model evaluation
Our models were fully evaluated at each step of the construction of fire-related PM2.5 concentrations (Extended Data Fig. 1): the performance of the GEOS-Chem PM2.5 simulations, the performance of the PM2.5 retrieval model and the performance of PM2.5 retrievals during fire events. Models based on the three fire emissions (GFED, QFED and GFAS) were compared in all the evaluations to understand the impacts of fire emissions on model performance. Total and fire-related PM2.5 concentrations estimated using the three fire emissions were compared with each other. We also quantified the wildfire-related PM2.5 estimates in previous years (2017 and 2021) and compared them with other studies as an additional evaluation of our methods.
Evaluation of GEOS-Chem PM2.5 simulations
The annual average GEOS-Chem PM2.5 concentrations driven by three fire emission inventories (GFED, QFED and GFAS) were evaluated against PM2.5 ground observations in Canada and the USA (Supplementary Fig. 15). In Canada, the modelled PM2.5 correlated reasonably well with ground observations, with R ranging between 0.41 and 0.69 and normalized model bias (NMB) between 0.01 and 0.72 for the three fire emissions. We noticed one outlier with unrealistically high PM2.5 simulations in Canada resulted from high fire emission estimates; therefore, we also reported evaluation statistics without this data point to avoid the effects of an outlier (Supplementary Fig. 15). After removing the outlier, R between modelled and observed PM2.5 increased from 0.41–0.69 to 0.66–0.82 and NMB was reduced from 0.01–0.72 to −0.04–0.30 in Canada. In the USA, simulations based on the three fire emissions had comparable performance against ground observations, with R between 0.46 and 0.47. As the only differences in these simulations are the underlying fire emission inventories, the differences in model performance can be solely attributed to the inventories. Compared with the GEOS-Chem simulation with the GFED inventory, using the QFED inventory helps to reduce both the overestimates of PM2.5 concentration in Canada (NMB from 0.30 to 0.16) and the underestimates in the USA (NMB from −0.22 to −0.06), whereas using the GFAS inventory reduces the overestimates of PM2.5 concentration in Canada (NMB from 0.30 to −0.04) but gives similar results to GFED in the USA (NMB from −0.22 to −0.27). By comparing estimates of fire-related PM2.5 exposure from different fire emission inventories in this way, we can evaluate the effects of GEOS-Chem model performance on the results (‘Comparisons of fire-related PM2.5 among three fire emissions and with previous studies’).
Evaluation of PM2.5 retrieval model
The PM2.5 retrieval model was evaluated by both station-based and sample-based 20-fold cross-validation. The model training dataset was randomly divided into 20 equal folds according to air-quality monitoring stations and station-day observations, separately. Then the model was trained on 19 of these folds and tested on the remaining fold. This process was repeated 20 times until each fold of the data was used for testing once. The model performance was quantified by the comparisons between cross-validation predictions and ground measurements at the daily, monthly and yearly levels to reflect model uncertainties at different temporal scales (Extended Data Fig. 8 and Supplementary Fig. 13). In addition, to highlight the model’s ability in retrieving daily variations in PM2.5 when controlling the seasonal and spatial variations, we included month-intercept (fix-month R2) as well as month and station intercepts (fix-month-and-station R2) when computing R2 (within R2) of the 20-fold cross-validation predictions following a previous study10 (Supplementary Table 5). We also calculated the station-specific R2 to show the variations in model performance in space. Evaluation results of the retrieval models trained with GEOS-Chem simulations using the three fire emissions in 2023 are listed in Extended Data Fig. 8, Supplementary Fig. 13 and Supplementary Table 5. As the HMS data are available for only North America, we assessed the impact of including HMS data on global PM2.5 retrievals. As shown in Supplementary Fig. 17, inclusion of the HMS data in the model led to substantial differences in PM2.5 retrievals over North America whereas it had minor impacts in other regions. Therefore, we incorporated the HMS data as a predictor to improve the accuracy of wildfire-related PM2.5 estimates in North America.
Globally, the retrieval model characterizes variations in PM2.5 well, and models based on the GFED, the QFED and the GFAS performed similarly well under both the station-based and sample-based 20-fold cross-validation. The station-based 20-fold cross-validation R2 ranged between 0.84 and 0.85 (RMSE between 8.55 μg m−3 and 8.62 μg m−3) at the daily scale, all equal to 0.88 (RMSE between 5.80 μg m−3 and 5.95 μg m−3) at the monthly scale, and ranged between 0.87 and 0.88 (RMSE between 4.80 μg m−3 and 5.00 μg m−3), in the year 2023. The sample-based cross-validation results are comparable to the station-based cross-validation results, indicating robust model performance in regions with limited observations. All models trained with the three fire emission inventories showed only a slight decrease in R2 when controlling the seasonal and spatial variations (that is, fix-month R2 between 0.86 and 0.88, and fix-month-and-station R2 between 0.80 and 0.82), indicating the model’s ability to capture daily variations in PM2.5 (Supplementary Table 5). Spatially, the global median station-specific R2 was 0.79, 0.79 and 0.79 with 90% of station-specific R2 above 0.37, 0.36 and 0.37 for models based on the GFED, the QFED and the GFAS, respectively. Our model performances are comparable to previous studies developing a global PM2.5 retrieval model, in that the cross-validation R2 of daily retrievals is around 0.91 and the RMSE ranges between 8.4 μg m−3 and 9.2 μg m−3 (refs. 10,41).
Regionally we found a lower R2 in Canada, which is mainly caused by several outliers in Canada resulting from unrealistically high GEOS-Chem simulations mentioned in ‘Evaluation of GEOS-Chem PM2.5 simulations’. To avoid the influence of these occasional outliers on model evaluation, we reported the model performance after removing data points outside 2σ of the Cook’s distance of linear regression (Supplementary Table 5). Models based on the three fire emissions performed comparably well after removing the outliers. As these outliers of the retrieved PM2.5 concentration occurred in a few days and in remote regions with sparse population, they do not considerably affect our exposure assessment and health burden quantification.
Evaluation of PM2.5 retrievals during fire events
To assess the model’s performance in capturing fire-related PM2.5 variations, we further evaluated the model performance during fire events (Extended Data Fig. 9 and Supplementary Table 4)—when the PM2.5 is dominated by wildfires. We first identified several major fire events and then compared the station-day observations during these fire events with the sample-based 20-fold cross-validation predictions. Thus, the station-day observations selected for evaluations were excluded from the model training to reveal the model performance in regions and periods without observations. We identified fire events with a similar protocol reported by an earlier work10, mainly according to the variations in PM2.5 observations. A station-day is labelled as affected by wildfires when (1) it is during one of the manually identified fire events that showed a substantial increase in national daily average PM2.5 monitoring time series data in Canada and the USA, separately. (2) The daily average PM2.5 concentration is higher than the median of all the station-day records during this fire event. Here we used the median as cut-off number rather than selecting one station with the largest increase because one station’s data were not sufficient to support the validation. (3) The PM2.5 concentration is higher than twice the background PM2.5 concentration before the fire season (first 2 months in 2023) at the corresponding station. (4) The PM2.5 concentration is higher than 15 μg m−3, the WHO air-quality guidelines level. In total, 11 and 9 fire events were identified in Canada and the USA, respectively. The events lasted between 3 days and 23 days. The national fire-related PM2.5 exposure during these fire events accounted for 83% and 81% of the annual fire-related PM2.5 exposure in Canada and the USA, respectively, indicating that most significant fire events were identified by this method.
In North America, the three retrieval models with different fire emissions correctly reflected PM2.5 variations during the fire events, with similar cross-validation R2 ranging between 0.78 and 0.80 in the USA, but moderately different R2 between 0.59 and 0.75 in Canada (Supplementary Table 4). We also noticed that as previously reported41,99, our model still slightly underestimated PM2.5 levels during extreme fire events with NMB ranging between −0.14 and −0.09 and between −0.13 and −0.10 for Canada and the USA, respectively.
Comparisons of fire-related PM2.5 among the three fire emissions and with previous studies
To further understand the impacts of different fire emissions on total and fire-related PM2.5 estimates, we compared the spatial and temporal distributions of our results among the three fire emissions (Supplementary Table 6). The spatial correlation was calculated using 3-year averaged gridded PM2.5 data among models based on the 3 fire emissions, whereas the temporal correlation was calculated using daily population-weighted mean PM2.5 concentrations for 3 years among models based on the 3 fire emissions. Some previous studies9,10,11,44 have investigated the impact of wildfires on PM2.5 exposure other than the 2023 Canadian extreme wildfires. We also compared our fire-related PM2.5 estimates with those previous studies as additional evaluations.
For spatial comparisons, the all-source PM2.5 estimates based on the three fire emission inventories showed similar spatial distributions to high correlations globally (all pairwise Pearson correlation coefficients r were 0.99). The spatial distributions of fire-related PM2.5 were also highly correlated between different emissions at the global scale, with pairwise r ranging between 0.88 and 0.95. The spatial patterns between the QFED-based and the GFAS-based estimates were more similar in specific regions than the GFED-base results, as fire emissions from the QFED and the GFAS were estimated using a similar approach based on FRP from the MODIS instrument. For temporal comparisons, the population-weighted daily mean all-source PM2.5 estimates showed high correlations at both the global scale (all pairwise Pearson correlation coefficients r were 0.99) and the regional scale. The fire-related PM2.5 estimates had some temporal differences globally among the 3 inventories, but showed even higher correlations in Canada and the USA, with pairwise r of fire-related PM2.5 ranging between 0.94 and 0.99 in Canada and between 0.94 and 0.97 in the USA.
Supplementary Fig. 14 shows the comparisons of estimated population-weighted fire-related PM2.5 in Canada and the USA based on the three inventories, as well as their comparisons with previous studies9,10,11,44. The fire-related PM2.5 in Canada and the USA based on the three inventories showed consistent increasing trends in 2017, 2021 and 2023; however, the magnitudes of the estimated fire-related PM2.5 varied between inventories, with the GFAS showing the lowest and the QFED showing the highest. In Canada, fire-related annual PM2.5 exposure estimated with the 3 inventories showed good agreement (3.75–4.49 μg m−3, relative difference within 20%). Given their similar fire-related PM2.5 exposure estimates in Canada, the choice of fire inventory has a relatively minor effect on our overall conclusions. In the USA, the GFED-based and the GFAS-based estimates on fire-related annual PM2.5 exposure in 2023 are also similar (1.96 μg m−3 and 1.83 μg m−3, respectively), whereas the QFED-based estimates (3.11 μg m−3) are 59–70% higher, consistent with several other recent studies (Supplementary Fig. 14). Because the GEOS-Chem PM2.5 simulations with the QFED inventory show better agreement with surface observations in the USA, our estimates of fire-related contributions to PM2.5 may be underestimated in the USA when using the GFED inventory.
Our population-weighted estimates of fire-related PM2.5 based on GFED emissions in 2017 (1.62 μg m−3 in Canada and 1.16 μg m−3 in the USA) are nearly identical to estimates by a previous study that used GFED emissions from 2000–201910 (1.50 μg m−3 in Canada and 1.21 μg m−3 in the USA). Another study44 estimated the source contribution to ambient PM2.5 in 2017 globally using GFED fire emissions and reported that the population-weighted fire-related PM2.5 in Canada and the USA was 1.35 μg m−3 and 0.90 μg m−3, respectively, again similar to the GFED-based and the GFAS-based estimates. In addition, ref. 11 developed a data fusion model to divide the fire-source and other-source PM2.5 in the USA and estimated that the mean fire-source PM2.5 concentrations inside and outside the vicinity of an EPA air-quality monitoring station (defined by a 5-km radius) in 2017 are 0.97 μg m−3 and 0.92 μg m−3, respectively, also close to our GFED-based and GFAS-based estimates.
In summary, although the GFED, the QFED and the GFAS inventories show varied air-pollutant emissions and led to different fire-related PM2.5 concentrations in GEOS-Chem simulations, the total PM2.5 retrieval models with different fire emission inventories showed generally consistent performance. For the fire-related PM2.5 concentrations, the GFED-based and the GFAS-based estimates are similar and more comparable to previous studies9,10,11,44. We understand that each fire inventory has its own advantages and disadvantages; therefore, we cannot justify which one is the best. Given that the GFED has frequently been used in recent analyses8,10,12,43, to put our results in the context with the literature, we choose to present the GFED-based estimates in the main text.
Fire source attribution
To further quantify the contributions of fires from different regions in North America, we divided Canada and the USA into four regions, Eastern Canada (CE), Western Canada (CW), Eastern USA (UE) and Western USA (UW) (Extended Data Fig. 3b). According to the EPA’s delineation of North American ecological regions (Supplementary Fig. 2; https://www.epa.gov/eco-research/ecoregions-north-america/), the Level I ecoregions divide North America into 15 broads. It can be seen that the majority of Canada is covered by forests, whereas the eastern and western parts of the USA are separated by plains and deserts. Therefore, we defined the above four regions according to administrative divisions in combination with the ecological regions in this study.
We then implemented simulations with regional wildfire emissions turned off using the GFED emissions (Supplementary Table 2). In addition, given that the zero-out approach may lead to additional bias owing to the nonlinear relationship between emissions and modelled PM2.5 concentrations, the case of ‘offNA’ is used to evaluate this impact and constrain the results. The regional contributions can be calculated using equations (7)–(11):
$$\begin{array}{l}{\rm{Scale}}=({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offNA}},{\rm{GFED}}})/[({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offCE}},{\rm{GFED}}})\\ \,\,+\,({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offCW}},{\rm{GFED}}})\\ \,\,+\,({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offUE}},{\rm{GFED}}})\\ \,\,+\,({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offUW}},{\rm{GFED}}})]\end{array}$$
(7)
$${F}_{{\rm{fireCE}},{\rm{GFED}}}={\rm{Scale}}\times ({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offCE}},{\rm{GFED}}})/{{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}$$
(8)
$${F}_{{\rm{fireCW}},{\rm{GFED}}}={\rm{Scale}}\times ({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offCW}},{\rm{GFED}}})/{{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}$$
(9)
$${F}_{{\rm{fireUE}},{\rm{GFED}}}={\rm{Scale}}\times ({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offUE}},{\rm{GFED}}})/{{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}$$
(10)
$${F}_{{\rm{fireUW}},{\rm{GFED}}}={\rm{Scale}}\times ({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offUW}},{\rm{GFED}}})/{{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}$$
(11)
where GCoffCE,GFED, GCoffCW,GFED, GCoffUE,GFED, GCoffUW,GFED and GCoffNA,GFED are the GEOS-Chem-simulated PM2.5 concentrations using the GFED emissions under the offCE, offCW, offUE, offUW and offNA scenarios, respectively. FfireCE,GFED, FfireCW,GFED, FfireUE,GFED and FfireUW,GFED are the fractional contributions to the PM2.5 concentrations from fire emissions in CE, CW, UE and UW, respectively.
The zero-out approach used here may introduce additional bias due to the nonlinear relationship between emissions and modelled PM2.5 concentrations. The bias can be calculated as equation (12):
$$\begin{array}{l}{\rm{Bias}}=({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offNA}},{\rm{GFED}}})-[({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offCE}},{\rm{GFED}}})\\ \,\,\,+\,({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offCW}},{\rm{GFED}}})\\ \,\,\,+\,({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offUE}},{\rm{GFED}}})\\ \,\,\,+\,({{\rm{GC}}}_{{\rm{base}},{\rm{GFED}}}-{{\rm{GC}}}_{{\rm{offUW}},{\rm{GFED}}})]\end{array}$$
(12)
The absolute and relative biases due to nonlinearity are presented in Supplementary Fig. 16. Over the four regions, the absolute biases range from −0.009 μg m−3 to 0.012 μg m−3, and the relative biases range from −1.35% to −0.09%, indicating that the biases related to nonlinear effects are relatively small.
To investigate the significance of the Canadian wildfire impact on the interannual variability in PM2.5 concentrations in downwind regions, we compared the estimated annual mean PM2.5 concentrations from the 2023 Canadian wildfires with 2000–2023 satellite-derived annual PM2.5 concentrations (https://sites.wustl.edu/acag/datasets/surface-pm2-5/#V6.GL.02.03/)100 and 2014–2023 annual mean PM2.5 observations over the USA and Europe (Supplementary Fig. 18). When comparing with the interannual variability in PM2.5 concentrations, the contribution of the 2023 Canadian wildfires was statistically significant in North America and Europe whereas it was statistically insignificant over other downwind regions. In the USA, the 2023 Canadian wildfires contributed to an average of 1.50 μg m−3 annual mean PM2.5 over the locations of the EPA sites, larger than the differences in observed the PM2.5 concentration between 2022 and 2023 (0.99 μg m−3) as well as the mean interannual variabilities in the observed PM2.5 concentration during 2014–2023 (0.59 μg m−3, P < 0.01). In Europe, the 2023 Canadian wildfires contributed to an average of 0.43 μg m−3 annual mean PM2.5 over the locations of European Environment Agency (EEA) sites, accounting for 20% of the differences in observed PM2.5 concentration between 2022 and 2023 (2.10 μg m−3) and 14% of the mean interannual variabilities in observed PM2.5 concentration during 2014–2023 (3.02 μg m−3, P < 0.05).
Long-range transport of the 2023 Canadian wildfire plumes to Europe
The Hybrid Single-Particle Lagrangian Integrated Trajectory (HYSPLIT) model has been widely used to track the back trajectories for air parcels arriving at the receptor101,102,103,104. To track the transports and sources of PM2.5 pollution in Europe related to the 2023 Canadian wildfires, we used the HYSPLIT model to calculate the 7-day back trajectories for the entire European region (excluding Russia) on a 0.1° × 0.1° grid basis from May to September in 2023 at an hourly resolution. The HYSPLIT model version 5.3.0 (https://www.ready.noaa.gov/ready2-bin/getlinuxtrial.pl/) was used for the analysis. Each backward trajectory was run for 7 days with 1-hour time steps, initialized at 0:00, 6:00, 12:00 and 18:00 coordinated universal time (UTC) daily. The arrival height was 500 m above ground, approximately within the planetary boundary layer. The HYSPLIT model was driven by three-dimensional meteorological fields from the Global Data Assimilation System of National Centers for Environmental Prediction (GDAS NCEP), with a time resolution of 3 hours, a horizontal resolution of 1° × 1° and a vertical resolution of 23 levels.
We then defined the transport trajectory density (TTD) to represent the capability of pollutant transport from source regions to the receptor region (Europe), which is the total number of trajectories passing through the 0.1° × 0.1° grid box of source regions during the study period. The TTD of each grid box during the study period can be calculated by equation (13):
$${{\rm{TTD}}}_{i}=\mathop{\sum }\limits_{j=1}^{m}{F}_{i,j}$$
(13)
where i is the index of each grid box and m is the total number of trajectories that passed through all receptor grid boxes (81,052 in this study). For each trajectory j (ranging from 1 to m), Fi,j is defined as 1 if trajectory j passes through grid box i; otherwise, Fi,j is defined as 0. Therefore, TTDi is the total number of trajectories that passed through the grid box i during the study period. Supplementary Fig. 9 shows the spatial distribution of TTD sums from May to September in 2023. High TTD values were observed over the majority of Canada between May and September 2023, indicating the frequent trans-Atlantic plumes that prompt pollution transported from the wildfire source regions to downwind regions, and finally reached Europe.
Supplementary Fig. 10 further illustrates the vertical transport process that brought the Canadian wildfire plumes to the surface of Europe during the late-June trans-Atlantic episode. High PM2.5 concentrations are observed on 29 June and 2 July over the large areas in Northern France and Belgium (Supplementary Fig. 8 and Supplementary Fig. 10a). During the pollution episode, the GEOS-Chem simulation shows that the enhancement of surface PM2.5 concentration from Canadian wildfires was accompanied by a high PM2.5 concentration at high altitude (Supplementary Fig. 10b). The modelled PM2.5 enhancement from Canadian wildfires corresponds well with the observed peak of PM2.5 concentration, demonstrating the vertical transport process during the episode. Meanwhile, 7-day back trajectories with 500-m arrival height indicate that the airflows originating from the Canadian wildfire source regions were transported into the boundary layer above the site location on 29 June and 2 July (Supplementary Fig. 10c), providing compelling evidence of long-range transport of wildfire-related PM2.5 from Canada to Europe.
Health impacts
Acute mortality attributable to exposure to Canadian wildfires
The acute and chronic mortality attributable to 2023-Canadian-wildfires-related PM2.5 exposure were estimated separately. The acute mortality was estimated for all grids ‘Canada smoke days’ in which both (1) grid daily mean PM2.5 concentrations exceeded 15 μg m−3 (the recommended 24-hour average guideline levels of the WHO) and (2) Canadian-wildfires-related PM2.5 accounted for at least half of the total daily PM2.5 (ref. 10). Details on the estimation of fire-specific PM2.5 concentrations are described above. Following previous studies50, the acute mortality attributable to Canadian wildfires PM2.5 exposure was assessed using equation (14):
$${D}_{i,j}=\mathop{\sum }\limits_{j=1}^{365}\{[({\rm{RR}}({C}_{i,j})-1)/{\rm{RR}}({C}_{i,j})]\times {P}_{i}\times ({I}_{i}({{\rm{Country}}}_{a})/365)\}$$
(14)
where Di,j represents the all-cause acute premature mortality attributable to Canadian-wildfires-related PM2.5 exposure in grid i on day j. RR(Ci,j) represents the relative risk at exposure level C in grid i on day j and the exposure level was assessed as the daily average Canadian-wildfires-related PM2.5 concentrations. A global RR estimates of 1.021 (95% CI, 1.018, 1.024)15 per 10 μg m−3 increase of wildfire PM2.5 exposure was used for all regions. Pi represents the population in grid i that was constructed as described in ‘Estimation of fire-specific PM2.5 exposure’. Ii represents the baseline all-cause death rate in grid i that belong to Country a, which was collected from the Global Burden of Disease (GBD) 2019 study (https://ghdx.healthdata.org/gbd-2019)45.
Compelling evidence for the increased toxicity of wildfire-related PM2.5 relative to all-source PM2.5 have been reported13,58, and various exposure–response functions for acute exposure to wildfire PM2.5 have been developed by recent studies15,46,48,49. Here we use a widely used global pooled relative risk of wildfire PM2.5 exposure15 to estimate the acute premature mortality considering the global nature of this study and the comparability across different regions. To investigate the impacts of the choice of relative risks on the acute premature mortality estimates, we further assess the acute premature mortality attributable to the 2023 Canadian wildfires by using the relative risks for wildfire PM2.5 from newly developed global meta-analysis46, two studies in North America48,49 and relative risk from a meta-analysis on all-source PM2.5 (ref. 50). The results of comparison are presented in Supplementary Table 7. The estimated global acute premature mortality attributable to Canada smoke day exposure varied by a factor of four when different RR estimates were used. Among the different functions derived for wildfire PM2.5, mortality estimates using global pooled relative risks15,46 (ranging from 2,800 to 5,400) were higher than those using regional relative risks48,49 (1,300–2,600). Estimates using relative risks derived for wildfire PM2.5 exposure15,46,48,49 generally yield higher mortality estimates (1,300–5,400) than those using all-source relative risk50 (1,800), implying increased toxicity of wildfire-related PM2.5 compared with all-source PM2.5.
Chronic premature mortality attributable to exposure to Canadian wildfires
All-cause premature mortality attributable to chronic smoke exposure from the 2023 Canadian wildfires was estimated with a meta-analysis relative risk estimate of 1.08 (95% CI, 1.06, 1.09) per 10 μg m−3 increase in PM2.5 exposure40.
The health burden attributable to chronic PM2.5 exposure was assessed using equation (15):
$${D}_{i}=[({\rm{RR}}({C}_{i})-1)/{\rm{RR}}({C}_{i})]\times {P}_{i}\times {I}_{i}({{\rm{Country}}}_{a})$$
(15)
where Dy,i,n represents the chronic premature mortality attributable to Canadian-wildfires-related PM2.5 exposure in grid i. RR(Ci) represents the relative risk at exposure level C in grid i. Ci represents the annual average PM2.5 concentration in grid i. Pi represents the population in grid i, and Ii represented the baseline all-cause death rate in grid i in Country a of year 2019, which was collected from the GBD 2019 study (https://ghdx.healthdata.org/gbd-2019)45. The theoretical minimum risk exposure level (TMREL) for the chronic health effects attributable to PM2.5 ranged between 2.4 μg m−3 and 5.9 μg m−3, as reported in the GBD 2019 study.
The chronic exposure mortality attributable to Canadian wildfires was assessed with the direct proportion approach105,106, which assumes that the increase in mortality is in proportion to the increases in PM2.5 exposure. Thus, the chronic exposure mortality attributable to all-source ambient PM2.5 exposure was assessed first and the Canadian-wildfires-associated chronic mortality was quantified by calculating the proportion of Canadian-wildfire-derived PM2.5 within all-source ambient PM2.5.
It should be noted that the RR used here is derived for all-source PM2.5 rather than wildfire PM2.5, owing to the limited epidemiological evidence of chronic health effects from wildfire-related PM2.5 exposure. We use the exposure–response function for all-cause mortality rather than the cause-specific exposure–response function (that is, the widely used GBD approach45) as this study aims to estimate the total mortality burden whereas the cause-specific may underestimate the total chronic mortality of ambient PM2.5 (ref. 107). Relative risk derived from regional meta-analysis47,49 may differ from those derived from global pooled analysis40. Given the global nature of this study, we choose the all-cause global pooled relative risk in our analysis40. We further conducted a sensitivity analysis to evaluate the impact of exposure–response functions on the chronic premature mortality45,47,49, as shown in Supplementary Table 8. Using cause-specific exposure–response function45 yields 31,000 global chronic premature mortality, lower than estimates using all-cause functions40,47,49 (ranging from 82,100 to 152,000). For all-cause premature mortality estimates with 3 different relative risks, mortalities estimated by global relative risk40 (82,100) are remarkably lower than those estimated by regional relative risks47,49 (that is, the USA, 117,500–152,000), indicating large variation of relative risks across different global regions.
We also reviewed the approaches of estimating chronic health burden from wildfire-related PM2.5 exposure used in different global and regional studies12,31,44,46,108,109,110,111,112,113,114,115,116,117,118,119,120,121,122,123, as shown in Supplementary Table 9. Although annual average exposure is widely used in those studies, it may not reflect the population experience of sporadic wildfire PM2.5 exposure owing to the different nature of wildfire smoke exposure and the urban and background pollution exposure. Substantial differences in chronic health burden assessment approaches were observed, which varied in mortality endpoints (cause-specific versus all-cause), the relative risk and the definition of exposed population. Many of those studies quantified PM2.5-related health impacts of landscape fires on populations that are annually impacted by fire-related air pollution from local or nearby fires, whereas our study investigated the transboundary air pollution and the health burden of a single extreme wildfire event. Specifically, ref. 31 quantified the transboundary PM2.5 health impact of wildfires in the Arctic Council, which is most comparable to the purpose of our study. In their analysis, areas where the increase in carbonaceous PM2.5 from wildfires was statistically insignificant were excluded from the health impact assessment. A previous study46 estimated global, regional, and national mortality burden attributable to fire-related PM2.5 exposure, which is estimated for all population with the same exposure–response function used in our study40. We then followed the widely accepted approach in our analysis and conducted significance tests for the estimated contribution of the annual mean PM2.5 concentration from the 2023 Canadian wildfires.
Uncertainty analysis
Our results are subject to a number of uncertainties and limitations. The uncertainty ranges (95% CI) in different steps of our analysis are discussed below.
First, the emission inventories used in this study bear large uncertainties. For example, the uncertainties in the GFED emissions mainly come from the inadequate representation of the natural variation of the emission factors during fire events, and the high uncertainties in the amount of fuel burned estimated from burned area. It is reported that a best-guess uncertainty assessment for GFEDv4.1s at regional scales could be a 1σ of about 50% in general but higher in areas where small fires burned area is important or where there is notable fuel consumption in organic soils37.
Second, the PM2.5 concentrations simulated by the global chemical transport model are affected by errors in emission inventories and the model’s representation of physical and chemical processes such as vertical transport51 and secondary organic aerosols124. Specifically, the smoke-injection height is not considered in the GEOS-Chem simulation, which may lead to overestimates and underestimates of fire-related PM2.5 concentrations in fire source regions and downwind regions, respectively51. Given the huge computational cost for model sensitivity simulations considering the uncertainties in emissions, we use the normalized root-mean-square deviation (NRMSD) between the modelled and the observed PM2.5 concentrations to represent the overall model errors in total PM2.5, and use the NRMSD between the modelled and the observed PM2.5 concentrations over fire events to represent the model errors in fire-related PM2.5. The NRMSD for GEOS-Chem-based total PM2.5 and fire-related PM2.5 ranged between 42.8% and 62.0% and between 44.3% and 53.0%, respectively, among Canada, the USA and Europe, but a bit higher in other regions globally (for example, close to 100.0%). Although the absolute errors in GEOS-Chem-based simulations are large, some errors are common between the total and fire-related PM2.5 and have limited impacts on their ratios.
Third, the multi-source fused PM2.5 data obtained from our retrieval model are influenced by errors in all the input data and the multilayer machine-learning model itself. We have fully evaluated the model’s performance using a cross-validation approach and the performance was comparable to previous studies10,41. We use the NRMSD between PM2.5 retrieval and observed PM2.5 concentrations to represent its uncertainties (2.0–7.1% among different regions).
As presented in equation (4), the overall uncertainties involved in fire-related PM2.5 exposures are determined by uncertainties in GEOS-Chem simulated fractional contributions of fire emissions (GCfire/GCbase) and in total PM2.5 exposures based on the retrieval model (CPM). The errors in GCbase, GCfire and CPM are defined above. The errors in the ratio (GCfire/GCbase) were then quantified by 10,000 trials of Monte Carlo simulation. Finally, the overall uncertainties of fire-related PM2.5 (Cfire) were derived from the aggregations of errors above.
The overall uncertainties (presented as 95% CI) of acute and chronic mortality attributable to Canadian wildfires were then assessed by Monte Carlo simulations with 1,000 iterations45. Uncertainties embedded in all input parameters of the risk assessment model were considered. The uncertainties in exposure levels are described above. The uncertainties in baseline mortality, exposure–response functions and TMREL were collected from the GBD 2019 study. The uncertainty in national total population was provided by United Nations data (high-fertility and low-fertility scenarios). All the parameters, except TMREL, which was simulated by a uniform distribution, were simulated by normal distributions.
Comparison with other relevant studies
The impacts of the 2023 Canadian wildfires on surface PM2.5 air quality have been reported in a few recent studies. By proposing a multidimensional air pollution correlation network framework, ref. 28 argued that the 2023 Canadian wildfires significantly impact the air pollution behaviour in the Northeastern USA region. Reference 29 estimated the PM2.5 concentration in the Northeastern USA in June 2023 by combining chemical transport model results and surface PM2.5 measurements. They identified two ‘smoke wake’ events in June 2023 (6–8 June and 28–30 June) with significant PM2.5 enhancement caused by the transport of Canadian wildfire plumes. Our results also capture these two events in the same region, although our estimates of PM2.5 concentrations are lower than those of ref. 29 during the first ‘smoke wave’ event (6–8 June), which might be attributed to the coarse model resolution in our analysis (Supplementary Table 3). By using a global chemical transport model, ref. 2 found that the Canadian wildfires significantly impacted air quality in the Northern Hemisphere, which was consistent with our findings. They identified six widespread air pollution episodes due to the Canadian wildfires from May to August, which are also captured in our GESO-Chem simulation (Supplementary Fig. 11) and retrieved PM2.5 concentration (Supplementary Fig. 12).
Canadian wildfire-related chronic PM2.5 exposure is associated with approximately 22, 10 and 3 deaths per 100,000 people in 2023 in Canada, the USA and Europe, respectively. For comparison, chronic mortalities from all-source PM2.5 exposure are 55, 57 and 80 per 100,000 people in 2023 in the three regions, respectively, underscoring the non-negligible contribution of wildfire smoke to public health burdens. Our estimates on Canadian wildfire-related per capita mortality rates in the USA and Canada are notably higher than wildfire-related per capita mortality rates reported in refs. 46,125, owing to high PM2.5 exposure levels attributable to the record-breaking Canadian wildfires in 2023 as well as the all-cause exposure–response function used in our analysis.
We estimated that the 2023 Canadian fires accounted for 3.82 μg m−3 (3.00–4.64 μg m−3) of annual mean PM2.5 exposure in Canada in 2023, which is lower than 2000–2019 annual mean fire-related PM2.5 exposure in typical wildfire hotspot regions10 such as sub-Saharan Africa (6.99 μg m−3), mainland Southeast Asia (5.77 μg m−3), Indonesia (6.28 μg m−3) and Brazil (5.68 μg m−3). In another study44, annual mean fire-related PM2.5 exposure in typical wildfire hotspot regions were estimated to be 0.72 μg m−3 in Indonesia, 1.26 μg m−3 in Brazil and 1.26 μg m−3 in Southeast Asia in 2017, indicating the large interannual variabilities in wildfire activities. Reference 109 estimated that global landscape fires alone result in 44 million and 4 million people annually being exposed to air quality considered ‘unhealthy’ and ‘hazardous’, respectively. In comparison, we estimated that the 2023 Canadian wildfires caused 139.3 million and 0.25 million people to be exposed to ‘unhealthy’ (PM2.5 > 55 μg m−3) and ‘hazardous’ (PM2.5 > 250.5 μg m−3) air quality. Large health impacts from wildfire-related air pollution have been reported in these hotspot regions46,109,110,111,112,114,115,116,118. For instance, wildfires-induced chronic mortality was estimated to be 160,200 in Africa in 2017110, 59,000 in Southeast Asia in 2014118, 13,700–44,000 in equatorial Asia during 2004–2015116, and 16,800 in South America in 2012111. In recent study46, it was estimated that 384,600, 144,300 and 79,300 people died annually in sub-Saharan Africa, Southeast Asia, and Latin American and the Caribbean, respectively, owing to chronic wildfire smoke exposure during 2000–2019. In contrast, we estimated that the 2023 Canadian fires resulted in 8,300 (95% CI, 5,800–10,800) PM2.5-attributable chronic premature deaths in Canada given the low population density close to fire regions. Globally, we estimated 82,100 (95% CI, 47,700–116,500) PM2.5-attributable chronic premature deaths owing to smoke exposure from the 2023 Canadian wildfires, indicating the large health impacts from long-range transported PM2.5 pollution.
Although notable impacts of the 2023 Canadian wildfires on surface PM2.5 concentration in Europe are observed, those impacts are remarkably low compared with the impact of smoke from local fires that are not diluted by long-range transport. For instance, ref. 126 estimated that 2.1 million people were exposed to concentrations above 36 μg m−3 for at least 1 day between 23 and 30 June owing to the Saddleworth Moor and Winter Hill fires in northern England. In comparison, our estimates shows that 0.81 million people in Europe were exposed to the same level of PM2.5 pollution for at least 1 day between 26 June and 7 July 2023 owing to the trans-Atlantic pollution of Canadian wildfires.
Dust has been recognized as another important natural source of air pollution44,127. A study127 estimated that Sahara dust contributed 5–20 μg m−3 of surface PM10 concentration in South Europe and 0.5–1.0 μg m−3 in North Europe during 2016–2017. Another study44 estimated that windblown dust increased the annual mean PM2.5 exposure in 2017 by 1.72 μg m−3, 1.50 μg m−3, 1.18 μg m−3, 0.07 μg m−3 and 0.19 μg m−3 in Central Europe, Eastern Europe, Western Europe, Canada and the USA, respectively. PM2.5 exposure was lower than dust-related PM2.5 exposure in Europe but higher than that in Canada and the USA. In turn, windblown dust contributed to 34,972, 33 and 1,126 annual chronic premature deaths in Europe, Canada and the USA, respectively44. Our estimates for Canadian-wildfire-related deaths are substantially lower than dust-related deaths in Europe, but much higher than that in Canada and the USA. It should be noted that both dust and fire activities have large interannual variabilities so the comparison could be different for other years.