incomplete_cholesky: Incomplete Cholesky Decomposition

View source: R/incomplete_cholesky.R

incomplete_choleskyR Documentation

Incomplete Cholesky Decomposition

Description

The incomplete Cholesky decomposition, which computes approximative low rank decompositions for either Gauss 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.

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 returns 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:

  • Gauss Kernel (fix v1.0.1: corrected sign in exponent): 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

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 indices representing the permutation matrix.

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
\Sexpr[results=rd]{tools:::Rd_expr_doi("10.1145/1102351.1102356")}.

Francis R. Bach, Michael I. Jordan
Kernel independent component analysis
Journal of Machine Learning Research 2002
\Sexpr[results=rd]{tools:::Rd_expr_doi("10.1162/153244303768966085")}

Examples

# approximation of a Gauss 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))

KernelICA documentation built on Aug. 4, 2026, 9:07 a.m.