Estimating PM2.5-related premature mortality and morbidity associated with future wildfire emissions in the western US

Wildfire activity in the western United States (US) has been increasing, a trend that has been correlated with changing patterns of temperature and precipitation associated with climate change. Health effects associated with exposure to wildfire smoke and fine particulate matter (PM2.5) include short- and long-term premature mortality, hospital admissions, emergency department visits, and other respiratory and cardiovascular incidents. We estimate PM2.5 exposure and health impacts for the entire continental US from current and future western US wildfire activity projected for a range of future climate scenarios through the 21st century. We use a simulation approach to estimate wildfire activity, area burned, fine particulate emissions, air quality concentrations, health effects, and economic valuation of health effects, using established and novel methodologies. We find that climatic factors increase wildfire pollutant emissions by an average of 0.40% per year over the 2006–2100 period under Representative Concentration Pathway (RCP) 4.5 (lower emissions scenarios) and 0.71% per year for RCP8.5. As a consequence, spatially weighted wildfire PM2.5 concentrations more than double for some climate model projections by the end of the 21st century. PM2.5 exposure changes, combined with population projections, result in a wildfire PM2.5-related premature mortality excess burden in the 2090 RCP8.5 scenario that is roughly 3.5 times larger than in the baseline period. The combined effect of increased wildfire activity, population growth, and increase in the valuation of avoided risk of premature mortality over time results in a large increase in total economic impact of wildfire-related PM2.5 mortality and morbidity in the continental US, from roughly $7 billion per year in the baseline period to roughly $36 billion per year in 2090 for RCP4.5, and $43 billion per year in RCP8.5. The climate effect alone accounts for a roughly 60% increase in wildfire PM2.5-related premature mortality in the RCP8.5 scenario, relative to baseline conditions.

Ecoregion delineations were adapted from Bailey et al. (1994) to ensure consistency with LOCA data -see Figure  S1 below and accompanying text below.

Data Type Description Data Documentation and Availability
Atmospheric modeling Open access GEOS-Chem model version 12.5.0 (September 2019). As the GCMs have very coarse spatial resolution, we use spatially downscaled climatic variables from the LOcalized Constructed Analogs (LOCA) dataset for each RCP and GCM combination.
The LOCA dataset provides daily projections through 2100 at a 1/16 degree resolution (~6.25 km) for three variables: daily maximum temperature (tmax), daily minimum temperature (tmin), and daily precipitation (see U.S. EPA 2017a for more details).
It is important to place wildfire activity that occurs in the western U.S. (31°-49°N, 101°-125°W) in the context of total wildfire activity across the contiguous U.S. As part of the Monitoring Trends in Burn Severity (MTBS) project, the USGS Center for Earth Resources Observation and Science (EROS) and the USDA-FS Remote Sensing Applications Center compiled fire data from federal agency databases (ICS 209) and state databases, then used Landsat remotely-sensed imagery data to map burn severity and perimeters of all large fires in the U.S. from 1984 to 2016. This dataset is limited by its exclusion of small fires (fires smaller than 500 acres in the east and smaller than 1,000 acres in the west). However, even with the exclusion of smaller fires, the dataset captured over 95% of area burned when compared to ground-based data in 2004. Hence, the MTBS data is used to estimate the fraction of contiguous U.S. wildfire activity that has historically occurred in the western U.S., both with respect to number of fires and area burned.
Prescribed wildfires and unknown wildfire types were excluded from the dataset before calculations spanning the baseline years (1986)(1987)(1988)(1989)(1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005) were performed. With respect to the number of fires, on average, there were 2.4 times as many fires in the West than in the East from 1986-2005. Put differently, 70% of fires in the contiguous U.S. occurred in the western U.S. on an average annual basis. The differential size thresholds for inclusion in the MTBS dataset result in an underestimate of contributions of Western fires, given that Western fires must be double the size of Eastern fires to be counted. Thus, not only does the West have a higher number of fires than the East, but the West also has larger fires. Accordingly, there were 5.6 times as many acres burned in the West than in the East on average, and 81% of area burned was located in the Western U.S. on an average annual basis from 1986-2005. In effect, this analysis captures a large fraction of contiguous U.S. wildfire activity in examining the western U.S.
The present analysis builds upon existing research relating climatic conditions to annual area burned to estimate the future health burden associated with wildfire emissions in the western U.S. The flow chart in main text Figure 1 summarizes our approach and highlights the components that are unique to the present analysis and those that have been adopted from earlier The spatial and temporal disaggregation described in step 3 of main text Figure 1 was adopted from Yue et al., (2013). As part of that process, the annual area burned in each ecoregion is distributed across the landscape and each month of the wildfire season while maintaining historic patterns of wildfire activity. Table S3 from Yue et al. (2013) (found here: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3763857/#FN3 ) provides the historic burn fraction for each ecosystem that was used as part of the spatial disaggregation process.
Yue et al. (2013) presents two alternative methodologies for relating climate variables to wildfire area burned in the western U.S. The first approach expands on earlier work conducted by Spracklen et al. (2009) and utilizes a stepwise regression analysis to relate observed annual area burned in each ecoregion to several meteorological predictor variables and fire indexes. The second approach is a physical parameterization of area burned based on monthly mean temperature, relative humidity and precipitation. The two approaches differ not only with respect to the specific climate variables utilized but also by their temporal resolutionthe regression models predict annual area burned, while the parameterization estimates monthly area burned.
Here, we use the six ecoregion-specific regression models because their predictions of present-day area burned are more consistent with observations than the parameterized model in every ecoregion other than the Great Plains. We modified and re-estimated the regression model for Rocky Mountain Forest ecoregion, as described in the section 3 belowthe result is summarized in Table S1. The regressions and R 2 statistics for five of the ecoregions (Pacific Northwest, California Coastal Shrub, Desert Southwest, Nevada Mnt./Semi-Desert, and Eastern Rocky Mtn./Great Plains)_are used directly from Yue et al., 2013, and are available in Table 1 here: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3763857/#FN3 ). Of particular importance is that the regression models perform better than prior work in the Pacific Northwest and California Coastal Shrub ecoregions. These two ecoregions contain a large percentage of total population in the western United States and are therefore expected to comprise a significant proportion of the future health impacts associated with wildfire smoke. Overall, our regressions capture a substantial amount of interannual variability in area burned that can be attributed to climatic variation, though they generally underestimate peak years across all ecoregions (Main Text Figure 2 and  (2013) were used (Fig. S1), These meteorological observations include the Global Surface Summary of the Day (GSOD) sites or United States Historical Climatology Network (USHCN) sites, depending on the type of climate variabletemperature and precipitation variables were calculated from USHCN sites and relative humidity variables and fire indexes were calculated from GSOD sites.

Relative Humidity and Rocky Mountain Ecoregion Regression
A complete estimate of area burned by wildfires requires the input of three additional climate variables: relative humidity (RH), wind speed and air pressure. We expand the LOCA climate dataset to include these variables using a binning approach described in the Technical Appendix to U.S EPA (2017a). Briefly, the binning approach effectively associate historical time series of climate variables to future projections. For example, to project future RH, we associated RH with temperature and precipitation based on historical daily data, additional climate variables, including relative humidity, wind speed and air pressure, using a binning approach. Relative humidity (RH), for example, is required to estimate key inputs for the area burned regressions, in particular the Duff Moisture Code (DMC) and Build-Up Index (BUI) variables, but RH is not available from the LOCA downscaled dataset. The historical pattern of relative humidity cannot be simply repeated to fill in the missing LOCA values, as the arrival times of the LOCA tmax, tmin, and precipitation outputs are drawn from the GCM projections rather than the historical time series. To fill in the missing values, we use a binning approach, which starts with the assumption that there is some relationship between temperature/precipitation and relative humidity. We effectively associate historical daily humidity with daily historical temperature and precipitation, and then used those binned associations to assign projected RH in LOCA future daily time series. Historical humidity for our study area is available from the Princeton Land Surface Hydrology Group (Sheffield et al. 2006).
Yue et al. originally calculated relative humidity (RH) as a function of daily mean temperature and dew point. This approach is different than the one utilized in the LOCA dataset, and resulted in minor differences in these regression variables that were amplified by the area burned equations.
To correct for this bias, a set of regression equations were developed to relate each RH predictor variable with temperature and precipitation variables. Then the RH variables were recalculated using those models and the LOCA temperature and precipitation data. While the fire indexes also rely on estimates of daily RH, the effect of these differences was not as significant for those variables, and they were not recalculated. The exception to this was the Rocky Mountain ecoregion, where a threshold effect related to the annual Duff Moisture Code index resulted in consistent underestimation of area burned. As a result, the regression model for that ecoregion was re-derived using the same stepwise regression approach utilized by Yue et al., but using only temperature and precipitation predictor variables.
Additional analysis was performed to improve the performance of the Rocky Mountain ecoregion area burned regressions using LOCA inputs.

Additional Details on Emissions Calculations
To calculate the quantity and spatial distribution of wildfire-related black carbon (BC) and organic carbon (OC) emissions, we followed the process described by Yue et al. (2013). Briefly, we first partitioned annual area burned into monthly totals using the observed historic seasonal pattern within each ecoregion, assuming that this intra-annual temporal pattern will not change significantly over the century. For each month, total area burned was distributed across a 0.5ºx 0.5º grid. To reflect observed spatial patterns in area burned, 70% of area burned was randomly distributed to only 10% of the grid cells while maintaining the historic burn fractions across ecosystems from Spracklen et al. (2009). The remaining 30% of area burned was distributed evenly to the remaining 90% of grid cells, also while maintaining historic burn fractions.
We then calculated the total fuel load for each grid cell using the U.S. Forest Service 1 km x1 km fuelbed map used by Yue et al. (2013)

GEOS-Chem Air Quality Results
Due to the computational expense of conducting 300 model simulation years (30 years, two climate scenarios, and five GCMs), we ran GEOS-Chem with offline chemistry (aerosol-only), which used monthly mean oxidants archived from a previous full-chemistry simulation, including OH, NO3, O3, total nitrate, and production and loss rates for H2O2. The lack of full chemistry feedback in our simulation could affect the concentrations of secondary aerosols (organic and inorganic) that are formed through chemical processes. Therefore, we used secondary aerosol concentrations from a full-chemistry simulation for 2013 conducted by Jun et al. (2019) using the same meteorological data and emissions as this study to correct for the potential bias in the offline simulation. We kept these secondary aerosol concentrations constant to 2013 to isolate the impacts of wildfire BC and OC emissions on PM2.5. Simulated PM2.5 concentrations were highly consistent with observations in both time and space ( Figure S5).       elasticity of VSL, we also perform a sensitivity test using an elasticity of 1.0, reflecting proportional growth in VSL with gross domestic product per capita. Using an elasticity of 1.0 compared with 0.4, our estimated costs associated with the mortality burdens are 87% higher in 2090. The sensitivity analysis suggests that our valuation of the mortality endpoints may be conservative.

BenMAP Modeling and Detailed Health Burden Results
wildfires burning 100+ acres or involving Incident Management Teams of type 1 or 2. We used geographic coordinate information in the web portal database to compile response costs at the state level, and assigned average response costs by state to our area burned estimates described above to generate total wildfire response costs for our projection period, as well as annual estimates, for each of the five GCMs and two RCPs. Response costs used in Mills et al. (2014) reflected data available for 2002 -2012; we updated the average regional response costs using newly available data for 2002 -2018, which led to a slight increase in the average response costs for most regions, and larger increases for the California regions.
This analysis used data from the National Wildfire Coordinating Group (NWCG) to monetize the projected changes in acreage burned by wildfires. Specifically, the analysis developed spatially resolved wildfire response costs from NWCG data on the size (i.e., acres burned), origin, and total response costs for distinct wildfires in the contiguous United States (U.S.) from 2002 to 2018. Figure S8 provides the Geographic Area Coordination Center (GACC) boundaries used by the NWCG to coordinate wildfire responses and to collect and report wildfire-related data (GACC 2011). Duplication Adjustments: Duplicate fires, with the same incident number, incident name, and start date were adjusted. The cost value that was the larger of the two was assumed to be cumulative and used for our calculation purposes. A total of 191 duplicates existed, of which only 23 had differing cost values.

Great Basin region:
In the years 2015-2018, Western Great Basin and Eastern Great Basin were combined. We manually altered the coding to reflect the two, WB and EB. According to the GACC, Nevada made up the Western Great Basin whereas the other states of the Great Basin region were part of the Eastern Great Basin.
Wildfire Mapping: Wildfires are mapped in ArcGIS based on latitude and longitude information provided in ICS-209 forms. Less than one percent of all fires do not have latitude and longitude. Those fires are excluded from the analysis. Fires with locations that do not fall within the contiguous US were also excluded. Some of these fires likely had latitudes and longitudes with errors (such as those placed within the Arctic Circle), but there was no clear pattern to the errors and thus the location data for these fires could not be reliably corrected. After excluding fires outside the contiguous US or without spatial information, 20,812 fires remain in the analysis. Each of these fires is assigned to a LOCA 0.5 o x 0.5 o grid cell.