Rmixmod: Classification with Mixture Modelling

###################################################################################
##                          GaussianModel.R                                      ##
###################################################################################

###################################################################################
##' @include global.R
##' @include Model.R
NULL
###################################################################################

###################################################################################
##' Constructor of [\code{\linkS4class{GaussianModel}}] class
##'
##' This class defines a gaussian Model. Inherits the [\code{\linkS4class{Model}}] class.
##' 
##' \describe{
##'   \item{family}{character defining a family of models.}
##' }
##'
##' @examples
##'   new("GaussianModel")
##'   new("GaussianModel",family="general")
##'
##'   getSlots("GaussianModel")
##' 
##' @name GaussianModel-class
##' @rdname GaussianModel-class
##' @exportClass GaussianModel
##'
setClass(
    Class="GaussianModel",
    representation=representation(
        family = "character"
    ),
    contains=c("Model"),
    prototype=prototype(
        family = character(0)
    ),
    validity=function(object){
      # spherical models
      spherical.free=c("Gaussian_pk_L_I", "Gaussian_pk_Lk_I")
      spherical.equal=c("Gaussian_p_L_I", "Gaussian_p_Lk_I")
      # all spherical
      spherical=c(spherical.free,spherical.equal)
      # diagonal models
      diagonal.free=c("Gaussian_pk_L_B", "Gaussian_pk_Lk_B", "Gaussian_pk_L_Bk", "Gaussian_pk_Lk_Bk")
      diagonal.equal=c("Gaussian_p_L_B", "Gaussian_p_Lk_B", "Gaussian_p_L_Bk", "Gaussian_p_Lk_Bk")
      # all diagonal
      diagonal=c(diagonal.free,diagonal.equal)
      #  general models
      general.free=c("Gaussian_pk_L_C", "Gaussian_pk_Lk_C", "Gaussian_pk_L_D_Ak_D", "Gaussian_pk_Lk_D_Ak_D", "Gaussian_pk_L_Dk_A_Dk", "Gaussian_pk_Lk_Dk_A_Dk", "Gaussian_pk_L_Ck", "Gaussian_pk_Lk_Ck")
      general.equal=c( "Gaussian_p_L_C", "Gaussian_p_Lk_C", "Gaussian_p_L_D_Ak_D", "Gaussian_p_Lk_D_Ak_D", "Gaussian_p_L_Dk_A_Dk", "Gaussian_p_Lk_Dk_A_Dk", "Gaussian_p_L_Ck", "Gaussian_p_Lk_Ck")
      # all general
      general=c(general.free,general.equal)
      # all models
      all.free=c(spherical.free,diagonal.free,general.free)
      all.equal=c(spherical.equal,diagonal.equal,general.equal)
      all=c(spherical,diagonal,general)
          
      # check listModels validity
      if ( sum(object@listModels %in% all) != length(object@listModels) )
          stop("At least one model is not a valid model. See ?mixmodGaussianModel for the list of all gaussian models.")
          
      # check proportions parameters validity
      if ( !object@equal.proportions & !object@free.proportions )
        stop("equal.proportions and free.porportions cannot be both as FALSE !")
      
      # check whether some family names are not valid
      if ( sum( !(object@family %in% c("all","general","diagonal","spherical")) ) ){
        warning( object@family[which(!(object@family %in% c("all","general","diagonal","spherical")))], ": unknown family name !")
      }
      
      # check proportions
      if ( !object@free.proportions & (sum(object@listModels %in% all.free)>0) )
        stop("At least one model has a free proportions but free.proportions is set as FALSE. See ?mixmodGaussianModel for the list of models with equal proportions.")
      if ( !object@equal.proportions & (sum(object@listModels %in% all.equal)>0) )
        stop("At least one model has an equal proportions but equal.proportions is set as FALSE. See ?mixmodGaussianModel for the list of models with free proportions.")
      
      return(TRUE)
    }
)
###################################################################################

###################################################################################
##' Create an instance of the [\code{\linkS4class{GaussianModel}}] class using new/initialize.
##' 
##' Initialization method. Used internally in the `Rmixmod' package.
##' 
##' @seealso \code{\link{initialize}}
##'
##' @keywords internal
##'
##' @rdname initialize-methods
##'
setMethod(
  f="initialize",
  signature=c("GaussianModel"),
  definition=function(.Object, listModels, family, free.proportions, equal.proportions){
          
    # spherical models
    spherical.free=c("Gaussian_pk_L_I", "Gaussian_pk_Lk_I")
    spherical.equal=c("Gaussian_p_L_I", "Gaussian_p_Lk_I")
    # all spherical
    spherical=c(spherical.free,spherical.equal)
    
    # diagonal models
    diagonal.free=c("Gaussian_pk_L_B", "Gaussian_pk_Lk_B", "Gaussian_pk_L_Bk", "Gaussian_pk_Lk_Bk")
    diagonal.equal=c("Gaussian_p_L_B", "Gaussian_p_Lk_B", "Gaussian_p_L_Bk", "Gaussian_p_Lk_Bk")
    # all diagonal
    diagonal=c(diagonal.free,diagonal.equal)
    
    #  general models
    general.free=c("Gaussian_pk_L_C", "Gaussian_pk_Lk_C", "Gaussian_pk_L_D_Ak_D", "Gaussian_pk_Lk_D_Ak_D", "Gaussian_pk_L_Dk_A_Dk", "Gaussian_pk_Lk_Dk_A_Dk", "Gaussian_pk_L_Ck", "Gaussian_pk_Lk_Ck")
    general.equal=c( "Gaussian_p_L_C", "Gaussian_p_Lk_C", "Gaussian_p_L_D_Ak_D", "Gaussian_p_Lk_D_Ak_D", "Gaussian_p_L_Dk_A_Dk", "Gaussian_p_Lk_Dk_A_Dk", "Gaussian_p_L_Ck", "Gaussian_p_Lk_Ck")
    # all general
    general=c(general.free,general.equal)
    
    # all models
    all.free=c(spherical.free,diagonal.free,general.free)
    all.equal=c(spherical.equal,diagonal.equal,general.equal)
    all=c(spherical,diagonal,general)
    
    if ( !missing(listModels) ){
      # save the list of models
      .Object@listModels <- listModels
      # set family 
      if ( missing(family) ){
        if ( sum(listModels %in% spherical) & sum(listModels %in% diagonal) & sum(listModels %in% general) ){
          .Object@family<-"all"
        }
        else{
          family<-character(0)
          if ( sum(listModels %in% spherical) ){ family<-c(family,"spherical") }
          if ( sum(listModels %in% diagonal) ){ family<-c(family,"diagonal") }
          if ( sum(listModels %in% general) ){ family<-c(family,"general") }
          .Object@family<-family
        }
      }
      else { .Object@family<-family }
      
      # set free.proportions
      if ( missing(free.proportions) ){
        if ( sum(listModels %in% all.free) ){ .Object@free.proportions<-TRUE }
        else{ .Object@free.proportions<-FALSE }
      }
      else{ .Object@free.proportions<-free.proportions }
      # set equal.proportions
      if ( missing(equal.proportions) ){
        if ( sum(listModels %in% all.equal) ){ .Object@equal.proportions<-TRUE }
        else{ .Object@equal.proportions<-FALSE }
      }
      else{ .Object@equal.proportions<-equal.proportions }
    }
    else{
      # check free.proportions option
      if ( missing(free.proportions) ){ .Object@free.proportions<-TRUE }
      else{ .Object@free.proportions<-free.proportions }
      
      # check equal.proportions option
      if ( missing(equal.proportions) ){ .Object@equal.proportions<-TRUE }
      else{ .Object@equal.proportions<-equal.proportions }
      
      # define an empty list of models
      list<-character(0)
      
      # set family as "all" if missing
      if ( missing(family) ){
        .Object@family<-"all"
        if ( .Object@free.proportions ){ list<-c(list,all.free) }
        if ( .Object@equal.proportions ){ list<-c(list,all.equal) }
      }
      else{
        # all gaussian models
        if (sum(family=="all")){
          # set family label in case of multiple entries
          .Object@family = "all"
          if ( .Object@free.proportions ){ list<-c(list,all.free) }
          if ( .Object@equal.proportions ){ list<-c(list,all.equal) }
        }
        else{
          # all spherical models
          if ( sum(family=="spherical") ){
            if ( .Object@free.proportions ){ list<-c(list,spherical.free) }
            if ( .Object@equal.proportions ){ list<-c(list,spherical.equal) }
          }
          # all diagonal models
          if ( sum(family=="diagonal") ){
            .Object@family = "diagonal"
            if ( .Object@free.proportions ){ list<-c(list,diagonal.free) }
            if ( .Object@equal.proportions ){ list<-c(list,diagonal.equal) }
          }
          # all general models
          if ( sum(family=="general") ){
            .Object@family = "general"
            if ( .Object@free.proportions ){ list<-c(list,general.free) }
            if ( .Object@equal.proportions ){ list<-c(list,general.equal) }
          }
          # set family label in case of multiple entries
          .Object@family = family
        }
      }
      # create the list of models depending on the proportions option
      .Object@listModels<-list
    }
    
    validObject(.Object)        
    return(.Object)
  }
)
###################################################################################


###################################################################################
##' Create an instance of the [\code{\linkS4class{GaussianModel}}] class
##'
##' Define a list of Gaussian model to test in MIXMOD.
##' 
##' In the Gaussian mixture model, following Banfield and Raftery (1993) and Celeux and Govaert (1995), we consider a parameterization of the variance matrices of the mixture components consisting of expressing the variance matrix \eqn{\Sigma_{k}} in terms of its eigenvalue decomposition \deqn{ \Sigma_{k}= \lambda_{k} D_{k} A_{k}D'_{k}} where \eqn{\lambda_{k}=|\Sigma_{k}|^{1/d}, D_{k}} is the matrix of eigenvectors of \eqn{\Sigma_{k}} and \eqn{A_{k}} is a diagonal matrix, such that \eqn{| A_{k} |=1}, with the normalized eigenvalues of \eqn{\Sigma_{k}} on the diagonal in a decreasing order. The parameter \eqn{\lambda_{k}} determines the \emph{volume} of the \eqn{k}th cluster, \eqn{D_{k}} its \emph{orientation} and \eqn{A_{k}} its \emph{shape}.  By allowing some but not all of these quantities to vary between clusters, we obtain parsimonious and easily interpreted models which are appropriate to describe various clustering situations.
##' 
##' In general family, we can allow the volumes, the shapes and the orientations of clusters to vary or to be equal between clusters. Variations on assumptions on the parameters \eqn{\lambda_{k}, D_{k}} and \eqn{A_{k}} \eqn{(1 \leq k \leq K)} lead to 8 general models of interest. For instance, we can assume different volumes and keep the shapes and orientations equal by requiring that \eqn{A_{k}=A} (\eqn{A} unknown) and \eqn{D_{k}=D} (\eqn{D} unknown) for \eqn{k=1,\ldots,K}. We denote this model \eqn{[\lambda_{k}DAD']}. With this convention, writing \eqn{[\lambda D_{k}AD'_{k}]} means that we consider the mixture model with equal volumes, equal shapes and different orientations.
##' In diagonal family, we assume that the variance matrices \eqn{\Sigma_{k}} are diagonal. In the parameterization, it means that the orientation matrices \eqn{D_{k}} are permutation matrices. We write \eqn{\Sigma_{k}=\lambda_{k}B_{k}} where \eqn{B_{k}} is a diagonal matrix with \eqn{| B_{k}|=1}.  This particular parameterization gives rise to 4 models: \eqn{[\lambda B]}, \eqn{[\lambda_{k}B]}, \eqn{[\lambda B_{k}]} and \eqn{[\lambda_{k}B_{k}]}.
##' 
##' In spherical family, we assume spherical shapes, namely \eqn{A_{k}=I}, \eqn{I} denoting the identity matrix. In such a case, two parsimonious models are in competition: \eqn{[\lambda I]} and \eqn{[\lambda_{k}I]}.
##'
##' @param family character defining a family of models. "general" for the general family, "diagonal" for the diagonal family, "spherical" for the spherical family and "all" for all families. Default is "general".
##' @param listModels a list of characters containing a list of models. It is optional.
##' @param free.proportions logical to include models with free proportions. Default is TRUE.
##' @param equal.proportions logical to include models with equal proportions. Default is TRUE.
##'
##' @return an object of [\code{\linkS4class{GaussianModel}}] which contains some of the 28 Gaussian Models:
##' \tabular{rlllll}{
##'     Model  \tab Family \tab Prop. \tab Volume \tab Shape \tab Orient. \cr
##'     Gaussian_p_L_C         \tab General \tab Equal \tab Equal \tab Equal  \tab Equal \cr
##'     Gaussian_p_Lk_C        \tab \tab \tab Free \tab Equal \tab Equal \cr
##'     Gaussian_p_L_D_Ak_D    \tab  \tab \tab Equal \tab Free \tab Equal \cr
##'     Gaussian_p_Lk_D_Ak_D   \tab  \tab \tab Free \tab Free \tab Equal \cr
##'     Gaussian_p_L_Dk_A_Dk   \tab  \tab \tab Equal \tab Equal \tab Free \cr
##'     Gaussian_p_Lk_Dk_A_Dk  \tab  \tab \tab Free \tab Equal \tab Free \cr
##'     Gaussian_p_L_Ck        \tab  \tab \tab Equal \tab Free \tab Free \cr
##'     Gaussian_p_Lk_Ck       \tab  \tab \tab Free \tab Free \tab Free \cr
##'     Gaussian_p_L_B         \tab Diagonal  \tab Equal \tab Equal \tab Equal \tab Axes \cr
##'     Gaussian_p_Lk_B        \tab  \tab \tab Free \tab Equal \tab Axes \cr
##'     Gaussian_p_L_Bk        \tab  \tab \tab Equal \tab Free \tab Axes \cr
##'     Gaussian_p_Lk_Bk       \tab  \tab \tab Free \tab Free \tab Axes \cr
##'     Gaussian_p_L_I         \tab Spherical \tab Equal \tab Equal \tab Equal \tab NA \cr
##'     Gaussian_p_Lk_I        \tab \tab \tab Free \tab Equal \tab NA \cr
##'     Gaussian_pk_L_C        \tab General  \tab Free \tab Equal \tab Equal \tab Equal \cr
##'     Gaussian_pk_Lk_C       \tab \tab \tab Free \tab Equal \tab Equal \cr
##'     Gaussian_pk_L_D_Ak_D   \tab \tab \tab Equal \tab Free \tab Equal \cr
##'     Gaussian_pk_Lk_D_Ak_D  \tab \tab \tab Free \tab Free \tab Equal \cr
##'     Gaussian_pk_L_Dk_A_Dk  \tab \tab \tab Equal \tab Equal \tab Free \cr
##'     Gaussian_pk_Lk_Dk_A_Dk \tab \tab \tab Free \tab Equal \tab Free \cr
##'     Gaussian_pk_L_Ck       \tab \tab \tab Equal \tab Free \tab Free \cr
##'     Gaussian_pk_Lk_Ck      \tab \tab \tab Free \tab Free \tab Free \cr
##'     Gaussian_pk_L_B        \tab  Diagonal  \tab Free \tab Equal \tab Equal \tab Axes \cr
##'     Gaussian_pk_Lk_B       \tab \tab \tab Free \tab Equal \tab Axes \cr
##'     Gaussian_pk_L_Bk       \tab \tab \tab Equal \tab Free \tab Axes \cr
##'     Gaussian_pk_Lk_Bk      \tab \tab \tab Free \tab Free \tab Axes \cr
##'     Gaussian_pk_L_I        \tab Spherical  \tab Free \tab Equal \tab Equal \tab NA \cr
##'     Gaussian_pk_Lk_I       \tab \tab \tab Free \tab Equal \tab NA \cr 
##' }
##'
##' @references C. Biernacki, G. Celeux, G. Govaert, F. Langrognet. "Model-Based Cluster and Discriminant Analysis with the MIXMOD Software". Computational Statistics and Data Analysis, vol. 51/2, pp. 587-600. (2006)
##' @examples
##'   mixmodGaussianModel()
##'   # all Gaussian models with equal proportions
##'   mixmodGaussianModel(family="all",free.proportions=FALSE)
##'   # Diagonal and Spherical Gaussian models
##'   mixmodGaussianModel(family=c("diagonal","spherical"))
##'   # Gaussian models with a pre-defined list
##'   mixmodGaussianModel(listModels=c("Gaussian_p_L_C","Gaussian_p_L_Ck","Gaussian_pk_L_I"))
##'
##' @author  Florent Langrognet and Remi Lebret and Christian Poli ans Serge Iovleff, with contributions from C. Biernacki and G. Celeux and G. Govaert \email{contact@@mixmod.org}
##' @export
##'
mixmodGaussianModel<- function( family="all", listModels=NULL, free.proportions=TRUE, equal.proportions=TRUE ){
    if ( is.null(listModels) ){
      new("GaussianModel", family=family, free.proportions=free.proportions, equal.proportions=equal.proportions)
    }else{
      new("GaussianModel", listModels=listModels)
    }
}
###################################################################################



###################################################################################
##' @rdname extract-methods
##' @aliases [,GaussianModel-method
##'
setMethod(
  f="[", 
  signature(x = "GaussianModel"),
  definition=function(x,i,j,drop){
    if ( missing(j) ){
      switch(EXPR=i,
        "listModels"={return(x@listModels)},
        "free.proportions"={return(x@free.proportions)},
        "equal.proportions"={return(x@equal.proportions)},
        "family"={return(x@family)},
        stop("This attribute doesn't exist !")
      )
    }else{
      switch(EXPR=i,
        "listModels"={return(x@listModels[j])},
        stop("This attribute doesn't exist !")
      )
    }
  }
)
###################################################################################



###################################################################################
##' 
##'
##' @name [
##' @rdname extract-methods
##' @aliases [<-,GaussianModel-method
##'
setReplaceMethod(
  f="[", 
  signature(x = "GaussianModel"), 
  definition=function(x,i,j,value){
    if ( missing(j) ){
      switch(EXPR=i,
        "listModels"={x@listModels<-value},
        "free.proportions"={x@free.proportions<-value},
        "equal.proportions"={x@equal.proportions<-value},
        "family"={return(x@family)},
        stop("This attribute doesn't exist !")
      )
    }else{
      switch(EXPR=i,
        "listModels"={x@listModels[j]<-value},
        stop("This attribute doesn't exist !")
      )
    }
    validObject(x)
    return(x)
  }
)
###################################################################################