Modeling and Mapping Forest Fire Occurrence from Aboveground Carbon Density in Mexico

: Understanding the spatial patterns of ﬁre occurrence is key for improved forest ﬁres management, particularly under global change scenarios. Very few studies have attempted to relate satellite-based aboveground biomass maps of moderate spatial resolution to spatial ﬁre occurrence under a variety of climatic and vegetation conditions. This study focuses on modeling and mapping ﬁre occurrence based on ﬁre suppression data from 2005–2015 from aboveground biomass—expressed as aboveground carbon density (AGCD)—for the main ecoregions in Mexico. Our results showed that at each ecoregion, unimodal or humped relationships were found between AGCD and ﬁre occurrence, which might be explained by varying constraints of fuel and climate limitation to ﬁre activity. Weibull equations successfully ﬁtted the ﬁre occurrence distributions from AGCD, with the lowest ﬁt for the desert shrub-dominated north region that had the lowest number of observed ﬁres. The models for predicting ﬁre occurrence from AGCD were signiﬁcantly di ﬀ erent by region, with the exception of the temperate forest in the northwest and northeast regions that could be modeled with a single Weibull model. Our results suggest that AGCD could be used to estimate spatial ﬁre occurrence maps; those estimates could be integrated into operational GIS tools for assistance in ﬁre danger mapping and ﬁre and fuel management decision-making. Further investigation of anthropogenic drivers of ﬁre occurrence and fuel characteristics should be considered for improving the operational spatial planning of ﬁre management. The modeling strategy presented here could be replicated in other countries or regions, based on remote-sensed measurements of aboveground biomass and ﬁre activity or ﬁre suppression records.

Globally, the country of Mexico ranks highly among areas of the largest fire activity [50]. Mexico has experienced a marked increase in the fire weather season due to climatic change [7]. Furthermore, fire activity in the North American region is projected to increase under climate warnings [6]. Therefore, a detailed knowledge of the spatial distribution of fire occurrence is increasingly necessary in the country to inform strategic fuel and fire management planning decisions. Moreover, from the scientific perspective, Mexico provides an opportunity to analyze the variability of wildfire occurrence and explore potential relationships to biomass by ecoregion, because it hosts a wide range of vegetation types and climate, from arid deserts dominated by scrubby vegetation, temperate pine-oak forests and tropical forests [51], that results in a variety of fire regimes in the country [52][53][54].
Previous studies in Mexico have related fire occurrence with vegetation type, mainly at a local scale [55][56][57]. Some studies have documented differences in fire occurrence by ecoregion at a national level [58], but, in spite of medium resolution aboveground biomass or carbon density maps from satellite images being available at national [59,60] and at local scales [61][62][63][64], no study has attempted to relate the spatial patterns of fire occurrence with aboveground biomass in Mexico.
The lack of an operational system for fire risk and danger mapping in Mexico, unlike countries like US or Canada [65][66][67][68][69], led to a national research project aiming at developing a fire danger forecast system for Mexico [70,71]. Previous studies have resulted in the development of temporal and spatial models for predicting the number and spatial distribution of wildfires from daily fuel dryness indices at 1 km resolution by vegetation types and regions in the country [72][73][74], but there is still no information on how these spatial patterns of fire occurrence might be affected by aboveground biomass within different vegetation types and ecoregions, limiting our knowledge of fire ecology in these variety of ecosystems and preventing the spatially explicit decision-making for fuel treatment and fire suppression strategic planning in Mexico.
The current study aimed at exploring relationships between aboveground biomass and fire occurrence for the main ecoregions of Mexico. The specific objectives of the study are: (i) To develop models for predicting the percentage of fire suppression records from aboveground biomass-expressed as carbon density-for the main ecoregions of Mexico; (ii) to map the expected percentage of fires by ecoregion from aboveground biomass for Mexico, at 1 km resolution.

Study Area and Forest Regions
The regionalization of the country was based on the regions proposed by [72], which was based on the North American ecoregions map (North American Level 3 Ecoregions Map, EPA, https://www.epa.gov/eco-research/ecoregions-north-america), together with previous analysis of the spatial patterns of fire occurrence in the country [70,75]. In addition, in this study we included two additional regions, north scrubland and chaparral, based on the land use map of the National Institute of Geography and Statistics (INEGI in Spanish, Land Use Map Series V, 1:250,000 http: //www.inegi.org.mx/geo/contenidos/recnat/usosuelo/) [51].
This resulted in a total of seven regions for analysis: CHAP: Chaparral forests; N: North scrubland forests; NW: Northwest forests; NE: Northeast forests; C: Center forests; SC: South-center forests and SE: Southeast tropical forests ( Figure 1). Note that here the term "forest" is used in a broad sense. In fact, the chaparral and north forest comprise a mix of shrubby and wooded vegetation, the former being by far the dominant cover in both regions.

Aboveground Carbon Density Data (AGCD)
The aboveground carbon density map (AGCD) from the study of Cartus et al. [59] was utilized, because of its public availability. This map was obtained from Landsat images, radar data from ALOS PALSAR, calibrated with ground data of 16,906 plots from the Mexican National Forest and Soil inventory [62]. Species-specific biomass allometries and biomass to carbon conversion factors-ranging from 0.440 to 0.516-were utilized [59].
The AGCD map was rescaled from 30 × 30 m to 1000 × 1000 m wall to wall for scale matching with the layers used in the ongoing Mexican Fire Danger System [70,72]. A low pass filter with a moving window of 3 × 3 pixels, using the package included in Spatial module of ERDAS IMAGINE 2014 software (ERDAS Inc., Atlanta, GA, USA), [76] was applied to the AGCD map.

Forest Fire Data
The fire suppression records database from the Mexican Government Forest Agency (CONAFOR) from the period 2005 to 2015 was utilized. The forest fires records in water bodies, cities and agriculture were excluded for the analysis based in the land use map from INEGI (2011) [51]. The database comprised a total of 53,000 fire suppression records in the period of study.

Aboveground Carbon Density Data (AGCD)
The aboveground carbon density map (AGCD) from the study of Cartus et al. [59] was utilized, because of its public availability. This map was obtained from Landsat images, radar data from ALOS PALSAR, calibrated with ground data of 16,906 plots from the Mexican National Forest and Soil inventory [62]. Species-specific biomass allometries and biomass to carbon conversion factors-ranging from 0.440 to 0.516-were utilized [59].
The AGCD map was rescaled from 30 × 30 m to 1000 × 1000 m wall to wall for scale matching with the layers used in the ongoing Mexican Fire Danger System [70,72]. A low pass filter with a moving window of 3 × 3 pixels, using the package included in Spatial module of ERDAS IMAGINE 2014 software (ERDAS Inc., Atlanta, GA, USA), [76] was applied to the AGCD map.

Forest Fire Data
The fire suppression records database from the Mexican Government Forest Agency (CONAFOR) from the period 2005 to 2015 was utilized. The forest fires records in water bodies, cities and agriculture were excluded for the analysis based in the land use map from INEGI (2011) [51]. The database comprised a total of 53,000 fire suppression records in the period of study.

Observed Fire Occurrence by AGCD Calculation
We firstly calculated the total number of forest fires, termed NFTi, for every region i. The aboveground carbon density (AGCD), in Mg C ha −1 , was retrieved for every fire suppression record, using the extract multiple values to points tool in ArcGIS 10.1 software (ESRI: Redlands, CA, USA) [77] software.
For every forest region (i), the fraction of fire occurrence recorded by AGCD value (j) was calculated following Equation (1): where: FOB: Fire occurrence index by biomass level (fraction of the number of fires recorded at a forest region i for every AGCD value j); NF ij = number of forest fires recorded a forest region i for every AGCD value j; NF Ti = total number of forest fires recorded for the i region.

Modeling and Mapping Fire Occurrence from AGCD
The two-parameters Weibull probability density function (PDF) model (Weibull, [78]) was tested to predict FOB from AGCD for every region following Equation (2): where: FOB: Fire occurrence index by biomass level at every forest region i for every AGCD value j (Equation (1)); AGCD is the observed aboveground carbon density value j for every i region; and b and c are the scale and the shape parameters of the model coefficients to be estimated, respectively. The Weibull model was selected based on visual inspection of the data distributions. This function has been widely applied (e.g., [79][80][81][82][83]) because of its flexibility to capture different shapes of distributions together with its relative simplicity. The models to predict observed FOB from AGCD (Equation (2)) were fitted with non-linear regression (NLS) using the package "nls" of the software R (R Core Team, Vienna, Austria [84]).
In addition, the nonlinear extra sum of squares method (Bates and Watts, [85], pp. 103-104) was used to assess whether the models were significantly different among different regions. This method has often been utilized to assess whether separate models are necessary for different species or different geographic regions (e.g., [73,[86][87][88]). It involves the fitting of full and reduced models. The reduced model corresponds to the same set of global parameters in the function (Equation (2)) for all regions. In the full models, different sets of global parameters are tested for different regions. This is tested by expanding each global parameter in Equation (2) by including an associated parameter and a dummy variable to differentiate the region.
For example, the expansion of a global parameter c of Equation (2) can be written as Equation (3): where: c1-c7 are the associated parameters of the full model (Equation (2)) for each of the seven regions, and d1-d7 are the dummy categorical variables for considering the regions, which are defined as follows: d1 = 1 if region = region 1, otherwise d1 = 0; . . . ; d7 = 1 if region = region 7, otherwise d7 = 0; The non-linear extra sum of squares follows an F distribution. We selected candidate models where the grouped coefficients were significantly different as detected by the F-test (p < 0.05).
Since the goal of our study was to evaluate the goodness-of-fit of Weibull distributions utilized as regression models, standard tests often used to evaluate the performance of probability distribution functions (e.g., Kolmogorov-Smirnov, Cramear-von Mises and Andersen-Darling tests) were not used. Such tests were not designed to answer questions about the error associated with predictions obtained from a regression model (e.g., [89][90][91]), which was the goal of our analysis. Candidate models were evaluated using the coefficient of determination for nonlinear regression (R 2 ) as a statistical criterion, together with root mean square error (RMSE) and model bias. R 2 is defined as the square correlation coefficient between the measured and estimated values and the root mean square error (RMSE) [92]. In addition, we considered the simultaneous F-test for intercept = 0 and slope = 1 [93], to evaluate the goodness-of-fit of the Weibull models.
The predicted FOB values were normalized from 0 to 100, scaling by the maximum predicted value for every region. Normalized predicted FOB (FOBN) was mapped for the whole country using ArcGIS 10.1software (ESRI: Redlands, CA, USA) [77].

Results
For the period of study (2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015), the Center region had the largest number of fires, with an NFT i of 13,220 fire suppression records, equivalent to a fire density of 395 fires per 10,000 km 2 of forest land. The lowest fire density was observed in the arid North region, with only 91 fires in the 10 years of the period of study, equivalent to a fire density of eight fires per 10,000 km 2 .
The estimates and the confidence interval of the Weibull PDF parameters and the values of the goodness-of-fit statistics for each region are shown in Table 1. The results of the nonlinear extra sum of squares method showed that all models were significantly different between regions, with the exception of NW and NE regions, which could be explained with a common set of parameters because of the absence of significant differences between these two temperate-forest dominated regions. Therefore, these two regions (NW and NE) were considered as a unique unit. The scale parameter of the two-parameters Weibull PDF (b) which influences the AGCD value at which the maximum FOB is predicted, varied from 17 to 33. The lowest predicted values of b coefficient were obtained for North scrubland (N) and chaparral (CHAP) regions, with values of 17 and 21, followed by South-center (SC), Northeast (NE) and Northwest (NW), with b values ranging from 25 to 29. The highest b values were obtained for SE and C regions, with values of 31 and 33, respectively. The parameter c is related to the shape of the two-parameters Weibull PDF, with higher values meaning steeper distributions result in FOB values focused on a narrower range of AGCD values. The highest values of c coefficient (steeper distributions) in the range 3.6-3.8, were observed for the regions of the NW-NE and N forests. Lower c values (flatter distributions) of around 2.5, were obtained for the SC and SE forest and the lowest c coefficient values (1.9) were obtained for the CHAP region.
The best goodness-of-fit statistics were obtained for the C and NW-NE regions with R 2 values of 0.95 and 0.98, respectively and RMSE values of 0.002 and 0.003, respectively. The lowest values of the goodness-of-fit statistics were obtained for the N region with a R 2 and RMSE of 0.48 and 0.029, respectively (Table 1). These goodness-of-fits are similar to or higher than the range of R 2 of 0.5-0.7 reported in the literature for spatial analysis of fuel availability and fire activity (e.g., [35,41]). Nevertheless, the fuel-related variables utilized in those studies-leaf area index, tree cover, among others, the size of the area of study, and the pixel size in the case of global analysis, are not directly comparable to the approach and spatial resolution utilized in the current study.
The curve of observed fire occurrence by biomass level FOB for each forest region is shown in Figure 2, together with the predicted values, against AGCD. In the chaparral and North scrubland forests (Figure 2a,b) the highest values of FOB occurred at low values of AGCD (from 7 to 25, 10 to 22 Mg C ha −1 , respectively). The curve of observed fire occurrence by biomass level FOB for each forest region is shown in Figure 2, together with the predicted values, against AGCD. In the chaparral and North scrubland forests (Figure 2a,b) the highest values of FOB occurred at low values of AGCD (from 7 to 25, 10 to 22 Mg C ha −1 , respectively).
For the regions of the combined Northwest, Northeast and South-centre (Figure 2c,e), highest values of FOB were related with medium values of AGCD (from 15 to 35 Mg C ha −1 ), with a steeper distribution for the NW and NE compared to SC. In the Center and Southeast regions, (Figure 2d,f), the highest FOB values were observed in slightly higher AGCD values, ranging from 20-22 to 35 Mg C ha −1 , with a flatter distribution for SE. Lower fire occurrence was observed at the highest biomass levels for all regions. Predicted against observed FOB values and the related goodness-of-fits statistics are shown in Figure 3. The results of the simultaneous F-test (Figure 3) show that there are no reasons to reject the null-hypothesis of intercept = 0 and slope = 1 meaning that there are no systematic over or underestimates in any ecoregion. The map of predicted FOB for Mexico based on the fitted models is shown in Figure 4. A good agreement can be seen between predicted FOB values mapped and the historical fire suppression records spatial distribution in the country.
Selected regional detailed examples of the mapped FOBN are shown in Figure 5. In the case of the chaparral forests, higher FOB values were observed in sites with low levels of AGCD (5 to 20, Mg C ha −1 ) (Figure 5a), for which maximum fire occurrence was predicted (Figure 5a, right figures). In the Northwest and Northeast forests (Figure 5b  For the regions of the combined Northwest, Northeast and South-centre (Figure 2c,e), highest values of FOB were related with medium values of AGCD (from 15 to 35 Mg C ha −1 ), with a steeper distribution for the NW and NE compared to SC.
In the Center and Southeast regions, (Figure 2d,f), the highest FOB values were observed in slightly higher AGCD values, ranging from 20-22 to 35 Mg C ha −1 , with a flatter distribution for SE. Lower fire occurrence was observed at the highest biomass levels for all regions.
Predicted against observed FOB values and the related goodness-of-fits statistics are shown in Figure 3. The results of the simultaneous F-test (Figure 3) show that there are no reasons to reject the null-hypothesis of intercept = 0 and slope = 1 meaning that there are no systematic over or underestimates in any ecoregion. The map of predicted FOB for Mexico based on the fitted models is shown in Figure 4. A good agreement can be seen between predicted FOB values mapped and the historical fire suppression records spatial distribution in the country.     ; where low, medium, high, very high and extreme represent FOB N values of <1%, 1%-25%, 25%-50%, 50%-75% and >75% respectively. (b) Location of suppressed fires in the study period. The three regional detailed windows included in Figure 4b are analyzed in detail in Figure 5.

Discussion
The varying relationships of AGCD with fire occurrence by region agree with studies that have found different effects of fuel availability on fire activity by climatic regions (e.g., [41][42][43]). At each ecoregion, biomass effects on fire activity might vary because of the interactions of climate and fuel availability on fire occurrence (e.g., [35,[44][45][46][47]49]).
For the chaparral region, (Figures 2a and 5a), a large percentage of fire suppression records occurred at AGCD values in the range 5-25 Mg C ha −1 , with the highest values of fire occurrence found in the range 15-25 Mg C ha −1 , approximately corresponding to fuel loads from 30 to 50 Mg ha −1 . The latter figure is likely the upper limit for fuel load in this type of vegetation in Baja California, where due to lower fire rotation periods than in a similar vegetation type in Southern California [94], fuel load is probably lower than in the latter region [95,96]. This suggests that fire suppression in this region seems to have mainly focused on areas with the highest chaparral fuel load and continuity, where fire spread is expected to be more frequent and dangerous [39,[97][98][99].  (1) show observed AGCD (Mg C ha −1 ), together with observed forest fire suppression records in red (2). Right panels show predicted FOB N (3) together with observed forest fire suppression records in blue (4). Figure 5. In the case of the chaparral forests, higher FOB values were observed in sites with low levels of AGCD (5 to 20, Mg C ha −1 ) (Figure 5a), for which maximum fire occurrence was predicted (Figure 5a, right figures). In the Northwest and Northeast forests (Figure 5b), highest fire occurrence was predicted for intermediate levels of AGCD (from 15 to 35 Mg C ha −1 ), with lower values of FOB for the highest AGCD values (>35 Mg C ha −1 ). In the C and Southeast forests (Figure 5c,d), highest FOB values were observed in intermediate levels of 20 to 35 Mg C ha −1 and lower fire occurrence was also observed at the highest biomass (>35 Mg C ha −1 ).

Discussion
The varying relationships of AGCD with fire occurrence by region agree with studies that have found different effects of fuel availability on fire activity by climatic regions (e.g., [41][42][43]). At each ecoregion, biomass effects on fire activity might vary because of the interactions of climate and fuel availability on fire occurrence (e.g., [35,[44][45][46][47]49]).
For the chaparral region, (Figures 2a and 5a), a large percentage of fire suppression records occurred at AGCD values in the range 5-25 Mg C ha −1 , with the highest values of fire occurrence found in the range 15-25 Mg C ha −1 , approximately corresponding to fuel loads from 30 to 50 Mg ha −1 . The latter figure is likely the upper limit for fuel load in this type of vegetation in Baja California, where due to lower fire rotation periods than in a similar vegetation type in Southern California [94], fuel load is probably lower than in the latter region [95,96]. This suggests that fire suppression in this region seems to have mainly focused on areas with the highest chaparral fuel load and continuity, where fire spread is expected to be more frequent and dangerous [39,[97][98][99].
In contrast, in areas of low biomass, in spite of dry climate and the associated low dead fuel moisture, fire spread is probably limited by low fuel availability like in other Mediterranean shrubby areas. Reduction in fire occurrence in the latter areas has been attributed to changes in fuel load, structure and connectivity (e.g., [100]).
The relatively lower number of recorded fires in the highest AGCD values in this chaparral ecoregion seems to be associated to the pine-oak "sky islands" forests, which accumulate comparatively higher biomass values than the genuine shrubby chaparral-dominated areas. This is consistent with observations in Mediterranean landscapes where fire occurrence is also higher on scrublands than in wooded areas [29,30,101].
There are further factors that should also be considered in future research to explain the spatial patterns of fire occurrence in chaparral regions. For example, Vega-Nieva et al. [71] and Monjarás-Vega et al. [102] suggest that the spatial pattern of ignitions in the chaparral area is very strongly influenced by anthropogenic factors such as by vicinity of roads and urban areas, which are related to accessibility for human ignition and also for allowing access for fire suppression. This vicinity of ignitions to human settlements in the Baja California chaparral area might be linked to recent urban growth, and its associated wildland-urban interface, similarly to that reported by Syphard et al. [103] on the adjacent Southern California.
Unlike the chaparral region, temperate pine and oak forests in the NW, NE (Figures 2c and 5b) and C regions (Figures 2e and 5c) did not experience a large fire occurrence for AGCD <15 Mg C ha −1 . In these regions, low biomass areas are mainly covered by xerophytic scrublands [51] which, unlike chaparral, are characterized by lower fuel continuity and generally lower fire propagation hazard [53,54]. The lower occurrence of fires observed in these low biomass semi-desert and desert shrublands, compared to forested areas, has been previously reported in these regions in Mexico [55][56][57][58]. This also is in line with observations in other countries where lower fire occurrence in desert scrubland areas has been related to lower fuel load and continuity, limited by the extremely dry climate [47,100].
Several reasons might explain the highest fire occurrence found at intermediate levels of AGCD (15-35 Mg C ha −1 ) in the above-mentioned temperate forest ecosystems. Intermediate biomass is associated to areas dominated by woody secondary vegetation. This primarily consists of young trees stands with low crown base height and bearing high loads of surface fuels that result in potentially more intense fires, including torching [53,54]. This agrees with very high occurrence of MODIS active fires observed in the secondary vegetation, classified as a very hazardous fire-prone vegetation type by Vega-Nieva et al. [73]. Likewise, with observed increases in fuel availability and fire spread potential related to decreasing tree cover in other countries [36,42,104]. On the other extreme within those temperate forests, the highest biomass areas, with lower fire occurrence, correspond to dense old-growth forests, mainly located from middle to high elevations. These stands, with moderate to high basal area, consist of tall dominant mature trees with high crown base height. This fuel structure limits wind penetration and sunlight incidence, maintaining low surface fuel loads and microclimatic and fuel structural conditions that hamper the occurrence of high severity fires [105]. This would explain the higher resilience to fire of these ecosystems, adapted to low severity natural fires, thus enabling the maintenance of ecosystem´s health and self-regulating surface loads at levels that prevent catastrophic large severe fires [16,106,107].
Finally, the response of fire occurrence in tropical SE forests to forest biomass might be partially explained by taking into account differences in forest structure, microclimatic conditions and human influences on ignition. Thus, the highest fire occurrence observed in the intermediate biomass levels corresponds to the dry and seasonally dry tropical forest, as opposed to lower fire occurrence in the high biomass areas associated to the evergreen wet tropical forests, where high fuel moisture throughout the year limits fire occurrence [52,53,71,72]. Areas of highest biomass and tree cover tend to be fire-limited by fuel moisture with wet fuels preventing ignition, together with a lower fire spread potential [36,42,43,108].
At the same time, in tropical forests, lower biomass areas are more connected to active forest degradation at landscape level, including slash and burn activities, which are important causes of fire occurrence in the SE region [55,56,109,110]. An example of the evident role of human accessibility in fire occurrence in these wet tropical ecosystems can be found in Ressl et al. [111]. In that study, the majority of ignitions in the Kalakmul reserve (south of the Yucatan Peninsula) were located in the vicinity of roads, whereas most isolated areas had an extremely small number of ignitions, very possibly because of the difficulties in accessing inside this very dense type of tropical forest. Posterior observations of human factors of fire occurrence for the whole Yucatan Peninsula [71,102] confirm the same pattern. This strongly suggests the need for considering anthropogenic factors in combination with biomass and fuel dryness [73] to gain insight into the multiple interactions behind the spatial distribution of fire occurrence and improve our ability of fire occurrence prediction.
Overall, our results suggest a gradient in the AGCD values for the respective maximum fire occurrence in the observed FOB distributions from the driest (N and CHAP, with peaks of FOB at <20 Mg C ha −1 ) to the wettest ecoregions (C and SE, with maximum FOB at >30 Mg C ha −1 ), with intermediate FOB peaks (between 25 and 30 Mg C ha −1 ) for the temperate forests of the NE and NW regions. These variations in the shape and maximum values of FOB across climatic gradients of aridity resemble the general tendency hypothesized by Dohn et al. [112] and Kahiu and Hanan [36] (Figure 2b), among others, relating tree cover and precipitation gradients with fire activity. Unlike the latter study [36], which utilized a MODIS leaf area index for explaining fire activity, our study related the variations in the relationships between aboveground biomass and fire occurrence at 1 km for contrasting climatic ecoregions, suggesting potential of biomass maps as explanatory variables to aid in fire occurrence modeling at operational scales for fire management planning.
Whereas these differences between the peaks and distribution shapes of fire occurrence between the arid N, temperate NW and NE and more humid C and SE regions seem to confirm a gradient in the AGCD values of maximum fire occurrence, one exception would be the very rainy SC region, where the highest fire occurrence would be expected at higher AGCD values, similar to the rainy tropical SE region, or to the C region. We speculate that other factors in addition to climate, such as lower road accessibility and urban areas density, might result in these areas being less accessible for ignition in the SC. That contrasts with higher road and urban areas density of the C region that might result in a higher percentage of fires occurring at higher biomass values in this densely populated region. Further research in the combined role of biomass and human factors in explaining fire occurrence spatial patterns by region might clarify these potential interactions.
The number of ecoregions utilized offered a good compromise between having contrasted climatic and vegetation conditions and containing enough internal variability in aboveground biomass to allow the finding of significant relationships with a representative number of fire suppression records. One exception would be the desert shrubby N region, which is characterized by rare fire activity. Longer periods of observation would be required for improving these initial models in this region.
In addition, future work might explore a finer regionalization to account for further internal variability within ecoregions.
Combining these data about the most likely locations where fire may occur based on region-specific AGCD relationships, with other information about road densities and urban locations-which also affect spatial patterns of fire ignitions [71,102], better maps of wildfire ignition likelihood could be produced. These ignition likelihood maps could be combined with weather-driven metrics to create improved spatial estimates of wildfire risk across Mexico that would provide better decision support information for forest and fire managers.

Conclusions
The current study presents the first available map for fire occurrence based on aboveground carbon density in Mexico. The results support the hypothesis of fuel limitation-which explains a substantial part of fire occurrence in the most arid regions-and climate limitation for the wetter regions.
The maximum values of fire occurrence occurred at different aboveground carbon density values by region, with the most arid regions (chaparral and N desert shrubs) at values of <20 Mg C ha −1 , and the wetter regions C and SE at >30 Mg C ha −1 . Whereas these first models help improve our understanding of the spatial pattern of fire occurrence in Mexico, the complexity of the phenomenon of fire occurrence demands to explore additional variables likely affecting these patterns.
Variations in fuel complex load and other structural attributes, associated to climate and tree cover differences, combined with other additional factors, not still quantified, such as anthropogenic drivers are probably influencing fire occurrence and they deserve further study. The role of humans and infrastructure-related variables such as distance to roads, human settlements, level of agricultural activities and the presence of respective interfaces with forest on the spatial patterns of fire occurrence needs to be explored in future studies.
Although the study period is not very large, it comprised dry and wet years and appears sufficiently representative to have a good picture of the situation for most of the ecoregions. Nevertheless, given the very low fire occurrence in the 10-year study period for the shrubby vegetation in the N region, a longer period of observation will be required to improve these first models in this region.
The successful use of the Weibull model to predict the observed fire occurrence distributions as a function of the total aboveground biomass was likely due to the latter acting as a proxy of fuel load and other fuel structural features [113,114] connected with fire flammability, a point that needs to be further investigated. Still, it would be desirable for future studies to utilize a fuel characteristics map, once such map becomes available nationally. Additionally, the latter would assist to better understand the connection between fire occurrence and fuel structure and connectivity at landscape scale.
Given the reasonably good predictive ability of most of the eco-regional models to estimate fire occurrence from AGCD, the quantitative approach developed can be included in operational GIS tools for assistance in fuel treatment prioritizing and suppression resources allocation decision-making at a national scale. These first models could be improved in the future with larger datasets, and the inclusion of fuel characteristics and human factors.
The modeling strategy presented here could be replicated in any other country or region, based on satellite measures of aboveground biomass and satellite measurements of fire activity or fire suppression records, for analyzing patterns of fire occurrence under contrasting climates, vegetation types, and human patterns of ignition. Funding: Funding for this work was provided by CONAFOR/CONACYT Projects "CO2-2014-3-252620" and "CO-2018-2-A3-S-131553" for the development and enhancement of a Forest Fire Danger Prediction System for Mexico. The scholarships from CONACYT for Carlos Briones to study a Masters program in Forest Sciences at the Universidad Juarez del Estado de Durango and to fund a research stay at the University of Santiago de Compostela (Spain) is also acknowledged.