Remote Sensing of Environment 92 (2004) 363 – 369 www.elsevier.com/locate/rse
Mapping wildfire occurrence at regional scaleB Juan de la Rivaa,*, Fernando Pe´rez-Cabelloa, Noemı´ Lana-Renaulta, Nikos Koutsiasb,1 a
Facultad de Filosofia y Letras, Department of Geography, University of Zaragoza, Calle Pedro Cerbuna 12, Zaragoza E-50009, Spain b Department of Geography, University of Zurich, Winterthurerstrasse 190, Zurich CH-8057, Switzerland Received 19 July 2003; received in revised form 31 May 2004; accepted 8 June 2004
Abstract When assessing fire danger, interpolation of the dependent variable—historic fire occurrence—is required in order to statistically compare and analyze it with human factors, environmental parameters and census statistics. To confirm the compatibility between the distinct data types, occasionally, for this kind of spatial analysis, historical observations of the primary wildland fire (given as x and y coordinates) must be transformed either to continuous surfaces or to area data. The simple overlay approach converts single point observations to area data. However, this procedure assumes lack of spatial uncertainties that would otherwise result in serious errors caused by the positional inaccuracies of the original point observations. Here, we used kernel density interpolation to convert the original data on wildland fire ignition into an expression of areal units, defined by a raster grid and, subsequently, by the administrative borders of the municipalities in two study areas in Spain. By overlaying a normal bivariate probability density function (kernel) over each point observation, each ignition point was considered an uncertain point location rather than an exact one. D 2004 Elsevier Inc. All rights reserved. Keywords: Wildfire; Kernel density interpolation; Ignition point
1. Introduction Within the framework of the Firerisk project,2 the spatialisation of fire occurrence as a dependent variable has been a necessary requirement in the fire risk modeling. Traditional methods based on occurrence indexes (e.g. number of fires related to wildland area) have been shown to be inadequate when exploring statistical relations with causal factors. Occurrence indexes are often calculated for vectorial units (e.g. municipalities) while casual factors can have a continuous behaviour (e.g. climatic variables), high B This article was previously published in Remote Sensing of Environment 92(2). Doi of original article: 10.1016/j.rse.2004.06.013. * Corresponding author. Tel.: +34 976 76 10 00; fax: +34 976 76 15 06. E-mail addresses:
[email protected] (J. de la Riva)8
[email protected] (F. Pe´rez-Cabello)8
[email protected] (N. Lana-Renault)8
[email protected] (N. Koutsias). 1 Fax: +4116356848. 2 http://www.geogra.uah.es/proyectos/firerisk/. 0034-4257/$ - see front matter D 2004 Elsevier Inc. All rights reserved. doi:10.1016/j.rse.2004.06.022
spatial variability (e.g. land use) or punctual or lineal representation (e.g. roads). Therefore, new solutions of fire pattern spatialisation must be investigated. Fire occurrence data in Spain were recorded before 1998 both on a UTM 1010-km grid (i.e. x and y coordinates for each fire at 10-km resolution) and at administrative level (i.e. number of fires per municipality). Because of coarse grid resolution, these records introduce an enhanced degree of uncertainty in fire location, which may be further propagated during modeling. The success to explain the spatial distribution of fire patterns and to propose and test hypotheses about underlying causal factors depends on the quality of fire ignition observations. Particularly, the use of other geo-referenced data for extracting complementary information (i.e., overlay point observations on other spatial data, proximity distances to other features) may introduce significant errors into analysis. The lack of appropriate fire occurrence data, in terms of their content and accuracy, has a significant
364
J. de la Riva et al. / Remote Sensing of Environment 92 (2004) 363–369
impact on the theoretical and applied research on wildland fires and on their management. Even countries heavily affected by wildland fires often do not have proper data on fire incidence (Martı´n et al., 1994). Multiple types of spatial data are available. For instance, point data are used to represent wildland fire ignition locations, while area data are used to represent census statistics that describe socio-economic and demographic characteristics. Many complications may arise in relating one type of data to the other when performing a statistical comparison and analysis (Flowerdew & Pearce, 2001). Moreover, the location uncertainty of these data raises further questions about the integrity of comparisons. The simple overlay approach converts single point observations to area data (i.e. number of points falling inside an areal unit). However, this procedure can produce serious artifacts caused by the positional inaccuracies of the original point observations. Koutsias et al. (in press) found that wildland fire occurrence patterns shown by the overlaying approach (i.e. a regular grid of quadrats superimposed over positional inaccurate ignition points) can be highly inconsistent depending on the magnitude of the positional errors and the resolution of the grid. In the same study, it was proposed that an increment in cell size eliminates these problems by simultaneously reducing the level of detail, which results in loss of spatial variability. Interpolation, as a method to predict attribute values at unsampled locations from observations sampled inside the study area, can be used to convert data from point observations to continuous fields (Burrough & McDonnel, 1998). Several interpolation techniques are available for this purpose. However, most require a variable to be estimated as a function of location. In contrast, kernel density estimation can be used as an interpolation technique for individual point observations (Levine, 2002). Originally, this approach was developed as an alternative method to obtain a smooth probability density function—univariate or multivariate—from a sample of observations, i.e. histogram (Bailey & Gatrell, 1995; Levine, 2002). Since the estimation of the intensity of point observations (given in x and y coordinates) is very similar to the bivariate probability density one, the kernel approach can be adapted for this purpose (Bailey & Gatrell, 1995). Kernel density estimation, a non-parametric statistical method for estimating probability densities, has been extensively used for home range estimation in wildlife ecology (Seaman & Powell, 1996; Tufto et al., 1996; Worton, 1989). By converting original wildland fire ignition locations to continuous density surfaces and simultaneously addressing some of their inherent positional inaccuracies, this technique is also very useful in defining spatial fire occurrence patterns at landscape level (Koutsias et al., in press). A kernel (i.e. normal bivariate probability density) is placed over each point observation and the intensity at each intersection of a superimposed grid is estimated (Seaman & Powell, 1996). The method
is similar to the bmoving windowQ concept, where a window of fixed-size is moved over the point observations, except in this case the window is replaced by a 3D function (Gatrell et al., 1996). Mathematically, the kernel density estimator is defined as (Seaman & Powell, 1996; Silverman, 1986; Worton, 1989): n 1 X x Xi fˆ ð xÞ ¼ 2 K nh i1 h where n is the number of point observations, h is the bandwidth, K is the kernel, x is a vector of coordinates that represent the location where the function is being estimated, and X i are vectors of coordinates that represent each point observation. The bandwidth expresses the size of the kernel and controls the interpolation results. Depending on whether a fixed value or multiple adaptive values are used for the bandwidth (smoothing parameter), the kernel is distinguished into the fixed and adaptive mode, respectively. Regardless of the mode, a kernel must be selected from a variety of functions, for instance normal distribution, triangular function, quartic function, etc., although the normal distribution is the most commonly used (Levine, 2002). Finally, the smoothing parameter must be taken in accordance with the rule that narrow kernels allow nearby point observations to have a greater effect on the density estimation than wide kernels (Seaman & Powell, 1996). In this regard, the size of the bandwidth controls the degree of smoothing of density estimations. The goal of the present study is to spatialize fire occurrence data as an input for fire risk modeling by using a kernel approach to interpolate historic fire observations. The analysis has been applied in two distinct mountain areas in Spain to prove that such a technique can be applied in the field of forest fires. In fire risk modeling, fire occurrence can be considered a dependent variable; this analysis implies continuous data and a more accurate location of the ignition points in order to obtain improved interpolation. When applying kernel interpolation, each fire ignition point was considered not as an exact point location but rather an uncertain one that defines a broader surrounding area within which the true point observation lies. A bivariate
Fig. 1. Density estimation is calculated by placing a kernel (i.e. bivariate normal probability density) over each wildland fire ignition point and estimating the intensity at each intersection of a superimposed grid. The mean density value was then obtained superimposing each administrative areal unit to the resulting kernel density surfaces.
J. de la Riva et al. / Remote Sensing of Environment 92 (2004) 363–369
365
Table 1 Analysis of the parameters of bandwidths Parameter
Pre-Pyrenees a
Mean polygon size (s) Diagonal of a theoretical square (D) Length of the theoretical radius (r) Mean number of ignitions points per polygon (N) Mean random distance (RDmean)2 Total acreage (including surrounding area) Number of ignition points (including surrounding area) Global mean random distance2 a
Fig. 2. Study areas.
probability density function was overlaid on each point observation (Fig. 1). Since fire managers frequently work with data that refers to administrative units (i.e. municipalities), fire occurrence is also represented at municipality level, by superimposing these units to the resulting kernel density surfaces and considering the mean density value.
2. Materials and methods 2.1. Study area Two study areas with similar physical characteristics but distinct administrative organizations and fire patterns were
2
Iberian range
28.6 km 7562.8 m 3781.4 m 3.8
18.58 km2 6097.2 m 3048.6 m 2.3
2757.4 m 9301.2 km2
2842 m 9111.7 km2
1220
1134
2761.1 m
2834.6 m
Polygons b5 km2 were not considered.
selected: the Central Spanish Pre-Pyrenees and the Eastcentral Iberian range (Fig. 2). These two areas are located in Mediterranean mountain environments, and they can be classified as high-risk areas for wildfires (Pe´rez-Cabello & de la Riva, 2001). The Central Spanish Pre-Pyrenees comprises an area of 4192 km2 with complex topography and altitudes that range from 500 to 1700 m. Vegetation is dominated by Pinus sylvestris, Pinus nigra (most afforested), Quercus faginea, Buxus sempervirens (indicative of some oceanic influences), Aphyllantes monspeliensis, etc. The size of the municipalities is highly heterogeneous since some cover more than 600 km2, while others are smaller than 15 km2. Socio-economic changes in the mid-20th century led to the abandonment of farming activities and intense emigration. Nowadays, recreational activities have increased in specific zones. Most fires are caused by humans: between 1983 and 2001 there were 616 fires, of
Fig. 3. Spatial reference units: polygons produced from overlaying the UTM grid (1010 km) and the municipality boundaries. (a) Pre-Pyrenees, and (b) Iberian range.
366
J. de la Riva et al. / Remote Sensing of Environment 92 (2004) 363–369
which 55% were caused by humans and 45% by lightning. The East-central Iberian range area occupies 4060 km2 with elevations ranging from 400 to 1300 m. The vegetation consists predominantly of Pinus pinaster, P. nigra, Quercus ilex rotundifoliae, Quercus coccifera and Brachypodium ramosum. The size of the municipalities is fairly homogeneous with a mean size of 39.8 km2. Similar to the Pre-Pyrenees, the area has suffered a drastic decline in population and agricultural activities over the last century but no recreational activities have been developed. Between 1983 and 2001, there were 572 fires, of which 56% were caused by humans and 44% by lightning.
2.2. The kernel approach in transforming point data to area data Because fire occurrence data, obtained from the official Spanish wildfire census, were provided on a UTM 1010km grid and at municipality level, there was no information on the exact x/y UTM position of the ignition points. To improve the accuracy of fire location, a new spatial reference system was designed. Data were referenced by randomly sampling within each polygon created after overlaying the UTM grid (1010 km) and the municipality boundaries (Fig. 3). Within each bnew polygonQ, where the number of fires is known, points were randomly positioned throughout the wildland area only (forest, shrub and grass
Fig. 4. Fire densities using the kernel density approach at various bandwidths to randomly distributed points in Pre-Pyrenees.
J. de la Riva et al. / Remote Sensing of Environment 92 (2004) 363–369
367
Fig. 5. Fire densities using the kernel density approach at various bandwidths to randomly distributed points in Iberian range.
areas). Using this random sampling, we established fire ignitions points at a finer spatial resolution. Fire data from a wider area were included to preserve the effect of the external points and to minimize problems associated with edge effect. Including the surrounding area, 1220 and 1134 random points were introduced for the Pre-Pyrennes and Iberian range areas, respectively. Kernel density interpolation was then applied to these fire ignition points using the fixed mode approach (i.e. constant value for the smoothing parameter) and a bivariate normal probability density function. We did not use the adaptive mode since the point observations were treated in a distinct way according to their concentration in space (Worton, 1989). Fire densities were estimated at a grid resolution of 100 by 100 m. CrimeStatR,3 a spatial statistics program for the analysis of crime incident locations, was used to perform kernel density interpolation (Levine, 2002). The size of bandwidth (i.e. standard deviation of the normal distribution) is critical because it determines the degree of smoothing in the density output surfaces. Bandwidth value depends on the scale adopted and the specific characteristics of the study case related to the spatial fire pattern. This implies knowledge of the mean polygon size and the mean number of ignition points within each. Several methods were tested to define the appropriate size of the smoothing parameter of the kernel. 3 CrimeStatR V. 2.0 is available on http://www.icpsr.umich.edu/ NACJD/crimestat.html.
–The first method was based solely on the mean polygon size, assuming the polygon as a theoretical square with the same size. In this case, a theoretical distance was estimated on the basis of the length of the theoretical radius (r): r ¼ D=2 where D is the diagonal of a theoretical square. –The second considered the mean random distance calculations (RDmean), on the basis of a local approach (i.e. mean polygon size and mean number of ignition points per polygon) and on a global one (i.e. total size of the study area and total number of ignition points). RDmean is mathematically defined as: 1 RDmean ¼ 2
rffiffiffiffiffi A N
where A is the mean size polygon and N is the mean number of ignitions points falling inside the polygons. On the basis of previous experience, the double of the RDmean value was decided to be used for bandwidth definition (Koutsias et al., in press). –In the third method, the effect of the randomly distributed points on kernel density outputs at certain bandwidths was evaluated. Random sampling was performed using a specific script of ArcView 3.2; each time the script was applied, a distinct sampling distribution was obtained. To test the sensitivity of the bandwidth to the randomness of the ignition points distribution, a correlation
368
J. de la Riva et al. / Remote Sensing of Environment 92 (2004) 363–369
analysis between the results obtained in the three random sampling for the each bandwidth was applied. The Pearson coefficient shows the bandwidth which is less affected by the randomness of the ignition points distribution. –Finally, a visual-subjective approach was used. These estimators define a range of values used as indicators for selecting the bandwidth. However, only after the analysis of the results was the final bandwidth chosen.
3. Results and discussion 3.1. Mapping fire densities The bandwidth parameters, estimated from the methods described in the previous section, are summarized in Table 1. Although the total area and the total number of ignition points in the two study areas were almost the same, the mean polygon size differed considerably because of the greater number of municipalities, and therefore polygons, in the Iberian range. This accounted for similar results: in the Pre-Pyrenees, the theoretical radius was 3781 m, the RDmean 2757 m and the global mean random distance 2761 m; while in the Iberian range these values were 3049, 2842 and 2835 m, respectively. To perform a visual-subjective evaluation, distinct bandwidths were tested: – –
2500, 3250, 4000, 5000 and 7500 m in the Pre-Pyrenees (Fig. 4), 2500, 3000 and 5000 m in the Iberian range area (Fig. 5),
and the best results were obtained with bandwidths of 3250 and 4000 m in the Pre-Pyrenees and 3000 m in the Iberian range. A narrower bandwidth allowed a high effect of the localization of the established ignition points while a wider one introduced excessive smoothing. The effect of the sampling method to establish the fire ignition points was also considered. The correlation analyses applied between the three kernel density outputs resulting from the three random samplings (Table 2) show that the density results for the Pre-Pyrenees using a 2500-m bandwidth are affected more by the method used to locate the ignition points (mean Pearson correlation coeffi-
Table 3 Correlation analysis to evaluate the effect of three random distribution points (1, 2, 3) in the Iberian range Bandwidth
1–2
1–3
2–3
Mean
2500 3000 5000
0.87 0.91 0.97
0.84 0.90 0.96
0.86 0.90 0.97
0.86 0.90 0.97
Table shows the Pearson correlation coefficients for each random pattern and bandwidth; mean value is also included.
cient=0.89) than for a 3250-m bandwidth (mean Pearson correlation coefficient=0.93). Differences between bandwidths of 4000, 5000 or 7500 m were not as significant (mean Pearson correlation coefficient is 0.95 to 0.99). For the Iberian range, the same analysis showed that results were less affected using a 3000 m bandwidth (Table 3, mean Pearson correlation coefficient=0.90). According to the previous calculations, the appropriate bandwidth in the Pre-Pyrenees ranges between 2750 and 3800 m and we chose a width of 3250 m. For the Iberian range area, the appropriate bandwidth ranges between 2800 and 3100 m and the bandwidth selected was 3000 m. 3.2. Summarizing fire densities at administrative level Application of the data on fire densities to administrative units involves homogenizing fire occurrence to a single value for each municipality and consequently the loss of local spatial distribution. However, the use of these units, at regional scale, is usually a requirement for fire management. The value densities obtained for each grid cell in the interpolation applied maintains the sample size. Therefore, these densities sum the total number of fires considered in the random sampling process, and express the probability of fire occurrence for each cell in relation with the total number of fires. The final result for each administrative unit
Table 2 Correlation analysis to evaluate the effect of three random distribution points (1, 2, 3) in the Pre-Pyrenees area Bandwidth
1–2
1–3
2–3
Mean
2500 3250 4000 5000 7500
0.90 0.94 0.96 0.97 0.99
0.88 0.92 0.95 0.97 0.99
0.90 0.93 0.95 0.97 0.99
0.89 0.93 0.95 0.97 0.99
Table shows the Pearson correlation coefficients for each random pattern and bandwidth; mean value is also included.
Fig. 6. Fire density at municipality level using the mean kernel density value (3250 m bandwidth) in Pre-Pyrenees.
J. de la Riva et al. / Remote Sensing of Environment 92 (2004) 363–369
369
several methods. The geometric estimator (RDmean) and the analysis of the effect of the sampling method provide the most appropriate bandwidth. However, these estimators define a range of values rather than a single one. Data on the spatial distribution of fire occurrence in administrative areas is useful for fire risk analysis and fire management, even if they homogenize the risk to a single value. However, representation as a continuous surface preserves a more realistic pattern of fire occurrence according to the considered scale, and thus allows the spatial analysis of the causal factors.
Acknowledgements This research was supported by the Spanish Ministry of Science and Technology (contract AGL2000-0842): FIRERISK project (Remote Sensing and Geographic Information Systems for forest fire risk estimation: an integrated analysis of natural and human factors). Fig. 7. Fire density at municipality level using the mean kernel density value (3000 m bandwidth) in Iberian range.
References (municipality) is the mean of values inside the study area. Fig. 6 (Pre-Pyrenees) and Fig. 7 (Iberian range) show the probability of occurrence in each municipality, expressed in five categories: very high, high, medium, low, and very low.
4. Concluding remarks Data on the spatial distribution of fire occurrence is one of the most common requirements for wildfire danger assessment. This data is essential in explaining wildfire causal factors. Although the current positioning system allows accurate location of ignition points, there is still a substantial lack of information, particularly for historic fire data. In Spain, x/y UTM coordinates to track fires have been used only since 1998; before, occurrence was recorded both on a UTM 1010-km grid and at municipality level. Here, we used kernel density interpolation to spatially define historic fire occurrence. In contrast to the overlay approach, where the locations of wildland fire ignition are considered as exact points, in the kernel approach they are taken as spatially uncertain points, achieved by placing a normal bivariate probability density over each event. Our results show that bandwidth is critical since it determines the degree of smoothing in fire density results. A procedure including several methods to define the bandwidth size was followed. Bandwidth value depends on both the scale adopted and the specific characteristics of the study area, especially those related to the spatial fire pattern. Therefore, two distinct study areas were chosen to provide rigorous results. The analysis reveals that there is no single method, as the best results are obtained by combining
Bailey, T. C., & Gatrell, A. C. (1995). Interactive spatial data analysis (pp. 84 – 88). England7 Longman. Burrough, P. A., & McDonnel, R. A. (1998). Principles of geographical information systems (pp. 98 – 99). Oxford7 Oxford Univ. Press. Flowerdew, R., & Pearce, J. (2001). Linking point and area data to model primary school performance indicators. Geographical and Environmental Modelling, 5, 23 – 41. Gatrell, A. C., Bailey, T. C., Diggle, P. J., & Rowlingsont, B. S. (1996). Spatial point pattern analysis and its application in geographical epidemiology. Transactions of the Institute of British Geographers, 21, 256 – 274. Koutsias, N., Kalabokidis, K. D., & Allgfwer, B. (in press). Fire occurrence patterns at landscape level: beyond positional accuracy of ignition points with kernel density estimation methods. Natural Resource Modeling, (in press). Levine, N., 2002. CrimeStat II: A Spatial Statistics Program for the Analysis of Crime Incident Locations (version 2.0). Ned Levine and Associates Annandale, VA and The National Institute of Justice, Washington, DC. Martı´n, M. P., Viedma, D., & Chuvieco, E. (1994). High versus low resolution satellite images to estimate burned areas in large forest fires. In D. X. Viegas (Ed.), 2nd International Conference of Forest Fire Research (pp. 653 – 663). University of Coimbra, Coimbra, Portugal7 ADAI. Pe´rez-Cabello, F., & de la Riva, J., 2001. Forest fires and land degradation in Spain. The Huesca Western Pre-Pyrenees case study. Keynote in the workshop bLandnutzungswandel und Landdegradation in SpanienQ, Frankfurt am Main, Germany. Seaman, D. E., & Powell, R. A. (1996). An evaluation of the accuracy of kernel density estimators for home range analysis. Ecology, 77, 2075 – 2085. Silverman, B. W. (1986). Density estimation for statistics and data analysis (pp. 7 – 94). London, England7 Chapman & Hall. Tufto, J., Andersen, R., & Linnell, J. (1996). Habitat use and ecological correlates of home range size in a small cervid: the roe deer. Journal of Animal Ecology, 65, 715 – 724. Worton, B. J. (1989). Kernel methods for estimating the utilization distribution in home-range studies. Ecology, 70, 164 – 168.