Package 'jointDiag' reference manual

Title:	Joint Approximate Diagonalization of a Set of Square Matrices
Description:	Different algorithms to perform approximate joint diagonalization of a finite set of square matrices. Depending on the algorithm, orthogonal or non-orthogonal diagonalizer is found. These algorithms are particularly useful in the context of blind source separation. Original publications of the algorithms can be found in Ziehe et al. (2004), Pham and Cardoso (2001) <doi:10.1109/78.942614>, Souloumiac (2009) <doi:10.1109/TSP.2009.2016997>, Vollgraff and Obermayer <doi:10.1109/TSP.2006.877673>. An example of application in the context of Brain-Computer Interfaces EEG denoising can be found in Gouy-Pailler et al (2010) <doi:10.1109/TBME.2009.2032162>.
Authors:	Cedric Gouy-Pailler <[email protected]>
Maintainer:	Cedric Gouy-Pailler <[email protected]>
License:	GPL (>= 2)
Version:	0.4
Built:	2024-11-04 19:54:09 UTC
Source:	CRAN

Wrapper: Joint approximate diagonalization of a set of matrices

Description

This function is mainly a wrapper to the different algorithms provided in the package. So see the help of the different algorithms for the details.

Usage

ajd(M, A0 = NULL, B0 = NULL, eps = .Machine$double.eps, 
	itermax = 200, keepTrace = FALSE, methods = c("jedi"))
ajd(M, A0 = NULL, B0 = NULL, eps = .Machine$double.eps, 
	itermax = 200, keepTrace = FALSE, methods = c("jedi"))

Arguments

`M`	DOUBLE ARRAY (KxKxN). Three-dimensional array with dimensions KxKxN representing the set of square and real-valued matrices to be jointly diagonalized. N is the number of matrices. Matrices are KxK square matrices.
`A0`	DOUBLE MATRIX (KxK). The initial guess of the inverse of a joint diagonalizer. If NULL, an initial guess is automatically generated by the algorithm.
`B0`	DOUBLE MATRIX (KxK). The initial guess of a joint diagonalizer. If NULL, an initial guess is automatically generated by the algorithm.
`eps`	DOUBLE. The algorithm stops when the criterion difference between two iterations is less than eps.
`itermax`	INTEGER. Alternatively, the algorithm stops when itermax sweeps have been performed without reaching convergence. If the maximum number of iteration is performed, a warning appears.
`keepTrace`	BOOLEAN. Do we want to keep the successive estimations of the joint diagonalizer.
`methods`	STRING. One or more methods, choosen among the set of available algorithms. Possible values are: jedi, ffdiag, jadiag, uwedge, qdiag

Details

This function is mainly a wrapper to use the different algorithms provided in the package (see help of the different functions).

Value

If the number of methods is one, the result is the structure provided by the algorithm used.

If the number of methods is more than one, a list of results provided by each algorithm is given. Names of the list correspond to methods.

Author(s)

Cedric Gouy-Pailler ([email protected])

Examples

# generating diagonal matrices
D <- replicate(30, diag(rchisq(df=1,n=10)), simplify=FALSE)
# Mixing and demixing matrices
B <- matrix(rnorm(100),10,10)
A <- solve(B)
C <- array(NA,dim=c(10,10,30))
for (i in 1:30) C[,,i] <- A %*% D[[i]] %*% t(A)
ajd(C,method=c("jedi","ffdiag"))
# generating diagonal matrices
D <- replicate(30, diag(rchisq(df=1,n=10)), simplify=FALSE)
# Mixing and demixing matrices
B <- matrix(rnorm(100),10,10)
A <- solve(B)
C <- array(NA,dim=c(10,10,30))
for (i in 1:30) C[,,i] <- A %*% D[[i]] %*% t(A)
ajd(C,method=c("jedi","ffdiag"))

Joint Approximate Diagonalization of a set of square, symmetric and real-valued matrices

Description

This function performs a Joint Approximate Diagonalization of a set of square and real-valued matrices.

Usage

ffdiag(C0, V0 = NULL, eps = .Machine$double.eps, itermax = 200, 
		keepTrace = FALSE)
ffdiag(C0, V0 = NULL, eps = .Machine$double.eps, itermax = 200, 
		keepTrace = FALSE)

Arguments

`C0`	DOUBLE ARRAY (KxKxN). Three-dimensional array with dimensions KxKxN representing the set of square and real-valued matrices to be jointly diagonalized. N is the number of matrices. Matrices are KxK square matrices.
`V0`	DOUBLE MATRIX (KxK). The initial guess of a joint diagonalizer. If NULL, an initial guess is automatically generated by the algorithm.
`eps`	DOUBLE. The algorithm stops when the criterium difference between two iterations is less than eps.
`itermax`	INTEGER. Alternatively, the algorithm stops when itermax sweeps have been performed without reaching convergence. If the maximum number of iteration is performed, a warning appears.
`keepTrace`	BOOLEAN. Do we want to keep the successive estimations of the joint diagonalizer.

Details

Given a set $C_i$ of N KxK real-valued matrices, the algorithm is looking for a matrix $B$ such that $\forall i \in [1,N]$ , $B C_i B^T$ is as close as possible of a diagonal matrix.

Value

`B`	Estimation of the Joint Diagonalizer.
`criter`	Successive estimates of the cost function across sweeps.
`B_trace`	Array of the successive estimates of B across iterations.

Author(s)

Cedric Gouy-Pailler ([email protected]), from the initial matlab code by A. Ziehe.

References

A. Ziehe, P. Laskov, G. Nolte and K.-R. Mueller; A Fast Algorithm for Joint Diagonalization with Non-orthogonal Transformations and its Application to Blind Source Separation; Journal of Machine Learning Research vol 5, pages 777-800, 2004

Examples

# generating diagonal matrices
D <- replicate(30, diag(rchisq(df=1,n=10)), simplify=FALSE)
# Mixing and demixing matrices
B <- matrix(rnorm(100),10,10)
A <- solve(B)
C <- array(NA,dim=c(10,10,30))
for (i in 1:30) C[,,i] <- A %*% D[[i]] %*% t(A)
B_est <- ffdiag(C)$B
# B_est should be an approximate of B=solve(A)
B_est %*% A
# close to a permutation matrix (with random scales)
# generating diagonal matrices
D <- replicate(30, diag(rchisq(df=1,n=10)), simplify=FALSE)
# Mixing and demixing matrices
B <- matrix(rnorm(100),10,10)
A <- solve(B)
C <- array(NA,dim=c(10,10,30))
for (i in 1:30) C[,,i] <- A %*% D[[i]] %*% t(A)
B_est <- ffdiag(C)$B
# B_est should be an approximate of B=solve(A)
B_est %*% A
# close to a permutation matrix (with random scales)

Joint Approximate Diagonalization of a set of square, symmetric and real-valued matrices

Description

This function performs a Joint Approximate Diagonalization of a set of square, symmetric and real-valued matrices.

Usage

jadiag(M, W_est0 = NULL, eps = .Machine$double.eps, itermax = 200, 
		keepTrace = FALSE)
jadiag(M, W_est0 = NULL, eps = .Machine$double.eps, itermax = 200, 
		keepTrace = FALSE)

Arguments

`M`	DOUBLE ARRAY (KxKxN). Three-dimensional array with dimensions KxKxN representing the set of square, symmetric and real-valued matrices to be jointly diagonalized. N is the number of matrices. Matrices are KxK square matrices.
`W_est0`	DOUBLE MATRIX (KxK). The initial guess of a joint diagonalizer. If NULL, an initial guess is automatically generated by the algorithm.
`eps`	DOUBLE. The algorithm stops when the criterium difference between two iterations is less than eps.
`itermax`	INTEGER. Alternatively, the algorithm stops when itermax sweeps have been performed without reaching convergence. If the maximumu number of iteration is performed, a warning appears.
`keepTrace`	BOOLEAN. Do we want to keep the successive estimations of the joint diagonalizer.

Details

Given a set $C_i$ of N KxK symmetric and real-valued matrices, the algorithm is looking for a matrix $B$ such that $\forall i \in [1,N]$ , $B C_i B^T$ is as close as possible of a diagonal matrix.

Value

`B`	Estimation of the Joint Diagonalizer.
`criter`	Successive estimates of the cost function across sweeps.
`B_trace`	Array of the successive estimates of B across iterations.

Author(s)

Cedric Gouy-Pailler ([email protected]), from the initial C code by Dinh-Tuan Pham.

References

Pham, D. & Cardoso, J.; Blind separation of instantaneous mixtures of nonstationary sources; IEEE Trans. Signal Process., 2001, 49, 1837-1848

Examples

# generating diagonal matrices
D <- replicate(30, diag(rchisq(df=1,n=10)), simplify=FALSE)
# Mixing and demixing matrices
B <- matrix(rnorm(100),10,10)
A <- solve(B)
C <- array(NA,dim=c(10,10,30))
for (i in 1:30) C[,,i] <- A %*% D[[i]] %*% t(A)
B_est <- jadiag(C)$B
# B_est should be an approximate of B=solve(A)
B_est %*% A
# close to a permutation matrix (with random scales)
# generating diagonal matrices
D <- replicate(30, diag(rchisq(df=1,n=10)), simplify=FALSE)
# Mixing and demixing matrices
B <- matrix(rnorm(100),10,10)
A <- solve(B)
C <- array(NA,dim=c(10,10,30))
for (i in 1:30) C[,,i] <- A %*% D[[i]] %*% t(A)
B_est <- jadiag(C)$B
# B_est should be an approximate of B=solve(A)
B_est %*% A
# close to a permutation matrix (with random scales)

Approximate non-orthogonal joint diagonalization of a set of square real-valued matrices

Description

This function performs a Joint Approximate Diagonalization of a set of square and real-valued matrices (not necessarily symmetric). The algorithm seeks the inverse of the joint diagonalizer (the mixing matrix in terms of source separation).

The algorithm uses Givens and hyperbolic rotations to find the inverse of a non-orthogonal joint diagonalizer. It is an extension of the JADE method (orthogonal joint diagonalization).

Usage

jedi(M, A0 = NULL, eps = .Machine$double.eps, itermax = 200, 
		keepTrace = FALSE)
jedi(M, A0 = NULL, eps = .Machine$double.eps, itermax = 200, 
		keepTrace = FALSE)

Arguments

`M`	DOUBLE ARRAY (KxKxN). Three-dimensional array with dimensions KxKxN representing the set of square and real-valued matrices to be jointly diagonalized. N is the number of matrices. Matrices are KxK square matrices.
`A0`	DOUBLE MATRIX (KxK). The initial guess of the inverse of a joint diagonalizer. If NULL, an initial guess is automatically generated by the algorithm.
`eps`	DOUBLE. The algorithm stops when the criterium difference between two iterations is less than eps.
`itermax`	INTEGER. Alternatively, the algorithm stops when itermax sweeps have been performed without reaching convergence. If the maximum number of iteration is performed, a warning appears.
`keepTrace`	BOOLEAN. Do we want to keep the successive estimations of the joint diagonalizer.

Details

Given a set $M_i$ of N K $\times$ K square and real-valued matrices, the algorithm is looking for a matrix $A$ such that $\forall i \in [1,N]$ , $A^{-1} C_i A^{-T}$ is as close as possible of a diagonal matrix.

Value

`A`	Estimation of the Joint Diagonalizer.
`criter`	Successive estimates of the cost function across sweeps.
`A_trace`	Array of the successive estimates of A across iterations.

Warning

This algorithm based on a combination of givens and hyperbolic rotations is covered by a patent (see A. Souloumiac, CEA Saclay).

Author(s)

Cedric Gouy-Pailler ([email protected]), with help from Antoine Souloumiac.

References

Souloumiac, A.; Non-Orthogonal Joint Diagonalization by Combining Givens and Hyperbolic Rotations; IEEE Trans. Signal Process., 2009

Examples

# generating diagonal matrices
D <- replicate(30, diag(rchisq(df=1,n=10)), simplify=FALSE)
# Mixing and demixing matrices
B <- matrix(rnorm(100),10,10)
A <- solve(B)
C <- array(NA,dim=c(10,10,30))
for (i in 1:30) C[,,i] <- A %*% D[[i]] %*% t(A)
A_est <- jedi(C)$A
# A_est should be an approximate of A
B %*% A_est
# close to a permutation matrix (with random scales)
# generating diagonal matrices
D <- replicate(30, diag(rchisq(df=1,n=10)), simplify=FALSE)
# Mixing and demixing matrices
B <- matrix(rnorm(100),10,10)
A <- solve(B)
C <- array(NA,dim=c(10,10,30))
for (i in 1:30) C[,,i] <- A %*% D[[i]] %*% t(A)
A_est <- jedi(C)$A
# A_est should be an approximate of A
B %*% A_est
# close to a permutation matrix (with random scales)

Joint Approximate Diagonalization of a set of square, symmetric and real-valued matrices

Description

This function performs a Joint Approximate Diagonalization of a set of square, symmetric and real-valued matrices.

Usage

qdiag(C, W0 = NULL, eps = .Machine$double.eps, itermax = 200, 
		keepTrace = FALSE)
qdiag(C, W0 = NULL, eps = .Machine$double.eps, itermax = 200, 
		keepTrace = FALSE)

Arguments

`C`	DOUBLE ARRAY (KxKxN). Three-dimensional array with dimensions KxKxN representing the set of square, symmetric and real-valued matrices to be jointly diagonalized. N is the number of matrices. Matrices are KxK square matrices.
`W0`	DOUBLE MATRIX (KxK). The initial guess of a joint diagonalizer. If NULL, an initial guess is automatically generated by the algorithm.
`eps`	DOUBLE. The algorithm stops when the criterium difference between two iterations is less than eps.
`itermax`	INTEGER. Alternatively, the algorithm stops when itermax sweeps have been performed without reaching convergence. If the maximum number of iteration is performed, a warning appears.
`keepTrace`	BOOLEAN. Do we want to keep the successive estimations of the joint diagonalizer.

Details

Given a set $C_i$ of N KxK symmetric and real-valued matrices, the algorithm is looking for a matrix $B$ such that $\forall i \in [1,N]$ , $B C_i B^T$ is as close as possible of a diagonal matrix.

Value

`B`	Estimation of the Joint Diagonalizer.
`criter`	Successive estimates of the cost function across sweeps.
`B_trace`	Array of the successive estimates of B across iterations.

Note

Two versions of the quadratic optimization are present in the paper referenced below. These two versions have different complexities, O(N K^3) and O(K^5). Currently only the version with O(N K^3) is implemented.

Author(s)

Cedric Gouy-Pailler ([email protected]), from the initial matlab code by R. Vollgraf.

References

R. Vollgraf and K. Obermayer; Quadratic Optimization for Approximate Matrix Diagonalization; IEEE Transaction on Signal Processing, 2006

Examples

# generating diagonal matrices
D <- replicate(30, diag(rchisq(df=1,n=10)), simplify=FALSE)
# Mixing and demixing matrices
B <- matrix(rnorm(100),10,10)
A <- solve(B)
C <- array(NA,dim=c(10,10,30))
for (i in 1:30) C[,,i] <- A %*% D[[i]] %*% t(A)
B_est <- qdiag(C)$B
# B_est should be an approximate of B=solve(A)
B_est %*% A
# close to a permutation matrix (with random scales)
# generating diagonal matrices
D <- replicate(30, diag(rchisq(df=1,n=10)), simplify=FALSE)
# Mixing and demixing matrices
B <- matrix(rnorm(100),10,10)
A <- solve(B)
C <- array(NA,dim=c(10,10,30))
for (i in 1:30) C[,,i] <- A %*% D[[i]] %*% t(A)
B_est <- qdiag(C)$B
# B_est should be an approximate of B=solve(A)
B_est %*% A
# close to a permutation matrix (with random scales)

Joint Approximate Diagonalization of a set of square, symmetric and real-valued matrices

Description

This function performs a Joint Approximate Diagonalization of a set of square, symmetric and real-valued matrices.

Usage

uwedge(M, W_est0 = NULL, eps = .Machine$double.eps, itermax = 200, 
		keepTrace = FALSE)
uwedge(M, W_est0 = NULL, eps = .Machine$double.eps, itermax = 200, 
		keepTrace = FALSE)

Arguments

`M`	DOUBLE ARRAY (KxKxN). Three-dimensional array with dimensions KxKxN representing the set of square, symmetric and real-valued matrices to be jointly diagonalized. N is the number of matrices. Matrices are KxK square matrices.
`W_est0`	DOUBLE MATRIX (KxK). The initial guess of a joint diagonalizer. If NULL, an initial guess is automatically generated by the algorithm.
`eps`	DOUBLE. The algorithm stops when the criterium difference between two iterations is less than eps.
`itermax`	INTEGER. Alternatively, the algorithm stops when itermax sweeps have been performed without reaching convergence. If the maximum number of iteration is performed, a warning appears.
`keepTrace`	BOOLEAN. Do we want to keep the successive estimations of the joint diagonalizer.

Details

Given a set $C_i$ of N KxK symmetric and real-valued matrices, the algorithm is looking for a matrix $B$ such that $\forall i \in [1,N]$ , $B C_i B^T$ is as close as possible of a diagonal matrix.

Value

`B`	Estimation of the Joint Diagonalizer.
`criter`	Successive estimates of the cost function across sweeps.
`B_trace`	Array of the successive estimates of B across iterations.

Author(s)

Cedric Gouy-Pailler ([email protected]), from the initial matlab code by P. Tichavsky.

References

Tichavsky, P. & Yeredor, A.; Fast Approximate Joint Diagonalization Incorporating Weight Matrices; IEEE Trans. Signal Process., 2009, 57, 878-891

Examples

# generating diagonal matrices
D <- replicate(30, diag(rchisq(df=1,n=10)), simplify=FALSE)
# Mixing and demixing matrices
B <- matrix(rnorm(100),10,10)
A <- solve(B)
C <- array(NA,dim=c(10,10,30))
for (i in 1:30) C[,,i] <- A %*% D[[i]] %*% t(A)
B_est <- uwedge(C)$B
# B_est should be an approximate of B=solve(A)
B_est %*% A
# close to a permutation matrix (with random scales)
# generating diagonal matrices
D <- replicate(30, diag(rchisq(df=1,n=10)), simplify=FALSE)
# Mixing and demixing matrices
B <- matrix(rnorm(100),10,10)
A <- solve(B)
C <- array(NA,dim=c(10,10,30))
for (i in 1:30) C[,,i] <- A %*% D[[i]] %*% t(A)
B_est <- uwedge(C)$B
# B_est should be an approximate of B=solve(A)
B_est %*% A
# close to a permutation matrix (with random scales)

Package 'jointDiag'

Help Index

Wrapper: Joint approximate diagonalization of a set of matrices

Description

Usage

Arguments

Details

Value

Author(s)

Examples

Joint Approximate Diagonalization of a set of square, symmetric and real-valued matrices

Description

Usage

Arguments

Details

Value

Author(s)

References

Examples

Joint Approximate Diagonalization of a set of square, symmetric and real-valued matrices

Description

Usage

Arguments

Details

Value

Author(s)

References

Examples

Approximate non-orthogonal joint diagonalization of a set of square real-valued matrices

Description

Usage

Arguments

Details

Value

Warning

Author(s)

References

Examples

Joint Approximate Diagonalization of a set of square, symmetric and real-valued matrices

Description

Usage

Arguments

Details

Value

Note

Author(s)

References

Examples

Joint Approximate Diagonalization of a set of square, symmetric and real-valued matrices

Description

Usage

Arguments

Details

Value

Author(s)

References

Examples