################################################################################## # Load the Old Faithful data data(faithful) # Load the teigen package library(teigen) # Fit all of the models in the teigen family for G=1 to 10 # Use the default starting value strategy (k-means) # The data is not scaled (scaling is the default in teigen) tfaith <- teigen(faithful, models="all", Gs=1:10,scale=FALSE,parallel.cores=TRUE) # Plot the BIC and ICL values for each model matplot(t(tfaith$allbic),ylab="BIC",type="p") matplot(t(tfaith$allbic),type="l",add=TRUE) # View the output tfaith # View the model parameters tfaith$parameters # Plot the resulting density estimate and uncertainty plot (dimensions 1 and 2) plot(tfaith,what="contour") plot(tfaith,what="uncertainty") ################################################################################## # View the highest ICL model parameters tfaith$iclresults$parameters # Plot the model with the highest ICL value # We create a new object called tfaithICL and overwrite # the highest BIC model with the corresponding highest # ICL model. # This is then used for plotting tfaithICL<-tfaith tfaithICL$parameters<-tfaithICL$iclresults$parameters tfaithICL$fuzzy<-tfaithICL$iclresults$fuzzy tfaithICL$G<-tfaithICL$iclresults$G tfaithICL$bestmodel<-tfaithICL$iclresults$bestmodel tfaithICL$classification<-tfaithICL$iclresults$classification # Plot the highest ICL model plot(tfaithICL,what="contour") plot(tfaithICL,what="uncertainty") ################################################################################## # Load the Old Faithful Data data(faithful) # Load the EMMIXskew package library(EMMIXskew) # Scale the data for ease of comparison with other models # The results are not invariant to scaling of the data X <- scale(faithful) # Fit the Skew Normal Model with G=1,2,...,10 # Fit five different scale matrix structures (EEE,EEI,VVV,VVI,VII) # Store the BIC and ICL values for each model BICmat <- matrix(NA,10,5) ICLmat <- matrix(NA,10,5) for (g in 1:10){ for (ncov in 1:5){ fit <- EmSkew(X,g=g,ncov=ncov,distr="msn",nrandom=500,nkmeans=500) BICmat[g,ncov] <- fit$bic ICLmat[g,ncov] <- fit$ICL } } # Plot the resulting BIC values (minus outputed BIC) matplot(-BICmat,ylab="BIC",type="l") matplot(-BICmat,type="p",add=TRUE) # Plot the resulting ICL values matplot(ICLmat,ylab="ICL",type="l") matplot(ICLmat,type="p",add=TRUE) # Refit the model with the highest BIC rmsnfaith <- EmSkew(X,g=2,ncov=2,distr="msn",nrandom=500,nkmeans=500) # Examine the fitted model rmsnfaith # Plot the resulting fit EmSkew.contours(faithful,rmsnfaith) ################################################################################## # Load the Old Faithful Data data(faithful) # Load the EMMIXskew package library(EMMIXskew) # Fit the Skew t Model with G=1,2,...,4 # Fit five different scale matrix structures (EEE,EEI,VVV,VVI,VII) BICmat<-matrix(NA,4,5) ICLmat<-matrix(NA,4,5) for (g in 1:4){ for (ncov in 1:5){ fit<-EmSkew(faithful,g=g,ncov=ncov,distr="mst",nrandom=100,nkmeans=100) BICmat[g,ncov]<-fit$bic ICLmat[g,ncov]<-fit$ICL } } # Plot the resulting BIC values (minus outputed BIC) matplot(-BICmat,ylab="BIC",type="l") matplot(-BICmat,type="p",add=TRUE) # Plot the resulting BIC values (minus outputed BIC) matplot(ICLmat,ylab="ICL",type="l") matplot(ICLmat,type="p",add=TRUE) # Refit the model with the highest BIC rmstfaith<-EmSkew(faithful,g=2,ncov=2,distr="mst",nrandom=100,nkmeans=100) # Look at fitted model rmstfaith # Plot the resulting fit EmSkew.contours(faithful,rmstfaith) ################################################################################## # Load the Old Faithful Data data(faithful) # Load the EMMIXuskew package library(EMMIXuskew) # Fit the unrestricted Skew t Model with G=1,...,5 BICvec<-rep(NA,5) for (g in 1:5){ BICvec[g]<-fmmst(g=g,faithful)$bic } # Plot the resulting BIC values (minus outputed BIC) plot(-BICvec,type="l",ylab="BIC",xlab="G") points(-BICvec) # Refit the model with the highest BIC umstfaith<-fmmst(g=2,faithful) # Examine the fitted model summary(umstfaith) # Plot the resulting fit fmmst.contour.2d(faithful,umstfaith) ################################################################################## # Load the Old Faithful data data(faithful) # Load the flowClust package library(flowClust) # Fit the Box-Cox transformed multivariate t with G=1,2,...,6 # Fit cases where degrees of freedom fixed, estimated, cluster specific (nu.est=0,1,2) # Fit cases where transformation is not used, estimated, cluster specific (trans=0,1,2) nuest<-c(0,1,2,0,1,2,0,1,2) lambdaest<-c(0,0,0,1,1,1,2,2,2) BICmat<-matrix(NA,6,9) for (j in 1:9){ BICmat[,j] <-criterion(flowClust(faithful,K=1:6,nu.est=nuest[j],trans=lambdaest[j],randomStart=100,level=1)) } # Fit the best model again BCfaith<-flowClust(faithful,expName="Old Faithful",K=2,nu.est=1,trans=0,randomStart=100,level=1) # Plot the results dens<-density(BCfaith, faithful) plot(dens) points(faithful,pch=3,col="darkgray") ################################################################################## # Load the Old faithful data data(faithful) # Load the MixGHD package library(MixGHD) # Fit the MGHD models for G=1,2, ..., 10. MGHDfaith<-MGHD(faithful,G=1:10,max.iter=1000,scale=FALSE,modelSel="BIC")