chol {base} | R Documentation |
Compute the Cholesky factorization of a real symmetric positive-definite square matrix.
chol(x, ...)
## Default S3 method:
chol(x, pivot = FALSE, LINPACK = FALSE, tol = -1, ...)
x |
an object for which a method exists. The default method applies to numeric (or logical) symmetric, positive-definite matrices. |
... |
arguments to be passed to or from methods. |
pivot |
logical: should pivoting be used? |
LINPACK |
logical. Defunct and gives an error. |
tol |
a numeric tolerance for use with |
chol
is generic: the description here applies to the default
method.
Note that only the upper triangular part of x
is used, so
that R'R = x
when x
is symmetric.
If pivot = FALSE
and x
is not non-negative definite an
error occurs. If x
is positive semi-definite (i.e., some zero
eigenvalues) an error will also occur as a numerical tolerance is used.
If pivot = TRUE
, then the Cholesky decomposition of a positive
semi-definite x
can be computed. The rank of x
is
returned as attr(Q, "rank")
, subject to numerical errors.
The pivot is returned as attr(Q, "pivot")
. It is no longer
the case that t(Q) %*% Q
equals x
. However, setting
pivot <- attr(Q, "pivot")
and oo <- order(pivot)
, it
is true that t(Q[, oo]) %*% Q[, oo]
equals x
,
or, alternatively, t(Q) %*% Q
equals x[pivot,
pivot]
. See the examples.
The value of tol
is passed to LAPACK, with negative values
selecting the default tolerance of (usually) nrow(x) *
.Machine$double.neg.eps * max(diag(x))
. The algorithm terminates once
the pivot is less than tol
.
Unsuccessful results from the underlying LAPACK code will result in an error giving a positive error code: these can only be interpreted by detailed study of the FORTRAN code.
The upper triangular factor of the Cholesky decomposition, i.e., the
matrix R
such that R'R = x
(see example).
If pivoting is used, then two additional attributes
"pivot"
and "rank"
are also returned.
The code does not check for symmetry.
If pivot = TRUE
and x
is not non-negative definite then
there will be a warning message but a meaningless result will occur.
So only use pivot = TRUE
when x
is non-negative definite
by construction.
This is an interface to the LAPACK routines DPOTRF
and
DPSTRF
,
LAPACK is from https://netlib.org/lapack/ and its guide is listed in the references.
Anderson. E. and ten others (1999)
LAPACK Users' Guide. Third Edition. SIAM.
Available on-line at
https://netlib.org/lapack/lug/lapack_lug.html.
Becker, R. A., Chambers, J. M. and Wilks, A. R. (1988) The New S Language. Wadsworth & Brooks/Cole.
chol2inv
for its inverse (without pivoting),
backsolve
for solving linear systems with upper
triangular left sides.
qr
, svd
for related matrix factorizations.
( m <- matrix(c(5,1,1,3),2,2) )
( cm <- chol(m) )
t(cm) %*% cm #-- = 'm'
crossprod(cm) #-- = 'm'
# now for something positive semi-definite
x <- matrix(c(1:5, (1:5)^2), 5, 2)
x <- cbind(x, x[, 1] + 3*x[, 2])
colnames(x) <- letters[20:22]
m <- crossprod(x)
qr(m)$rank # is 2, as it should be
# chol() may fail, depending on numerical rounding:
# chol() unlike qr() does not use a tolerance.
try(chol(m))
(Q <- chol(m, pivot = TRUE))
## we can use this by
pivot <- attr(Q, "pivot")
crossprod(Q[, order(pivot)]) # recover m
## now for a non-positive-definite matrix
( m <- matrix(c(5,-5,-5,3), 2, 2) )
try(chol(m)) # fails
(Q <- chol(m, pivot = TRUE)) # warning
crossprod(Q) # not equal to m