Nothing
optPenaltyVAR2 <- function (Y,
lambdaMin,
lambdaMax,
lambdaInit=(lambdaMin+lambdaMax)/2,
optimizer="nlm",
...){
########################################################################
#
# DESCRIPTION:
# Automatic selection of the ridge penalties for the estimation of the
# parameters of the VAR(2) model. The penalties are selected through
# maximization of the LOOCV log-likelihood.
#
# ARGUMENTS:
# -> Y : Three-dimensional array containing the data. The
# first, second and third dimensions correspond to
# covariates, time and samples, respectively. The
# data are assumed to centered covariate-wise.
# -> lambdaMin : Numeric of length three, containing the minimum
# values of ridge penalty parameters to be
# considered. The first element is the ridge
# parameter corresponding to the penalty on A1,
# the matrix with lag one auto-regression
# coefficients, the second to matrix A2 containing
# the lag two auto-regression coefficients, while the
# third parameter relates to the penalty on Omega,
# the precision matrix of the errors.
# -> lambdaMax : Numeric of length three, containing the maximum
# values of ridge penalty parameters to be
# considered. The first element is the ridge
# parameter corresponding to the penalty on A1, the
# matrix with lag one auto-regression coefficients,
# the second to matrix A2 containing the lag two
# auto-regression coefficients, while the third
# parameter relates to the penalty on Omega, the
# precision matrix of the errors.
# -> lambdaInit : Numeric of length three, containing the initial
# values of ridge penalty parameters to be considered.
# The first element is the ridge parameter
# corresponding to the penalty on A1, the matrix with
# lag one auto-regression coefficients, the second to
# matrix A2 containing the lag two auto-regression
# coefficients, while the third parameter relates to
# the penalty on Omega, the precision matrix of the
# errors.
# -> optimizer : A character : which optimization function should be
# used: "nlm" (default) or "optim"?
# -> ... : Additional arguments passed on to loglikLOOCVVAR2
#
# DEPENDENCIES:
# library(base) # functions: nlminb, constrOptim
# library(ragt2ridges) # functions: loglikLOOCVVAR2 and its dependencies.
#
# NOTES:
# ....
#
########################################################################
# input checks
if (!is(Y, "array")){
stop("Input (Y) is of wrong class.")
}
if (length(dim(Y)) != 3){
stop("Input (Y) is of wrong dimensions: either covariate, time or sample dimension is missing.")
}
if (!is(lambdaMin, "numeric")){
stop("Input (lambdaMin) is of wrong class.")
}
if (length(lambdaMin) != 3){
stop("Input (lambdaMin) is of wrong length.")
}
if (any(is.na(lambdaMin))){
stop("Input (lambdaMin) does not comprise of positive numbers.")
}
if (any(lambdaMin < 0)){
stop("Input (lambdaMin) does not comprise of positive numbers.")
}
if (!is(lambdaMax, "numeric")){
stop("Input (lambdaMax) is of wrong class.")
}
if (length(lambdaMax) != 3){
stop("Input (lambdaMax) is of wrong length.")
}
if (any(is.na(lambdaMax))){
stop("Input (lambdaMax) does not comprise of positive numbers.")
}
if (any(lambdaMax < 0)){
stop("Input (lambdaMax) does not comprise of positive numbers.")
}
if (any(lambdaMax <= lambdaMin)){
stop("Input (lambdaMax) must be larger (element-wise) than lambdaMin")
}
if (!is(lambdaInit, "numeric")){
stop("Input (lambdaInit) is of wrong class.")
}
if (length(lambdaInit) != 3){
stop("Input (lambdaInit) is of wrong length.")
}
if (any(is.na(lambdaInit))){
stop("Input (lambdaInit) does not comprise of positive numbers.")
}
if (any(lambdaInit < 0)){
stop("Input (lambdaInit) does not comprise of positive numbers.")
}
if (any(lambdaInit <= lambdaMin)){
stop("Input (lambdaInit) must be larger (element-wise) than lambdaMin")
}
if (any(lambdaInit >= lambdaMax)){
stop("Input (lambdaInit) must be smaller (element-wise) than lambdaMax")
}
# optimize LOOCV log-likelihood w.r.t. the penalty parameters
if (optimizer=="optim"){
optLambdas <- constrOptim(lambdaInit,
loglikLOOCVVAR2,
grad=NULL,
ui=rbind(diag(rep(1, 3)), diag(rep(-1, 3))),
ci=c(lambdaMin, -lambdaMax),
Y=Y,
control=list(reltol=0.001),
...)$par
}
if (optimizer=="nlm"){
optLambdas <- nlminb(lambdaInit,
loglikLOOCVVAR2,
gradient=NULL,
lower=lambdaMin,
upper=lambdaMax,
Y=Y,
control=list(rel.tol=0.001),
...)$par
}
return(optLambdas)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.