robregcc (version 1.0)

robregcc_sp: Robust model estimation approach for regression with compositional covariates.

Description

Fit regression model with compositional covariates for a range of tuning parameter lambda. Model parameters is assumed to be sparse.

Usage

robregcc_sp(X, y, C, beta.init = NULL, gamma.init = NULL,
  beta.wt = NULL, gamma.wt = NULL, control = list(),
  penalty.index = 3, alpha = 1, verbose = TRUE)

Arguments

X

CLR transformed predictor matrix. Centered X

y

model response vector, centered y

C

sub-compositional matrix

beta.init

initial value of model parameter beta

gamma.init

inital value of shift parameter gamma

beta.wt

specify weight of model parameter

gamma.wt

initial value of shift parameter for weight construction/initialization

control

a list of internal parameters controlling the model fitting

penalty.index

a vector of length 2 specifying type of penalty for model parameter and shift parameter respectively. 1, 2, 3 corresponding to adaptive, soft and hard penalty

alpha

elastic net penalty

verbose

TRUE/FALSE for showing progress of the cross validation

Value

Method

Type of penalty used

betapath

model parameter estimate along solution path

gammapath

shift parameter estimate along solution path

lampath

sequence of fitted lambda)

k0

scaling factor

References

Mishra, A., Mueller, C.,(2019) Robust regression with compositional covariates. In prepration. arXiv:1909.04990.

Examples

# NOT RUN {
library(robregcc)
library(magrittr)

data(simulate_robregcc_sp)
X <- simulate_robregcc_sp$X;
y <- simulate_robregcc_sp$y
C <- simulate_robregcc_sp$C
n <- nrow(X); p <- ncol(X); k <-  nrow(C)

# Predictor transformation due to compositional constraint:
# Equivalent to performing centered log-ratio transform 
Xt <- svd(t(C))$u %>% tcrossprod() %>% subtract(diag(p),.) %>% crossprod(t(X),.)
#
Xm <- colMeans(Xt)
Xt <- scale(Xt,Xm,FALSE)                  # centering of predictors 
mean.y <- mean(y)
y <- y - mean.y                           # centering of response 
Xt <- cbind(1,Xt)                         # accounting for intercept in predictor
C <- cbind(0,C)                           # accounting for intercept in constraint
bw <- c(0,rep(1,p))                       # weight matrix to not penalize intercept 

example_seed <- 2*p+1               
set.seed(example_seed) 

# Breakdown point for tukey Bisquare loss function 
b1 = 0.5                    # 50% breakdown point
cc1 =  1.567                # corresponding model parameter
# b1 = 0.25; cc1 =  2.937   

# }
# NOT RUN {
# Initialization [PSC analysis for compositional data]
control <- robregcc_option(maxiter=1000,tol = 1e-4,lminfac = 1e-7)
fit.init <- cpsc_sp(Xt, y,alp=0.4, cfac=2, b1=b1,cc1=cc1,C,bw,1,control) 

# Robust procedure
# control parameters
control <- robregcc_option()
beta.wt <- fit.init$betaR           # Set weight for model parameter beta
beta.wt[1] <- 0
control$gamma = 2                   # gamma for constructing  weighted penalty
control$spb = 40/p                  # fraction of maximum non-zero model parameter beta
control$outMiter = 1000             # Outer loop iteration
control$inMiter = 3000              # Inner loop iteration
control$nlam = 50                   # Number of tuning parameter lambda to be explored
control$lmaxfac = 1                 # Parameter for constructing sequence of lambda
control$lminfac = 1e-8              # Parameter for constructing sequence of lambda 
control$tol = 1e-20;                # tolrence parameter for converging [inner  loop]
control$out.tol = 1e-16             # tolerence parameter for convergence [outer loop]
control$kfold = 5                   # number of fold of crossvalidation


# Robust regression using adaptive elastic net penalty [case III, Table 1]
fit.ada <- robregcc_sp(Xt,y,C, beta.init=fit.init$betaR, 
                       gamma.init = fit.init$residualR,
                       beta.wt=abs(beta.wt), 
                       gamma.wt = abs(fit.init$residualR),
                       control = control, 
                       penalty.index = 1, alpha = 0.95) 
                       
# Robust regression using lasso penalty [Huber equivalent]   [case II, Table 1]
fit.soft <- robregcc_sp(Xt,y,C, beta.init=NULL, gamma.init = NULL,
                        beta.wt=bw, gamma.wt = NULL,
                        control = control, penalty.index = 2, 
                        alpha = 0.95)


# Robust regression using hard thresholding penalty [case I, Table 1]
control$lmaxfac = 1e2        # Parameter for constructing sequence of lambda
control$lminfac = 1e-3       # Parameter for constructing sequence of lambda
fit.hard <- robregcc_sp(Xt,y,C, beta.init=fit.init$betaf, 
                        gamma.init = fit.init$residuals,
                        beta.wt=bw, gamma.wt = NULL,
                        control = control, penalty.index = 3, 
                        alpha = 0.95)
                        
                        
 
# }