-
Notifications
You must be signed in to change notification settings - Fork 4
/
scfm.Rmd
172 lines (151 loc) · 9.92 KB
/
scfm.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
---
title: "scfm"
author: "Ian"
date: "December 18, 2018"
output: pdf_document
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```
## Overview
Fire is modelled over a grid of raster cells (pixels) and is treated as a three stage stochastic process of ignition, escape from the cell of origin, and subsequent spread. The original fire model is described by Cumming et al. (1998), and more accessibly in Armstrong and Cumming (2003), and was recently implemented as a collection of SpaDES modules by Cumming, McIntire, and Eddy (in prep). Each fire module corresponds to one of the processes being modelled). The Ignition stage determines if a fire actually starts in a particular grid cell in a given model time step. The Escape stage models the effect of fire suppression, or other ecological or sampling effects that alter the lower end of the fire size distribution. This is simulated in the model by determining a probability that a fire stays within the cell of origin, and accounts for lakes and other non-flammable geographic features. The Spread stage simulates fire growth, using a spread parameter that determines the probability of a fire spreading across a cell boundary, from a burning cell into an unburnt neighbour.
In addition to these core fire modules, the model also makes use of several support modules involved in parameter estimation and model calibration. Fire regime parameters are estimated from historic fire data using the Canadian National Fire Database (Canadian Forest Service, n.d.) in order to reproduce realistic fire patterns. Parameterization and calibration of the model can be done separately for each geographic region used in the model (e.g., by ecoregion), allowing finer scale parameterization across large geographic areas. Individual spread parameters are determined in each geographic region following a four-part process. First, the region is buffered by a set distance and a flammability map is generated for the buffered region using the 2005 Land Cover Map of Canada (Latifovic et al., 2005). Next, fires are repeatedly simulated on the landscape across a range of spread probabilities. These fires cannot start in the buffered area, but may spread to and from the buffer. This reduces the influence of edge effects in determining mean fire size. Then the spread probabilities and resulting fire sizes are fit with a shape-constrained additive model (SCAM) (Pya and Wood, 2015). The SCAM is monotonic to ensure fire size increases for any incremental increase in spread probability. Lastly, an optimization algorithm is used to predict the spread probability that will generate the estimated mean fire size for the region.
This version of the model considers lightning-caused fires only. Fires are started by applying a random number sequence to each cell independently. Fires are extinguished when no further spread events (jumps) occur. After a fire, forested cells are reset to age 1. Thus, fires are presumed to be stand initiating. These processes result in irregular patches of variously sized burns. The locations of fires are not constant across simulations because of stochastic patterns of fire ignition and growth. The model tracks the sizes and compositions of each burn. Fire sizes follow a probability distribution largely determined by the spread parameter(s). The current implementation does not account for differences in vegetation type with respect to fire spread (e.g., flammability), nor is it climate sensitive. The more experimental fireSense model (Marchal, Cumming, and McIntire 2017), recently implemented in SpaDES, could be used in future updates to explicitly account for climate and vegetation impacts on fires (in prep).
## References
Armstrong, Glen W, and Steven G Cumming. 2003. “Estimating the Cost of Land Base Changes Due to Wildfire Using Shadow Prices.” Forest Science 49 (5): 719–30. https://doi.org/10.1093/forestscience/49.5.719.
Cumming, Steve G, D Demarchi, and C Walters. 1998. “A Grid-Based Spatial Model of Forest Dynamics Applied to the Boreal Mixedwood Region.” Working Paper 1998–8. Sustainable Forest Management Network. https://doi.org/10.7939/R35N33.
Latifovic, Rasim and Darren Pouliot, 2005. "Multi-temporal land cover mapping for Canada: Methodology and Products," Canadian Journal of Remote Sensing, Vol. 31, no. 5, pp. 347-363.
Marchal, Jean, Steve G Cumming, and Eliot J B McIntire. 2017. “Land Cover, More than Monthly Fire Weather, Drives Fire-Size Distribution in Southern Québec Forests: Implications for Fire Risk Management.” PLoS ONE 12 (6): 1–17. https://doi.org/10.1371/journal.pone.0179294.
Pya, Natalya, and Simon N Wood, 2015. "Shape constrained additive models," Statistics and Computing, Vol. 25, Iss. 3, pp 543-559.
```{r load-SpaDES, eval=TRUE}
library(magrittr)
library(raster)
library(SpaDES)
```
## Usage example
```{r module usage example, eval = FALSE}
cloudFolderID <- "https://drive.google.com/open?id=1PoEkOkg_ixnAdDqqTQcun77nUvkEHDc0"
#Parameters
timeunit <- "year"
times <- list(start = 0, end = 20)
defaultInterval <- 1.0
defaultPlotInterval <- 1.0
defaultInitialSaveTime <- NA #don't be saving nuffink
parameters <- list(
.progress = list(type = "text", interval = 1),
scfmLandcoverInit = list(
.plotInitialTime = times$start,
.plotInterval = defaultPlotInterval,
.saveInitialTime = defaultInitialSaveTime,
.saveInterval = defaultInterval),
ageModule = list( #age is really time since fire, regardless of landcover
initialAge = 100,
maxAge = 200,
returnInterval = defaultInterval,
startTime = times$start,
.plotInitialTime = times$start,
.plotInterval = defaultPlotInterval,
.saveInitialTime = defaultInitialSaveTime,
.saveInterval = defaultInterval),
scfmIgnition = list(
pIgnition = 0.0001,
returnInterval = defaultInterval,
startTime = times$start,
.plotInitialTime = NA,
.plotInterval = defaultPlotInterval,
.saveInitialTime = defaultInitialSaveTime,
.saveInterval = defaultInterval),
scfmEscape = list(
p0 = 0.05,
returnInterval = defaultInterval,
startTime = times$start,
.plotInitialTime = NA,
.plotInterval = defaultPlotInterval,
.saveInitialTime = defaultInitialSaveTime,
.saveInterval = defaultInterval),
scfmSpread = list(
pSpread = 0.235,
returnInterval = defaultInterval,
startTime = times$start,
.plotInitialTime = times$start,
.plotInterval = defaultPlotInterval,
.saveInitialTime = defaultInitialSaveTime,
.saveInterval = defaultInterval),
scfmRegime = list(fireCause=c("L", "H")),
scfmDriver = list(targetN = 1000) #this is an important parameter.
#increase targetN for more robust estimates, longer run-time
)
modules <- list("group_scfm")
#Paths
paths <- list(
cachePath = file.path("cache"),
modulePath = file.path("modules"),
inputPath = file.path("inputs"),
outputPath = file.path("outputs")
)
#RasterToMatch:
#if you supply studyArea you should supply rtm to make sure the crs are identical.
#if run with no studyArea, the default is a small area in southwest Alberta
objects <- list(
"cloudFolderID" = cloudFolderID
# studyArea = tempStudyArea,
# rasterToMatch = rtm
)
#Run module
#The calibration process may take hours (it's cached)
options("reproducible.cachePath" = paths$cachePath)
options("reproducible.useMemoise" = FALSE)
mySim <- simInit(times = times, params = parameters, modules = modules,
objects = objects, paths = paths)
dev()
clearPlot()
set.seed(23657)
outSim <- SpaDES.core::spades(mySim, progress = FALSE, debug = TRUE)
```
```{r module evaluate module output, eval = FALSE}
#for comparing parameters of interest
comparePredictions <- function(polyList, simList) {
out <- lapply(polyList, FUN = function(x, sim = simList) {
regime <- sim$scfmRegimePars[[x]]
driver <- sim$scfmDriverPars[[x]]
landscapeAttr <- sim$landscapeAttr[[x]]
firePoints <- outSim$firePoints[outSim$firePoints$PolyID == x,]
hist_MAAB <- sum(firePoints$SIZE_HA[firePoints$SIZE_HA > landscapeAttr$cellSize])/
(landscapeAttr$burnyArea*(sim@params$scfmRegime$fireEpoch[2] - sim@params$scfmRegime$fireEpoch[1] + 1)) * 100
#This is a long way of saying, sum of fires/ (flammable landscape * fire epoch )
#hist_mfs will be NaN if there were no fires larger than one pixel
pSpread <- driver$pSpread
burnSum <- sim$burnSummary[sim$burnSummary$polyID == x,]
burnSum$N <- as.numeric(burnSum$N)
burnSum$areaBurned <- as.numeric(burnSum$areaBurned)
burnSum <- burnSum[burnSum$N > 1]
mod_MAAB <- sum(burnSum$areaBurned)/(landscapeAttr$burnyArea * (times(sim)$end - times(sim)$start)) * 100
pred <- data.frame("PolyId" = x, #Polygon ID
"histMeanSize" = regime$xBar, #The predicted (empirical) mean size of fires
"modMeanSize" = mean(burnSum$areaBurned), #The modeled mean size of fires
"pSpread" = pSpread, # The spread probability estimated from the SCAM model
"hist_MAAB" = hist_MAAB,#The empirical mean annual area burned (from NFDB 1970-2000)
"mod_MAAB" = mod_MAAB) #The modelled mean annual area burned
return(pred)
})
return(out)
}
df <- comparePredictions(names(outSim$scfmDriverPars), outSim) %>%
rbindlist(.)
#Some useful plots
breaks = c(0.20, 0.21, 0.22, 0.23, 0.24, 0.25)
ggplot(df, aes(x = histMeanSize, y = modMeanSize, color = pSpread)) +
geom_point(aes(histMeanSize, modMeanSize)) +
scale_color_gradient(low = "yellow", high = "red")+
labs(y = "modeled mean fire size size", x = "historical mean size fire size") +
theme_minimal() +
#geom_text(aes(label = PolyId), hjust = 0, vjust = 0)
geom_abline(slope = 1)
ggplot(df, aes(x = hist_MAAB, y = mod_MAAB, col = pSpread)) +
geom_point(aes(hist_MAAB, mod_MAAB)) +
labs(y = "model mean annual area burned (%)", x = "empirical mean annual area burned (%)") +
scale_color_gradient(low = "yellow", high = "red")+
theme_minimal() + #+
geom_abline(slope = 1)
# xlim(c(0, 8000)) +
# ylim(c(0, 12000)) # +