Fits a ridge regularised weighted least squares regression using a chosen solver subject to optional nonnegativity or box constraints on the coefficients.

Fit Constrained Ridge Regularised Weighted Least Squares Regression using chosen Solver

Usage

lmridge(
  X,
  y,
  lambdas = rep(0, ncol(X)),
  weights = NULL,
  nonnegative = FALSE,
  lower = if (nonnegative) {
     rep(0, ncol(X))
 } else {
    
    rep(-.Machine$double.xmax, ncol(X))
 },
  upper = rep(.Machine$double.xmax, ncol(X)),
  scaling = rep(1, ncol(X)),
  prior.mean = NULL,
  solver = "eigen",
  method = 1L,
  start = rep(0, ncol(X)),
  iteration = 1,
  control.solver = do.call("control.solver.gen", list(solver = solver, method = method)),
  verbose = FALSE
)

Arguments

X: A numeric matrix or sparse matrix of predictors (design matrix).
y: A numeric vector of responses.
lambdas: Numeric value or vector specifying the ridge penalty parameter(s).
weights: Optional numeric vector of observation weights. Default is NULL.
nonnegative: Logical; if TRUE, nonnegativity constraints are imposed on the coefficients (equivalent to setting lower = 0). Default is FALSE.
lower: Numeric vector specifying lower bounds on the coefficients. Default is a vector of -Inf.
upper: Numeric vector specifying upper bounds on the coefficients. Default is a vector of Inf.
scaling: Numeric vector specifying scaling factors for the columns of X. Default is a vector of ones (no scaling). This allows scaling of columns of X, e.g., by their L2 norm.
prior.mean: Numeric vector specifying prior means for Gaussian priors implied by the ridge penalties. Default is a vector of zeros (i.e. zero centered means) to shrink to zero. Nonzero values will result in a targeted nonzero-centered ridge regression being used instead.
solver: Character string specifying which solver to use. Options are "eigen", "bigGlm" (alias "glmnet"), "osqp". Default is "eigen".
method: Integer code specifying the method to be used within the solver. Default is 1L, i.e. Cholesky LLT.
start: Numeric vector of starting values for the coefficients. Default is a vector of zeros.
iteration: Integer specifying the iteration number, used internally for iterative algorithms. Default is 1.
control.solver: List of control parameters for the solver. Default is generated by control.solver.gen(solver = solver, method = method).
verbose: Logical; if TRUE, additional information is printed during computation. Default is FALSE.

Value

A numeric vector of estimated coefficients.

Details

The lmridge function fits a ridge regression model of the form: $$\beta = \arg \min_{\beta} \left\{ \frac{1}{2} \| W^{1/2}(y - X \beta) \|_2^2 + \frac{1}{2}\sum_j \lambda_j (\beta_j - \mu_j)^2 \right\}$$ where $W$ contains the observation weights, $\lambda_j$ are the ridge penalty parameters, and $\mu$ is the prior mean of the coefficients (specified by prior.mean).

The function supports optional bounds on the coefficients via the lower and upper parameters and can impose nonnegativity constraints if nonnegative = TRUE.

The scaling parameter allows scaling of the predictors. If scaling is applied, the coefficients are rescaled back to the original scale before being returned.

Different solvers can be specified via the solver parameter:

"eigen": Uses the Eigen package to solve the ridge regression problem using one of several solvers (see method option).
"bigGlm": Uses the bigGlm function in the glmnet package to solve the ridge regression problem. "glmnet" is an alias.
"osqp": Uses the osqp package (Operator Splitting Quadratic Program solver) to solve the problem with constraints. Most appropriate when box or nonnegativity constraints are applied.

The method parameter specifies the specific method to be used when the "eigen" solver is used (1=Cholesky LLT, 2=Cholesky LDLT, 3=CG with incomplete Cholesky preconditioning for sparse systems or diagonal preconditioning for dense systems, 4=LSCG with diagonal least-squares preconditioner, 5=column-pivoted QR, 6=BDCSVD for dense systems or Jacobi SVD for sparse systems).

LSCG uses Eigen's diagonal least-squares preconditioner without additional external column equilibration. A finite approximate LSCG solution is retained with a warning if the requested tolerance is not reached; a missing or nonfinite solution triggers the Cholesky fallback. OSQP uses the original coefficient coordinates and reports a warning when its requested accuracy is not reached. On a failed OSQP setup or solve, an infeasible/nonconvex status, or an invalid returned coefficient vector, a warning is issued and the previous IRLS coefficients (clipped to the bounds) are retained. Such a warning does not establish convergence. Finite-tolerance solutions can change the subsequent adaptive ridge selection path.

Eligible wide Eigen ridge systems use the Woodbury identity to solve an observation-space system of dimension n, rather than a predictor-space system of dimension p. The requested Eigen method is retained. For LLT, LDLT and CG, the default forms the Gram matrix directly for dense predictors. Sparse predictors, or dense-stored working matrices with at least one million entries and fewer than 10 percent nonzeros, use sparse row augmentation. Set control.solver$rowaug explicitly to override this automatic choice. Dense and sparse designs retain their storage type; optional row augmentation constructs a sparse augmented matrix. Eigen's SVD requires a dense factorization internally even when its input is sparse. OSQP requires sparse arrays at its solver boundary.

Author

Tom Wenseleers tom.wenseleers@kuleuven.be \ Laboratory of Socioecology and Social Evolution, KU Leuven \ https://bio.kuleuven.be/eeb/tw/research \ ORCID: 0000-0002-1434-861X

Examples

if (FALSE) { # \dontrun{
# Simulate some data
set.seed(123)
n <- 100
p <- 10
X <- matrix(rnorm(n * p), n, p)
beta_true <- rnorm(p)
y <- X %*% beta_true + rnorm(n)

# Fit ridge regression model with default settings
beta_est <- lmridge(X, y, lambdas = 1)

# Fit ridge regression with nonnegativity constraints using the osqp solver
beta_est_nn <- lmridge(X, y, lambdas = 1, nonnegative = TRUE, solver = "osqp")

# Compare estimated coefficients to true coefficients
plot(beta_true, beta_est, xlab = "True Coefficients", ylab = "Estimated Coefficients",
     main = "Ridge Regression Coefficients")
abline(0, 1, col = "red")
} # }