enetLTS-package: Robust and Sparse Methods for High Dimensional Linear and...

In enetLTS: Robust and Sparse Methods for High Dimensional Linear and Binary and Multinomial Regression

Robust and Sparse Methods for High Dimensional Linear and Binary and Multinomial Regression

Description

Compute fully robust versions of the elastic net estimator, which allows for sparse model estimates, for linear regression and binary and multinomial logistic regression.

Usage

enetLTS(
    xx,
    yy,
    family=c("gaussian","binomial","multinomial"),
    alphas=seq(0,1,length=41),
    lambdas=NULL,
    lambdaw=NULL,
    intercept=TRUE,
    scal=TRUE,
    hsize=0.75,
    nsamp=c(500,10),
    nCsteps=20,
    nfold=5,
    repl=1,
    ncores=1,
    tol=-1e6,
    seed=NULL,
    del=0.0125,
    crit.plot=FALSE,
    typegrouped=FALSE,
    type.response=c("link","response","class")
  )
  

Arguments

xx	a numeric matrix containing the predictor variables.
yy	response variable. Quantitative for family="gaussian". For family="binomial" should be a factor with two levels which is coded as 0 and 1. For family="multinomial" should be a factor with the number of categories (NC) which is coded as 1,2,...,NC.
family	a description of the error distribution and link function to be used in the model. "gaussian", "binomial" and family="multinomial" options are available.
alphas	a user supplied alpha sequence for the elastic net penalty, which is the mixing proportion of the ridge and lasso penalties and takes value in [0,1]. α=1 is the lasso penalty, and α=0 the ridge penalty. If not provided a sequence, default is 41 equally spaced values.
lambdas	a user supplied lambda sequence for the strength of the elastic net penalty. If not provided a sequence, default is chosen with steps of size -0.025 lambda0 with 0≤λ≤ lambda0 for linear regression and -0.025 lambda00 with 0≤λ≤ lambda00 for binary logistic regression. lambda0 is determined based on the Pearson correlation between y and the jth predictor variable x_j on winsorized data for linear regression. In lambda00 for logistic regression, the Pearson correlation is replaced by a robustified point-biserial correlation. Default is chosen with steps of size -0.05 from 0.95 to 0.05 for multinomial logistic regression.
lambdaw	a user supplied lambda sequence for reweighting step. If not provided, default is computed by using k-fold cross-validation via cv.glmnet function.
intercept	a logical indicating whether a constant term should be included in the model (the default is TRUE).
scal	a logical value indicating whether scale the predictors by their arithmetic means and standard deviations. For family="gaussian", it also indicates if mean-center the response variable or not. The default is TRUE. Note that scaling is performed on the subsamples rather than the full data set.
hsize	a user supplied numeric value giving the percentage of the residuals for which the elastic net penalized sum of squares for linear regression or for which the elastic net penalized sum of deviances for binary and multinomial logistic regression should be minimized. The default is 0.75.
nsamp	a numeric vector giving the number of subsamples to be used in the beginning of the algorithm, which gives the number of initial subsamples to be used. The default is to first perform C-steps on 500 initial subsamples, and then to keep the s1 subsamples with the lowest value (or highest value based on which model is used - family="gaussian" or family="binomial" or family="multinomial") of the objective function for additional C-steps until convergence.
nCsteps	a positive integer giving the number of C-steps to perform on determined s1 subsamples. The default is 20.
nfold	a user supplied numeric value for fold number of k-fold cross-validation which used in varied functions of the algorithm. The default is 5-fold cross-validation.
repl	a user supplied positive number for more stable results, repeat the k-fold CV repl times and take the average of the corresponding evaluation measure. The default is 1.
ncores	a positive integer giving the number of processor cores to be used for parallel computing (the default is 1 for no parallelization). If this is set to NA, all available processor cores are used. For prediction error estimation, parallel computing is implemented on the R level using package parallel.
tol	a small numeric value for convergence. The default is -1e6.
seed	optional initial seed for the random number generator (see.Random.seed) when determine initial subsets at thebeginning of the algorithm. The default is NULL.
del	The default is 0.0125.
crit.plot	a logical value indicating if produces a plot for k-fold cross-validation based on alpha and lambda combinations. The default is TRUE.
typegrouped	This argument is available for only family="multinomial" in the last fit based on the best subset. TRUE means "grouped" and FALSE means "ungrouped". If "TRUE" then a grouped lasso penalty is used on the multinomial coefficients for a variable. The default is FALSE.
type.response	type of prediction required. type="link" gives the linear predictors. type="response" gives the fitted probabilities for family="multinomial" and family="binomial" and gives the fitted values for family="gaussian". type="class" is available only for family="binomial" and family="multinomial", and produces the class label corresponding to the maximum probability.

Details

The idea of repeatedly applying the non-robust classical elastic net estimators to data subsets only is used for linear and logistic regression. The algorithm starts with 500 elemental subsets only for one combination of α and λ, and takes the warm start strategy for subsequent combinations. This idea saves the computation time. To choose the elastic net penalties, k-fold cross-validation is used and the replication option is provided for more stable results. Robustness has been achieved by using trimming idea, therefore a reweighting step is introduced in order to improve the efficiency. The outliers are identified according to current model. For family="gaussian", standardized residuals are used. For family="binomial", the Pearson residuals which are approximately standard normally distributed is used. Then the weights are defined by the binary weight function using del=0.0125, which allows to be flagged as outliers of the 2.5% of the observations in the normal model. For family="multinomial", group-wise scaled robust distances are used. The the binary weights defined using the constant $c_2=5$. Therefore, binary weight function produces a clear distinction between the "good observations" and "outliers".

Value

objective	a numeric vector giving the respective values of the enetLTS objective function, i.e., the elastic net penalized sums of the h smallest squared residuals from the raw fits for family="gaussian" and the elastic net penalized sums of the h deviances from the raw fits for family="binomial".
raw.rmse	root mean squared error for raw fit, which is available for only family="gaussian".
rmse	root mean squared error for reweighted fit, which is available for only family ="gaussian".
raw.mae	mean absolute error for raw fit.
mae	mean absolute error for reweighted fit.
best	an integer vector containing the respective best subsets of h observations found and used for computing the raw estimates.
raw.wt	an integer vector containing binary weights that indicate outliers from the respective raw fits, i.e., the weights used for the reweighted fits.
wt	an integer vector containing binary weights that indicate outliers from the respective reweighted fits, i.e., the weights are 1 for observations with reasonably small reweighted residuals and 0 for observations with large reweighted residuals.
raw.coefficients	a numeric vector containing the respective coefficient estimates from the raw fit.
coefficients	a numeric vector containing the respective coefficient estimates from the reweighted fit.
raw.fitted.values	a numeric vector containing the respective fitted values of the response from the raw fits.
fitted.values	a numeric vector containing the respective fitted values of the response from the reweighted fits.
raw.residuals	a numeric vector containing the respective residuals for family="gaussian" and respective deviances for family="binomial" and family="multinomial" from the raw fits.
residuals	a numeric vector containing the respective residuals for family="gaussian" and respective deviances for family="binomial" and family="multinomial" from the reweighted fits.
alpha	an optimal elastic net mixing parameter value obtained with k-fold cross-validation.
lambda	an optimal value for the strength of the elastic net penalty obtained with k-fold cross-validation.
lambdaw	an optimal value for the strength of the elastic net penalty re-obtained with k-fold cross-validation for reweighted fit.
num.nonzerocoef	the number of the nonzero coefficients in the model.
n	the number of observations.
p	the number of variables.
h	the number of observations used to compute the raw estimates.
classnames	class names for logistic model, which is available for only family="binomial" and family="multinomial".
classize	class sizes for logisitic model, which is available for only family="binomial" and family="multinomial".
inputs	all inputs used in the function enetLTS.R.
call	the matched function call.

Author(s)

Fatma Sevinc KURNAZ, Irene HOFFMANN, Peter FILZMOSER

References

Kurnaz, F.S., Hoffmann, I. and Filzmoser, P. (2017) Robust and sparse estimation methods for high dimensional linear and logistic regression. Chemometrics and Intelligent Laboratory Systems.

See Also

print, predict, coef, nonzeroCoef.enetLTS, plot, plotCoef.enetLTS, plotResid.enetLTS, plotDiagnostic.enetLTS, residuals, fitted, weights

Examples

## for gaussian

set.seed(86)
n <- 100; p <- 25                             # number of observations and variables
beta <- rep(0,p); beta[1:6] <- 1              # 10% nonzero coefficients
sigma <- 0.5                                  # controls signal-to-noise ratio
x <- matrix(rnorm(n*p, sigma),nrow=n)
e <- rnorm(n,0,1)                             # error terms
eps <- 0.1                                    # contamination level
m <- ceiling(eps*n)                           # observations to be contaminated
eout <- e; eout[1:m] <- eout[1:m] + 10        # vertical outliers
yout <- c(x %*% beta + sigma * eout)        # response
xout <- x; xout[1:m,] <- xout[1:m,] + 10      # bad leverage points


# determine user supplied alpha and lambda sequences
# alphas=seq(0,1,length=11)
# l0 <- robustHD::lambda0(xout,yout)          # use lambda0 function from robustHD package
# lambdas <- seq(l0,0,by=-0.1*l0)
# fit <- enetLTS(xout,yout,alphas=alphas,lambdas=lambdas)


## for binomial

eps <-0.05                                     # %10 contamination to only class 0
m <- ceiling(eps*n)
y <- sample(0:1,n,replace=TRUE)
xout <- x
xout[y==0,][1:m,] <- xout[1:m,] + 10;          # class 0
yout <- y                                      # wrong classification for vertical outliers


# determine user supplied alpha and lambda sequences
# alphas=seq(0,1,length=11)
# l00 <- lambda00(xout,yout,normalize=TRUE,intercept=TRUE)
# lambdas <-  seq(l00,0,by=-0.01*l00)
# fit <- enetLTS(xout,yout,family="binomial",alphas=alphas,lambdas=lambdas)


## for multinomial

n <- 120; p <- 15
NC <- 3           # number of groups
X <- matrix(rnorm(n * p), n, p)
betas <- matrix(1:NC, ncol=NC, nrow=p, byrow=TRUE)
betas[(p-5):p,]=0; betas <- rbind(rep(0,NC),betas)
lv <- cbind(1,X) %*% betas
probs <- exp(lv)/apply(exp(lv),1,sum)
y <- apply(probs,1,function(prob){sample(1:NC, 1, TRUE, prob)})
xout <- X
eps <-0.05                          # %10 contamination to only class 0
m <- ceiling(eps*n)
xout[1:m,] <- xout[1:m,] + 10       # bad leverage points
yout <- y


# determine user supplied alpha and lambda sequences
alphas=seq(0,1,length=11)
lambdas <- seq(from=0.95,to=0.05,by=-0.05)
fit <- enetLTS(xout,yout,family="multinomial",alphas=alphas,lambdas=lambdas)

enetLTS documentation built on May 22, 2022, 1:05 a.m.