# Base R diagnostic teaching benchmark. No external packages or network calls.
# Rscript residual-diagnostics.R mmm-practice-synthetic.csv [output.pdf]
# Data is synthetic/noiseless: near-zero error is expected, NOT real MMM performance.
args <- commandArgs(trailingOnly = TRUE)
if (length(args) < 1) stop('Supply the path to mmm-practice-synthetic.csv')
out <- if (length(args) > 1) args[2] else 'residual-diagnostics.pdf'
d <- read.csv(args[1])
required <- c('week', 'revenue', 'search_spend', 'social_spend', 'video_spend', 'promotion')
stopifnot(all(required %in% names(d)), nrow(d) >= 20)
d$week <- as.Date(d$week, format = '%Y-%m-%d')
stopifnot(!anyNA(d[, required]), !anyDuplicated(d$week), all(diff(d$week) == 7))
stopifnot(all(vapply(d[setdiff(required, 'week')], function(x) is.numeric(x) && all(is.finite(x)), logical(1))))
d$t <- seq_len(nrow(d)) - 1
# Known synthetic seasonal term. Do not assume this specification fits real data.
d$season_sin <- sin(2*pi*d$t/52)
n_train <- floor(nrow(d)*0.8)
train <- seq_len(n_train)
test <- seq.int(n_train+1, nrow(d))
fit <- lm(revenue ~ search_spend + social_spend + video_spend + promotion + season_sin,
          data = d[train, ])
pred <- predict(fit, newdata = d[test, ])
stopifnot(all(is.finite(pred)))
res <- residuals(fit)
rmse <- function(x) sqrt(mean(x^2))
cat('Training weeks:', length(train), 'Holdout weeks:', length(test), '\n')
cat('Train RMSE:', rmse(res), 'Holdout RMSE:', rmse(d$revenue[test]-pred), '\n')
cat('Synthetic noiseless benchmark, not evidence of causal identification.\n')
pdf(out, width=10, height=8)
par(mfrow=c(2,2), mar=c(4,4,3,1))
plot(d$week[train], res, type='b', xlab='Week', ylab='Residual', main='Training residuals over time')
abline(h=0, col='grey')
plot(fitted(fit), res, xlab='Fitted revenue', ylab='Residual', main='Residual scale')
abline(h=0, col='grey')
qqnorm(res, main='Q-Q: near-zero synthetic residuals'); qqline(res)
plot(d$week[test], d$revenue[test], type='l', xlab='Week', ylab='Revenue', main='Future-period prediction')
lines(d$week[test], pred, col='blue', lty=2)
legend('topleft', c('Observed', 'Predicted'), col=c('black','blue'), lty=c(1,2), bty='n')
invisible(dev.off())
cat('Wrote', out, '\n')