version 1.2.4

DESCRIPTION

@@ -2,7 +2,7 @@ Package: BNPMIXcluster
 Type: Package
 Title: Bayesian Nonparametric Model for Clustering with Mixed Scale
         Variables
-Version: 1.2.2
+Version: 1.2.4
 Authors@R: c(person("Christian", "Carmona", role = c("aut", "cre"), email = "carmona@stats.ox.ac.uk"),
              person("Luis", "Nieto-Barajas", role = "aut", email = " lnieto@itam.mx"),
              person("Antonio", "Canale", role = "ctb", email = "antonio.canale@unito.it")
@@ -11,15 +11,16 @@ Description: Bayesian nonparametric approach for clustering that is capable to c
 License: GPL (>= 2)
 LazyData: TRUE
 Depends: R (>= 2.10),
-Imports: Rcpp, matrixcalc, truncnorm, mvtnorm, plyr, MASS, compiler
-LinkingTo: Rcpp, RcppArmadillo,
+Imports: compiler, gplots, MASS, matrixcalc, mvtnorm, plyr, Rcpp,
+        truncnorm
+LinkingTo: Rcpp, RcppArmadillo
 Suggests: scatterplot3d
 RoxygenNote: 6.0.1
 NeedsCompilation: yes
-Packaged: 2017-09-19 01:09:42 UTC; carmona
+Packaged: 2017-09-26 19:24:38 UTC; Chris
 Author: Christian Carmona [aut, cre],
   Luis Nieto-Barajas [aut],
   Antonio Canale [ctb]
 Maintainer: Christian Carmona <carmona@stats.ox.ac.uk>
 Repository: CRAN
-Date/Publication: 2017-09-19 07:43:13 UTC
+Date/Publication: 2017-09-26 19:46:36 UTC

MD5

@@ -1,6 +1,6 @@
-5394fba36bc0ff277bdbd6d065a55e05 *DESCRIPTION
-bb5f146f968bfa35e763205c0fbc6015 *NAMESPACE
-accdb60cc818d276e859a50e775c17c4 *R/MIXclustering.R
+a1aca3d7720752c4918f41e438313cfc *DESCRIPTION
+ee8dc04c07cf29bf971aded815ffa1dd *NAMESPACE
+e5f2255fedc354da33a2be5c2cbcd75c *R/MIXclustering.R
 29757b2aa5891398a21a9ae3b2624ffb *R/RcppExports.R
 f03e938d710427a1eb94702927098660 *R/Y_ex_5_1.R
 eb32ef2a7789a076f4d2874dd9a9b351 *R/Z_latent_ex_5_1.R
@@ -11,29 +11,31 @@ d394a878ad63b74f4059f2e18cd880e2 *R/log_f_post_Omega.R
 8732b15933a3369213a6ab87b3760c32 *R/log_f_post_a.R
 195ccc0b2811e15eb218b51049ecbf5c *R/log_f_post_b.R
 384a68a0ae5a8f9afe81f3fd2670b54c *R/meta_param_ex.R
-f6194db2220a60696f73dd6af1533141 *R/plot.MIXcluster.R
+32a050d2e63112c3fc67f498d26d0fca *R/plot.MIXcluster.R
 69be32fb98977f5fdca6425606c70577 *R/poverty.data.R
 acc85741b57f9b59eb93c9270b5332fb *R/print.MIXcluster.R
-4be068a71f5261a7bb7a3e3fab1de5bc *R/sampling_Lambda_jj.R
-686751d117afde7f59329f29218e6ae7 *R/sampling_Omega_ij.R
-da3441f8baff0c7f76cab85c5f25d5d3 *R/sampling_a.R
-9a720f7ae035aa6390a2faa11f190eb0 *R/sampling_b.R
+e0d5b6ad9d258ce1d53580b4d6724c1e *R/sampling_Lambda_jj.R
+5700bd37e4246b9b7ebc2dd8fbfc5028 *R/sampling_Omega_ij.R
+1b1d886146de50df10cbf36505d0a858 *R/sampling_a.R
+cb61793c325136318eb796401ea46bf7 *R/sampling_b.R
 cf1d4c55ebf03e9a5905fc28c5dd1a86 *R/summary.MIXcluster.R
 e74971d952c037671b5c05e5ee7f7692 *data/Y_ex_5_1.rda
 834cac29dc53fe3031328b93c508b31f *data/Z_latent_ex_5_1.rda
 6fc61e48a9f16f1a1c671d96fbddad9f *data/meta_param_ex.rda
 5b4971a336869c5c5160fc0fc2e39c69 *data/poverty.data.rda
 de98b5f45e6bba3970ee1b3fba520a1a *inst/include/BNPMIXcluster.h
 1cac0c598598395619c503dde0b410a5 *inst/include/BNPMIXcluster_RcppExports.h
-2a2e7b03194fa719a2399e04d549fae1 *man/MIXclustering.Rd
+aaca22bf91d115a60a060f65756dae8a *man/MIXclustering.Rd
 2105272bd48fddfe899d1c4967bba3a4 *man/Y_ex_5_1.Rd
 48ae820b9ba3277af28cd2eafa053995 *man/Z_latent_ex_5_1.Rd
 cc5626248654855a2b86be9a5b30913c *man/get_latents.Rd
 56fce3947c463cce12478971b57daeb4 *man/meta_param_ex.Rd
 82ae95703ed254011e5f844ee9a76f86 *man/plot.MIXcluster.Rd
 454e8a1e27e92faf0c683dd6f91ce0ec *man/poverty.data.Rd
-0672ba05c9f92b356c9f5e0013db945b *man/sampling_Omega_ij.Rd
-aca607affca259c0e21d120ceb29240a *man/sampling_a.Rd
+57f6472d997d6a9dff8a617241f2bff6 *man/sampling_Lambda_jj.Rd
+6f4f8652ebd2c039213a9f3b0c2711aa *man/sampling_Omega_ij.Rd
+e9f43ecfef04b3387a73daefa075e61a *man/sampling_a.Rd
+f436ffef64f3d9705b8ca3e605e15e4b *man/sampling_b.Rd
 b599ae226c33b831bd2f909146b2c644 *man/summary.MIXcluster.Rd
 4d751a816a417548d9de308ba62c3ba7 *src/Makevars
 4d751a816a417548d9de308ba62c3ba7 *src/Makevars.win

NAMESPACE

@@ -7,6 +7,10 @@ export(MIXclustering)
 import(plyr)
 importFrom(MASS,ginv)
 importFrom(Rcpp,sourceCpp)
+importFrom(compiler,cmpfun)
+importFrom(gplots,heatmap.2)
+importFrom(grDevices,colorRamp)
+importFrom(grDevices,rgb)
 importFrom(graphics,abline)
 importFrom(graphics,barplot)
 importFrom(graphics,hist)
@@ -19,7 +23,6 @@ importFrom(stats,cor)
 importFrom(stats,dbeta)
 importFrom(stats,dgamma)
 importFrom(stats,diffinv)
-importFrom(stats,heatmap)
 importFrom(stats,model.matrix)
 importFrom(stats,pnorm)
 importFrom(stats,qnorm)

R/MIXclustering.R

@@ -40,6 +40,7 @@
 #' @param d_1_mu Hyperparameter in the prior distribution of the variance of the location in each cluster. See \code{details}.
 #'
 #' @param log_file Specifies a file to save the details with the execution time and the parameters used.
+#' @param keep_param_chains Indicates wheter the simulations of parameters \code{a}, \code{b}, \code{lambda} and \code{omega} should be returned as output.
 #'
 #' @details
 #'
@@ -290,7 +291,8 @@ MIXclustering <- function( Y,
                            sampling_prob=NULL,
                            expansion_f=NULL,
 
-                           log_file=NULL ) {
+                           log_file=NULL,
+                           keep_param_chains=FALSE ) {
 
   if( !any( c(is.matrix(Y),is.data.frame(Y)) ) ) { stop("Y has to be a Matrix or data frame, with each column representing a variable, and each row representing an individual") }
 
@@ -446,6 +448,12 @@ MIXclustering <- function( Y,
       cat('\nError: There is an inconsistency between "sampling_prob" and the number of rows in "Y"\n')
       stop('There is an inconsistency between "sampling_prob" and the number of rows in "Y"')
     }
+  } else {
+    if(is.null(expansion_f)) {
+    expansion_f <- 1/sampling_prob
+    } else {
+      stop('Only one of "sampling_prob" or "expansion_f" should be specified')
+    }
   }
 
   #sampling_prob <- sampling_prob/sum(sampling_prob) # should add up 1
@@ -550,18 +558,27 @@ MIXclustering <- function( Y,
 
   ##### Monitoring acceptance rate for MH #####
 
-
-
+  if(keep_param_chains) {
     Lambda_sim <- matrix(as.numeric(NA),nrow=n_iter,ncol=ncol(sigma_Z) )
     Omega_sim <- array(as.numeric(NA),dim=c(nrow(sigma_Z),ncol(sigma_Z),n_iter))
     a_sim <- as.numeric(NULL)
     b_sim <- as.numeric(NULL)
 
-    Lambda_accept <- matrix(as.numeric(NA),nrow=n_iter,ncol=ncol(sigma_Z) )
-    Omega_accept <- array(as.numeric(NA),dim=c(nrow(sigma_Z),ncol(sigma_Z),n_iter))
-    a_accept <- as.numeric(NULL)
-    b_accept <- as.numeric(NULL)
+    #Lambda_accept <- matrix(as.numeric(NA),nrow=n_iter,ncol=ncol(sigma_Z) )
+    #Omega_accept <- array(as.numeric(NA),dim=c(nrow(sigma_Z),ncol(sigma_Z),n_iter))
+    #a_accept <- as.numeric(NULL)
+    #b_accept <- as.numeric(NULL)
+  } else {
+    Lambda_sim <- NULL
+    Omega_sim <- NULL
+    a_sim <- NULL
+    b_sim <- NULL
 
+    #Lambda_accept <- NULL
+    #Omega_accept <- NULL
+    #a_accept <- NULL
+    #b_accept <- NULL
+  }
 
 
   #if(dev_verbose) {
@@ -763,7 +780,7 @@ MIXclustering <- function( Y,
                                          USING_CPP=USING_CPP )
 
       Lambda_new[j_sigma,j_sigma] <- sqrt(aux_Lambda[[1]])
-      Lambda_accept[iter_i,j_sigma] <- aux_Lambda[[2]]
+      #Lambda_accept[iter_i,j_sigma] <- aux_Lambda[[2]]
 
       # Element with unitary variance in Lambda_new
       # diag(Lambda_new)[aux_var1_Z] <- 1
@@ -819,7 +836,7 @@ MIXclustering <- function( Y,
 
 
           omega_ij_new <- aux_omega_ij_new[[1]]
-          Omega_accept[i_omega,j_omega,iter_i] <- Omega_accept[j_omega,i_omega,iter_i] <- aux_omega_ij_new[[2]]
+          #Omega_accept[i_omega,j_omega,iter_i] <- Omega_accept[j_omega,i_omega,iter_i] <- aux_omega_ij_new[[2]]
 
           Omega_new[i_omega,j_omega] <- Omega_new[j_omega,i_omega] <- omega_ij_new
           if(dev_verbose) {
@@ -885,11 +902,11 @@ MIXclustering <- function( Y,
                                USING_CPP=USING_CPP)
 
       a_new <- aux_a_new$a.chain
-      a_accept <- c( a_accept , aux_a_new$accept.indic )
+      #a_accept <- c( a_accept , aux_a_new$accept.indic )
 
     } else {
       a_new <- a_fix
-      a_accept <- NULL
+      #a_accept <- NULL
     }
 
     a <- a_new
@@ -913,13 +930,13 @@ MIXclustering <- function( Y,
                                USING_CPP=USING_CPP )
 
       b_new <- aux_b_new[[1]]
-      b_accept <- c(b_accept, aux_b_new[[2]])
+      #b_accept <- c(b_accept, aux_b_new[[2]])
 
       if(!all(b_new>-a)){stop('There is a problem sampling from "b", it should be >-a\nb=',b,"\n-a=",-a,sep="")}
 
     } else {
       b_new <- b_fix
-      b_accept <- NULL
+      #b_accept <- NULL
     }
     b <- b_new
     if(dev_verbose) {
@@ -958,10 +975,12 @@ MIXclustering <- function( Y,
     mu_star_map_sim[,paste("iter_",iter_i,sep="")] <- mu_star_map
     #mu_star_n_r_sim[[iter_i]] <- mu_star_n_r
 
-    Lambda_sim[iter_i,] <- diag(Lambda)
-    Omega_sim[,,iter_i] <- Omega
-    a_sim[iter_i] <- a
-    b_sim[iter_i] <- b
+    if(keep_param_chains) {
+      Lambda_sim[iter_i,] <- diag(Lambda)
+      Omega_sim[,,iter_i] <- Omega
+      a_sim[iter_i] <- a
+      b_sim[iter_i] <- b
+    }
 
     time_sim[iter_i+1] <- Sys.time()
 
@@ -1074,6 +1093,7 @@ MIXclustering <- function( Y,
     nr<-table(cluster)[1:nc]
 
     hmm<-matrix(NA,nc,q)
+
     for (i in 1:nc){
       # i <- 1
       for (j in 1:q){

R/plot.MIXcluster.R

@@ -27,22 +27,35 @@
 #' \code{\link{MIXclustering}}
 #'
 #' @importFrom graphics abline barplot hist layout par plot
-#' @importFrom stats heatmap
-#'
+#' @importFrom grDevices colorRamp rgb
+#' @importFrom gplots heatmap.2
+#' 
 #' @export
 
-plot.MIXcluster <- function(x,
-                            type=c("heatmap","chain")[1],
-                            chain.obj=c("n.cluster","a","b","Lambda","Omega","all")[1],
-                            ...
-                            ) {
-
+plot.MIXcluster <- function( x,
+                             type=c("heatmap","chain")[1],
+                             chain.obj=c("n.cluster","a","b","Lambda","Omega","all")[1],
+                             ... ) {
+
   if(is.element("heatmap",type)) {
-    heatmap(1-x$cluster.matrix.avg)
+    ramp <- grDevices::colorRamp(c("white","blue"))
+    blue <- grDevices::rgb( ramp(seq(0, 1, length = 25)), max = 255)
+    gplots::heatmap.2( x$cluster.matrix.avg,
+                       dendrogram='none',
+                       key=FALSE,
+                       #key.title=NA,
+                       keysize=0.1,
+                       Rowv=TRUE,
+                       Colv=TRUE,
+                       trace='none',
+                       density.info="none",
+                       margins=c(1,1),
+                       labRow="", labCol="",
+                       col=blue )
   }
-
+  
   if( is.element("chain",type) ) {
-
+    
     plot_mcmc <- function(x,...){
       op <- par()
       hist_n <- hist(x,plot=F)
@@ -53,7 +66,7 @@ plot.MIXcluster <- function(x,
       barplot( hist_n$counts,axes=F,space=0,xlim=c(0,max(hist_n$counts)),horiz=T)
       suppressWarnings(par(op))
     }
-
+    
     if( (is.element("all",chain.obj)|is.element("n.cluster",chain.obj)) ) {
       plot_mcmc(x$MC.values$n.clusters,main="Number of clusters",ylab="Number of clusters",xlab="iteration")
     }

R/sampling_Lambda_jj.R

@@ -1,4 +1,17 @@
+#' @title
+#'   MCMC sampling of parameter "\eqn{\Lambda_{j,j}}" in the \emph{mixdpcluster} model for bayesian clustering.
+#'
+#' @description
+#'   Generates a sample from the the posterior distribution of the j-th diagonal element \eqn{(i,j)} of the \eqn{\Lambda} matrix in the \emph{mixdpcluster} model for bayesian clustering.
+#'   The simulation is done via Metropolis-Hastings method.
+#' 
+#' 
 #' @importFrom stats dgamma
+#' @importFrom compiler cmpfun
+#' 
+#' @keywords internal
+#' 
+
 sampling_Lambda_jj <- function( n_sim_mh=1, sigma_jj_ini,j,
                                 d_0_z, d_1_z, kappa=1,
                                 Z, mu_Z, sigma_Z, sampling_prob,

R/sampling_Omega_ij.R

@@ -1,8 +1,9 @@
-#' Simulation of "\eqn{\Omega_{i,j}}" in the \emph{mixdpcluster} model for bayesian clustering.
+#' @title
+#'   MCMC sampling of parameter "\eqn{\Omega_{i,j}}" in the \emph{mixdpcluster} model for bayesian clustering.
 #'
 #' @description
-#' Generates simulation of the posterior distribution of the \eqn{(i,j)} element of the \eqn{\Omega} matrix in the \emph{mixdpcluster} model for bayesian clustering.
-#' The simulation is done via Metropolis-Hastings method.
+#'   Generates a sample from the the posterior distribution of the \eqn{(i,j)} element of the \eqn{\Omega} matrix in the \emph{mixdpcluster} model for bayesian clustering.
+#'   The simulation is done via Metropolis-Hastings method.
 #'
 #' @param n number of simulations to be generated
 #' @param Omega.ini matrix \eqn{\Omega} with an initialization value for \eqn{\Omega_{i,j}}.
@@ -26,8 +27,11 @@
 #'
 #' @references
 #' Carmona C., Nieto-Barajas L., Canale A. (2017). \emph{Model based approach for household clustering with mixed scale variables.}
-#'
+#' 
+#' @importFrom compiler cmpfun
+#' 
 #' @keywords internal
+#' 
 
 sampling_Omega_ij <- function( n=1,Omega.ini,i,j,delta=4,
                                Z, mu_Z, Lambda, sampling_prob,

R/sampling_a.R

@@ -1,8 +1,9 @@
-#' Simulation of "\eqn{a}" in the \emph{mixdpcluster} model for bayesian clustering.
+#' @title
+#'   MCMC sampling of parameter "\eqn{a}" in the \emph{mixdpcluster} model for bayesian clustering.
 #'
 #' @description
-#' Generates simulation of the posterior distribution of \eqn{a} in the \emph{mixdpcluster} model for bayesian clustering.
-#' The simulation is done via Metropolis-Hastings method.
+#'   Generates a sample from the posterior distribution of \eqn{a} in the \emph{mixdpcluster} model for bayesian clustering.
+#'   The simulation is done via Metropolis-Hastings method.
 #'
 #' @param n number of simulations to generate
 #' @param a.ini initialization value
@@ -27,7 +28,8 @@
 #' Carmona C., Nieto-Barajas L., Canale A. (2017). \emph{Model based approach for household clustering with mixed scale variables.}
 #'
 #' @importFrom stats rbinom
-#'
+#' @importFrom compiler cmpfun
+#' 
 #' @keywords internal
 
 sampling_a <- function( n=1, a.ini,

R/sampling_b.R

@@ -1,3 +1,14 @@
+#' @title
+#'   MCMC sampling of parameter "\eqn{b}" in the \emph{mixdpcluster} model for bayesian clustering.
+#'
+#' @description
+#'   Generates a sample from the posterior distribution of \eqn{b} in the \emph{mixdpcluster} model for bayesian clustering.
+#'   The simulation is done via Metropolis-Hastings method.
+#'   
+#' @importFrom compiler cmpfun
+#' 
+#' @keywords internal
+#' 
 
 sampling_b <- function( n_sim_mh=1, b_ini,
                         a, d_0_b, d_1_b, eta=1,