Natural Hazards and Earth System Sciences A stochastic Forest Fire Model for future land cover scenarios assessment

Land cover is affected by many factors including economic development, climate and natural disturbances such as wildfires. The ability to evaluate how fire regimes may alter future vegetation, and how future vegetation may alter fire regimes, would assist forest managers in planning management actions to be carried out in the face of anticipated socio-economic and climatic change. In this paper, we present a method for calibrating a cellular automata wildfire regime simulation model with actual data on land cover and wildfire size-frequency. The method is based on the observation that many forest fire regimes, in different forest types and regions, exhibit power law frequency-area distributions. The standard Drossel-Schwabl cellular automata Forest Fire Model (DS-FFM) produces simulations which reproduce this observed pattern. However, the standard model is simplistic in that it considers land cover to be binary – each cell either contains a tree or it is empty – and the model overestimates the frequency of large fires relative to actual landscapes. Our new model, the Modified Forest Fire Model (MFFM), addresses this limitation by incorporating information on actual land use and differentiating among various types of flammable vegetation. The MFFM simulation model was tested on forest types with Mediterranean and sub-tropical fire regimes. The results showed that the MFFM was able to reproduce structural fire regime parameters for these two regions. Further, the model was used to forecast future land cover. Future research will extend this model to refine the forecasts of future land cover and fire regime scenarios under climate, land use and socio-economic change. Correspondence to: P. Fiorucci (p.fiorucci@cimafoundation.org)


Introduction
Landscapes are dynamic systems that reflect the complex interplay of many factors including climate, natural disturbance, natural succession, economic development and public policy.Public planning efforts directed at providing the best mix of social benefits derived from complex landscapes requires an understanding of the key factors driving landscape change and how economic and ecological factors interact over time and space (Wear and Bolstad, 1998;Drechsler et al., 2009).Forecasts of anticipated future conditions required for landscape planning and public policy will improve as scientists come to understand the slow variables that constrain fast ecological and economic processes (Carpenter and Turner, 2000;Clark et al., 2001).Within the context of wildfires and land use planning, species change and fuel accumulation resulting from fire suppression (Dombeck et al., 2004), climate change (Westerling et al., 2006), and landscape homogenisation (Loepfe et al., 2010) are slow moving variables that influence the occurrence and severity of wildfires (a fast variable).Public planning efforts designed to protect the public from the undesirable consequences of wildfires would benefit from a better understanding of the feedback between slow and fast variables at a landscape scale (Holmes et al., 2008b).
Recent research indicates that several feedback pathways exist between landscape structure and human populations in anthropogenic fire regimes.Fuel connectivity is a key factor in explaining the area burned by wildfires (Miller and Urban, 2000).However, connectivity in anthropogenic fire regimes can be disrupted by fuel fragmentation (due to factors such as road building or conversion of forests to crop land) and fire suppression effort (Guyette et al., 2002).Alternatively, the connectivity of fuels can be increased by agricultural land abandonment and the succession to shrublands -a process referred to as landscape homogenisation -which is currently occurring in many regions of Western Europe (Loepfe et al., 2010).The homogenisation of landscapes tends to increase fire occurrence which, in turn, increases landscape heterogeneity as burned over lands, for example, are converted to pastures (Lloret et al., 2003).Finally, humans are the source of many, if not most, wildfires in densely populated areas and are an integral factor in the analysis of anthropogenic fire regimes (Cardille et al., 2001;Prestemon et al., 2002).
In a broad sense, the term "cultural landscape" has been used to describe geographic areas in which the interactions between human activity and the environment have created typical patterns of plant and animal communities (Farina, 2000).While many historical-cultural landscapes are valuable to society, contemporary socio-economic forces pose challenges to the sustainability of these areas (Vos and Meekes, 1999).It has been suggested that preserving and/or promoting the traditional rural mosaic in some regions of Western Europe may increase fuel fragmentation and reduce the risk of large fires (Badia et al., 2002;Loepfe et al., 2010).Similarly, the historical-cultural practice of burning by indigenous people to reduce vegetative cover has been adopted in prescribed forest burning practices which are widely used on private lands in the Southern United States (Cronon, 1983;Fowler and Knopik, 2007).
Because shifts in human cultural traditions cause substantive transitions in fire regimes, it is of great interest to develop methodologies that are capable of simulating the dynamic relationships between vegetation, wildfire, socio-economic factors and public policy.In this paper, we present a novel simulation model which is used to examine the dynamic relationships between land cover and wildfires.The goal of the modelling strategy is to provide a platform for investigating a suite of public policy issues related to land-use planning, fire management and climate change.

Modelling approach
Forest fire regimes are generally characterised by a large number of low-intensity events (which have little impact on society) and a small number of catastrophic events (which have severe impacts on society).For example, during the period 1980-2002, 94 percent of fire suppression costs on US Forest Service land resulted from a mere 1.4 percent of the total number of fires (Holmes et al., 2008a).Size-frequency models for these types of disturbances can be represented by heavy-tailed distributions, that is, distributions for which unusual events occur more often than generally anticipated.Currently there is great interest in, and considerable debate over, the exact form of the probability distributions that best describe this behaviour (Millington et al., 2006).Among the heavy-tailed distributions, the power-law distribution is one of the most frequently used for the statistical characterisation of wildfire regimes (Malamud and Turcotte, 1998;Ricotta et al., 1999;Song et al., 2001;Turcotte and Malamud, 2004;Telesca et al., 2005;Malamud et al., 2005;Fiorucci et al., 2008).
Power-laws are revealed both by real world distributions of fire size-frequencies as well as being produced by cellular automata simulation models.The standard cellular automata forest-fire model was first described by Drossel and Schwabl (1992) -which we refer to as DS-FFM -and is defined on a square lattice with L 2 regularly spaced cells.Each cell is either occupied by a tree (or, more generally, flammable vegetation), or is empty.At each discrete time step, all sites are updated in parallel according to the following rules: (i) an empty site becomes occupied by a tree with probability p, (ii) a tree is ignited with probability f .The ignited tree and the entire cluster of trees connected to it (by nearest-neighbour coupling) burn and become empty sites.These seemingly simple rules result, when iterated over many time steps, in a fire size-frequency density characterised by a power-law: where A is the area burned, f (A) is the frequency associated with the area burned, α and β are parameters and β is called the scaling exponent and is typically negative (β < 0).The DS-FFM has been linked with phenomena known as self-organized criticality (Bak, 1991;Drossel and Schwabl, 1992).Criticality, a term drawn from statistical mechanics, describes the state of a system which is poised at a phase transition.The fact that the system is referred to as selforganized derives from the fact that although there is no local heterogeneity governing behaviour of the system, large scale patterns are created from purely local interactions.As suggested by Pascual and Guichard (2005), the power-law spatial pattern that results from many iterations of this model over many time steps depends on a double separation of time scales.Fire spread is a fast variable, forest growth is a slow variable, and the rate of ignition (ignitions per unit area) is a very slow variable.
It should be noted that processes other than described by the cellular automata model can generate power laws, and that the observation of power law phenomena does not necessarily imply the system is self-organized.For example, Reed (2001) demonstrated that phenomena which evolve following a proportional growth rate, and are observed at time periods which are random draws from an exponential distribution, should follow distributions with power-law behaviour (in the tail of the distribution).Since actual wildfires may grow proportionally until the area burned is stopped at random by rain or snow, the distribution of actual wildfires may follow a power-law (at least in part of its range) and yet not be generated by self-organized critical behaviour.A similar explanation for power laws in evolution and ecology was put forth by Solow (2005) who argued that the observation of power laws should not necessarily be taken as evidence of self-organized criticality.The traditional DS-FFM is clearly an oversimplification of forest dynamics.In addition, the power-law exponents, obtained using the DS-FFM, are typically larger than the exponents estimated for actual wildfire regimes.In response to this concern, Krenn and Hergarten (2009) sought to modify the rules governing a cellular automaton in an attempt to make the model more realistic and correspond more closely with real world observations.In particular, recognizing that many real world wildfires are ignited by people, they reasoned that anthropogenic ignitions should occur at the border of a forest (such as along roads or trails) and not within a forest as specified by the DS-FFM.By combining a new rule for anthropogenic ignitions, and linking this rule with lightning induced fires, they found that the scaling exponent increased by roughly 1/3 which is more similar to real world observations.
In keeping with this line of thinking that seeks to make the cellular automata forest fire model more realistic, the innovation of this paper is to calibrate the standard cellular automata forest-fire model for actual landscapes.In the standard cellular automata forest-fire model, the rules governing the simulation are held constant across a homogeneous square lattice.In this paper, we modify the standard model by: (1) populating the lattice with vegetation types observed on actual landscapes, (2) linking unique probability values p i for vegetation recovery with each vegetation type (i), and (3) creating  new probability values g i that describe the probability of fire spread within each vegetation type.By allowing the recovery (regeneration) probability p i to vary according to vegetation type, recovery dynamics are local.Finally, by allowing the probability of fire spread to vary within each vegetation type, the model can simulate changes in flammability which may result from management actions, such as prescribed fires.
We refer to this new model as the modified forest fire model (MFFM).The model is initiated by a random number of ignitions determined by the probability f for each 100 m 2 cell, which is constant over the entire grid.Each cell that is ignited burns with probability g i .Fire spreads from a burning cell to contiguous cells with probability g j where (j ) is the class of vegetation cover in the neighbouring cell.Each iteration terminates when no fires are burning and empty sites become occupied by vegetation class i with probability p i .At that point, the frequency of fires of all sizes is counted where size is determined as the sum of the number of burned cells in each fire.
The initial values of the probabilities g i are obtained using GIS analysis of the actual fires data series.The values of p i are defined arbitrarily considering the number of years necessary to have the same degree of coverage existing before the fire.The parameter f is initially undefined and it is obtained as a result of the tuning phase.It is necessary to "tune" the values of the initial probabilities p i , and g i to obtain a frequency-area distribution of simulated fires that closely approximates the power-law distribution associated with an actual fire regime.This is accomplished by using data on the sizes and frequencies of historical fires occurring within a study landscape.Taking the logarithms of both sides of Eq. ( 1), it can be seen that a linear relationship exists between frequencies and fire sizes.Thus, it has been common practice to estimate the power-law parameters using Ordinary Least Squares regression.However, it has been recently recognized that this procedure can produce seriously biased estimates of power-law parameters (Bauke, 2007;White et al., 2008;Clauset et al., 2009).Maximum likelihood fitting methods are preferred.Using Monte Carlo methods, White et al. (2008) demonstrated that maximum likelihood estimates are superior to other methods (specifically OLS methods using linear binning, logarithmic binning, normalized logarithmic binning or fitting the cumulative distribution function).Consequently, we used the following maximum likelihood estimator for the power-law parameter β: Equation ( 2) was used to estimate values for β from the actual data.Similarly, this equation was used to estimate values for β using simulated data on fire sizes and frequencies provided by the MFFM.A bootstrap procedure, based on 1000 iterations, is used to obtain an empirical distribution of the β parameter based on actual fires.Likewise, the MFFM can be iterated many (1000) times, and the β values from each iteration of the simulation model can be recorded.The p i , f , and g i values are tuned until the total number of fires and the total burned area from the actual and simulated fires are closest.Then, the β values from the actual and simulated fires are compared.
Given the tuned parameter values, it is possible to use the MFFM to simulate future scenarios based on relatively small perturbations of the tuned model parameters.A wide variety of future scenarios could be evaluated using this approach.To demonstrate this modelling approach, we present a longterm simulation to evaluate possible future trends in land cover change.This is accomplished by adding new rules to the MFFM.In particular, these rules (which are described in more detail below) allow the vegetation in each cell to randomly change during the recovery phase after a fire.

Data
The proposed model has been tested in two different geographic areas characterised by different fire-adapted ecosystems and wildfire regimes.The first test area is the administrative Region of Liguria (Italy) and the second is Alachua County (Florida, USA).Liguria is a region of 5400 km 2 lying between the Cote d'Azur (France) and Tuscany (Italy) on the northwest coast of the Tyrrhenian Sea (Fig. 1).For this Mediterranean region, characterised by complex topography and dense vegetation, wildfires are recurrent phenomena both in summer and winter and an average of 365 wildfires of size ≥0.01 km 2 burn an area of 55 km 2 per year.Alachua County covers 2500 km 2 in north-central Florida, United States, and is located in a sub-tropical ecosystem (Fig. 1).
About 50 wildfires of size ≥0.01 km 2 burn an area of about 9 km 2 per year in this County.In Liguria, most of the fires are caused by humans.Less than 1% of the fires are naturally caused in this region.On the contrary, although many of the fires in Alachua County are human caused, lightning is also an important ignition source in this ecosystem.The population density in Liguria is about 300 km 2 whereas in Alachua County is less than 100 km 2 .
In both cases, a grid of 100 m 2 cells is superimposed over the total area.Land cover types are classified into 6 main categories: no fuel, shrub, grass, oak, conifer and other broadleaves.The land cover map for the two regions is shown in Fig. 1.As can be seen, Alachua County land cover is mainly composed of coniferous species whereas Liguria land cover is mainly composed of other (non-oak) broadleaf species (Castanea sativa, Fagus sylvatica).

Analysis and results
A preliminary analysis of the fire regimes in the two regions has been carried out.In the first case study, we used the actual occurrence of fires in Liguria between 1986 and 2006.1000 different trials were conducted randomly selecting 500 fires, for every trial, from the complete time series of about 8000 fires.Thus, 1000 different β power-law parameter values were obtained.The same approach was used in the second case study, randomly selecting again 500 fires, for every trial, from the complete time series of about 1400 fires that occurred between 1981 and 2006 in Alachua County.Successively, 1000 simulations corresponding to the time series length were generated for the two case studies.Each simulation was made by n iterations, where n is the time series length in years, and each iteration represents one year of fire occurrences.
As discussed above, the initial values of the probabilities g i were obtained by a GIS analysis of the actual fires data series.For the time window 1997-2009, the perimeters of the area burned by the actual fires in Liguria were available.Superimposing the burned areas to the actual vegetation cover classified in the 5 classes, the total burned area of each vegetation class was normalized to the total area covered by each  1 for the two case studies area.
The tuning phase was based on the comparison between the actual number of fires and the total burned area during the considered time window and the number of simulated fires on the same time window, respectively, for the two case study areas.A close correspondence between the actual and simulated mean number of fires and total burned area was obtained for the two case studies (Table 2).
A second analysis was conducted to compare the distribution of the power-law parameter β as previously described.The parameters of DS-FFM were tuned to obtain the closest value of the total number of fires and the total burned area as well.The value of the parameters, after the tuning phase, is reported in Table 3 for the two case studies area.
In Table 4, the observed number of fires and total burned area are compared with their counterpart values obtained using the DS-FFM simulation model.
The distribution of β obtained using the time series of simulated fires using the MFFM (5 classes of vegetation) and the DS-FFM (1 class of vegetation) can be compared with the distribution of β obtained using the time series of the actual fires (Fig. 2).As can be seen, in both regions the powerlaw parameter estimates from the MFFM simulation model are more similar to the parameter estimates from the historical fire data than are the parameters computed using the DS-FFM.Thus, the MFFM appears to be very adequate for simulating historical fire regimes in these two case studies.In Fig. 2, it is evident that the use of 5 classes of vegetation in Liguria is essential because the vegetation cover is very heterogeneous.This suggests that the power-law behaviour observed in actual fire data may be directly related to spatial arrangement of vegetation clusters.On the contrary, Alachua is quite homogeneous in vegetation cover.In this case, the values of parameter β obtained by DS-FFM simulation are www.nat-hazards-earth-syst-sci.net/10/2161/2010/Nat.Hazards Earth Syst.Sci., 10, 2161-2167, 2010 not so far from the values obtained by the actual fires.However, in this case, the introduction of 5 classes of vegetation allows us to obtain results that are again more similar with the actual fire data.Finally, a long term simulation was performed to evaluate possible future trends in land cover change introducing new rules in MFFM.The new rules were: (1) after 3 years without fires in herbaceous areas, the vegetation cover changes randomly to one of the neighbourhood vegetation classes; (2) in the shrub areas, the same process occurs with a time lapse of 10 years; (3) a burned area previously not vegetated by shrub or grass becomes occupied by grass with probability 0.7, shrub with probability 0.2 or by the previous vegetation class with probability 0.1.In Figs. 3 and 4, the obtained results are reported for Liguria and Alachua County, respectively.
The long term trend analysis in the first case study shows a decrease of grass and conifer.The area of shrub vegetation, after an initial period of expansion, slowly decreases after 60 years of simulation.On the contrary, oaks, after an initial decrease, start to slowly increase after 60 years of simulation.Other broadleaved land cover show a positive trend increasing to about 20% after 100 years.
Different trends resulted from the simulations conducted for Alachua County.In this case, only grass vegetation decreases towards an asymptotic value corresponding to 20% of the initial value.Coniferous vegetation increases toward an asymptotic value of about 115% of the initial value.Oaks, shrub and the other broadleaf vegetation also show a positive trend.

Conclusions
The MFFM is able to reproduce the variability of the powerlaw parameters of actual fires in both the considered regions.As the case studies represent widely divergent, fire-adapted ecosystems (Mediterranean and Sub-tropical), these preliminary results suggest that the methodology developed for the MFFM is robust and promising for future analyses.
Some preliminary conclusions are suggested concerning the values of the tuned parameters, f , p i and g i in the 2 regions.First, the ignition probability parameter, f , is very different between Liguria and Alachua County.This difference may be attributable to differences in population density.In fact, Liguria has a population density about three times greater than in Alachua County leading to a greater likelihood of human-caused fires.Second, the vegetation regrowth parameters, p i , have the same values in the two regions for the various classes of vegetation except for shrubs.This result requires further analysis but may be due to the different species of shrubs and different climatic conditions found in the two regions.Third, the probability of fire spread, g i , for grass and coniferous vegetation are greater for Liguria than for Alachua County.This difference might be explained by noting that Liguria has a rougher topography than Alachua County, which may affect the ability to suppress wild fires.
Finally, the MFFM is able to provide scenarios describing possible long-term trends in land cover change.Although this type of analysis requires further calibration beyond the preliminary models presented here, we conclude that this modelling approach appears to be promising.Future research will refine the preliminary models described here and be directed at issues including the impact of climate change, fuel reduction strategies and economic development on successional dynamics and wildfire regimes.

Figure 1 .Fig. 1 .
Figure 1.Location map of the study areas.The land cover of Alachua County and Liguria (in red within the box) are reported on the left and on the right, respectively.
For example, using the DS-FFM, Malamud et al. (1998) estimated power-law exponents falling within the range −1.16≤ β ≤−1.02.However, for actual wildfires occurring on US Forest Service land across the contiguous United States, Malamud et al. (2005) found that power laws ranged from −1.30 (in the Mediterranean ecosystem) to −1.81 (in the subtropical ecosystem).Apparently the DS-FFM overestimates the frequency of large fires.

Figure 2 .
Figure 2. Distribution of parameter β respectively simulated by MFFM (blue), DS-FFM (green) and actual fires (red).The results obtained for the Liguria and Alachua case studies are reported above and below, respectively.

Fig. 2 .
Fig. 2. Distribution of parameter β respectively simulated by MFFM (blue), DS-FFM (green) and actual fires (red).The results obtained for the Liguria and Alachua case studies are reported above and below, respectively.

Figure 3 .Fig. 3 .
Figure 3.Long term simulation obtained using the vegetation succession model for Liguria case study.

Figure 4 .Fig. 4 .
Figure 4. Long term simulation obtained using the vegetation succession model for Alachua County case study.

Table 1 .
Parameters of the Modified Forest Fire Model; f is the probability of ignition, g is the probability of spread, p is the probability of re-growth.

Table 2 .
Comparison between observed number of fires, total burned area, and simulated distribution of the same parameters obtained by MFFM (5 classes).The mean values and the standard deviation of the number of fires and the total burned area obtained simulating 1000 different trials are reported for the 2 considered regions.5 classes # of fires Total burned area [km 2 ]

Table 3 .
Parameters of the Drossel-Schwabl Forest Fire Model; f is the probability of ignition, g is the probability of spread, p is the probability of re-growth.

Table 4 .
Comparison between observed number of fires, total burned area and simulated distribution of the same parameters by using DS-FFM (1 class).The mean values and the standard deviation of the number of fires and the total burned area obtained simulating 1000 different trials are reported for the 2 considered regions.the initial values of g i .The initial values obtained for Liguria were used both for Liguria and Alachua.The value of the parameters, after the tuning phase, is reported in Table