mvplnMCMCclus: Clustering Using mixtures of MVPLN via MCMC-EM

In anjalisilva/mixMVPLN: Mixtures of Matrix Variate Poisson-log Normal Distributions

Clustering Using mixtures of MVPLN via MCMC-EM

Description

Performs clustering using mixtures of matrix variate Poisson-log normal (MVPLN) via MCMC-EM with an option for parallelization. Model selection can be done using AIC, AIC3, BIC and ICL.

Usage

mvplnMCMCclus(
  dataset,
  membership = "none",
  gmin,
  gmax,
  nChains = 3,
  nIterations = NA,
  initMethod = "kmeans",
  nInitIterations = 2,
  normalize = "Yes",
  numNodes = NA
)

Arguments

dataset	A list of length nUnits, containing Y_j matrices. A matrix Y_j has size r x p, and the dataset will have 'j' such matrices with j = 1,...,n. If a Y_j has all zeros, such Y_j will be removed prior to cluster analysis.
membership	A numeric vector of size length(dataset) containing the cluster membership of each Y_j matrix. If not available, leave as "none".
gmin	A positive integer specifying the minimum number of components to be considered in the clustering run.
gmax	A positive integer, greater than or equal to gmin, specifying the maximum number of components to be considered in the clustering run.
nChains	A positive integer specifying the number of Markov chains. Default is 3, the recommended minimum number.
nIterations	A positive integer specifying the number of iterations for each MCMC chain (including warmup). The value should be greater than 40. The upper limit will depend on size of dataset.
initMethod	A method for initialization. Current options are "kmeans", "random", "medoids", "clara", or "fanny". If nInitIterations is set to >= 1, then this method will be used for initialization. Default is "kmeans".
nInitIterations	A positive integer or zero, specifying the number of initialization runs prior to clustering. The run with the highest loglikelihood is used to initialize values, if nInitIterations > 1. Default value is set to 0, in which case no initialization will be performed.
normalize	A string with options "Yes" or "No" specifying if normalization should be performed. Currently, normalization factors are calculated using TMM method of edgeR package. Default is "Yes".
numNodes	A positive integer indicating the number of nodes to be used from the local machine to run the clustering algorithm. Else leave as NA, so default will be detected as parallel::makeCluster(parallel::detectCores() - 1).

Value

Returns an S3 object of class mvplnParallel with results.

• dataset - The input dataset on which clustering is performed.

• nUnits - Number of units in the input dataset.

• nVariables - Number of variables in the input dataset.

• nOccassions - Number of occassions in the input dataset.

• normFactors - A vector of normalization factors used for input dataset.

• gmin - Minimum number of components considered in the clustering run.

• gmax - Maximum number of components considered in the clustering run.

• initalizationMethod - Method used for initialization.

• allResults - A list with all results.

• loglikelihood - A vector with value of final log-likelihoods for each cluster size.

• nParameters - A vector with number of parameters for each cluster size.

• trueLabels - The vector of true labels, if provided by user.

• ICL_all - A list with all ICL model selection results.

• BIC_all - A list with all BIC model selection results.

• AIC_all - A list with all AIC model selection results.

• AIC3_all - A list with all AIC3 model selection results.

• totalTime - Total time used for clustering and model selection.

Author(s)

Anjali Silva, [e-mail], Sanjeena Dang, [e-mail].

References

Aitchison, J. and C. H. Ho (1989). The multivariate Poisson-log normal distribution. Biometrika 76.

Akaike, H. (1973). Information theory and an extension of the maximum likelihood principle. In Second International Symposium on Information Theory, New York, NY, USA, pp. 267–281. Springer Verlag.

Arlot, S., Brault, V., Baudry, J., Maugis, C., and Michel, B. (2016). capushe: CAlibrating Penalities Using Slope HEuristics. R package version 1.1.1.

Biernacki, C., G. Celeux, and G. Govaert (2000). Assessing a mixture model for clustering with the integrated classification likelihood. IEEE Transactions on Pattern Analysis and Machine Intelligence 22.

Bozdogan, H. (1994). Mixture-model cluster analysis using model selection criteria and a new informational measure of complexity. In Proceedings of the First US/Japan Conference on the Frontiers of Statistical Modeling: An Informational Approach: Volume 2 Multivariate Statistical Modeling, pp. 69–113. Dordrecht: Springer Netherlands.

Robinson, M.D., and Oshlack, A. (2010). A scaling normalization method for differential expression analysis of RNA-seq data. Genome Biology 11, R25.

Schwarz, G. (1978). Estimating the dimension of a model. The Annals of Statistics 6.

Silva, A. et al. (2019). A multivariate Poisson-log normal mixture model for clustering transcriptome sequencing data. BMC Bioinformatics 20. Link

Silva, A. et al. (2018). Finite Mixtures of Matrix Variate Poisson-Log Normal Distributions for Three-Way Count Data. arXiv preprint arXiv:1807.08380.

Examples

## Not run: 
# Generating simulated matrix variate count data
set.seed(1234)
trueG <- 2 # number of total G
truer <- 2 # number of total occasions
truep <- 3 # number of total responses
trueN <- 100 # number of total units

# Mu is a r x p matrix
trueM1 <- matrix(rep(6, (truer * truep)),
                 ncol = truep,
                 nrow = truer, byrow = TRUE)

trueM2 <- matrix(rep(1, (truer * truep)),
                 ncol = truep,
                 nrow = truer,
                 byrow = TRUE)

trueMall <- rbind(trueM1, trueM2)

# Phi is a r x r matrix
# Loading needed packages for generating data
# if (!require(clusterGeneration)) install.packages("clusterGeneration")
# library("clusterGeneration")

# Covariance matrix containing variances and covariances between r occasions
# truePhi1 <- clusterGeneration::genPositiveDefMat("unifcorrmat",
#                                                   dim = truer,
#                                                   rangeVar = c(1, 1.7))$Sigma
truePhi1 <- matrix(c(1.075551, -0.488301, -0.488301, 1.362777), nrow = 2)
truePhi1[1, 1] <- 1 # For identifiability issues

# truePhi2 <- clusterGeneration::genPositiveDefMat("unifcorrmat",
#                                                   dim = truer,
#                                                   rangeVar = c(0.7, 0.7))$Sigma
truePhi2 <- matrix(c(0.7000000, 0.6585887, 0.6585887, 0.7000000), nrow = 2)
truePhi2[1, 1] <- 1 # For identifiability issues
truePhiall <- rbind(truePhi1, truePhi2)

# Omega is a p x p matrix
# Covariance matrix containing variances and covariances between p responses
# trueOmega1 <- clusterGeneration::genPositiveDefMat("unifcorrmat", dim = truep,
#                                    rangeVar = c(1, 1.7))$Sigma
trueOmega1 <- matrix(c(1.0526554, 1.0841910, -0.7976842,
                       1.0841910,  1.1518811, -0.8068102,
                       -0.7976842, -0.8068102,  1.4090578),
                       nrow = 3)
# trueOmega2 <- clusterGeneration::genPositiveDefMat("unifcorrmat", dim = truep,
#                                    rangeVar = c(0.7, 0.7))$Sigma
trueOmega2 <- matrix(c(0.7000000, 0.5513744, 0.4441598,
                       0.5513744, 0.7000000, 0.4726577,
                       0.4441598, 0.4726577, 0.7000000),
                       nrow = 3)
trueOmegaAll <- rbind(trueOmega1, trueOmega2)

# Generated simulated data
sampleData <- mixMVPLN::mvplnDataGenerator(nOccasions = truer,
                                           nResponses = truep,
                                           nUnits = trueN,
                                           mixingProportions = c(0.79, 0.21),
                                           matrixMean = trueMall,
                                           phi = truePhiall,
                                           omega = trueOmegaAll)

# Clustering simulated matrix variate count data
clusteringResults <- mixMVPLN::mvplnMCMCclus(dataset = sampleData$dataset,
                                      membership = sampleData$truemembership,
                                      gmin = 1,
                                      gmax = 2,
                                      nChains = 3,
                                      nIterations = 300,
                                      initMethod = "kmeans",
                                      nInitIterations = 1,
                                      normalize = "Yes")

## End(Not run)

anjalisilva/mixMVPLN documentation built on Jan. 15, 2024, 1:10 a.m.