Nothing
#' @title Main function for CUBE models
#' @description Main function to estimate and validate a CUBE model for given ratings,
#' explaining uncertainty, feeling and overdispersion.
#' @aliases CUBE
#' @usage CUBE(Formula,data,...)
#' @param Formula Object of class Formula.
#' @param data Data frame from which model matrices and response variables are taken.
#' @param ... Additional arguments to be passed for the specification of the model, Including Y, W, Z for
#' explanatory variables for uncertainty, feeling and overdispersion. Set expinform=TRUE if inference should
#' be based on expected information matrix for model with no covariate. Set starting = ... to pass initial
#' values for EM iterations.
#' @return An object of the class "GEM"-"CUBE" is a list containing the following results:
#' \item{estimates}{Maximum likelihood estimates: \eqn{(\pi, \xi, \phi)}}
#' \item{loglik}{Log-likelihood function at the final estimates}
#' \item{varmat}{Variance-covariance matrix of final estimates}
#' \item{niter}{Number of executed iterations}
#' \item{BIC}{BIC index for the estimated model}
#' @details It is the main function for CUBE models, calling for the corresponding functions whenever
#' covariates are specified: it is possible to select covariates for explaining all the three parameters
#' or only the feeling component. \cr
#' The program also checks if the estimated variance-covariance matrix is positive definite: if not,
#' it prints a warning message and returns a matrix and related results with NA entries.
#' The optimization procedure is run via "optim". If covariates are included only for feeling,
#' the variance-covariance matrix is computed as the inverse of the returned numerically differentiated
#' Hessian matrix (option: hessian=TRUE as argument for "optim"), and the estimation procedure is not
#' iterative, so a NULL result for $niter is produced.
#' If the estimated variance-covariance matrix is not positive definite, the function returns a
#' warning message and produces a matrix with NA entries.
#' @references
#' Iannario M. (2014). Modelling Uncertainty and Overdispersion in Ordinal Data,
#' \emph{Communications in Statistics - Theory and Methods}, \bold{43}, 771--786 \cr
#' Piccolo D. (2015). Inferential issues for CUBE models with covariates,
#' \emph{Communications in Statistics. Theory and Methods}, \bold{44}(23), 771--786. \cr
#' Iannario M. (2015). Detecting latent components in ordinal data with overdispersion by means
#' of a mixture distribution, \emph{Quality & Quantity}, \bold{49}, 977--987 \cr
#' Iannario M. (2016). Testing the overdispersion parameter in CUBE models.
#' \emph{Communications in Statistics: Simulation and Computation}, \bold{45}(5), 1621--1635.\cr
#' @seealso \code{\link{probcube}}, \code{\link{loglikCUBE}}, \code{\link{loglikcuben}}, \code{\link{inibestcube}},
#' \code{\link{inibestcubecsi}}, \code{\link{inibestcubecov}},
#' \code{\link{varmatCUBE}}
#' @keywords internal #models
CUBE<-function(Formula,data,...){
ellipsis.arg<-list(...)
mf<-model.frame(Formula,data=data,na.action=na.omit)
ordinal<-as.numeric(model.response(mf))
#formula given as Formula=ordinal~covpai|covcsi|covphi
covpai<-model.matrix(Formula,data=mf,rhs=1)
covcsi<-model.matrix(Formula,data=mf,rhs=2)
covphi<-model.matrix(Formula,data=mf,rhs=3)
if (ncol(covpai)==0){
Y<-NULL
} else {
Y<-covpai[,-1]
}
if (ncol(covcsi)==0){
W<-NULL
} else {
W<-covcsi[,-1]
}
if (ncol(covphi)==0){
Z<-NULL
} else {
Z<-covphi[,-1]
}
lista<-ellipsis.arg[[1]]
m<-lista[['m']]
maxiter<-lista[['maxiter']]
toler<-lista$toler
starting<-lista$starting
expinform<-lista$expinform
# if (is.null(maxiter)){
# maxiter <- 1000
# }
# if (is.null(toler)){
# toler <- 1e-6
# }
if (is.null(expinform)){
expinform <- FALSE
}
if(is.null(starting)){
if(is.null(W) & is.null(Y) & is.null(Z)) {
starting<-inibestcube(m,ordinal)
}else if(is.null(Y) & is.null(Z) & !is.null(W)){
W<-as.matrix(W)
initial<-inibestcube(m,ordinal)
starting<-inibestcubecsi(m,ordinal,W,initial,maxiter=500,toler=1e-6)
} else {
W<-as.matrix(W); Y<-as.matrix(Y); Z<-as.matrix(Z);
starting<-inibestcubecov(m,ordinal,Y,W,Z)
}
}
if(is.null(W) & is.null(Y) & is.null(Z)) {
mod<-cube000(m,ordinal,starting,maxiter,toler,expinform)
} else if (is.null(Y) & is.null(Z) & !is.null(W)){
W<-as.matrix(W)
mod<-cubecsi(m,ordinal,W,starting,maxiter,toler)
} else if(!is.null(Y) & !is.null(W) & !is.null(Z)){
W<-as.matrix(W); Y<-as.matrix(Y); Z<-as.matrix(Z);
mod<-cubecov(m,ordinal,Y,W,Z,starting,maxiter,toler=1e-2)
} else {
cat("CUBE models not available for this variables specification")
}
stime<-mod$estimates
durata<-mod$time
loglik<-as.numeric(mod$loglik)
niter<-mod$niter
varmat<-mod$varmat
BIC<-as.numeric(mod$BIC)
ordinal<-factor(mod$ordinal,ordered=TRUE)
results<-list('estimates'=stime,'ordinal'=ordinal,'time'=durata,
'loglik'=loglik,'niter'=niter,'varmat'=varmat,
'BIC'=BIC)
#class(results)<-"cube"
return(results)
}
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.