bounded.reg {quadrupen} | R Documentation |
Fit a linear model with infinity-norm plus ridge-like regularization
Description
Adjust a linear model penalized by a mixture of a (possibly
weighted) \ell_\infty
-norm (bounding the
magnitude of the parameters) and a (possibly structured)
\ell_2
-norm (ridge-like). The solution path is computed
at a grid of values for the infinity-penalty, fixing the amount of
\ell_2
regularization. See details for the criterion
optimized.
Usage
bounded.reg(
x,
y,
lambda1 = NULL,
lambda2 = 0.01,
penscale = rep(1, p),
struct = NULL,
intercept = TRUE,
normalize = TRUE,
naive = FALSE,
nlambda1 = ifelse(is.null(lambda1), 100, length(lambda1)),
min.ratio = ifelse(n <= p, 0.01, 0.001),
max.feat = ifelse(lambda2 < 0.01, min(n, p), min(4 * n, p)),
control = list(),
checkargs = TRUE
)
Arguments
x |
matrix of features, possibly sparsely encoded
(experimental). Do NOT include intercept. When normalized os
|
y |
response vector. |
lambda1 |
sequence of decreasing |
lambda2 |
real scalar; tunes the |
penscale |
vector with real positive values that weight the infinity norm of each feature. Default set all weights to 1. See details below. |
struct |
matrix structuring the coefficients. Must be at
least positive semidefinite (this is checked internally if the
|
intercept |
logical; indicates if an intercept should be
included in the model. Default is |
normalize |
logical; indicates if variables should be
normalized to have unit L2 norm before fitting. Default is
|
naive |
logical; Compute either 'naive' of 'classic' bounded
regression: mimicking the Elastic-net, the vector of parameters is
rescaled by a coefficient |
nlambda1 |
integer that indicates the number of values to put
in the |
min.ratio |
minimal value of infinity-part of the penalty
that will be tried, as a fraction of the maximal |
max.feat |
integer; limits the number of features ever to
enter the model: in our implementation of the bounded regression,
it corresponds to the variables which have left the boundary along
the path. The algorithm stops if this number is exceeded and
|
control |
list of argument controlling low level options of the algorithm –use with care and at your own risk– :
|
checkargs |
logical; should arguments be checked to
(hopefully) avoid internal crashes? Default is
|
Details
The optimized criterion is
βhatλ1,λ2 = argminβ 1/2 RSS(β) + λ1 | D β |∞ + λ/2 2 βT S β,
where
D
is a diagonal matrix, whose diagonal terms are provided
as a vector by the penscale
argument. The \ell_2
structuring matrix S
is provided via the struct
argument, a positive semidefinite matrix (possibly of class
Matrix
).
Note that the quadratic algorithm for the bounded regression may
become unstable along the path because of singularity of the
underlying problem, e.g. when there are too much correlation or
when the size of the problem is close to or smaller than the
sample size. In such cases, it might be a good idea to switch to
the proximal solver, slower yet more robust. This is the strategy
adopted by the 'bulletproof'
mode, that will send a warning
while switching the method to 'fista'
and keep on
optimizing on the remainder of the path. When bulletproof
is set to FALSE
, the algorithm stops at an early stage of
the path of solutions. Hence, users should be careful when
manipulating the resulting 'quadrupen'
object, as it will
not have the size expected regarding the dimension of the
lambda1
argument.
Singularity of the system can also be avoided with a larger
\ell_2
-regularization, via lambda2
, or a
"not-too-small" \ell_\infty
regularization, via
a larger 'min.ratio'
argument.
Value
an object with class quadrupen
, see the
documentation page quadrupen
for details.
See Also
See also quadrupen
,
plot,quadrupen-method
and crossval
.
Examples
## Simulating multivariate Gaussian with blockwise correlation
## and piecewise constant vector of parameters
beta <- rep(c(0,1,0,-1,0), c(25,10,25,10,25))
cor <- 0.75
Soo <- toeplitz(cor^(0:(25-1))) ## Toeplitz correlation for irrelevant variables
Sww <- matrix(cor,10,10) ## bloc correlation between active variables
Sigma <- bdiag(Soo,Sww,Soo,Sww,Soo)
diag(Sigma) <- 1
n <- 50
x <- as.matrix(matrix(rnorm(95*n),n,95) %*% chol(Sigma))
y <- 10 + x %*% beta + rnorm(n,0,10)
## Infinity norm without/with an additional l2 regularization term
## and with structuring prior
labels <- rep("irrelevant", length(beta))
labels[beta != 0] <- "relevant"
plot(bounded.reg(x,y,lambda2=0) , label=labels) ## a mess
plot(bounded.reg(x,y,lambda2=10), label=labels) ## good guys are at the boundaries