library(agriPriceRisk)
checks <- 0L
assert <- function(x) {stopifnot(isTRUE(x)); checks <<- checks + 1L}
errors <- function(expr) inherits(tryCatch({force(expr); NULL}, error = function(e) e), "error")
d <- simulate_commodity_data()
p <- prepare_price_data(d)
assert(nrow(p) == 120L)
assert(all(diff(as.integer(p$date)) > 0))
assert(errors(prepare_price_data(rbind(d, d[1, ]))))
assert(errors(prepare_price_data(d[-20, ])))
bad <- d; bad$wholesale_price[1] <- 0
assert(errors(prepare_price_data(bad)))
bad <- d; bad$wholesale_price[1] <- NA
assert(errors(prepare_price_data(bad)))
assert(is.na(prepare_price_data(bad, missing = "keep")$price[1]))
s <- classify_price_state(p$price, p$month)
assert(all(is.na(s$state[1:60])))
assert(all(s$status[1:60] == "insufficient_history"))
assert(all(!is.na(s$state[61:120])))
# Perturbing future data cannot change earlier classifications.
changed <- p$price; changed[100:120] <- changed[100:120] * 4
s2 <- classify_price_state(changed, p$month)
assert(identical(s[1:99, ], s2[1:99, ]))
short <- classify_price_state(p$price[1:96], p$month[1:96])
assert(identical(s[1:96, ], short))
base <- data.frame(month = 1:12, center = rep(100, 12), scale = rep(10, 12), n = rep(10L, 12))
bands <- classify_price_state(c(109.9, 110, 115, 125, 85, 75), rep(1, 6), baseline = base)
assert(identical(bands$state, c("Normal", "Warning", "High Risk", "Extreme Shock", "High Risk", "Extreme Shock")))
assert(identical(bands$direction, c(rep("Upward", 4), rep("Downward", 2))))
base$scale <- 0
assert(all(is.na(classify_price_state(100, 1, baseline = base)$state)))
assert(errors(classify_price_state(100, 1, thresholds = c(1, 1, 2))))
assert(errors(classify_price_state(100, 13)))
assert(errors(seasonal_baseline(p$price, p$month, min_history = 1)))
robust <- seasonal_baseline(p$price, p$month, method = "median_mad", window = 5)
assert(all(robust$n == 5L))
fc <- forecast_price_risk(p, model = "seasonal_naive", h = 15)
assert(nrow(fc) == 15L)
assert(fc$date[1] == as.Date("2020-01-01"))
assert(fc$date[15] == as.Date("2021-03-01"))
assert(all(fc$lower <= fc$forecast & fc$forecast <= fc$upper))
assert(all(fc$p_warning >= fc$p_high_risk & fc$p_high_risk >= fc$p_extreme))
assert(all(fc$p_warning >= 0 & fc$p_warning <= 1))
assert(all(abs(fc$p_up_high + fc$p_down_high - fc$p_high_risk) < 1e-12))
assert(all(abs(fc$forecast[1:12] - tail(p$price, 12)) < 1e-12))
assert(abs(fc$forecast_se[13] / fc$forecast_se[1] - sqrt(2)) < 1e-12)
assert(errors(forecast_price_risk(p[1:24, ])))
assert(errors(forecast_price_risk(p, h = 0)))
assert(errors(forecast_price_risk(p, level = 1)))
# Optimized or documented fixed-coefficient Holt-Winters always reports fitting mode.
hw <- forecast_price_risk(p)
assert(nzchar(attr(hw, "fitting_note")))
assert(all(is.finite(hw$forecast)))
# Constant series does not manufacture an alert from zero seasonal scale.
constant <- d; constant$wholesale_price <- 100
cf <- forecast_price_risk(prepare_price_data(constant), model = "seasonal_naive")
assert(all(is.na(cf$state)))
bt <- backtest_ews(p, initial = 108, h = 3, model = "seasonal_naive")
assert(identical(bt$metrics$n_forecasts, c(12L, 11L, 10L)))
assert(all(bt$predictions$origin < bt$predictions$date))
q <- bt$predictions[bt$predictions$horizon == 1, ]
assert(abs(bt$metrics$MAE[1] - mean(abs(q$actual - q$forecast))) < 1e-10)
assert(abs(bt$metrics$Brier[1] - mean((q$p_high_risk - q$actual_event)^2)) < 1e-12)
# Later perturbations cannot alter a backtest forecast from an earlier origin.
p2 <- p; p2$price[117:120] <- p2$price[117:120] * 2
bt2 <- backtest_ews(p2, initial = 108, h = 1, model = "seasonal_naive")
assert(identical(bt$predictions$forecast[bt$predictions$horizon == 1][1:9], bt2$predictions$forecast[1:9]))
cal <- calibrate_alerts(bt)
assert(cal$best_cutoff %in% cal$comparison$cutoff)
assert(cal$comparison$mean_cost[cal$comparison$cutoff == cal$best_cutoff] == min(cal$comparison$mean_cost))
assert(errors(calibrate_alerts(bt, cutoffs = 0)))
assert(errors(calibrate_alerts(bt, false_alarm_cost = 0, missed_shock_cost = 0)))
ch <- estimate_transmission_chain(d, farm_price_col = "farm_gate_price")
assert(ch$coefficients$estimate[1] > 0)
assert(ch$coefficients$estimate[3] < 0)
assert(ch$rainfall_lag == 1L)
assert(isFALSE(ch$income_price_proxy))
assert(nrow(ch$coefficients) == 6L)
assert(all(ch$coefficients$lower <= ch$coefficients$estimate & ch$coefficients$upper >= ch$coefficients$estimate))
# Negative rainfall is supported, but non-positive logged outcomes are rejected.
neg <- d; neg$production[2] <- 0
assert(errors(estimate_transmission_chain(neg)))
assert(errors(estimate_transmission_chain(d, rainfall_lag = 120)))
sc <- simulate_shock_propagation(ch, nsim = 200)
assert(nrow(sc) == 5L)
assert(sc$percent_change[1] < 0 && sc$percent_change[3] > 0)
zero <- simulate_shock_propagation(ch, rainfall_change = 0, nsim = 0)
assert(all(zero$percent_change == 0))
channels <- attr(sc, "income_channels")
assert(abs(channels$log_change[3] - sum(channels$log_change[1:2])) < 1e-12)
assert(abs(sc$log_change[sc$stage == "income"] - channels$log_change[3]) < 1e-12)
assert(identical(sc, simulate_shock_propagation(ch, nsim = 200)))
set.seed(989); prior <- .Random.seed
invisible(simulate_commodity_data(seed = 123)); assert(identical(prior, .Random.seed))
invisible(simulate_shock_propagation(ch, seed = 123, nsim = 20)); assert(identical(prior, .Random.seed))
base <- seasonal_baseline(p$price, p$month)
sc2 <- simulate_shock_propagation(ch, baseline_price = tail(p$price, 1), baseline = base, month = 1, nsim = 0)
assert(is.finite(attr(sc2, "projected_wholesale_price")))
assert(errors(simulate_shock_propagation(ch, baseline_price = 100)))
assert(grepl("Conditional", interpret_ews_state(s, hw, sc)))
# Groups are isolated: a second commodity cannot affect the first.
d2 <- d; d2$commodity <- "Second commodity"; d2$wholesale_price <- d2$wholesale_price * 2
panel <- run_price_ews(rbind(d, d2), commodity_col = "commodity", market_col = "market", model = "seasonal_naive")
assert(length(panel) == 2L)
assert(abs(panel[[2]]$forecast$forecast[1] / panel[[1]]$forecast$forecast[1] - 2) < 1e-12)
# Exercise every plotting method without touching the user's directories.
f <- tempfile(fileext = ".pdf")
grDevices::pdf(f)
plot(s, dates = p$date); plot(hw); plot(bt); plot(ch); plot(sc)
grDevices::dev.off(); assert(file.info(f)$size > 0); unlink(f)
cat("All", checks, "invariant and numerical checks passed.\n")