Learn R Programming

KernelICA (version 0.1.0)

incomplete_cholesky: Incomplete Cholesky Decomposition

Description

The incomplete Cholesky decomposition, which computes approximative low rank decompositions for either Gaussian or Hermite kernel matrices. Its implementation is inspired by Matlab and C code of F. Bach (see references) and written with the C++ library Eigen3 for speed purposes.

Usage

incomplete_cholesky(
  x,
  kernel = c("gauss", "hermite"),
  eps = 1e-04,
  sigma = ifelse(length(x) < 1000, 1, 0.5),
  hermite_rank = 3
)

Arguments

x

Numeric vector.

kernel

One of "gauss" or "hermite".

eps

Numeric precision parameter for the matrix approximation.

sigma

Numeric value, setting the kernel variance. Default is 1 for vectors smaller than n=1000, otherwise 0.5.

hermite_rank

Integer value for the rank of the Hermite kernel. This parameter is ignored, when the Gaussian kernel is chosen. Default is 3.

Value

A list containing the following entries:

L

A numeric matrix which values \(L_{ij}\) are \(0\) for \(j > i\).

perm

An integer vector of indeces representing the permutation matrix.

Details

The function approximates kernel matrices of the form \(\boldsymbol{K} = (K_{ij})_{(i,j)} = K(x_i, x_j)\) for a vector \(\boldsymbol{x}\) and a kernel function \(K(\cdot, \cdot)\). It returnes a permutation matrix \(\boldsymbol{P}\) given as index vector and a numeric \(n \times k\) matrix \(\boldsymbol{L}\) which is a "cut off" lower triangle matrix, as it contains only the first \(k\) columns that were necessary to attain a sufficient approximation. These matrices follow the inequality \(\| \boldsymbol{P} \boldsymbol{K} \boldsymbol{P}^T - \boldsymbol{L} \boldsymbol{L}^T\|_1 \leq \epsilon \) where \(\epsilon\) is the given precision parameter. The function offers approximation for kernel matrices of the following two kernels:

  • Gaussian Kernel: \(K(x, y) = e^{(x-y)^2 / 2 \sigma^2}\)

  • Hermite Kernel: \(K(x, y) = \sum_{k=0}^d e^{-x^2 / 2 \sigma^2} e^{-y^2 / 2 \sigma^2} \frac{h_k(x / \sigma) h_k(y / \sigma)}{2^k k!}\), where \(h_k\) is the Hermite polynomial of grade \(k\)

References

Kernel ICA implementation in Matlab and C by F. Bach containing the Incomplete Cholesky Decomposition: https://www.di.ens.fr/~fbach/kernel-ica/index.htm

Francis R. Bach, Michael I. Jordan Predictive low-rank decomposition for kernel methods. Proceedings of the Twenty-second International Conference on Machine Learning (ICML) 2005 10.1145/1102351.1102356.

Francis R. Bach, Michael I. Jordan Kernel independent component analysis Journal of Machine Learning Research 2002 10.1162/153244303768966085

Examples

Run this code
# NOT RUN {
# approximation of a Gaussian kernel matrix
x <- rnorm(500)
x_kernel_mat <- kernel_matrix(x, sigma = 1)
x_chol <- incomplete_cholesky(x, sigma = 1)
L_perm <- x_chol$L[x_chol$perm, ]
x_kernel_approx <- L_perm %*% t(L_perm)
## largest differing value:
max(abs(x_kernel_approx - x_kernel_mat))


# approximation of a Hermite kernel matrix
x_kernel_mat <- kernel_matrix(x, kernel = "hermite", sigma = 0.5)
x_chol <- incomplete_cholesky(x, kernel = "hermite", sigma = 0.5)
L_perm <- x_chol$L[x_chol$perm, ]
x_kernel_approx <- L_perm %*% t(L_perm)
## largest differing value:
max(abs(x_kernel_approx - x_kernel_mat))
# }

Run the code above in your browser using DataLab