diff --git a/docker/depends/pecan_package_dependencies.csv b/docker/depends/pecan_package_dependencies.csv index 7ee8e35d12..de9868d889 100644 --- a/docker/depends/pecan_package_dependencies.csv +++ b/docker/depends/pecan_package_dependencies.csv @@ -252,6 +252,7 @@ "ncdf4","*","models/stics","Imports",FALSE "ncdf4","*","modules/assim.sequential","Imports",FALSE "ncdf4","*","modules/data.remote","Imports",FALSE +"ncdf4","*","modules/uncertainty","Suggests",FALSE "ncdf4",">= 1.15","base/utils","Imports",FALSE "ncdf4",">= 1.15","base/visualization","Imports",FALSE "ncdf4",">= 1.15","models/biocro","Imports",FALSE diff --git a/modules/uncertainty/DESCRIPTION b/modules/uncertainty/DESCRIPTION index 7acb7e7791..4b983d1f17 100644 --- a/modules/uncertainty/DESCRIPTION +++ b/modules/uncertainty/DESCRIPTION @@ -47,6 +47,7 @@ Imports: tidyr Suggests: mockery, + ncdf4, testthat (>= 1.0.2), withr License: BSD_3_clause + file LICENSE diff --git a/modules/uncertainty/NEWS.md b/modules/uncertainty/NEWS.md index ac0c0ee648..4a261963ac 100644 --- a/modules/uncertainty/NEWS.md +++ b/modules/uncertainty/NEWS.md @@ -10,7 +10,15 @@ sensitivity analysis. This function ensures non-parameter inputs (met, IC, soil) are held constant, which is required for valid variance decomposition in one-at-a-time (OAT) sensitivity analysis (#3729) -* `get.parameter.samples()` now consistently accepts PFTs with no `outdir` specified. These previously failed when no database connection was available. +* Fixed `read.sa.output()` dropping the median (q50) row when using the + manifest-based OAT sensitivity path. `write.sa.configs()` writes a single + shared median run with `pft_name = "NA"` and `trait = "NA"`, but the lookup + previously required an exact per-trait/per-PFT match, so the cell was always + `NA`. The fix adds a fallback: when `quantile == "50"` and no exact match is + found, the function now resolves the shared median entry before emitting a + warning. +* `get.parameter.samples()` now consistently accepts PFTs with no `outdir` specified. + These previously failed when no database connection was available. * Improved row alignment and text scaling in `plot_variance_decomposition()`. * `plot_variance_decomposition()` gains argument `order_by`, sorting highest variance to the top by default. To reproduce old behavior, use `order_by = "none"`. diff --git a/modules/uncertainty/R/sensitivity.R b/modules/uncertainty/R/sensitivity.R index 39a58105d4..93a22ec966 100644 --- a/modules/uncertainty/R/sensitivity.R +++ b/modules/uncertainty/R/sensitivity.R @@ -38,11 +38,15 @@ read.sa.output <- function(traits, quantiles, pecandir, outdir, pft.name = "", for (trait in traits) { for (quantile in quantiles) { # We look for the row that matches the current pft, trait, and quantile. + # Use which() to handle NA values in manifest columns (read.csv converts + # the string "NA" to actual R NA, and comparison with NA yields NA, not FALSE). subset_df <- manifest[ - manifest$type == "Sensitivity" & - manifest$pft_name == pft.name & - manifest$trait == trait & - as.character(manifest$quantile) == as.character(quantile), + which( + manifest$type == "Sensitivity" & + manifest$pft_name == pft.name & + manifest$trait == trait & + as.character(manifest$quantile) == as.character(quantile) + ), ] if (nrow(subset_df) == 1) { @@ -50,6 +54,26 @@ read.sa.output <- function(traits, quantiles, pecandir, outdir, pft.name = "", } else if (nrow(subset_df) > 1) { PEcAn.logger::logger.warn("Multiple runs found for", trait, quantile, "- using the last one.") run.id <- utils::tail(subset_df$run_id, 1) + } else if (as.character(quantile) == "50") { + # The median run is written as a single shared entry with pft_name = "NA" + # and trait = "NA". Fall back to that shared row when no exact match exists. + # Note: read.csv() converts the string "NA" to actual R NA values, + # so we must use is.na() rather than == "NA" for the comparison. + median_df <- manifest[ + manifest$type == "Sensitivity" & + is.na(manifest$pft_name) & + is.na(manifest$trait) & + as.character(manifest$quantile) == "50", + ] + if (nrow(median_df) >= 1) { + PEcAn.logger::logger.debug( + "Using shared median run for", trait, "quantile 50" + ) + run.id <- utils::tail(median_df$run_id, 1) + } else { + PEcAn.logger::logger.warn("No run found in manifest for", trait, quantile) + next # Skip this quantile + } } else { PEcAn.logger::logger.warn("No run found in manifest for", trait, quantile) next # Skip this quantile diff --git a/modules/uncertainty/tests/testthat/test-read_sa_output_median.R b/modules/uncertainty/tests/testthat/test-read_sa_output_median.R new file mode 100644 index 0000000000..a7eb8049aa --- /dev/null +++ b/modules/uncertainty/tests/testthat/test-read_sa_output_median.R @@ -0,0 +1,99 @@ +# Regression test for the shared-median manifest fallback in read.sa.output() +# +# Background: write.sa.configs() writes the median (q50) run as a single shared +# entry with pft_name = "NA" and trait = "NA". read.sa.output() previously did +# a strict per-trait/per-PFT lookup, so the median row always failed to match +# and the 50th-quantile cell was returned as NA. +# +# The fix adds a fallback: when quantile == "50" and no exact match is found, +# look for the shared median entry before giving up. +# +# Reference: https://github.com/pecanproject/pecan/issues/3882 + +test_that("read.sa.output resolves median (q50) via shared manifest fallback", { + withr::with_tempdir({ + pft <- "temperate.coniferous" + trait <- "growth_resp_factor" + yr <- 2004 + + # ---- Build a manifest that matches what write.sa.configs() actually writes ---- + # One shared median row (pft_name = "NA", trait = "NA") + no trait-specific row + # for quantile 50. This is the exact situation that triggered the original bug. + median_run_id <- "SA-median--1" + manifest <- data.frame( + type = "Sensitivity", + pft_name = "NA", + trait = "NA", + quantile = "50", + run_id = median_run_id, + stringsAsFactors = FALSE + ) + write.csv(manifest, "runs_manifest.csv", row.names = FALSE) + + # ---- Stub out the model output directory so read.output() finds a file ---- + run_outdir <- file.path(getwd(), median_run_id) + dir.create(run_outdir, recursive = TRUE) + + # Create a minimal NetCDF file with an NPP variable for the target year + nc_path <- file.path(run_outdir, paste0(yr, ".nc")) + nc_obj <- ncdf4::nc_create( + nc_path, + list(ncdf4::ncvar_def("NPP", "kg m-2 s-1", list(), missval = NA_real_)) + ) + ncdf4::ncvar_put(nc_obj, "NPP", 1.23) + ncdf4::nc_close(nc_obj) + + # ---- Run read.sa.output() ---- + # Should NOT warn "Run ID invalid or missing" for quantile 50 any more. + expect_no_warning( + out <- PEcAn.uncertainty::read.sa.output( + traits = trait, + quantiles = "50", + pecandir = getwd(), + outdir = getwd(), + pft.name = pft, + start.year = yr, + end.year = yr, + variable = PEcAn.utils::convert.expr("NPP")$variable.eqn + ), + regexp = "Run ID invalid or missing" + ) + + # The median cell must not be NA + expect_false( + is.na(out[["50", trait]]), + label = "median (q50) output should not be NA" + ) + }) +}) + + +test_that("read.sa.output still warns when no median fallback row exists", { + withr::with_tempdir({ + # Manifest with quantile 50 row that has the wrong pft_name (not "NA"), + # meaning neither exact match nor shared-median fallback can be found. + manifest <- data.frame( + type = "Sensitivity", + pft_name = "some.other.pft", + trait = "NA", + quantile = "50", + run_id = "SA-median--1", + stringsAsFactors = FALSE + ) + write.csv(manifest, "runs_manifest.csv", row.names = FALSE) + + expect_warning( + PEcAn.uncertainty::read.sa.output( + traits = "growth_resp_factor", + quantiles = "50", + pecandir = getwd(), + outdir = getwd(), + pft.name = "temperate.coniferous", + start.year = 2004, + end.year = 2004, + variable = PEcAn.utils::convert.expr("NPP")$variable.eqn + ), + regexp = "No run found in manifest" + ) + }) +})