diff --git a/NEWS.md b/NEWS.md index 5a13914..8ac6d13 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,5 +1,6 @@ # fireSense_dataPrepFit (development version) +- New parameter `spreadIntercept` (default `FALSE`): with `TRUE` the spread formula is `~ 1 + ...`, so `fireSense_spreadFit` fits a free intercept with centred covariates (`fireSenseUtils::spreadInterceptTxt`) and the coefficients describe variation about the centre. Covariates are rescaled to 0-1, so without an intercept the level of the linear predictor is set by the coefficients alone and the climate coefficient trades off against the fuel and non-forest ones (ELF 13.1: the other coefficients explain a median 73% of CMD's variance in a final population). `FALSE` keeps the formula string `~ 0 + ...` byte for byte. Needs the `fireSenseUtils` change in PredictiveEcology/fireSenseUtils#131 (floor to be set once it has a version). - New parameter `minCovariateProp` (default `0.05`): a land-cover class covering less than this share of the ELF's flammable pixels (most recent data year) gets no spread coefficient of its own. A rare non-forest class joins the non-forest group nearest its burn rate; with less treed wetland than this there is no `treedWetland_agb` and its tree AGB stays in `dom_agb_*`/`sec_agb_*`. Treed wetland is under 5% in 7 of the 10 fitted ELFs. Needs fireSenseUtils >= 0.2.3.9081 (PredictiveEcology/fireSenseUtils#129). Fits and cached `dataPrepBuild`/`prepSpreadFitData` results of ELFs whose covariates change are invalid. - Fixed: with one fitted ELF in the SpreadFit ledger, `Init()` now sets `nonForestedLCCGroups` and `missingLCCgroup` to that ELF's groups, so this run's ignition and escape fits use the same non-forest columns (`nfLCC_*`) that diff --git a/R/spreadFormula.R b/R/spreadFormula.R new file mode 100644 index 0000000..055af76 --- /dev/null +++ b/R/spreadFormula.R @@ -0,0 +1,12 @@ +#' The spread formula +#' +#' `~ 0 + RHS`, which has no intercept, or `~ 1 + RHS` with `intercept = TRUE` +#' (`fireSenseUtils::spreadInterceptTxt` is its coefficient in `fireSense_spreadFit`). Without an +#' intercept the string is what it always was, so cache keys that hold it do not move. +#' +#' @param RHS character; the covariate terms joined with `" + "`. +#' @param intercept logical; `TRUE` for a formula with an intercept. +#' @return a character string. +spreadFormula <- function(RHS, intercept = FALSE) { + paste0("~ ", if (isTRUE(intercept)) "1" else "0", " + ", RHS) +} diff --git a/fireSense_dataPrepFit.R b/fireSense_dataPrepFit.R index 56804cb..74e2f97 100644 --- a/fireSense_dataPrepFit.R +++ b/fireSense_dataPrepFit.R @@ -113,6 +113,14 @@ defineModule(sim, list( NA, NA, paste("Name of the ledger file in `spreadFitGoogleDriveFolder`: study area polygons with", "their fitted SpreadFit parameters. `\"latest\"` (the default) takes each polygon's", "fit from the most recent ledger file that has it (`fireSenseUtils::latestSpreadFits()`).")), + defineParameter("spreadIntercept", "logical", FALSE, NA, NA, + paste("If `TRUE`, the spread formula has an intercept (`~ 1 + ...`) and `fireSense_spreadFit` fits", + "it, with the covariates centred, so the coefficients describe variation and the intercept", + "the level. Covariates are rescaled to 0-1 and so are all >= 0: without an intercept the", + "level of the linear predictor is set by the coefficients alone, and a climate coefficient trades off", + "against the fuel and non-forest ones (ELF 13.1). `FALSE` (default) keeps `~ 0 + ...`.", + "Ignored if `sim$fireSense_spreadFormula` is supplied. Needs a `fireSense_spreadFit` and a", + "`fireSense_spreadPredict` that know the intercept.")), defineParameter("targetFuelClasses", "numeric", 5, 1, 7, "the target number of unique fuel classes when using semi-automated approach"), defineParameter("useRasterizedFireForSpread", "logical", FALSE, NA, NA, @@ -939,7 +947,7 @@ prepare_SpreadFit <- function(sim) { sim$fireSense_nonAnnualSpreadFitCovariates <- nonAnnuals if (is.null(sim$fireSense_spreadFormula)) { - sim$fireSense_spreadFormula <- paste0("~ 0 + ", RHS) + sim$fireSense_spreadFormula <- spreadFormula(RHS, Par$spreadIntercept) } return(invisible(sim)) diff --git a/tests/testthat/test-metadata.R b/tests/testthat/test-metadata.R index 2304290..49f4d86 100644 --- a/tests/testthat/test-metadata.R +++ b/tests/testthat/test-metadata.R @@ -97,7 +97,7 @@ test_that("parameters are the expected names", { "forestedLCC", "fuelClassCol", "fuelCovariates", "heldOutFold", "igAggFactor", "minBufferSize", "minCovariateProp", "nonflammableLCC", "nonForestCanBeYoungAge", "scanfiVersion", "sppEquivCol", "spreadFitFilename", - "spreadFitGoogleDriveFolder", "targetFuelClasses", + "spreadFitGoogleDriveFolder", "spreadIntercept", "targetFuelClasses", "useRasterizedFireForSpread", "whichModulesToPrepare")) ) }) diff --git a/tests/testthat/test-spreadFormulaYoungAgeOnce.R b/tests/testthat/test-spreadFormulaYoungAgeOnce.R index b831f02..a194add 100644 --- a/tests/testthat/test-spreadFormulaYoungAgeOnce.R +++ b/tests/testthat/test-spreadFormulaYoungAgeOnce.R @@ -23,6 +23,7 @@ test_that("the spread formula lists a veg-table youngAge column once, not twice" env$sim <- new.env() env$sim$climateVariablesForFire <- list(spread = "CMD") env$sim$fireSense_spreadFormula <- NULL + env$Par <- list(spreadIntercept = FALSE) # the module's `Par`, as it is with the default ## the shape fireSenseVegData has after joinFireBuffersToVeg() and setnames(..., "buffer", ## "burned"): a youngAge column that is not all-zero, as in ELFs 14.3/14.4 env$fireSenseVegData <- data.table::data.table( diff --git a/tests/testthat/test-spreadIntercept.R b/tests/testthat/test-spreadIntercept.R new file mode 100644 index 0000000..7a6c3bf --- /dev/null +++ b/tests/testthat/test-spreadIntercept.R @@ -0,0 +1,35 @@ +## `spreadIntercept` (default FALSE): the spread formula has an intercept, `~ 1 + RHS`, which +## fireSense_spreadFit fits with centred covariates. Off, the formula string is byte-identical to what it +## has always been, so every cache key that holds it is unchanged. + +test_that("spreadFormula() gives '~ 0 + RHS' by default and '~ 1 + RHS' with an intercept", { + rhs <- "CMD + youngAge + fuelA" + expect_identical(spreadFormula(rhs), "~ 0 + CMD + youngAge + fuelA") + expect_identical(spreadFormula(rhs, intercept = FALSE), "~ 0 + CMD + youngAge + fuelA") + expect_identical(spreadFormula(rhs, intercept = TRUE), "~ 1 + CMD + youngAge + fuelA") + ## the default string is the one the module built before the parameter existed + expect_identical(spreadFormula(rhs), paste0("~ 0 + ", rhs)) + ## and the intercept is what terms() sees + expect_equal(attr(terms(as.formula(spreadFormula(rhs, TRUE))), "intercept"), 1L) + expect_equal(attr(terms(as.formula(spreadFormula(rhs))), "intercept"), 0L) + expect_identical(fireSenseUtils::spreadDesignCols(spreadFormula(rhs, TRUE)), + c(fireSenseUtils::spreadInterceptTxt, "CMD", "youngAge", "fuelA")) +}) + +test_that("the module defines spreadIntercept, a logical that is FALSE by default", { + def <- Filter(function(x) is.call(x) && identical(x[[1]], as.name("defineModule")), + parse(testthat::test_path("..", "..", "fireSense_dataPrepFit.R"), keep.source = FALSE)) + params <- eval(as.list(def[[1]][[3]])[["parameters"]], envir = asNamespace("SpaDES.core")) + i <- match("spreadIntercept", params$paramName) + expect_false(is.na(i)) + expect_identical(params$paramClass[[i]], "logical") + expect_identical(params$default[[i]], FALSE) +}) + +test_that("prepare_SpreadFit() builds the formula with spreadFormula() from the spreadIntercept parameter", { + src <- paste(deparse(body(prepare_SpreadFit)), collapse = "\n") + expect_match(src, "spreadFormula\\(RHS, Par\\$spreadIntercept\\)") + expect_no_match(src, '"~ 0 \\+ "', fixed = FALSE) + ## the formula is only built when none was supplied + expect_match(src, "is.null\\(sim\\$fireSense_spreadFormula\\)") +})