Global Wildfire Susceptibility Mapping Based on Machine Learning Models

Abstract

Wildfires are a major natural hazard that lead to deforestation, carbon emissions, and loss of human and animal lives every year. Effective predictions of wildfire occurrence and burned areas are essential to forest management and firefighting. In this paper we apply various machine learning (ML) methods on a 0.25° monthly resolution global dataset of wildfires. We test the prediction accuracies of four different fire occurrence classifiers: random forest (RF), eXtreme Gradient Boosting (XGBoost), multilayer perceptron (MLP) neural network, and a logistic regression. Our best ML model predicts wildfire occurrence with over 90% accuracy, compared to approximately 70% using a logistic regression. We then train ML regression models to predict the size of burned areas and obtain an MAE score of 3.13 ${\mathrm{km}}^2$, compared to 7.48 ${\mathrm{km}}^2$ using a linear regression. To the best of our knowledge, this is the first study to be conducted in such resolution on a global dataset. We use the developed models to shed light on the influence of various factors on wildfire occurrence and burned areas. We suggest building upon these results to create ML-based fire weather indices.

1. Introduction

In August 2020 a series of lightning strikes ignited hundreds of wildfires across California, which resulted in the largest wildfire in California’s documented history [1]. This catastrophe occurred less than a year after Australia’s “Black Summer”, during which the largest bushfires in the continent’s history burned 11 million hectares [2]. It is believed that climate change is a major cause of these extreme events, as fire weather is becoming more frequent (e.g., [3,4]).

Wildfire danger estimation is an essential tool that stands at the base of forest management and firefighting strategies. However, predicting the occurrence and size of areas burned by wildfires is extremely challenging as fire characteristics depend on numerous factors: meteorological factors, such as temperature, relative humidity and wind speed; fuel characteristics, such as vegetation type and cover; anthropogenic factors; and so forth (e.g., [5,6,7,8]).

The challenge of predicting wildfire occurrence has been drawing the interest of both scholars and forest managers for decades. Attempts to understand the conditions that lead to the ignition and spread of wildfires can be dated at least 150 years back ([9] cited in [10]).

The common method of predicting wildfire occurrence is based on empirical and statistical models. Statistical models are data-driven, and their accuracy relies on the availability and quality of the studied data. Statistical models have been widely used in wildfire science. Research interest in the field includes (but is not limited to) predicting wildfire occurrence (reviewed in [11]), burned areas (e.g., [12]), lightning-related ignitions (e.g., [13]), effects of climate change (e.g., [14]), and more.

Although statistical analyses have made an extraordinary contribution to wildfire science and to forest management, the nature of wildfire danger estimation—which is influenced by nonlinear interactions between numerous factors—has limited the achievements of traditional statistical models [15]. In recent years, machine learning models have become increasingly popular and have been applied in almost every scientific field. This trend, combined with new datasets made available by satellite observations, has introduced machine learning (ML) methods to wildfire science. The nonlinear nature of ML algorithms has been recognized by scholars as potentially beneficial in the complex prediction of wildfire occurrence [15].

Previous scholars built upon state-of-the-art ML models for the purpose of predicting wildfire occurrence and burned areas. A recently published review of ML applications in wildfire science [15] presents a comprehensive summary of the work in the field. These include various advanced ML models such as random forest (RF), AdaBoost, artificial neural networks (ANNs), support-vector machines (SVMs), k-nearest-neighbor (KNN), and more. According to a review of almost 30 papers that compared ML model accuracies to logistic regression, almost all papers found ML models to be preferable [15]. However, most ML studies in the field have been conducted on specific regions and not on global datasets. (These include [16] on Alberta, Canada; [17] on the Galicia region of Spain; [18] on the Greek island of Lesvos; [19] on Lebanon; [20] on Turkey; [21,22] on Australia; [23] on Slovenia; [24] on the central plateau of the Iberian Peninsula; [25] and [26] on the Montesinho Natural Park in Portugal; [27] on Southern China; [28] on Cambodia; [29] on Puerto Rico; [30] on the Canton Ticino in Switzerland; [31] on the Daxinganling Mountains of Northeastern China; [32] on the Kroumirie Mountains in Northwestern Tunisia; [33] on South Central US; and [34] on Chile).

Although these regional wildfire danger estimation models are of great value, they are limited in several aspects. First, regional datasets have a relatively small number of observations. The performance of ML models improves significantly when they are trained on large datasets; thus, applying ML models on limited datasets might not realize the full potential of these models. Second, regional models are, by definition, only appropriate in specific locations. A model trained in one region does not necessarily perform well in other regions. Finally, focusing on one specific region does not provide information on the influence of region-dependent characteristics. For example, population density is known to have a strong impact on wildfire occurrence and burned areas; a model which is focused on a single region in which the population density is approximately constant would not identify the significance of this important factor.

In recent years, data from satellites combined with advanced image processing techniques have made global datasets with millions of wildfire observations available [35,36,37]. Such large datasets provide an opportunity to improve the predictions of current ML models and to accurately identify the most hazardous fires. Analyzing a global dataset also enables a location-dependent prediction that can improve current fire weather alerts. Recent global analyses include [38], which performed a statistical analysis on a global dataset with an annual resolution; [39], which applied ML regression models on a monthly 1.89° × 0.25° resolution global dataset to predict the size of regional burned areas; and [40], which applied neural network classification models on a monthly 0.25° × 0.25° resolution global dataset to create wildfire susceptibility maps.

In this paper, we apply various ML models to predict both fire occurrence and burned areas based on meteorological data, fuel characteristics, topography, anthropogenic factors, and regional fire history. In addition to the raw meteorological data, we use the commonly applied fire weather indices as predictors. The dataset we use includes a million spatiotemporal observations of monthly 0.25° × 0.25° regions [37].

The objectives of the study are as follows: (a) to develop and evaluate the prediction accuracies of ML models and linear/logistic regression (LR) models in the prediction of wildfire occurrence and burned areas; (b) to develop models which are not region-dependent and are valid in different climate zones; (c) to analyze the effect of various variables on wildfire occurrence and burned areas; (d) to compare and evaluate the performance of the developed models in different regions around the globe; and (e) to demonstrate the potential of ML models in wildfire danger estimation as a step towards ML-based fire weather indices.

Although [40] is similar to the current study in some respects, our study differs from it in several aspects: first, [40] focused on classification models, whereas in the current study we also develop regression models to estimate the size of burned areas; second, in the current study we include several additional factors which we believe are important to wildfire danger estimation (these include the traditional fire weather indices, topographic slope, population density, and long-term values of meteorological factors such as yearly precipitation); and finally, in the current study we include RF and XGBoost models in addition to neural networks, as these models have been shown to be advantageous for tabular data analysis (e.g., [41,42]).

The paper is organized as follows. We begin by providing a detailed summary of our data, followed by a description of the ML methods we apply in the paper. The next two sections present the models’ principal results, divided into prediction of wildfire occurrence and prediction of burned area size. In the two final sections, we discuss the contribution, implications, and limitations of the study, and propose several directions for future research.

2. Materials and Methods

2.1. Data

We used the Terra MODIS satellite data (FireCCI51) to create the dependent variable in the study. The data were obtained from the ECMWF website [37,43]. The dataset is comprised of monthly burned areas in 0.25° × 0.25° regions across the entire globe for the year 2015. The original MODIS product was at a 250 m resolution and was aggregated to a 0.25° grid by ECMWF by summing the burned areas in each pixel. We eliminated regions where no fires occurred in the decade before 2015, including oceans and seas. Including regions in which wildfires seldom occur would lead to a higher predictive performance, but would not provide a clear representation of the model’s performance in the regions in which wildfire danger estimation is of interest. The remaining data were composed of approximately a million (986, 606) spatiotemporal observations in 121,476 regions throughout the globe, including 128,618 observations with non-zero burned areas. Figure 1 and Figure 2 presents descriptive statistics of the data.

The features (independent variables) include meteorological factors, fuel characteristics, topography, anthropogenic factors, and regional fire history. Most meteorological data were taken from the ERA5 monthly averaged reanalysis dataset [44]. We used 2 m temperature, precipitation, and relative humidity and 10 m wind velocity and direction. For precipitation and relative humidity, we included both present and past value variables: mean value in the month of the observation, mean value in the previous month, and mean value in the previous year. We also included a variable for incoming short-wave solar radiation, obtained as a monthly mean in a 0.25° × 0.25° resolution [45].

For each month we include the mean monthly value and the mean of the highest seven values of various fire weather indices. These indices include three groups: (1) Canadian Forest Service’s Fire Weather Index Rating System; (2) the Australian McArthur Mark 5 Rating System; (3) U.S. Forest Service’s National Fire Danger Rating System. The variables in each group include: (1a) fire weather index; (1b) build-up index; (1c) danger index; (1d) drought code; (1e) duff moisture code; (1f) initial fire spread index; (1g) fine fuel moisture code; (1h) fire daily severity rating; (2a) Keetch–Byram drought index; (2b) fire danger index; (3a) spread component; (3b) energy release component; (3c) burning index; (3d) ignition component. All data are available in 0.25° resolution and were obtained from the Copernicus Climate Change Service [46].

Previous studies have shown that population density has a substantial effect on the area burned by wildfires (e.g., [39]). We included population density based on the dataset of the Center for International Earth Science Information Network [47]. Although the original dataset is provided with a resolution of ~1 km, we calculated the mean population density in a 0.25° resolution region.

Leaf area index (LAI) is a variable that describes the leaf material in a given location. LAI is a dimensionless variable that varies between 0 and approximately 10. LAI data at a 1/112° (~1 km) resolution were taken from [48] and aggregated to the relevant resolution. In addition, we used the normalized difference vegetation index (NDVI), a dimensionless parameter which is commonly used to estimate the density of live green vegetation. NDVI is calculated as the difference between near-infrared (NIR) and red reflectance, divided by their sum [49]: NDVI = (NIR − Red)/(NIR + Red). We include the NDVI value of the 0.25° region of the wildfire, obtained from the NASA Earth Observations website [49]. Finally, we also included the percent of burnable area in each pixel taken from [37].

We included a variable of soil moisture from the Copernicus Climate Change Service [50]. The variable indicates the content of liquid water in a surface soil layer at a depth of 2 to 5 cm, expressed as the percentage of total saturation derived from satellite-based observations [50]. The data were obtained in a 0.25° spatial resolution and monthly mean temporal resolution.

Topography is known to affect the growth rate of a wildfire, and consequently, the area burned by the fire (e.g., [51]). We included the mean slope in each region, based on the dataset in [52].

Regional fire history was obtained from the ECMWF dataset that was also used for the burned area target variable [37,43]. For each region we created features of the mean and median monthly burned area from 2003 up to the month prior of that observation. These variables did not include the month of the observation to prevent data leakage. As the wildfire history variables are not always used in wildfire danger estimation models, we presented both the models which include them as independent variables and models which do not. Examining the models which do not use the wildfire history as predictors could also shed light on the effect that different features have on the models—the wildfire history variables could potentially conceal the effect of certain regional characteristics.

As a categorical variable, month of the year was transformed by 1-of-C dummy encoding [53].

Table 1. Summary of data sources.

Variable	Abbreviation	Source
monthly burned area (target variable)	BA	[37]
2 m temperature	temp	[44]
relative humidity	RH
10 m wind velocity	wind_speed
precipitation	prec
mean relative humidity in previous month	RH_1_month
mean precipitation in previous month	prec_1_month
mean relative humidity in previous year	RH_12_months
mean precipitation in previous year	prec_12_months
percentage of burnable area	burnable	[37]
median burned area	median_burned
mean burned area	mean_burned
latitude	lat	-
longitude	lon	-
month (categorical)	month_1, month_2, etc.	-
leaf area index—low vegetation	LAI_low	[48]
leaf area index—high vegetation	LAI_high
total leaf area index	LAI_tot
normalized difference vegetation index	NDVI	[49]
incoming short-wave solar radiation	radiation	[45]
soil moisture	soil	[50]
mean slope	slope	[52]
population density	population	[47]
fire weather index—mean	FWI_mean	[46]
fire weather index—highest 7 in month	FWI_7
build up index—mean	BUI_mean
build up index—highest 7 in month	BUI_7
danger index—mean	danger_mean
danger index—highest 7 in month	danger_7
drought code—mean	drought_mean
drought code—highest 7 in month	drought_7
duff moisture code—mean	DM_mean
duff moisture code—highest 7 in month	DM_7
initial fire spread index—mean	ISI_mean
initial fire spread index—highest 7 in month	ISI_7
fine fuel moisture code—mean	FFMC_mean
fine fuel moisture code—highest 7 in month	FFMC_7
fire daily severity rating—mean	severity_mean
fire daily severity rating—highest 7 in month	severity_7
Keetch–Byram drought index—mean	KBDI_mean
Keetch–Byram drought index—highest 7 in month	KBDI_7
fire danger index—mean	FFDI_mean
fire danger index—highest 7 in month	FFDI_7
spread component—mean	SC_mean
spread component—highest 7 in month	SC_7
energy release component—mean	energy_mean
energy release component—highest 7 in month	energy_7
burning index—mean	BI_mean
burning index—highest 7 in month	BI_7
ignition component—mean	IC_mean
ignition component—highest 7 in month	IC_7

presents a summary of the variables used in the study.

2.2. Methodology

We developed two types of ML models. The first model differentiates observations to a burned/unburned binary classification. One important issue to address in this regard is that the full dataset is imbalanced: the number of negative (unburned) observations is approximately seven times higher than the number of positive (burned) observations. Not balancing the data would lead to higher prediction accuracies but would not reflect the actual model performance. For example, if 99% of the observations were unburned, simply guessing that all pixels are unburned would lead to a 99% accuracy. We applied two common methods to address this issue. The first method is to balance the data by including all of the burned observations and randomly sampling an equal number of unburned observations (as performed in [54], for example). The second method is to use the entire imbalanced dataset but to impose a penalty on misclassifications of the minority class—in this case, positive observations (e.g., [55]). It is most common to impose a penalty value equal to the ratio between the number of majority and minority observations. We examined this value as well as several lower and higher values.

Since different studies refer to datasets with different spatial and temporal resolutions, it is difficult to indicate a single benchmark accuracy score. It is more challenging to predict fire occurrence in datasets with a finer resolution, and lower accuracies are to be expected for their models. In the absence of a global accuracy benchmark, we followed [15], which found that a considerable number of papers in the field use the linear or logistic regression results as a benchmark for a specific dataset.

Different studies use different accuracy metrics for binary classification, two of which are the accuracy metric, which is the percentage of correctly classified observations, and the area under the curve (AUC) metric, which is the area under the relative operating characteristic (ROC) curve [56]. The AUC metric is preferable to accuracy in binary classifications [57]. A random classification would produce a 0.5 AUC score; AUC scores between 0.5 and 0.7 are considered poor prediction accuracies; scores between 0.7 and 0.9 are considered moderate; and AUC scores of 0.9 or above are considered excellent prediction accuracies [27,58]. It is recommended to present more than a single accuracy metric [59]. We presented several additional metrics: accuracy, true positive rate (TPR), and true negative rate (TNR). These metrics can be defined using the four components of the confusion matrix: true positives (TP), false positives (FP), true negatives (TN), and false negatives (FN). TPR (also known as sensitivity) describes the percentage of correctly classified positive observations, whereas TNR (also known as specificity) represents the percentage of correctly classified negative observations.

$$\mathrm{accuracy}=\frac{\mathrm{T}\mathrm{P}+\mathrm{T}\mathrm{N}}{+\mathrm{FP}+\mathrm{TN}+\mathrm{FN}}$$

(1)

$$\mathrm{TRP}=\frac{\mathrm{TP}}{\mathrm{TP}+\mathrm{FN}}$$

(2)

$$\mathrm{TNR}=\frac{\mathrm{TN}}{\mathrm{TN}+\mathrm{FP}}$$

(3)

The second group of ML models we developed is of regression models which predict the size of burned areas. We evaluated the accuracy of these models by measuring the mean absolute error (MAE), as is common in similar papers (e.g., [25,60,61]). We also presented the mean-square error (MSE) and the root-mean-square error (RMSE) for each model. It is recommended to present these metrics alongside MAE, as MAE is known to be less sensitive to outliers [62].

We included maps describing the prediction accuracies of the testing data for both the classification and the regression models. In the classification models, the values in these maps represent the mean prediction accuracy in each region. For the regression models, these maps describe the MAE scores in each region.

We applied four different classification and regression models: (i) random forest (RF) [63] (ii) eXtreme Gradient Boosting (XGBoost) [64] (iii) multilayer perceptron (MLP), which is a form of neural network [65], and (iv) logistic regression [66] for fire occurrence or linear regression for fire size. We performed a 10-fold cross-validation to reliably compare the accuracies of the different models. The analyses were performed using Python’s scikit-learn package [67], except for the eXtreme Gradient Boosting model which is based on the XGBoost package [64]. As mentioned in the previous section, we reran the models with or without the regional wildfire history variables.

The RF, XGBoost, and MLP models are tuned by various hyperparameters. We performed hyperparameter optimization for both models to achieve optimal predictions. The following hyperparameters were examined for the RF and XGBoost: number of estimators (“n_estimators”) between 100 and 550 and maximal tree depth (“max_depth”) between 8 and 10. The accuracies presented in the Results section are of the optimal hyperparameters. The following hyperparameters were examined for the MLP model: number and size of hidden layers (“hidden_layer_sizes”)—(1–3 hidden layers of 50–150 neurons).

To improve the model interpretability, we analyzed and present the SHAP (SHapley Additive exPlanations) values of the most important features in the highest performance model. For each feature, SHAP values were calculated by comparing the predictions without the feature (assigning it some baseline value, usually its mean value) with the predictions including the feature (for example, [68]). Each dot in each feature represents the feature’s effect on the prediction in a specific observation. The horizontal axis represents the effect on the dependent variable, whereas the colors represent the values of the different independent variables. The independent variables were ordered from the most significant feature to the least significant one. The SHAP graph representation enables a clear interpretation of the effect each variable had on the model. The colors in the graph provide information on the sign of each feature’s impact on the model’s output. For example, if the red dots (high values) of a specific feature are strongly concentrated around the negative side of the horizonal axis, this would indicate a negative correlation between this feature and the dependent variable, and vice versa.

3. Results

3.1. Wildfire Occurrence

In this section we present the results of wildfire occurrence prediction models. We first present the results of various classification models which predict whether a wildfire will occur in each observation. We show the results for the balanced dataset, as well as the analysis of the full weighted dataset. We present both the various metrics described in the previous section (AUC, accuracy, TPR, and FPR), and the ROC curves, and finally an analysis of feature importance.

Table 2. Summary of model performances—wildfire occurrence.

Including Regional Wildfire History
Model	AUC	Accuracy	TPR	TNR	Hyperparameters Tested	Best Parameters
RF	0.92	0.83	0.81	0.86	max depth: 8–10	10
					n_estimators: 100–550	550
XGBoost	0.97	0.92	0.90	0.93	max depth: 8–10	10
					n_estimators: 100–550	550
MLP	0.69	0.67	0.45	0.83	hidden layers: 1–3	1
					# neurons in layer: 50–150	100
LR	0.73	0.69	0.52	0.79	-	-
Excluding Regional Wildfire History
Model	AUC	Accuracy	TPR	TNR	Hyperparameters Tested	Best Parameters
RF	0.89	0.80	0.77	0.83	max depth: 8–10	10
					n_estimators: 100–550	550
XGBoost	0.94	0.86	0.86	0.86	max depth: 8–10	10
					n_estimators: 100–550	400
MLP	0.90	0.80	0.81	0.84	hidden layers: 1–3	2
					# neurons in layer: 50–150	150
LR	0.81	0.73	0.71	0.76	-	-

Optimized prediction scores of the different models. The hyperparameters tested for each model are presented, alongside the chosen values.

and Figure 3 present the different prediction scores and ROC curves for the different models. The purple line represents a random classifier with an AUC score of 0.5. The green line represents a logistic regression which serves as the benchmark accuracy as a non-ML model and achieves an AUC score of 0.73. The best performing models are the XGBoost model, followed by the RF classifier with excellent AUC scores of 0.97 and 0.92, respectively. The XGBoost model was also superior to the other ML and LR models in terms of accuracy, TPR, and FPR. The MLP model performed worse than the LR model. These results address objective (a) by comparing the performances of ML and LR models.

To address objective (b) and verify that the ML models we developed provide accurate wildfire danger estimations which are not region-dependent, we present the prediction accuracies by region for the XGBoost and LR models in Figure 4. The XGBoost model outperformed the LR model in all of the regions around the globe.

Table 3. Summary of model performances—weighted wildfire occurrence.

Including Regional Wildfire History
Model	AUC	Accuracy	TPR	TNR
XGBoost (w = 1)	0.94	0.93	0.59	0.98
XGBoost (w = 5)	0.94	0.89	0.81	0.91
XGBoost (w = 7)	0.94	0.88	0.85	0.88
XGBoost (w = 9)	0.94	0.86	0.87	0.86
XGBoost (w = 20)	0.94	0.80	0.92	0.78
XGBoost (w = 50)	0.93	0.72	0.95	0.69
LR (w = 1)	0.64	0.87	0.11	0.98
LR (w = 5)	0.77	0.75	0.71	0.75
LR (w = 7)	0.78	0.70	0.78	0.69
LR (w = 9)	0.79	0.66	0.83	0.63
LR (w = 20)	0.80	0.46	0.95	0.38
LR (w = 50)	0.69	0.14	0.99	0.01
Excluding Regional Wildfire History
Model	AUC	Accuracy	TPR	TNR
XGBoost (w = 1)	0.92	0.92	0.52	0.98
XGBoost (w = 5)	0.92	0.88	0.78	0.89
XGBoost (w = 7)	0.92	0.85	0.82	0.86
XGBoost (w = 9)	0.92	0.84	0.85	0.84
XGBoost (w = 20)	0.92	0.76	0.91	0.74
XGBoost (w = 50)	0.92	0.67	0.95	0.62
LR (w = 1)	0.81	0.89	0.26	0.98
LR (w = 5)	0.81	0.80	0.64	0.83
LR (w = 7)	0.81	0.74	0.72	0.74
LR (w = 9)	0.81	0.67	0.78	0.66
LR (w = 20)	0.81	0.43	0.93	0.36
LR (w = 50)	0.81	0.22	0.99	0.10

Prediction scores for the weighted classification models applied on the entire (unbalanced) dataset. The value marked as “w” represents the weighting ratio between the positive and negative prediction error penalties. The ratio between negative and positive observations was approximately 7 (rows marked in bold), but additional weighting values were examined.

and Figure 5 present the AUC scores and ROC curve in the full imbalanced dataset. We present the results for the highest performing ML model, XGBoost, compared to a logistic regression. As expected, in terms of accuracy some of the XGBoost models outperformed the previous analysis for the reasons mentioned in the Methods section; however, in terms of AUC the current models obtained a lower score in comparison to the balanced data (0.94 versus 0.97 when regional fire history is included, or 0.92 versus 0.94 otherwise).

The XGBoost model obtained the highest score on this dataset with a weighting ratio of 7, approximately equal to the ratio between negative and positive observations in the dataset. The XGBoost models outperformed the logistic regressions in all weighting values.

To address objective (c), we applied the SHAP-value analysis which provides a better understanding of the most influential features in the models. Figure 6 presents the SHAP values of the most important features in the XGBoost model of the balanced dataset. For each feature, SHAP values were calculated by comparing the predictions without the feature with the predictions including the feature (for example, [68]). Each dot in each feature represents the feature’s effect on the prediction in a specific observation. The color of each observation indicates the value of the variable in that observation, and its position on the horizontal axis indicates its effect on the model output. For clarity of viewing the figure, 5000 observations were randomly chosen from the entire dataset.

When available (left subplot), the historical mean burned area had the highest impact on the model. This result is in line with the prediction performances of these two models—the model which did not include regional fire history obtained a lower AUC score (0.94 versus 0.97). The most significant meteorological factor in both models was RH, both in the month of the observation and in the previous year. Precipitation was strongly and negatively correlated with the dependent variable in both models. Although the current precipitation was negatively correlated with the probability of burning, the 12 months of precipitation had an opposite effect, as expected. The NDVI index had a substantial effect on both models. Some fire weather indices were found to be influential; the variable describing the seven highest FWI values in each month (FWI_7) had the highest effect in the left subplot. Regional characteristics such as population and slope had a higher impact on the model that did not include the regional wildfire history variables (right subplot).

3.2. Size of Burned Areas

In this section we present the results of the burned area prediction models. We present the results of various ML regression models which estimate the size of burned areas for each observation, compared to a linear regression as a benchmark accuracy.

Table 4. Summary of model performances—size of burned areas.

Including Regional Wildfire History
Model	MAE	RMSE	MSE	Parameters Tested	Best Parameters
	$(k{m}^2$)	$(k{m}^2$)	$(k{m}^4$)
RF	3.44	16.21	262.87	max depth: 8–10	10
				n_estimators: 100–550	550
XGBoost	3.13	14.30	204.61	max depth: 8–10	10
				n_estimators: 100–550	100
MLP	4.78	21.34	455.63	Hidden layers: 1–3	2
				# neurons in layer: 50–150	150
Linear Regression	7.48	21.28	452.94	-	-
Excluding Regional Wildfire History
Model	MAE	RMSE	MSE	Parameters Tested	Best Parameters
	$(k{m}^2$)	$(k{m}^2$)	$(k{m}^4$)
RF	4.08	17.78	316.15	max depth: 8–10	10
				n_estimators: 100–550	400
XGBoost	3.75	15.67	245.50	max depth: 8–10	10
				n_estimators: 100–550	100
MLP	3.90	17.11	292.62	Hidden layers: 1–3	3
				# neurons in layer: 50–150	150
Linear Regression	7.51	22.03	485.52	-	-

Optimized prediction scores of the different regression models. The hyperparameters tested for each model are presented, alongside the chosen values.

presents the MAE values for the different models. Once again, the XGBoost model outperformed all the other models including the RF, MLP, and linear regression (MAE scores of 3.13, 3.44, 4.78, and 7.48, respectively). The MSE and RMSE scores were higher than the MAE scores in all models, but lower in the XGBoost model compared to the other models. In similarity to the fire occurrence analysis, removing the regional fire history variables reduced the accuracy of the highest performing model (XGBoost) and increased the MAE score from 3.13 to 3.75. The LR models had significantly lower accuracies compared to the various ML models in both types of analyses.

To compare and evaluate the performance of the developed models in different regions around the globe (objective (d)), we present the prediction accuracies by region in Figure 7. Unlike the fire occurrence analysis, the advantage of the XGBoost model over the LR model was not uniform in all continents. Its performance was only slightly better in Australia and North America, as opposed to a major advantage in Africa, Asia, and South America.

Figure 8 presents the SHAP values for the XGBoost regression model. Similar to the fire occurrence analysis, the historical mean regional wildfire variable was the most influential when included. RH and precipitation were once again the most important meteorological factors. In the current analysis, the fire weather indices were considerably more influential compared to the wildfire occurrence models. Specifically, the mean monthly FWI value was ranked high in both models. Additional indices such as ISI and FFMC were also highly influential in the model excluding wildfire history. Similar to the wildfire occurrence analysis, the NDVI was the most important vegetation index and had a higher impact than the LAI variables and the percentage of the burnable area.

4. Discussion

In this paper we applied multiple ML models to predict the occurrence and size of wildfires around the globe. The models included random forest, XGBoost, multilayer perceptron, and logistic and linear regressions. The models were trained using a large dataset which includes wildfire observations from around the globe over a full year. Once trained, the best models provided promising prediction accuracies and were able to predict wildfire occurrence with over 90% accuracy and burned area size with an MAE score of 3.13 ${\mathrm{km}}^2$. The XGBoost model obtained the most accurate results, followed by random forest. All models showed a substantial improvement in prediction accuracies after being trained on a large dataset with a million observations.

The prediction accuracies varied significantly between different regions of the globe. In the wildfire occurrence prediction, the XGBoost classification model outperformed the logistic regression in all the examined regions. As for burned area predictions, however, the improved performance of the XGBoost regression model was most significant in Africa, Asia, and South America. Its advantage over the linear regression in North America, Australia, and Europe was significantly smaller. One possible explanation is that North America and Australia are the central regions for which the traditional fire weather indices were planned, and thus the simple linear model which includes the traditional indices provides a relatively good prediction accuracy in these regions. The advantage of the ML model is most distinct in other regions, where fire weather indices were not originally planned to be used.

As expected, the most important predictor of wildfire occurrence and burned area size was the regional wildfire history, with the mean burned area variable providing higher predictive performance than the median burned area. These variables have not been used in most studies in the field [16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34]; although they provide little insight as to the factors which dominate wildfire occurrence, they improve the performance of the models, so we recommend including them in future studies and in models used for practical applications. Relative humidity, both in the month of the observation and in earlier times, was found to be the most important meteorological factor, followed by precipitation. This is in line with most previous studies in the field (e.g., [33]). As expected, the longer-term precipitation and RH had the opposite effect than the values of these factors in the month of the observations. Although precipitation substantially decreases the risk of wildfires, its occurrence in earlier times is correlated to increased vegetation growth and, therefore, increased wildfire risk (e.g., [69]).

NDVI was also an especially influential factor for both types of models and had a substantially larger effect on the models compared to LAI or the percent of burnable area. Whereas some previous studies (e.g., [40]) included both NDVI and LAI in their models and found NDVI to have a higher contribution, some papers only included LAI as an independent variable (e.g., [33]). We suggest including both in future studies. Other regional factors were substantially less significant than expected, especially population density and topographic slope. This is an unexpected result, as previous studies have found a strong negative correlation between population density and burned areas (e.g., [39]) and a positive correlation between topographic slope and wildfire propagation (e.g., [51]). The most probable explanation is that these factors are steady and do not change significantly over time, and their effect is therefore reflected in the mean regional burned area variable. For example, a steeper slope is known to increase the wildfire propagation rate, and thus regions in which the mean slope is high also suffer from higher historical burned areas.

The models which did not include the wildfire history variables mostly had somewhat lower prediction accuracies, both in the classification analysis and in the regression analysis. The ML models still significantly outperformed the logistic and linear regressions. When examining the most influential variables in these models, the relative importance of some regional variables increased in comparison to the models which included the wildfire history variables. Namely, topographic slope, elevation, and population density had a much larger effect on the models which did not include the regional wildfire history (Figure 6 and Figure 8). This is probably because the regional wildfire history variables already reflected certain regional characteristics with little or no temporal variability, and thus including them concealed the effects of the additional regional characteristics.

Different fire weather indices were important for both types of models. Among the various indices, the Canadian FWI index and the FFMC index had the most substantial effect on the models. Both the seven highest values in each month and the mean values of the entire month had an impact on the model, depending on the index and the specific model. Although some similar studies include fire weather indices as independent variables and find them to be beneficial (e.g., [26]), many other studies do not include them and instead only include the raw meteorological factors. As these indices clearly contribute to wildfire danger estimation, we suggest including them in future models to build on the accumulated knowledge in the field.

Accurate wildfire danger estimation by ML models presents a promising opportunity to improve wildfire alerts and provide forest managers with tools of assessing regional wildfire risk. We propose to build upon the results of this research and develop regional wildfire indices that take local fire history into account. Such indices would require training ML models on a dataset of daily resolution, as opposed to the monthly resolution of the current study. The models developed in this study could benefit from validation on additional wildfire datasets. The monthly resolution of the current study limits its applicability for fire weather indices, as these require a daily resolution. In addition, it is very likely that the results and the most influential features in a daily resolution model would differ from the current one. For example, in the current study mean monthly wind velocity is one of the least significant factors, whereas previous studies have demonstrated the importance of wind velocity on the day of the fire. Regional fire indices could vary based on local firefighting abilities, fuel characteristics, and climate. The disadvantage of developing regional fire weather indices is that the data in each region are more limited than data from the entire globe. However, as satellites continue to gather information and global datasets continue to grow, it is a matter of time until wildfire datasets will be large enough to train region-specific ML models which would probably provide improved prediction accuracies and could potentially replace the traditional fire weather indices.

5. Conclusions

Wildfire occurrence and burned areas are extremely difficult to predict as they depend on nonlinear interactions between numerous factors. In this study we compared the prediction performance of ML models to LR models on a global dataset. We demonstrated the overwhelming advantage of ML models over LR models in these tasks (objective (e)). Working with a global dataset provides vast amounts of data for model training and validation, and in addition, allows researching the effect of various factors in a global context, rather than region-dependent phenomena. We analyzed and presented the main contributing factors of each model and their relative importance. We found that LR models which include traditional fire weather indices perform well in North America and Australia, but less so in other continents. The ML models, however, obtain excellent predictive performance worldwide. The results of this study demonstrate the advantage of ML models over traditional fire weather indices in wildfire danger estimation. We propose to work towards the application of ML-based fire weather indices, which we believe will gradually replace the traditional indices.