Apply youngAge mutual exclusivity at prediction time, matching the fit

NEWS.md

@@ -1,5 +1,11 @@
 # fireSense_SpreadPredict (development version)
 
+- `spreadProbOneELF()` (fireSense_SpreadPredict.R:282 pre-fix) called
+  `fireSenseUtils::spreadProbFromIntegerCovs()` with `mutuallyExclusive = NULL`, so a young pixel's
+  fuel biomass and non-forest land-cover columns reached the logistic unchanged instead of being
+  zeroed with `youngAge`, as the fit requires. Prediction now derives the same
+  `youngAge`-exclusivity rule the fit uses, via the new `fireSenseUtils::youngAgeExclusiveCols()`.
+  Requires `fireSenseUtils@development (>= 0.2.3.9048)`. Version 1.0.0.9006.
 - New parameter `.studyAreaName` (default `NA`), the name PredictiveEcology modules use for the study area. This module does not use it yet.
 - The fitted per-year random effect (`yearSpreadSD`, fireSense_SpreadFit / fireSenseUtils >= 0.2.3.9041) is no
   longer read as a covariate coefficient: with it, a single ELF treated it as a fourth logistic parameter and several

fireSense_SpreadPredict.R

@@ -10,14 +10,14 @@ defineModule(sim, list(
     person("Alex M.", "Chubaty", email = "achubaty@for-cast.ca", role = "ctb")
   ),
   childModules = character(),
-  version = list(fireSense_SpreadPredict = "1.0.0.9005", SpaDES.core = "0.1.0"),
+  version = list(fireSense_SpreadPredict = "1.0.0.9006", SpaDES.core = "0.1.0"),
   timeframe = as.POSIXlt(c(NA, NA)),
   timeunit = "year",
   citation = list("citation.bib"),
   documentation = list("README.txt", "fireSense_SpreadPredict.Rmd"),
   reqdPkgs = list("magrittr", "Matrix", "methods", "terra", "SpaDES.core (>=3.0.4)", "stats",
                   "ggplot2", "viridis",
-                  "PredictiveEcology/fireSenseUtils@development (>= 0.2.3.9038)"),
+                  "PredictiveEcology/fireSenseUtils@development (>= 0.2.3.9048)"),
   parameters = bindrows(
     defineParameter(name = "lowerSpreadProb", class = "numeric", default = 0.13,
                     desc = "Lower asymptote of the 2- and 3-parameter logistic."),
@@ -248,9 +248,12 @@ spreadProbOneELF <- function(covs, params, covMinMax, formula, yr, maxFireSpread
   ## fireSenseUtils::fuelLinearRange, c(0, 1e4), as that covariate's covMinMax_spread, and its
   ## coefficients only mean anything for biomass / 1e4: undo the log with the function the fit used.
   ## A fit made on the log scale has the log range there, and its covariates are left as they are.
+  fuelCols <- character()
   for (cn in intersect(names(covMinMax), names(covs))) {
-    if (fireSenseUtils::isLinearFuelRange(covMinMax[[cn]]))
+    if (fireSenseUtils::isLinearFuelRange(covMinMax[[cn]])) {
       covs[[cn]] <- fireSenseUtils::fuelLogToLinear(covs[[cn]])
+      fuelCols <- c(fuelCols, cn)
+    }
   }
 
   ## the covariates this ELF was fitted with
@@ -275,11 +278,22 @@ spreadProbOneELF <- function(covs, params, covMinMax, formula, yr, maxFireSpread
   shortAnnDTx1000 <- toX1000(list(covs))[[1]] |> setDT()
   colsToUse <- setdiff(names(covs), "pixelID")
 
+  ## youngAge is mutually exclusive with every other non-climate covariate, exactly as in the fit
+  ## (fireSense_SpreadFit::spreadFitPrep()): wherever youngAge is non-zero, fuel biomass, nfLCC_*
+  ## and treedWetland are zero. fireSenseUtils::fireSenseCovariatesCreate() already applies this
+  ## when the covariates are built, but that must not be this module's only defence -- it is
+  ## re-applied here so a young pixel's covariates are correct regardless of how they arrived.
+  mutuallyExclusive <- if (fireSenseUtils::youngAgeTxt %in% colsToUse) {
+    fireSenseUtils::youngAgeExclusiveCols(colsToUse, fuelCols = fuelCols)
+  } else {
+    NULL
+  }
+
   shortAnnDT <-
     spreadProbFromIntegerCovs(shortAnnDTx1000 = shortAnnDTx1000,
                               yr = yr,
                               covMinMax = covMinMax,
-                              mutuallyExclusive = NULL, # alraedy done in dataPrepPredict
+                              mutuallyExclusive = mutuallyExclusive,
                               colsToUse = colsToUse,
                               doAssertions = FALSE,
                               logisticPars = params,

tests/testthat/test-spreadProb-values.R

@@ -176,8 +176,9 @@ test_that("a fit on linear fuel biomass is predicted with biomass / 1e4", {
   ## x = 1.5 * MDC/200 - 0.5 * youngAge + 1 * biomass/1e4
   ## cell 1: 0 + 0 + 0 = 0                (absent fuel is exactly 0, not the log floor)
   ## cell 3: 0.75 + 0.5 = 1.25            cell 4: 1.5 + 2 = 3.5   (biomass/1e4 = 2: above 1, not clamped)
-  ## cell 6: 0.75 + 0.25 = 1              cell 7: 1.125 + 0       cell 8: 0.75 - 0.5 + 2 = 2.25
-  expect_equal(v[c(1, 3, 4, 6, 7, 8)], handLogistic3(c(0, 1.25, 3.5, 1, 1.125, 2.25)), tolerance = 1e-7)
+  ## cell 6: 0.75 + 0.25 = 1              cell 7: 1.125 + 0
+  ## cell 8: youngAge = 1 zeroes fuelA (mutually exclusive, as in the fit): 0.75 - 0.5 + 0 = 0.25
+  expect_equal(v[c(1, 3, 4, 6, 7, 8)], handLogistic3(c(0, 1.25, 3.5, 1, 1.125, 0.25)), tolerance = 1e-7)
 })
 
 test_that("a fit made on the log scale is predicted on the log scale, as before", {

tests/testthat/test-youngAgeExclusivity.R

@@ -0,0 +1,34 @@
+## youngAge must be mutually exclusive with every other non-climate covariate at prediction time,
+## exactly as in the fit (fireSense_SpreadFit::spreadFitPrep()). Before the fix,
+## spreadProbOneELF() (fireSense_SpreadPredict.R:282 pre-fix) called
+## fireSenseUtils::spreadProbFromIntegerCovs() with `mutuallyExclusive = NULL`, so a young pixel's
+## fuel biomass and non-forest land-cover columns reached the link unchanged instead of being
+## zeroed alongside youngAge = 1.
+
+test_that("a young pixel's covariates entering the link are youngAge = 1 and all others 0", {
+  covs <- toyCovariates()
+  ## rows are pixelID 7, 1, 3, 2, 8, 4, 6 (see setup-toySpread.R)
+  covs$fuelA <- fireSenseUtils::logMinB(c(0, 0, 5000, 0, 20000, 20000, 2500))  # biomass, supplied logged
+  covs$nfLCC_40 <- c(0, 0, 0, 1, 1, 0, 0)                                     # land-cover indicator
+
+  formula <- "~ MDC + youngAge + fuelA + nfLCC_40 - 1"
+  p <- toyParams(nfLCC_40 = 2)
+  ins <- toyInputs(covs = covs, params = p, formula = formula)
+  ins$covMinMax_spread$fuelA <- c(0, 1e4)      # a fit on linear fuel biomass, as fireSense_SpreadFit makes
+  ins$covMinMax_spread$nfLCC_40 <- c(0, 1)
+
+  v <- predVals(toyRun(ins))
+
+  ## cell 8 (pixelID 8): MDC 100, youngAge 1, fuelA biomass 20000, nfLCC_40 = 1.
+  ## Without the fix: x = 0.75 - 0.5 + 2 + 2 = 4.25 (fuel and land cover both leak through).
+  ## With the fix, youngAge zeroes both: x = 0.75 - 0.5 + 0 + 0 = 0.25
+  expect_equal(v[8], handLogistic3(0.25), tolerance = 1e-7)
+
+  ## cell 2 (pixelID 2): MDC 50, youngAge 1, fuelA already 0, nfLCC_40 = 1.
+  ## Without the fix: x = 0.375 - 0.5 + 0 + 2 = 1.875. With the fix: x = 0.375 - 0.5 = -0.125
+  expect_equal(v[2], handLogistic3(-0.125), tolerance = 1e-7)
+
+  ## cell 3 (pixelID 3) is NOT young: its fuel and land cover are left alone
+  ## x = 0.75 + 0 + 0.5 + 0 = 1.25
+  expect_equal(v[3], handLogistic3(1.25), tolerance = 1e-7)
+})