R/ClusterVII.R

In GridOnClusters: Cluster-Preserving Multivariate Joint Grid Discretization

#' Fast VII Gaussian Mixture Model
#'
#' @param data A numeric matrix or data frame.
#' @param G An integer vector specifying the numbers of mixture components. 
#'        If NULL, defaults to 1:9.
#' @param threads Number of threads for OpenMP parallelization.
#' @param min_pts an integer to specify the minimum number of points inside 
#'        a cluster. Default: 5
#' @param verbose logical. if \code{TRUE}, print out additional messages. 
#'        Default: \code{FALSE} 
#' @return An object of class \code{ballclust}, mimicking \code{mclust}.
#' @export
ballclust <- function(data, G = NULL, threads = parallel::detectCores() - 1, 
                      min_pts = 5, verbose = FALSE) {
   X <- as.matrix(data)
   
   # Auto K-selection fallback
   if (is.null(G)) {
      G <- 1:9
   }
   
   # Call the Rcpp backend
   res <- fast_VII_mclust(X, G, threads = threads, min_pts = min_pts, verbose = verbose)
   
   # --- ADD THIS SAFETY CHECK ---
   if (is.null(res$z)) {
      stop(paste("Clustering failed for all requested K values (", paste(G, collapse=","), ").",
                 "This usually occurs because k-means failed to initialize or 
                 clusters collapsed below the 'min_pts' ", min_pts, " threshold."))
   }
   
   # Determine hard classifications
   classification <- max.col(res$z)
   
   # Format output to match mclust structure
   out <- list(
      call = match.call(),
      modelName = "VII",
      n = nrow(X),
      d = ncol(X),
      G = res$K,
      BIC = res$BIC,
      loglik = res$loglik,
      df = (res$K - 1) + (res$K * ncol(X)) + res$K,
      bic = res$BIC_all, # All tested BICs
      parameters = list(
         pro = as.numeric(res$pro),
         mean = t(res$mean), # mclust expects p x K
         variance = list(
            modelName = "VII",
            d = ncol(X),
            G = res$K,
            sigmasq = as.numeric(res$variance), # Vector of length K
            scale = as.numeric(res$variance)
         )
      ),
      z = res$z,
      classification = classification
   )
   
   class(out) <- c("mclustVII", "mclust")
   return(out)
}