Nothing
choose.MSMSNC <-
function(y, X=NULL, max.iter=1000, prec=1e-4, est.var=TRUE, criteria="AIC",cluster=FALSE)
{
se.est<-function(P,y,X,dist="MN", nu=3, gamma=0.5)
{
#mf <- match.call(expand.dots = FALSE)
#if("nu" %in% names(mf$...)) nu=mf$...$nu
#if("gamma" %in% names(mf$...)) gamma=mf$...$gamma
if (!any(dist == c("MN","MT","MSL","MCN","MSN","MSNC","MSTEC","MSTN","MSTT",
"MSSLEC","MSSL","MSSL2","MSCN","MSCN2","MSCEC")))
stop("distribution is not recognized")
y<-as.matrix(y)
if(!is.matrix(y))
stop("y must have at least one element")
if(is.null(X)){X<-array(c(diag(ncol(y))),c(ncol(y),ncol(y),nrow(y)))}
if(is.array(X)==FALSE & is.list(X)==FALSE)
stop("X must be an array or a list")
if(is.array(X))
{Xs<-list()
if(ncol(y)>1 | !is.matrix(X)){
for (i in 1:nrow(y)){
Xs[[i]]<- matrix(t(X[,,i]),nrow=ncol(y))}}
if(ncol(y)==1 & is.matrix(X)){
for (i in 1:nrow(y)){
Xs[[i]]<- matrix(t(X[i,]),nrow=1)}}
X<-Xs}
if (ncol(y) != nrow(X[[1]]))
stop("y does not have the same number of columns than X")
if (nrow(y) != length(X))
stop("y does not have the same number of observations than X")
FI.dist <- get(paste("FI.", dist, sep = ""), mode = "function")
if(dist=="MSTT" | dist=="MSSL2" | dist=="MSTEC" | dist=="MSSLEC") {
#if(!exists("nu")) nu=3
if(!is.numeric(nu)) stop("nu must be a number greater than 1")
if(as.numeric(nu)<1) stop("nu must be a number greater than 1")}
if(dist=="MSCN2" | dist=="MSCEC") {
#if(!exists("nu")) nu=0.1
#if(!exists("gamma")) gamma=0.5
if(!is.numeric(nu)) stop("nu must be a number between 0 and 1")
if((as.numeric(nu)<0 | as.numeric(nu)>1)) stop("nu must be a number between 0 and 1")
if(!is.numeric(gamma)) stop("gamma must be a number between 0 and 1")
if(as.numeric(gamma)<0 | as.numeric(gamma)>1) stop("gamma must be a number between 0 and 1")}
P<-matrix(P,ncol=1)
if (!any(dist == c("MSTEC","MSSLEC","MSCEC","MSTT","MSSL2","MSCN2"))) MI.obs<-FI.dist(P=P, y=y, X=X)
if (dist == "MSTEC" | dist=="MSSLEC" | dist=="MSTT" | dist=="MSSL2") MI.obs<-FI.dist(P=P, y=y, X=X, nu=nu)
if (dist == "MSCEC" | dist=="MSCN2") MI.obs<-FI.dist(P=P, y=y, X=X, nu=nu, gamma=gamma)
test=try(solve2(MI.obs),silent=TRUE)
se=c()
if(is.numeric(test) & max(diag(test))<0)
{
se=sqrt(-diag(test))
}
else stop("Standard errors can't be estimated: Numerical problems with the inversion of the information matrix")
se
}
if(cluster!=TRUE) cluster=FALSE
if(all(criteria!=c("AIC","BIC")))
stop("criteria must be AIC or BIC")
if (cluster==FALSE){
fit.msnc=try(estimate.MSNC(y=y,X=X,max.iter=max.iter,prec=prec,est.var=FALSE),silent=TRUE)
fit.mstec=try(estimate.MSTEC(y=y,X=X,max.iter=max.iter,prec=prec,est.var=FALSE, nu.fixed=FALSE),silent=TRUE)
fit.msslec=try(estimate.MSSLEC(y,X,max.iter=max.iter,prec=prec,est.var=FALSE, nu.fixed=FALSE),silent=TRUE)
fit.mscec=try(estimate.MSCEC(y,X,max.iter=max.iter,prec=prec,est.var=FALSE, nu.fixed=FALSE, gamma.fixed=FALSE),silent=TRUE)
}
else {
cl0 <- parallel::detectCores()-1
cl <- parallel::makeCluster(cl0)
doParallel::registerDoParallel(cl);j=1
kname=c("MSNC","MSTEC","MSSLEC","MSCEC")
fit4 <- foreach(j = 1:4) %dopar%{
#source('SMSNCfunctions.r')
require(skewMLRM)
XX=kname[j]
estimate.dist.XX <- get(paste("estimate.", XX, sep = ""), mode = "function")
if(XX=="MSNC") fit4= try(estimate.dist.XX(y=y,X=X,max.iter=max.iter,prec=prec,est.var=FALSE),silent=TRUE)
if(XX=="MSTEC" | XX=="MSSLEC") fit4= try(estimate.dist.XX(y=y,X=X,max.iter=max.iter,prec=prec,est.var=FALSE,nu.fixed=FALSE),silent=TRUE)
if(XX=="MSNC") fit4= try(estimate.dist.XX(y=y,X=X,max.iter=max.iter,prec=prec,est.var=FALSE,nu.fixed=FALSE,gamma.fixed=FALSE),silent=TRUE)
}
fit.msnc=fit4[[1]]
fit.mstec=fit4[[2]]
fit.msslec=fit4[[3]]
fit.mscec=fit4[[4]]
parallel::stopCluster(cl)
}
if(grepl("Error",fit.mstec)[1] & grepl("Error",fit.msslec)[1] &
grepl("Error",fit.mscec)[1] & grepl("Error",fit.msnc)[1])
{stop("estimation problem in all the models in the SMSNC class")}
index<-c(); maxi=c()
if(criteria=="AIC")
{
if(grepl("Error",fit.mstec)[1]==FALSE){ index<-c(index,1); maxi<-c(maxi, fit.mstec$AIC)}
if(grepl("Error",fit.msslec)[1]==FALSE){ index<-c(index,2); maxi<-c(maxi, fit.msslec$AIC)}
if(grepl("Error",fit.mscec)[1]==FALSE){ index<-c(index,3); maxi<-c(maxi, fit.mscec$AIC)}
if(grepl("Error",fit.msnc)[1]==FALSE){ index<-c(index,4); maxi<-c(maxi, fit.msnc$AIC)}
}
if(criteria=="BIC")
{
if(grepl("Error",fit.mstec)[1]==FALSE){ index<-c(index,1); maxi<-c(maxi, fit.mstec$BIC)}
if(grepl("Error",fit.msslec)[1]==FALSE){ index<-c(index,2); maxi<-c(maxi, fit.msslec$BIC)}
if(grepl("Error",fit.mscec)[1]==FALSE){ index<-c(index,3); maxi<-c(maxi, fit.mscec$BIC)}
if(grepl("Error",fit.msnc)[1]==FALSE){ index<-c(index,4); maxi<-c(maxi, fit.msnc$BIC)}
}
fitted<-c("MSTEC","MSSLEC","MSCEC","MSNC")[index]
index<-index[which.min(maxi)]
fit<-switch(index, "1"=fit.mstec, "2"=fit.msslec, "3"=fit.mscec, "4"=fit.msnc)
selected<-c("MSTEC","MSSLEC","MSCEC","MSNC")[index]
se<-"Error"
if(est.var)
{
if(selected=="MSCEC") se<-try(se.est(P=fit$coefficients,y=y,X=X,dist=selected,nu=fit$nu,gamma=fit$gamma),silent=TRUE)
else se<-try(se.est(P=fit$coefficients,y=y,X=X,dist=selected,nu=fit$nu),silent=TRUE)
}
if(!grepl("Error",se[1])){
names(se)<-names(fit$coefficients)
if(selected!="MSCEC")
{
if(selected!="MSNC") RVAL<- list(coefficients=fit$coefficients, se=se, nu=fit$nu, logLik=fit$logLik,
AIC=fit$AIC, BIC=fit$BIC, iterations=fit$iterations, conv=fit$conv, fitted.models=fitted,
class="MSMSNC", selected.model=selected)
if(selected=="MSNC") RVAL<- list(coefficients=fit$coefficients, se=se, logLik=fit$logLik,
AIC=fit$AIC, BIC=fit$BIC, iterations=fit$iterations, conv=fit$conv, fitted.models=fitted,
class="MSMSNC", selected.model=selected)
}
if(selected=="MSCEC") RVAL<- list(coefficients=fit$coefficients, se=se, nu=fit$nu, gamma=fit$gamma, logLik=fit$logLik,
AIC=fit$AIC, BIC=fit$BIC, iterations=fit$iterations, conv=fit$conv, fitted.models=fitted,
class="MSMSNC", selected.model=selected)
}
if(grepl("Error",se[1]))
{
if(selected!="MSCEC")
{
if(selected!="MSNC") RVAL<- list(coefficients=fit$coefficients, nu=fit$nu, logLik=fit$logLik,
AIC=fit$AIC, BIC=fit$BIC, iterations=fit$iterations, conv=fit$conv, fitted.models=fitted,
class="MSMSNC", selected.model=selected)
if(selected=="MSNC") RVAL<- list(coefficients=fit$coefficients, logLik=fit$logLik,
AIC=fit$AIC, BIC=fit$BIC, iterations=fit$iterations, conv=fit$conv, fitted.models=fitted,
class="MSMSNC", selected.model=selected)
}
if(selected=="MSCEC") RVAL<- list(coefficients=fit$coefficients, nu=fit$nu, gamma=fit$gamma, logLik=fit$logLik,
AIC=fit$AIC, BIC=fit$BIC, iterations=fit$iterations, conv=fit$conv, fitted.models=fitted,
class="MSMSNC", selected.model=selected)
}
class(RVAL)<- "skewMLRM"
RVAL$y<-y; RVAL$X<-X
RVAL$dist<-RVAL$selected.model
RVAL$choose.crit<-criteria
RVAL$"function"<-"choose.MSMSNC"
RVAL
}
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.