set.seed(7869670) knitr::opts_knit$set(self.contained = TRUE)
This is a quick tutorial for using the general GMCM for unsupervised clustering.
The GMCM^1^ package is loaded.
#install.packages("GMCM") # Uncomment to install the GMCM package library("GMCM")
If GMCM is not installed, please uncomment the above line and rerun the script to install from CRAN. The development version can be installed from GitHub following the instructions there.
First, we simulate some toy data. We wish to simulate, say, 1000 observations of 2-dimensions each of which stems from one of 3 components.
In order to do so, we construct a parameter object theta
and simulate from the GMCM.
true_theta <- rtheta(m = 3, d = 2) plot(true_theta) ds <- SimulateGMCMData(n = 1000, theta = true_theta) str(ds)
As can be seen, the SimulateGMCMData
function returns a list
containing the copula observations u
, the unobserved process z
, the latent groups K
, and the theta
used for simulation.
Plotting the true copula realizations and coloring by the true classes shows what we intend to estimate and recover:
plot(ds$u, col = ds$K)
A ranking of u
(or z
) corresponds to what would be the data in a real application.
uhat <- Uhat(ds$u) plot(uhat)
To fit a general GMCM, we must choose some initial parameters. The choose.theta
is a (sometimes) helpful default we here invoke explicitly:
init_theta <- choose.theta(uhat, m = 3) print(init_theta)
The function needs to know how many components we want to estimate though this number may very well be unknown in practice.
With the data loaded and defined initial parameters, the model is now fitted.
est_theta <- fit.full.GMCM(u = uhat, # Ranking function is applied automatically theta = init_theta, method = "NM", max.ite = 5000, verbose = FALSE) print(est_theta)
The fitting method is set to "NM"
with a maximum number of iterations of 100
.
The estimated parameters are used to calculated posterior component probabilities on which the classification is based:
membership <- classify(uhat, est_theta) str(membership)
If it is of interest, the posterior probability can be computed directly using
post_prob <- get.prob(uhat, est_theta) # Compute component probabilities
The number of observations in each class can be e.g. counted by
table(membership)
The results are also displayed by plotting
par(mfrow = c(1,2)) plot(uhat, col = membership, asp = 1) # Plot of estimated copula values z <- GMCM:::qgmm.marginal(uhat, theta = est_theta) # Estimate latent process plot(z, col = membership, asp = 1) # Plot of estimated latent process
The fitted theta
object can also be plotted directly:
plot(est_theta)
This report was generated using rmarkdown^2^ and knitr^3^ under the session given below.
sessionInfo()
Please cite the GMCM paper^1^ if you use the package or shiny app.
cites <- lapply(c("GMCM", "knitr", "rmarkdown"), citation) fmt_cites <- unlist(lapply(cites, format, style = "text")) cat(paste0(" ", seq_along(fmt_cites), ". ", fmt_cites, "\n"))
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.