GMCM (version 1.1.1)

GMCM-package: Fast optimization of Gaussian Mixture Copula Models

Description

Gaussian mixture copula models (GMCM) can be used for unsupervised clustering and meta analysis. In meta analysis, GMCMs are used to identify features which have been reproduced across multiple experiments. This package provides a fast and general implementation of GMCM cluster analysis and serves as an extension of the features available in the idr package.

Details

If the meta analysis of Li et al. (2011) is to be performed, the function fit.meta.GMCM is used to identify the maximum likelihood estimate of the special Gaussian mixture copula model (GMCM) defined by Li et al. (2011). The function get.IDR computes the local and adjusted Irreprodcible Discovery Rates defined by Li et al. (2011) to determine the level of reproduciblity.

Tewari et. al. (2011) proposed using GMCMs as an general unsupervised clustering tool. If such a general unsupervised clustering is needed, like above, the function fit.full.GMCM computes the maximum likelihood estimate of the general GMCM. The function get.prob is used to estimate the class membership probabilities of each observation.

SimulateGMCMData provide easy simulation from the GMCMs.

References

Li, Q., Brown, J. B. J. B., Huang, H., & Bickel, P. J. (2011). Measuring reproducibility of high-throughput experiments. The Annals of Applied Statistics, 5(3), 1752-1779. doi:10.1214/11-AOAS466

Tewari, A., Giering, M. J., & Raghunathan, A. (2011). Parametric Characterization of Multimodal Distributions with Non-gaussian Modes. 2011 IEEE 11th International Conference on Data Mining Workshops, 286-292. doi:10.1109/ICDMW.2011.135

See Also

Core user functions: fit.meta.GMCM, fit.full.GMCM, get.IDR, get.prob, SimulateGMCMData, SimulateGMMData, rtheta, Uhat, choose.theta, full2meta, meta2full

Package by Li et. al. (2011): idr.

Examples

# Loading data
data(u133VsExon)

# Subsetting data to reduce computation time
u133VsExon <- u133VsExon[1:5000, ]

# Ranking and scaling,
# Remember large values should be critical to the null!
uhat <- Uhat(1 - u133VsExon)

# Visualizing P-values and the ranked and scaled P-values
par(mfrow = c(1,2))
plot(u133VsExon, cex = 0.5, pch = 4, col = "tomato", main = "P-values"
     xlab = "P   (U133)", ylab = "P   (Exon)"
plot(uhat, cex = 0.5, pch = 4, col = "tomato", main = "Ranked P-values",
     xlab = "rank(1-P)   (U133)", ylab = "rank(1-P)   (Exon)"

# Fitting using BFGS
fit <- fit.meta.GMCM(uhat, init.par = c(0.5, 1, 1, 0.5), pgtol = 1e-2,
                     method = "L-BFGS", positive.rho = TRUE, verbose = TRUE)

# Compute IDR values and classify
idr <- get.IDR(uhat, par = fit)
table(idr$K) # 1 = irreproducible, 2 = reproducible

# See clustering results
par(mfrow = c(1,2))
plot(u133VsExon, cex = 0.5, pch = 4, main = "Classified genes",
     col = c("tomato", "steelblue")[idr$K],
     xlab = "P   (U133)", ylab = "P   (Exon)")
  plot(uhat, cex = 0.5, pch = 4, main = "Classified genes",
       col = c("tomato", "steelblue")[idr$K],
       xlab = "rank(1-P)   (U133)", ylab = "rank(1-P)   (Exon)")