An object to store estimated mixture model densities

Description

Instances of DensityModel store the estimated densities for each component and the overall (marginal) estimate of the density. The derived class DensityBatchModel additionally stores the density for each batch / component combination (i.e., if there are 3 components and 10 batches, there are 30 estimated densities). The intended use-case of the DensityModel class is to faciliate visualization of the estimated densities (see examples) as well as to provide an estimate of the number of modes in the overall density. If the number of estimated modes is smaller than the number of components of the best-fitting mixture model, post-hoc merging of components may be useful.

Slots

component

The component densities.

overall

The overall (marginal across batches and components) estimate of the density.

modes

A numeric vector providing the estimated modes in the overall density. The modes are defined by a crude estimate of the first derivative of the overall density (see findModes).

data

A numeric vector containing the data

clusters

A vector providing the k-means clustering of the component means using the modes as centers. If an object of class DensityModel is instantiated with merge=FALSE, this slot takes values 1, ..., K, where K is the number of components.

See Also

DensityModel

Examples

## marginal model
truth <- simulateData(N=2500, p=rep(1/3, 3),
                      theta=c(-1, 0, 1),
                      sds=rep(0.1, 3))
dm <- DensityModel(truth)
print(dm)
dm.merged <- DensityModel(truth, merge=TRUE)
print(dm.merged)
## here, because there are 3 distinct modes, specifying merge=TRUE
## does not change the resulting clusters
identical(clusters(dm), clusters(dm.merged))
## These objects can be plotted
plot(dm)
## Note that calling plot on a MixtureModel-derived object returns
## a density object as a side-effect of the plotting
dm2 <- CNPBayes::plot(truth)
identical(dm, dm2)
## batch model
k <- 3
nbatch <- 3
means <- matrix(c(-1.2, -1.0, -0.8,
                 -0.2, 0, 0.2,
                  0.8, 1, 1.2), nbatch, k, byrow=FALSE)
sds <- matrix(0.1, nbatch, k)
N <- 1500
truth <- simulateBatchData(N=N,
                           batch=rep(letters[1:3], length.out=N),
                           theta=means,
                           sds=sds,
                           p=c(1/5, 1/3, 1-1/3-1/5))
dm <- DensityModel(truth)
dm.merged <- DensityModel(truth, merge=TRUE)
print(dm)
dm2 <- CNPBayes::plot(truth)
identical(dm, dm2)
## suppress plotting of the batch-specific densities
CNPBayes::plot(dm2, show.batch=FALSE)

Results

R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(CNPBayes)
Loading required package: GenomicRanges
Loading required package: BiocGenerics
Loading required package: parallel

Attaching package: 'BiocGenerics'

The following objects are masked from 'package:parallel':

    clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
    clusterExport, clusterMap, parApply, parCapply, parLapply,
    parLapplyLB, parRapply, parSapply, parSapplyLB

The following objects are masked from 'package:stats':

    IQR, mad, xtabs

The following objects are masked from 'package:base':

    Filter, Find, Map, Position, Reduce, anyDuplicated, append,
    as.data.frame, cbind, colnames, do.call, duplicated, eval, evalq,
    get, grep, grepl, intersect, is.unsorted, lapply, lengths, mapply,
    match, mget, order, paste, pmax, pmax.int, pmin, pmin.int, rank,
    rbind, rownames, sapply, setdiff, sort, table, tapply, union,
    unique, unsplit

Loading required package: S4Vectors
Loading required package: stats4

Attaching package: 'S4Vectors'

The following objects are masked from 'package:base':

    colMeans, colSums, expand.grid, rowMeans, rowSums

Loading required package: IRanges
Loading required package: GenomeInfoDb

Attaching package: 'CNPBayes'

The following object is masked from 'package:stats':

    sigma

> png(filename="/home/ddbj/snapshot/RGM3/R_BC/result/CNPBayes/DensityModel-class.Rd_%03d_medium.png", width=480, height=480)
> ### Name: DensityModel-class
> ### Title: An object to store estimated mixture model densities
> ### Aliases: DensityBatchModel-class DensityModel-class
> 
> ### ** Examples
> 
> ## marginal model
> truth <- simulateData(N=2500, p=rep(1/3, 3),
+                       theta=c(-1, 0, 1),
+                       sds=rep(0.1, 3))
> dm <- DensityModel(truth)
> print(dm)
An object of class 'DensityModel'
   component densities:  list of 3 vectors 
   overall density:  vector of length 250 
   modes:  -1, 0, 1 
> dm.merged <- DensityModel(truth, merge=TRUE)
> print(dm.merged)
An object of class 'DensityModel'
   component densities:  list of 3 vectors 
   overall density:  vector of length 250 
   modes:  -1, 0, 1 
> ## here, because there are 3 distinct modes, specifying merge=TRUE
> ## does not change the resulting clusters
> identical(clusters(dm), clusters(dm.merged))
[1] TRUE
> ## These objects can be plotted
> plot(dm)
> ## Note that calling plot on a MixtureModel-derived object returns
> ## a density object as a side-effect of the plotting
> dm2 <- CNPBayes::plot(truth)
> identical(dm, dm2)
[1] TRUE
> ## batch model
> k <- 3
> nbatch <- 3
> means <- matrix(c(-1.2, -1.0, -0.8,
+                  -0.2, 0, 0.2,
+                   0.8, 1, 1.2), nbatch, k, byrow=FALSE)
> sds <- matrix(0.1, nbatch, k)
> N <- 1500
> truth <- simulateBatchData(N=N,
+                            batch=rep(letters[1:3], length.out=N),
+                            theta=means,
+                            sds=sds,
+                            p=c(1/5, 1/3, 1-1/3-1/5))
> dm <- DensityModel(truth)
> dm.merged <- DensityModel(truth, merge=TRUE)
> print(dm)
An object of class 'DensityBatchModel'
   batch: list of  3 matrices 
   component densities:  list of 3 vectors 
   overall density:  vector of length 250 
   modes:  -1.12, -0.85, -0.02, 1.03 
> dm2 <- CNPBayes::plot(truth)
> identical(dm, dm2)
[1] TRUE
> ## suppress plotting of the batch-specific densities
> CNPBayes::plot(dm2, show.batch=FALSE)
> 
> 
> 
> 
> 
> dev.off()
null device 
          1 
>