R Graphical Manual

Browse All

Last data update: 2014.03.03

R: Product Density Independent Component Analysis

ProDenICA

R Documentation

Product Density Independent Component Analysis

Description

Fits an ICA model by directly estimating the densities of the independent components using Poisson GAMs. The densities have the form of tilted Gaussians, and hense directly estimate the contrast functions that lead to negentropy measures. This function supports Section 14.7.4 of 'Elements of Statistical Learning (Hastie, Tibshirani and Friedman, 2009, 2nd Edition)'. Models include 'FastICA'.

Usage

ProDenICA(x, k = p, W0 = NULL, whiten = FALSE, maxit = 20, thresh = 1e-07,
restarts = 0, trace = FALSE, Gfunc = GPois, eps.rank = 1e-07, ...)

Arguments

`x`	input matrix
`k`	Number of components required, less than or equal to the number of columns of x
`W0`	Optional initial matrix (for comparing algorithms)
`whiten`	Logical variable - should x be whitened. If TRUE, the SVD of X=UDV' is computed, and U is used (up to rank(X) columns). Also k is reduced to min(k,rank(X)). If FALSE (default), it is assumed that the user has pre-whitened x (and if not, the function may not perform properly)
`maxit`	Maximum number of iterations; default is 20
`thresh`	Convergence threshold, in terms of relative change in Amari metric; dfault is 1e-7
`restarts`	Number of random restarts; default is 0
`trace`	Trace iterations; default is FALSE
`Gfunc`	Contrast functional which is basis for negentropy measure. Default is `'GPois'` which fits a tilted Gaussian density using a Poisson GAM. Other options are `'G1'` (cosh negentropy) and `'G0'` (kurtosis negentropy)
`eps.rank`	Threshold for deciding rank of x if option `whiten=TRUE`. Any singular value less than `eps.thresh` smaller than the largest is treated as zero
`...`	Additional arguments for `Gfunc` areguments

Details

See Section 14.7.4 of Elements of Statistical Learning (Hastie, Tibshirani and Friedman, 2009, 2nd Edition)

Value

An object of S3 class "ProDenICA" is returned, with the following components:

`W`	Orthonormal matrix that takes the whitened version of x to the independent components
`negentropy`	The total negentropy measure of this solution
`s`	the matrix of k independent components
`whitner`	if `whiten=TRUE`, the matrix that whitens `x`, else `NULL`
`call`	the call that produced this object
`density`	If `Gfunc=GPois`, an list of length `k` with the density estimates for each component

Author(s)

Trevor Hastie and Rob Tibshirani

References

Hastie, T. and Tibshirani, R. (2003) Independent Component Analysis through Product Density Estimation in Advances in Neural Information Processing Systems 15 (Becker, S. and Obermayer, K., eds), MIT Press, Cambridge, MA. pp 649-656
Hastie, T., Tibshirani, R. and Friedman, J. (2009) Elements of Statistical Learning (2nd edition), Springer.
http://www-stat.stanford.edu/~hastie/Papers/ESLII.pdf

Examples

p=2
### Can use letters a-r below for dist
dist="n" 
N=1024
A0<-mixmat(p)
s<-scale(cbind(rjordan(dist,N),rjordan(dist,N)))
x <- s %*% A0
###Whiten the data
x <- scale(x, TRUE, FALSE)
sx <- svd(x)	### orthogonalization function
x <- sqrt(N) * sx$u
target <- solve(A0)
target <- diag(sx$d) %*% t(sx$v) %*% target/sqrt(N)
W0 <- matrix(rnorm(2*2), 2, 2)
W0 <- ICAorthW(W0)
W1 <- ProDenICA(x, W0=W0,trace=TRUE,Gfunc=G1)$W
fit=ProDenICA(x, W0=W0,Gfunc=GPois,trace=TRUE, density=TRUE)
W2 <- fit$W
#distance of FastICA from target
amari(W1,target)
#distance of ProDenICA from target
amari(W2,target)
par(mfrow=c(2,1))
plot(fit)

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(ProDenICA)
Loading required package: gam
Loading required package: splines
Loading required package: foreach
Loaded gam 1.12

> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/ProDenICA/ProDenICA.Rd_%03d_medium.png", width=480, height=480)
> ### Name: ProDenICA
> ### Title: Product Density Independent Component Analysis
> ### Aliases: ProDenICA
> ### Keywords: multivariate distribution
> 
> ### ** Examples
> 
> p=2
> ### Can use letters a-r below for dist
> dist="n" 
> N=1024
> A0<-mixmat(p)
> s<-scale(cbind(rjordan(dist,N),rjordan(dist,N)))
> x <- s %*% A0
> ###Whiten the data
> x <- scale(x, TRUE, FALSE)
> sx <- svd(x)	### orthogonalization function
> x <- sqrt(N) * sx$u
> target <- solve(A0)
> target <- diag(sx$d) %*% t(sx$v) %*% target/sqrt(N)
> W0 <- matrix(rnorm(2*2), 2, 2)
> W0 <- ICAorthW(W0)
> W1 <- ProDenICA(x, W0=W0,trace=TRUE,Gfunc=G1)$W
Iter 1 G 0.3780068 crit 0.9205832 
Iter 2 G 0.3764168 crit 0.6160539 
Iter 3 G 0.3763918 crit 0.7774054 
Iter 4 G 0.3769515 crit 0.8009358 
Iter 5 G 0.3777908 crit 0.1089143 
Iter 6 G 0.3771324 crit 0.08831069 
Iter 7 G 0.3776729 crit 0.06994085 
Iter 8 G 0.3772463 crit 0.05617898 
Iter 9 G 0.3775914 crit 0.04446085 
Iter 10 G 0.3773191 crit 0.03553323 
Iter 11 G 0.3775377 crit 0.02813767 
Iter 12 G 0.377365 crit 0.02242656 
Iter 13 G 0.377503 crit 0.01777256 
Iter 14 G 0.3773938 crit 0.0141438 
Iter 15 G 0.3774808 crit 0.01121587 
Iter 16 G 0.3774119 crit 0.008918075 
Iter 17 G 0.3774668 crit 0.007075231 
Iter 18 G 0.3774233 crit 0.00562284 
Iter 19 G 0.3774578 crit 0.004462349 
Iter 20 G 0.3774304 crit 0.003545228 
Iter 21 G 0.3774522 crit 0.002814125 
> fit=ProDenICA(x, W0=W0,Gfunc=GPois,trace=TRUE, density=TRUE)
Iter 1 G 0.009631627 crit 0.2246009 
Iter 2 G 0.01956424 crit 0.03856108 
Iter 3 G 0.01988045 crit 0.009155474 
Iter 4 G 0.01987206 crit 0.002358653 
Iter 5 G 0.01984736 crit 0.0007049329 
Iter 6 G 0.01984961 crit 0.0002502941 
Iter 7 G 0.01985741 crit 0.000129817 
Iter 8 G 0.01985393 crit 4.851388e-05 
Iter 9 G 0.01985152 crit 8.057595e-06 
Iter 10 G 0.01985214 crit 4.091124e-06 
Iter 11 G 0.01985153 crit 4.410169e-06 
Iter 12 G 0.01985214 crit 4.384489e-06 
Iter 13 G 0.01985153 crit 4.386473e-06 
Iter 14 G 0.01985214 crit 4.386314e-06 
Iter 15 G 0.01985153 crit 4.386327e-06 
Iter 16 G 0.01985214 crit 4.386326e-06 
Iter 17 G 0.01985153 crit 4.386326e-06 
Iter 18 G 0.01985214 crit 4.386327e-06 
Iter 19 G 0.01985153 crit 4.386327e-06 
Iter 20 G 0.01985214 crit 4.386326e-06 
Iter 21 G 0.01985153 crit 4.386326e-06 
> W2 <- fit$W
> #distance of FastICA from target
> amari(W1,target)
[1] 0.1761545
> #distance of ProDenICA from target
> amari(W2,target)
[1] 0.01348925
> par(mfrow=c(2,1))
> plot(fit)
> 
> 
> 
> 
> 
> dev.off()
null device 
          1 
>