Variance-Guided Regression Improving Upon OLS and ANOVA

Fits variance-guided linear regression models that provide an alternative to ordinary least squares (OLS) for general linear-model design matrices, including ANOVA-style encodings. The methods use an iteratively reweighted least squares estimator or an iteratively reweighted lasso estimator and implement the global linear mean-variance model from the associated 2026 Statistics in Medicine article . Under the assumptions in that paper, the estimator matches the homoscedastic baseline in population predictive quasi-risk when variance is constant and improves on it when the variance depends on covariates. The grouping-based nonlinear prediction extension from Section 3 is available in the development version on GitHub.

varGuid: Variance-Guided Regression Improving Upon OLS and ANOVA

Authors

Sibei Liu ([email protected]) and Min Lu ([email protected])

Reference

Liu, S. and Lu, M. (2026). Variance-Guided Regression for Heteroscedastic Data with a Grouping-Based Extension for Nonlinear Prediction. Statistics in Medicine. 45(13-14):e70632. DOI: 10.1002/sim.70632

Description

The CRAN release of varGuid implements the global linear mean-variance model from Section 2 of the paper via iteratively reweighted least squares (IRLS) and iteratively reweighted lasso estimation. It supports general linear-model design matrices, including ANOVA-style encodings. For the grouping-based nonlinear prediction extension from Section 3 of the paper, please use the development version available at:

devtools::install_github("luminwin/varGuid")

Installation

install.packages("varGuid")
library(varGuid)

Examples

data(cobra2d, package = "varGuid")
dat <- cobra2d
set.seed(1)
tid <- sample(seq_len(nrow(dat)), 200)
train <- dat[-tid, ]
test <- dat[tid, ]
yid <- which(colnames(dat) == "y")

IRLS / weighted least squares fit

o <- lmv(X = train[, -yid], Y = train[, yid], lasso = FALSE)
summary(o$obj.varGuid)  # varGuid weighted fit
summary(o$obj.OLS)      # baseline OLS fit

head(prd(o, train[, -yid], model = "baseline"))
head(prd(o, train[, -yid], model = "varGuid"))

# The corresponding base R calls are:
head(stats::predict(o$obj.OLS, newdata = as.data.frame(train[, -yid])))
head(stats::predict(o$obj.varGuid, newdata = as.data.frame(train[, -yid])))

Iteratively reweighted lasso fit

o2 <- lmv(X = train[, -yid], Y = train[, yid], lasso = TRUE)
o2$beta            # varGuid-lasso coefficients
o2$obj.lasso$beta  # baseline lasso coefficients

head(prd(o2, train[, -yid], model = "baseline"))
head(prd(o2, train[, -yid], model = "varGuid"))

# The corresponding glmnet calls are:
head(drop(glmnet::predict.glmnet(o2$obj.lasso,
                                 newx = as.matrix(train[, -yid]))))
head(drop(glmnet::predict.glmnet(o2$obj.varGuid,
                                 newx = as.matrix(train[, -yid]))))