ELM2.1-XGBfire1.0: improving wildfire prediction by integrating a machine learning fire model in a land surface model

Wildfires have shown increasing trends in both frequency and severity across the contiguous United States (CONUS). However, process-based fire models have difficulties in accurately simulating the burned area over the CONUS due to a simplification of the physical process and cannot capture the interplay among fire, ignition, climate, and human activities. The deficiency of burned area simulation deteriorates the description of fire impact on energy balance, water budget, and carbon fluxes in the Earth system models (ESMs). Alternatively, fire models based on machine learning (ML), which capture statistical relationships between the burned area and environmental factors, have shown promising burned area predictions and corresponding fire impact simulation. We develop a hybrid framework (ELM2.1-XGBFire1.0) that integrates an eXtreme Gradient Boosting (XGBoost) wildfire model with the Energy Exascale Earth System Model (E3SM) land model (ELM) version 2.1. A Fortran–C–Python deep learning bridge is adapted to support online communication between ELM and the ML fire model. Specifically, the burned area predicted by the ML-based wildfire model is directly passed to ELM to adjust the carbon pool and vegetation dynamics after disturbance, which are then used as predictors in the ML-based fire model in the next time step. Evaluated against the historical burned area from Global Fire Emissions Database 5 from 2001–2019, the ELM2.1-XGBFire1.0 outperforms process-based fire models in terms of spatial distribution and seasonal variations. The ELM2.1-XGBFire1.0 has proven to be a new tool for studying vegetation–fire interactions and, more importantly, enables seamless exploration of climate–fire feedback, working as an active component of E3SM.

1 Introduction

Recent wildfire outbreaks worldwide have raised alarms due to wildfires burning longer and more intensely in many regions, posing significant threats to human livelihoods and biodiversity. In the past 2 decades, satellite-derived data suggest that the global total burned area has declined by over 20 %, which is primarily attributed to human influences (Jones et al., 2022; Andela et al., 2017). However, the contiguous United States (CONUS) has emerged as a hotspot for wildfires, where both climate change and human activities have fueled a 42 % increase in burned area (Jones et al., 2022). Such expansive burned areas release an average of 162×10⁶ t of CO₂ and 0.9×10⁶ t of PM_2.5 annually into the atmosphere, resulting in over USD 200 billion health costs due to exposure to wildfire smoke (Global Wildfire Information System, 2024; JEC, 2023). Accurate prediction of wildfire risks has become an urgent need.

Traditional fire models, predominantly process-based models, simulate the behavior of individual wildfires using theoretical equations for ignitions and fire spread (Hantson et al., 2016). These models explicitly simulate the number and size of individual fires by incorporating parameterizations and parameters derived from laboratory or field experiments and typically estimate the burned area by scaling up to the grid-cell level (Lasslop et al., 2014; Pfeiffer et al., 2013; Yue et al., 2014; Li et al., 2012; Thonicke et al., 2010; Huang et al., 2020, 2021; Arora and Boer, 2005; Burton et al., 2019). While process-based wildfire models are effective in simulating global burned area distribution (Hantson et al., 2020), they often fall short when accurately predicting the extent and temporal changes of wildfires over the CONUS (Forkel et al., 2019; Teckentrup et al., 2019). The climate and vegetation controls on the CONUS burned area and their relative importance are incorrectly represented, leading to failures in burned area predictions regarding both spatial distribution and temporal variations (Forkel et al., 2019). Human ignition and suppression are assumed to be linearly or log-linearly related to population density and the gross domestic product (GDP), respectively (Jones et al., 2023; Li et al., 2013). This assumption overlooks a more nuanced picture of human activities, such as road density, cultural differences, agricultural activities, and forest management policy (Jones et al., 2022; Villarreal et al., 2022; Hanan et al., 2021; Miller et al., 2009; Turco et al., 2023; Haas et al., 2022). Process-based fire models are often integrated with biogeochemical process-enabled land models (hereafter referred to as BGC models) within Earth system models (ESMs) to predict fire disturbances to carbon allocation, which is then used to update energy balance, water budget, and carbon fluxes in the land model. Incorrect simulation of burned areas over the CONUS induces large uncertainties in the assessment of fire impacts using ESMs.

Recent advances have explored the application of machine learning (ML) techniques in wildfire prediction (e.g., Buch et al., 2023; Li et al., 2023; Wang et al., 2021; Zhu et al., 2022). ML models offer the advantage of capturing nonlinear dependencies and complex interactions between driving factors and fire dynamics without the need for the explicit understanding of physical processes (Rodrigues and de la Riva, 2014). Zhu et al. (2022) presented a deep neural network (DNN) scheme that surrogated the process-based wildfire model with the Energy Exascale Earth System Model (E3SM) interface, demonstrating over 90 % higher accuracy in simulating the global burned area. Wang et al. (2021) combined the local predictors, large-scale meteorological patterns, and the eXtreme Gradient Boosting (XGBoost) algorithm to build an ML wildfire model, which improves the temporal correlations of burned areas in several regions over the CONUS by 14 %–44 %. Buch et al. (2023) developed a novel stochastic machine learning (SML) framework, SMLFire1.0, with a high spatial resolution of 12 km over the western US (WUS).

The newly developed ML fire models often focus on wildfire properties such as burned area, fire count, and fire emissions (Wang et al., 2021; Buch et al., 2023). Despite the improved fire predictions, fire impacts on the ecosystem, climate, and human community cannot be evaluated without integrating the wildfire process into the Earth system. In addition, climate change impacts on the burned area, either directly through fire weather conditions or indirectly through ecosystem productivity, vegetation type, fuel loads, and fuel moisture, cannot be fully understood without explicitly representing the complex interplays between climate, ecosystems, and fire. For instance, a warmer and drier climate has been shown to cause an 8-fold rise in the high-severity burned area from 1985 to 2017 over the WUS (Parks and Abatzoglou, 2020). The corresponding changes in fire dynamics may shift the vegetation species distribution from those originally low in resistance to wildfire to those in high resistance or even benefiting from regular fire occurrence (Rogers et al., 2015; Huang et al., 2024). The fire-adapted vegetation species, in turn, facilitate the frequent occurrence of wildfires. Taking this into consideration, a full coupling of fire, ecosystem, and climate is required to better predict fire changes and the corresponding impacts in a future climate.

Leveraging the accuracy of ML-based wildfire models and the representation of ecosystem–climate interactions in ESMs, in this study, we develop a novel hybrid framework to integrate a pretrained ML wildfire model with the E3SM land model (ELM) to study the full atmosphere–vegetation–wildfire feedback. This integration facilitates a dynamic feedback loop where outputs from the ML model (i.e., predicted burned areas) inform the land surface processes in ELM, which in turn update the inputs for the ML model for subsequent predictions. This approach leverages the detailed physical understanding of surface biogeophysical and biogeochemical processes provided by ELM and the predictive power of ML-based wildfire models to create a more accurate and robust framework for wildfire prediction and impact assessment. The remaining sections are arranged as follows: Sect. 2 introduces the ELM and ML wildfire model training method, coupling strategy, and datasets used in this study; Sect. 3 presents the simulated burned area compared with observations and several process-based fire models; and Sect. 4 contains the discussion and conclusion.

2 Materials and methods

2.1 Model description

2.1.1 Default wildfire model in ELM

The ELM is part of the E3SM project which started with a version of the Community Earth System Model (CESM1). The ELM default wildfire module originated from the Community Land Model (CLM4.5) (Li et al., 2012). This wildfire model calculates burned areas by multiplying the number of wildfires and burned area per fire on a grid-cell level. The number of wildfires (fire count) is derived using anthropogenic and natural ignition sources, fuel load and combustibility, surface meteorology, and anthropogenic suppression. The natural ignition source is derived from the number of cloud-to-ground lightning flashes multiplied by a constant ignition efficiency (Prentice and Mackerras, 1977). Anthropogenic ignitions are simply parameterized using a fixed number of potential anthropogenic ignitions by a person and population density (Venevsky et al., 2002). Humans also suppress wildfires. The capability of fire suppression is assumed to be a function of GDP and population density. The ignition efficiency is also altered by fuel conditions, including the fuel load (aboveground biomass) and fuel combustibility (approximated using relative humidity, temperature, and top or root zone soil moisture). The spread of each fire is approximated using an ellipse shape, with its length-to-breadth ratio determined by wind speed and fuel moisture (Rothermel, 1972). This simple concept captures the major constraints for predicting the global wildfire distribution and seasonal variations well (Rabin et al., 2017; Li et al., 2014; Huang et al., 2020).

Like many other process-based wildfire models, the default fire model in ELM benefits from the full ecosystem interactions from its hosting land model as well as the potential to be coupled with atmospheric models. With the BGC processes being turned on, ELM-BGC reallocates carbon and nitrogen in leaf, wood, root, litter, and soil pools after fire based on carbon combustion and mortality rate dependent on plant functional type (PFT). The biogeochemical changes subsequently influence biogeophysical properties such as leaf area index (LAI), vegetation canopy height, and albedo, disturbing the land–atmosphere exchanges of energy and water fluxes. The post-fire vegetation recovery in ELM-BGC depends on the plant photosynthesis processes and PFT competition strategy for soil resources. The interactions between wildfire and vegetation under historical climate have been thoroughly assessed in CLM long-term simulations (Li and Lawrence, 2017). The model framework is illustrated in Fig. 1. Hereafter the ELM coupled with the process-based fire model is referred to as ELM-BGC.

Figure 1 Schematic diagram of the hybrid model framework.

2.1.2 Machine learning wildfire model

The XGBoost-based wildfire model has proven to outperform process-based models in predicting burned areas over the CONUS (Wang et al., 2021). XGBoost is a highly efficient and scalable implementation of gradient boosting, designed for performance and speed (Chen and Guestrin, 2016). It builds sequential decision trees to correct errors from previous models, using techniques like regularization to prevent overfitting and parallel processing for faster computation. In this study, we adapt the XGBoost algorithm used in Wang et al. (2021) to develop an offline ML fire model using variables directly provided by ELM at each grid cell. Wang et al. (2021) integrated large-scale meteorological patterns alongside local weather, land surface properties, and socioeconomic data to enhance the prediction of burned areas. The large-scale patterns were identified using singular value decomposition (SVD) to capture influential atmospheric conditions that develop over days to weeks and cumulatively impact the monthly burned area. The feature importance analysis in their study noted that while large-scale patterns improved prediction, however, they played a secondary role. Therefore, we exclude the large-scale patterns from predictors without significantly affecting the model accuracy. Hereafter the uncoupled XGBoost fire model is referred to as offline-XGB.

2.1.3 Hybrid modeling framework

The offline-XGB model is integrated with the ELM using the ML4ESM coupling framework. The ML4ESM framework offers a robust and flexible solution for integrating ML parameterizations into ESMs through a Fortran–Python interface (Zhang et al., 2024). It supports popular ML libraries such as PyTorch, TensorFlow, and Scikit-learn, enabling the seamless incorporation of ML algorithms to represent complex climate processes like convection and wildfire dynamics. The interface leverages C language as an intermediary for efficient data transfer by accessing the same memory reference instead of the extra data copy or through files, minimizing memory overhead and computational inefficiencies. A C hub is then used to communicate variables from the Fortran-written ELM and the Python-written ML fire model. In our application, all surface meteorology, lightning, and socioeconomic data alongside the ELM simulated fuel conditions are passed to the ML-based fire model to predict the burned area. The burned area is returned to ELM to calculate fire impacts and update surface properties.

2.2 Datasets and processing

2.2.1 Burned area datasets

The primary dataset for training and validating the ML-based model is the Global Fire Emissions Database version 5 (GFED5) (Chen et al., 2023). The GFED5 is a succession of GFED4s (van der Werf et al., 2017), which we also use as an additional reference dataset. GFED4 is generated by fusing multiple streams of remote sensing data to create a 24-year (1997–2020) dataset of the monthly burned area at 0.25° spatial resolution. During 2001–2020, the GFED5 comprises the Moderate Resolution Imaging Spectroradiometer (MODIS) MCD64A1 burned area product (Hall et al., 2016; Giglio et al., 2016; Giglio et al., 2018) with adjustments for the errors of commission and omission. Adjustment factors are estimated based on region, land cover, and tree cover fraction using spatiotemporally aligned burned areas from Landsat or Sentinel-2 (Claverie et al., 2018). Because of a new fire detection method that significantly boosts the area of small fires, the CONUS annual burned area increases from 2.36 Mha in GFED4s to 6.04 in GFED5, which is primarily contributed to by the increase of crop fires from 0.83 to 3.09 Mha.

The FireCCI5.1 is obtained as another reference dataset (Chuvieco et al., 2019). FireCCI5.1 maps fires at 250 m resolution using the spectral information from MODIS in combination with the thermal anomalies. FireCCI5.1 has been reported to heavily underestimate the total burned area mainly due to the underrepresentation of small fires (Lizundia-Loiola et al., 2020).

Besides observations, we also obtain burned area from seven state-of-the-art process-based wildfire models participating in the Inter-Sectoral Impact Model Intercomparison Project (ISIMIP3a) (Burton et al., 2024b), including the Canadian Land Surface Scheme Including Biogeochemical Cycles (CLASSIC) (Melton et al., 2020), the Simplified Simple Biosphere model coupled with the Top-down Representation of Interactive Foliage and Flora Including Dynamics model (SSiB4-TRIFFID-Fire) (Huang et al., 2020, 2021), the SPread and InTensity of FIRE (SPITFIRE) coupled with the Organizing Carbon and Hydrology In Dynamic Ecosystems (ORCHIDEE) (Yue et al., 2014), the Joint UK Land Environment Simulator (JULES) coupled with the INFERNO fire model (Mathison et al., 2023; Mangeon et al., 2016), the LPJ-GUESS dynamic global vegetation model coupled with the SPITFIRE (LPJ-GUESS-SPITFIRE) and SIMple FIRE model (SIMFIRE) (Knorr et al., 2016) and BLAZe induced biosphere–atmosphere flux Estimator (BLAZE) (LPJ-GUESS-SIMFIRE-BLAZE) (Rabin et al., 2017), and the Vegetation Integrative Simulator for Trace gases (VISIT) (Ito, 2019). Driven by GSWP3-W5E5 historical climate forcing (Cucchi et al., 2020; Lange et al., 2021), these models provide the monthly burned area at 0.5° spatial resolution from 1901–2019. The multi-model output during 2001–2019 is used in this study. We also performed the benchmarking simulation using the built-in process model in ELM-BGC.

The process-based models differ from one another in not only their dynamic global vegetation models (DGVMs) but also the complexity of their fire models. ELM-BGC and SSiB4-TRIFFID utilized the same fire model from Li et al. (2012), LPJ-GUESS-SPITFIRE and ORCHIDEE, both coupled with SPITFIRE. Other models incorporate their own unique fire modules. The representation of fires over croplands and pastures varies across models (Burton et al., 2024a; Teckentrup et al., 2019). Most models, except for JULES, classify croplands as non-burnable. JULES treats croplands similarly to natural grasslands, while all other models exclude croplands from burning. Most models do not include pasture as a PFT and, therefore, do not distinguish pastures from grasslands in terms of both growth and fire behavior. In LPJ-GUESS-SIMFIRE-BLAZE, pastures are harvested, leading to reduced biomass and consequently a smaller burned area. The difference among process-based models is discussed in Sect. 4.

2.2.2 Surface meteorological, lightning, and socioeconomic datasets

Surface meteorological variables, including temperature, humidity, wind speed, downward shortwave radiation, downward longwave radiation, precipitation, and surface pressure, are obtained from NLDAS-2 (Phase 2 of the North American Land Data Assimilation System), forcing fields to both drive the ELM and construct the training set for the ML fire model. This dataset combines multiple sources of observations (such as precipitation gauge data, satellite data, and radar precipitation measurements) to produce estimates of climatological properties at or near the Earth's surface at hourly temporal resolution and $\mathrm{1}/\mathrm{8}$° grid spacing. We use the temperature, relative humidity, specific humidity, wind speed, and precipitation directly from NLDAS-2 to train the ML fire model. Additionally, we calculate the standardized precipitation evapotranspiration index (SPEI) following Beguería et al. (2014) and the vapor pressure deficit (VPD) based on NLDAS-2 dataset as additional input for the ML model (Table 1). We coarsen this dataset to 0.25° to align with burned area datasets.

In addition to surface meteorological forcing, while identical to those used by ISIMIP3a fire models, we acquire lightning and socioeconomic datasets from multiple sources. The 2 h climatology lightning flash data from NASA LIS/OTD v2.2 at 2.5° resolution are used to calculate the number of natural ignitions. Lightning data are aggregated by summing the 2 h data to derive monthly climatological means, and these monthly climatologies are repeated across all years, disregarding interannual variations. The annual gridded population density data are acquired from Klein Goldewijk et al. (2017), while the GDP per capita is from the World Bank (https://data.worldbank.org/, last access: 20 June 2024), and both are assigned constant values for all months within each corresponding year. All datasets are spatially resampled to a 0.25° × 0.25° grid using bilinear interpolation. To train the ML model, additional inputs, including top-layer soil moisture, LAI, and spatial fraction of each plant functional type (PFT), are simulated by ELM (explained further in Sect. 2.3).

Table 1 Meteorological forcing, land surface properties, and fire-specific inputs for driving the ELM-BGC and training the offline-XGB fire model.

2.3 Model configuration and offline-XGB training and coupling processes

In ELM-BGC, vegetation properties, including canopy height and LAI, vary with carbon allocation and distribution driven by climate variability and disturbances such as wildfires. To bring the model's carbon and nitrogen pools into equilibrium, we first conduct long-term spin-up simulations as suggested by Lawrence et al. (2011). We adopt a two-step approach consisting of a 400-year accelerated decomposition (AD) spin-up followed by a 400-year regular spin-up, driven by cycling NLDAS-2 meteorological forcing from 1981 through 2000. In the AD spin-up, acceleration factors are applied to accelerate decomposition in soil organic matter pools and for plant dead stem and coarse root mortality. The terrestrial carbon pools and vegetation distribution after spin-up simulations reach quasi-equilibrium states after the 800-year simulations.

With the quasi-equilibrium state from the spin-up simulation, we conduct transient simulations with the process-based fire model in the ELM-BGC driven by hourly NLDAS-2 meteorological forcings at a 0.25° resolution from 2001 to 2020. The process-based fire model operates on an hourly basis, matching the frequency of the meteorological inputs, while the ML fire model is trained and applied at a monthly interval consistent with GFED5 data intervals. For training the offline-XGB model, the ELM-BGC outputs, including LAI, surface soil moisture, and PFT fractions, are averaged to monthly intervals combined with monthly mean meteorological conditions, socioeconomic variables (GDP, population density), and lightning (as detailed in Table 1) to learn the relationship between predictors and burned area. To reduce overfitting, the 20-year dataset is split, with 80 % used for training and 20 % for validation. During training, grid cells with fewer than 30 months of non-zero burned area (∼ two-thirds of the total number of grid cells) are masked. This step is important to avoid feeding the ML model distinct predictor combinations that all correspond to zero burned areas, which could skew the model's learning process. Model performance was evaluated based on its accuracy in predicting the spatial distribution and temporal variation of burned areas. Validation metrics included root mean square error (RMSE) and the coefficient of determination (R²).

We then integrate the offline-XGB to ELM-BGC, forming the coupled model ELM2.1-XGBfire1.0. The coupled model runs at 0.25° and hourly resolutions, where the hourly model predictions are accumulated to calculate monthly means. At the end of each month, the ML fire model is called to predict the monthly burned area, updating the land surface properties (e.g., LAI and vegetation height), carbon cycling (biotic carbon in each pool), and ecohydrology processes (photosynthesis and soil moisture) in ELM-BGC.

2.4 Ecoregion

We evaluate the model simulated burned area for each ecoregion adopted from the US Environmental Protection Agency (EPA). Ecoregions are areas where ecosystems (and the type, quality, and quantity of environmental resources) are generally similar (Omernik and Griffith, 2014), and, generally, wildfire properties in each ecoregion are similar. A combination of level I and level II ecoregions is used, and some types have been combined to focus on the broad vegetation distribution. As shown in Fig. 2, the Western Forested Mountains include NW Forested Mountains, Marine West Coast Forests, and Mediterranean California from level I ecoregions. The North American (NA) Deserts include NA Deserts and small portions of Temperate Sierras and Southern Semi-Arid Highlands. The Northeast (NE) Temperate Forests include Mixed Wood Shield, Mixed Wood Plains, and Atlantic Highlands from level II ecoregions. The Southeast (SE) Temperate Forests include Southeastern US Plains Ozark, Ouachita–Appalachian Forests, and Mississippi Alluvial and Southeast US Coastal Plains level II ecoregions.

Figure 2 Ecoregions used in fire model evaluation. 1: Western Forested Mountains, 2: NA Desert, 3: Great Plains, 4: SE Temperate Forests, and 5: NE Temperate Forests.

3 Results

3.1 Evaluation of the burned area spatial distribution

The burned areas across the CONUS exhibit a strong spatial variation (Fig. 3a), which is primarily influenced by climate, vegetation, and human activities. According to the GFED5, the CONUS experiences an average burned area fraction (BAF) of 0.6 % yr⁻¹–0.9 % yr⁻¹ (4.8–7.1 Mha yr⁻¹). The BAF in the WUS (Western Forested Mountains and NA Desert) ranges between 0.4 % yr⁻¹–0.9 % yr⁻¹ (1.1–2.3 Mha yr⁻¹). States like California, Oregon, and Nevada, as well as the Rocky Mountain region, including parts of Colorado and Wyoming, experience large wildfires. The wildfires in the Pacific Northwest and northern California are generally lightning-caused and occur in boreal forests (Balch et al., 2017), whereas those in southern California are primarily caused by human ignition in dry forests and shrublands. The southwest, including Arizona and New Mexico, also sees significant burned areas in shrublands and dry forests. In the Great Plains, states such as Kansas and North Dakota also exhibit high burned areas alongside Texas and Oklahoma, with a BAF ranging between 0.7 % yr⁻¹ and 1.3 % yr⁻¹ (1.6–2.9 Mha yr⁻¹). These high burned areas are primarily contributed by agricultural fires, particularly for cleaning crop residues and managing pastures (Donovan et al., 2020). The southeastern US experiences 0.9 % yr⁻¹–1.5 % yr⁻¹ (1.5–2.6 Mha yr⁻¹) BAF annually, while the temperate forested areas covering Florida, Georgia, and the Carolinas, show lower burned areas compared to the West. The Midwest and northeast exhibit sparse burned areas, with BAF mostly less than 0.16 % yr⁻¹–0.25 % yr⁻¹ (0.2–0.3 Mha yr⁻¹). Burned areas in GFED4s and FireCCI5.1 are much smaller than GFED5 due to the underrepresentation of small fires. The overall spatial distributions are generally consistent across the three datasets, as shown by the high spatial correlation coefficients (R_p).

Figure 3 Observed burned area fraction (% yr⁻¹). (a) GFED5 (2001–2019), (b) GFED4s (2001–2016), and (c) FireCCI5.1 (2001–2019). The numbers indicate the mean (M) burned area fraction and burned area (in Mha) in brackets for each dataset. The pattern correlation (R) against GFED5 is also shown, with an asterisk (∗) denoting significance at the 0.01 level. Black contours outline the ecoregions.

The offline-XGB wildfire model reproduces the burned area distribution over the CONUS well (Fig. 3b), with a R_p of 0.98 (p<0.01) and a bias of −1.0 Mha yr⁻¹. While integrated with ELM, the performance degraded (R_p=0.59, p<0.01, bias = 1.9 Mha yr⁻¹) (Fig. 3d). This degradation is likely due to the fire–vegetation feedback. The aboveground biomass and fuel moisture from ELM-BGC have been used to train the offline-XGB prior to the coupled run within ELM. In the coupled simulation, ELM2.1-XGBfire1.0 updates the biotic carbon and fuel moisture based on the burned area simulated in the previous time step. Consequently, differences in the simulated burned area compared to the process-based models are reflected in the biotic carbon and fuel moisture, accumulating over the 20-year simulation period and influencing the burned area simulation in subsequent time steps.

Figure 4 Same as Fig. 3 but showing model outputs. The pattern correlation (R) and bias (B) against GFED5 are denoted.

In various ecoregions, the offline-XGB model demonstrates minimal biases, and the ELM2.1-XGBfire1.0 consistently outperforms all process-based fire models in predicting annual mean burned area (Fig. 4a–b). The accurate simulation of burned area over the Western Forest Mountains indicates that the ELM2.1-XGBfire1.0 framework generally captures the complex interplays between climate, vegetation, and human activities, with both climate forcings and predicted vegetation status acquired from ELM-BGC. Meanwhile, the ELM2.1-XGBfire1.0 shows superior performance over the Great Plains, indicating that the ML model effectively describes crop fire, thereby utilizing data on crop fraction and LAI.

The performance of the eight process-based fire models in simulating burned areas over the CONUS shows both strengths and weaknesses (Figs. 4c–j and 5). All models generally capture the high burned areas in key regions such as the WUS and southeast US, except for ORCHIDEE showing a concentrated burned area in the Great Plains and LPJ-GUESS-SIMFIRE-BLAZE model missing fires in SE US. However, all process-based models tend to overestimate burned areas in various regions across the CONUS. ELM-BGC has moderate overestimations over the CONUS, with 3.83 Mha yr⁻¹. The burned areas are doubled in CLASSIC, ORCHIDEE, JULES, and VISIT simulations, with values of up to 20.7 Mha yr⁻¹ (Fig. 4a).

Figure 5 Observed and simulated mean burned area fraction ( % yr⁻¹) over the CONUS and ecoregions. The red line in each panel indicates the observed burned area. Modeled burned areas greater than 4 % yr⁻¹ are truncated with the value denoted on the bar.

In the Western Forest Mountains, where fuel is abundant due to dense forest coverage, all process-based models except ORCHIDEE simulate 2 to 5 times of the GFED5 burned area. This overestimation can be related to many factors, including overestimation of fuel combustibility and underrepresentation of anthropogenic fire suppression (Balch et al., 2017). In contrast, wildfires in the NA Desert are primarily constrained by the fuel load. ELM-BGC and CLASSIC produce smaller overestimations, while SSiB4-TRIFFID-Fire, VISIT, JULES, and LPJ-GUESS models significantly overestimate the burned area (4–16 times of GFED5), likely due to overestimations of fuel load, which might be attributed to insufficient water stress on vegetation growth in the arid region (Liu and Xue, 2020; Zhang et al., 2015). Although none of the process-based models accurately capture the spatial distribution of burned area over the Great Plains (Fig. 1), ELM-BGC, SSiB4-TRIFFID-Fire, and VISIT produce comparable burned areas to observations, while CLASSIC and ORCHIDEE overpredict them (4–7 times of GFED5). The inaccurate description of the spatial pattern and large inter-model spread in the Great Plains may be caused by inaccurate treatments of cropland fires and pasture fires (Donovan et al., 2020). As noted by Teckentrup et al. (2019) and Burton et al. (2024a), none of the process-based models has activated the explicit cropland fire model. While LPJ-GUESS-SIMFIRE-BLAZE incorporates harvesting in pastures, reducing biomass and influencing fire dynamics, all other process-based vegetation models do not distinguish pastures from natural grasslands for both vegetation growth and fire processes. Therefore, information on how fuel properties, including the amount as well as physical (e.g., bulk density) and chemical characteristics, and fire ignitions differ between pastures and natural grasslands could help to improve burned area simulation in the process-based fire models (Rabin et al., 2017). Fuel management practices, such as prescribed burning and grazing, can significantly alter fire dynamics but are generally absent in current models. In the eastern US (EUS) forests (Southeast and Northeast Temperate Forests ecoregions), fires are more managed by prescribed burning, leading to fewer uncontrolled extreme wildfires. Although prescribed burning as an additional ignition source is not included in the process-based models, ignition is not a limiting factor in this region due to the abundance of lightning, which provides sufficient natural ignition sources. Consequently, the burned area is primarily controlled by fire spread, which is influenced by natural conditions such as fuel availability and wind, allowing the models to perform well in simulating fire dynamics.

3.2 Evaluation of the burned area temporal variability

We evaluate the model performance in simulating the monthly burned area and depicting fire seasons. Fire season is defined as a monthly burned area greater than $\mathrm{1}/\mathrm{12}$ of the annual total burned area. The CONUS has two fire seasons, i.e., March–April–May and August–September–October, affected by both climate and human activities (Fig. 6a). The WUS fire season spans from early summer to late fall and is primarily determined by the dry conditions and high temperature during these months (Safford et al., 2022; Schoennagel et al., 2017). Specifically, over the Western Forest Mountains, the fire season includes July-to-November (Fig. 6b). Most models capture the July to October fire season, except for ORCHIDEE (May–August). However, only offline-XGB, SSiB4-TRIFFID-Fire, and CLASSIC simulate the peak fire month in August, while others simulate a peak ∼ 1–2 months late. Similar fire season and model performance are observed over the NA Desert (Fig. 5c). In wildfire-dominant regions, the shift in fire peak months might be related to the representation of seasonality in vegetation production and fuel build-up in the BGC model (Hantson et al., 2020).

Figure 6 Monthly mean burned area fraction (% yr⁻¹) over each ecoregion. Vertical shadings indicate the fire seasons, monthly burned area greater than $\mathrm{1}/\mathrm{12}$ of the total burned area, and shadings along x axes indicate 1 standard deviation across the years.

Human activities can also change the timing of fire occurrences (Le Page et al., 2010). Over the Great Plains, pasture fires are conducted during late winter to early spring to control pests, recycle nutrients, and prepare fields for planting (Gates et al., 2017). During the late summer to early fall, crop fires are conducted to clear crop residues. However, these fires may become uncontrolled, leading to larger fires that significantly impact the region. The fire seasons due to pasture fires and crop fires are evident in observations and are captured in offline-XGB and ELM2.1-XGBfire1.0, despite ELM2.1-XGBfire1.0 slightly underestimating the peak in March. Except for LPJ-GUESS-SPITFIRE, none of the process-based models is able to simulate these periods, and instead, a summer fire season is predicted. LPJ-GUESS-SPITFIRE produces peaks in both spring and summer. In SE Temperate Forests, routinely prescribed burns reduce large fire occurrences across the year (Mitchell et al., 2014). The dry condition and/or fallen vegetation fuel larger burned areas in February–March and September–November. The ML-based models generally reproduce the fire seasons in March–April and September–November while none of the process-based models captures the bimodal seasonality. The results of NE Temperate Forests are similar to Great Plains, expect no peak burned area appears in November. The offline-XGB and SSiB4-TRIFFID-Fire models capture the spring peak. To the best of our knowledge, ELM-BGC is one of the few process-based models capable of explicitly simulating crop fires; however, this feature was not enabled in our study. None of the models used here includes explicit representations of pasture burning. Our evaluation suggests that including anthropogenic fires could help to improve model simulations in central and eastern US. However, this requires a better understanding of how fire is used for land management under different socioeconomic and cultural conditions (Pfeiffer et al., 2013; Li et al., 2013).

Over the CONUS, the observed interannual variability (IAV), measured using standard deviation, is 0.7 Mha yr⁻¹, representing 12 % of the annual total burned area in GFED5 (Fig. 7a). GFED4s and FireCCI5.1 suggest 1.1 Mha yr⁻¹ (45 %) and 0.9 Mha yr⁻¹ (30 %), respectively. Process-based models greatly overestimate the IAV, ranging from 2.5 (LPJ-GUESS-SIMFIRE-BLAZE) to 6.6 Mha yr⁻¹ (VISIT). The relative IAV regarding the modeled annual mean value, ranging from 12 % (JULES) to 41 % (ELM-BGC), generally within the range of observations. The machine learning models, offline-XGB and ELM2.1-XGBfire1.0 produce IAV of 0.6 Mha yr⁻¹ (11 %) and 0.8 Mha yr⁻¹ (10 %), respectively.

Despite the magnitude of IAV being amplified by process-based models, after extracting the mean values and dividing by standard deviation, the standardized time series correlated well with the observation (Fig. 7b). Since the modeled IAV is generally influenced by climate variability and the climate-driven fuel variability, both process-based and ML-based models capture the timing of the fluctuations.

Figure 7 Annual total burned area (Mha yr⁻¹). (a) Annual total value and (b) standardized by removing mean and standard deviation.

Monthly temporal variability in burned areas demonstrates significant regional differences across the ecoregions (Fig. 8). Over the entire simulation period, the ML-based models generally capture the timing of wildfires across the CONUS with a temporal correlation coefficient greater than 0.5 (p<0.01), whereas the process-based models exhibit a correlation of only 0.3 (p>0.01). The ML-based models also effectively capture the temporal variability across the ecoregions, although there is a slight decrease in the ELM2.1-XGBfire1.0 in the Great Plains and EUS. This decrease is likely related to the fire–vegetation feedback, which alters the fuel condition differently from the training set. In contrast, the process-based models show correlations comparable to the ML-based models in the WUS but fail to accurately predict burned area temporal variations in the Great Plains and EUS. Again, climatic factors play a dominant role in shaping the temporal variability of BAF in the WUS, while human activities largely influence the BAF in the Great Plains and EUS (Kupfer et al., 2020; Chen et al., 2023). Process-based models tend to better describe the responses of fuel load and combustibility to climate than the responses of fire ignition and suppression to human activities (Hantson et al., 2016).

Figure 8 Monthly correlation coefficient between simulations and GFED5 over each ecoregion.

4 Discussion and conclusion

4.1 Overview of the hybrid framework

This study introduces a hybrid framework integrating an XGBoost wildfire model into an Earth system model (ELM-BGC), resulting in ELM2.1-XGBfire1.0. Both offline and coupled versions of the ML model were evaluated against observations and compared to eight state-of-the-art process-based models. The offline-XGB model significantly reduces burned area biases, particularly in the WUS, while the ELM2.1-XGBfire1.0 model retains the spatial and temporal accuracies with slightly reduced performance. In regions such as the Great Plains and EUS, where human activities are major influences, offline-XGB and ELM2.1-XGBfire1.0 outperform all process-based models.

4.2 Challenges and insights for process-based models

We acknowledge that the simulation biases in process-based models may come from multiple sources. All ISIMIP3a fire models were driven by daily GSWP3-W5E5 forcings at a 0.5° spatial resolution. Differences in forcing data could lead to variations in burned area predictions. However, given that both ELM-BGC and ELM2.1-XGBfire1.0 are driven by the same set of forcings yet produce markedly different burned area predictions, we suggest that limitations in physical understanding may play a dominant role in hindering the performance of the process-based model. By contrast, the ML model incorporates the crop PFT fraction and is trained with data that include agricultural burning, allowing it to capture burning patterns that are often missing or underrepresented in process-based models. Meanwhile, all process-based fire models used in this study have used GFED4s or an earlier versions as a reference for calibration. GFED5 captures significantly more small fires compared to GFED4s, making the CONUS annual burned area increase by 156 %, with crop fire increasing by 240 % (Chen et al., 2023). The inclusion of crop fires is particularly impactful in the CONUS.

The process-based fire models used in this study differ in both fire models and DGVMs. VISIT, JULES, and LPJ-GUESS-SIMFIRE-BLAZE employ the semi-empirical fire models (Thonicke et al., 2001; Pechony and Shindell, 2009; Knorr et al., 2014), in which burned area is calculated without an explicit rate-of-spread model (Hantson et al., 2016). The CLM-Li fire model (Li et al., 2012), a fire model of intermediate complexity, is incorporated into both ELM-BGC and SSiB4-TRIFFID-Fire and partially used in CLASSIC (Melton and Arora, 2016). Consequently, similar performance is observed among these models, although CLASSIC tends to exhibit a larger overestimation. The highly complex SPITFIRE model (Thonicke et al., 2010) provides a more comprehensive description of fire behavior (e.g., fire duration and flame height) and is coupled with ORCHIDEE and LPJ-GUESS to describe the fire impact depending on plant traits (bark thickness and crown height). Although SPITFIRE provides more comprehensive description of fire, it does not outperform other fire models with regard to burned area simulation (Hantson et al., 2020).

With more sophisticated parameterization and fire parameters introduced, more observational analyses are required to understand the mechanism behind and to constrain the parametric uncertainty. The fire–vegetation feedback further complicates this problem, with more complex dynamic vegetation models being slow to reach equilibrium after disturbances. The choice of prescribed or dynamic vegetation could also play a role; note that among all the process-based models, CLASSIC, VISIT, and ELM used prescribed vegetation, while all others used dynamic vegetation. It is noteworthy that parameters involved in wildfire prediction are calibrated to align with the research interests of the institutes developing and managing these models. Advancing the physical understanding of wildfire processes for the CONUS and fine-tuning model parameters towards the new burned area dataset hold the potential to improve model performance (Huang et al., 2020).

4.3 Impact on carbon dynamics and broader application

Although ELM2.1-XGBfire1.0 significantly improves the simulation of burned areas, its impact on terrestrial carbon fluxes remains limited. Within the CONUS, fires primarily affect the terrestrial carbon cycle at localized scales due to the relatively small burned areas. ELM-BGC, for instance, underestimates gross primary production (GPP) by approximately 30 % (figure not shown). With more accurate fire predictions, ELM2.1-XGBfire1.0 helps to slightly reduce this negative bias (less than 1 %). Additionally, while ELM-BGC using prescribed PFT distributions can suppress the effects of fires on the ecosystem, it does not account for fire-induced shifts in vegetation species, where species with greater resistance or fire-adaptive traits may gradually dominate. Nonetheless, the coupling remains valuable, especially when the model is configured at higher resolutions. It is particularly important for evaluating fire-induced tree mortality, post-fire recovery, and fire emissions and their subsequent impacts on air quality, cloud formation, and surface meteorology, particularly when ELM is run as part of the E3SM.

The development and application of ML4Fire-XGB represent a significant step forward in our ability to model wildfire dynamics in regions with complicated interactions between fires, ecosystems, climate, and human activities, bypassing the explicit understanding of physical processes. By incorporating an ML wildfire model into a land surface model, we address the critical need for enhanced predictive capabilities at subseasonal to seasonal scales. Meanwhile, the predictability can adapt to the evolving nature of fire regimes under climate change. This research not only contributes to the scientific community's understanding of fire–ecosystem–climate interactions but also provides a practical tool for policymakers and resource managers engaged in wildfire preparedness and response.