JADE (version 2.0-4)

djd: Function for Joint Diagonalization of k Square Matrices in a Deflation Based Manner

Description

This function jointly diagonalizes k real-valued square matrices by searching an orthogonal matrix in a deflation based manner.

Usage

djd(X, G = "max", r = 2, eps = 1e-06, maxiter = 500)

Value

The matrix W

Arguments

X

an array containing the k p times p real valued matrices of dimension c(p, p, k).

G

criterion function used for the the algorithm. Options are max, pow and log. See details.

r

power value used if G="pow" or G="max". 0 is not meaningful for this value. See details.

eps

convergence tolerance.

maxiter

maximum number of iterations.

Author

Klaus Nordhausen, Jari Miettinen

Details

Denote the square matrices as \(A_i\), \(i=1,\ldots,k\). This algorithm searches then an orthogonal matrix W so that \(D_i=W'A_iW\) is diagonal for all \(i\). If the \(A_i\) commute then there is an exact solution. If not, the function will perform an approximate joint diagonalization by maximizing \(\sum G(w_j' A_i w_j)\) where \(w_j\) are the orthogonal vectors in W.

The function G can be choosen to be of the form \(G(x) = |x|^r\) or \(G(x) = log(x)\). If G="max" is chosen, the function G is of the form \(G(x) = |x|^r\), and the diagonalization criterion will be maximized globally at each stage by choosing an appropriate initial value from a set of random vectors. If G="pow" or G="log" are chosen, the initial values are the eigenvectors of \(A_1\) which plays hence a special role.

References

Nordhausen, K., Gutch, H. W., Oja, H. and Theis, F.J. (2012): Joint Diagonalization of Several Scatter Matrices for ICA, in LVA/ICA 2012, LNCS 7191, pp. 172--179.

Miettinen, J., Nordhausen, K., Oja, H. and Taskinen, S. (2014), Deflation-based Separation of Uncorrelated Stationary Time Series, Journal of Multivariate Analysis, 123, 214--227.

Miettinen, J., Nordhausen, K. and Taskinen, S. (2017), Blind Source Separation Based on Joint Diagonalization in R: The Packages JADE and BSSasymp, Journal of Statistical Software, 76, 1--31, <doi:10.18637/jss.v076.i02>.

Examples

Z <- matrix(runif(9), ncol = 3)
U <- eigen(Z %*% t(Z))$vectors
D1 <- diag(runif(3))
D2 <- diag(runif(3))
D3 <- diag(runif(3))
D4 <- diag(runif(3))

X.matrix <- array(0, dim=c(3, 3, 4))
X.matrix[,,1] <- t(U) %*% D1 %*% U
X.matrix[,,2] <- t(U) %*% D2 %*% U
X.matrix[,,3] <- t(U) %*% D3 %*% U
X.matrix[,,4] <- t(U) %*% D4 %*% U

W1 <- djd(X.matrix)
round(U %*% W1, 4) # should be a signed permutation 
                     # matrix if W1 is correct.

W2 <- djd(X.matrix, r=1)
round(U %*% W2, 4) # should be a signed permutation 
                     # matrix if W2 is correct.

W3 <- djd(X.matrix, G="l")
round(U %*% W3, 4) # should be a signed permutation 
                     # matrix if W3 is correct.