kappa {base}  R Documentation 
Compute or Estimate the Condition Number of a Matrix
Description
The condition number of a regular (square) matrix is the product of the norm of the matrix and the norm of its inverse (or pseudoinverse), and hence depends on the kind of matrixnorm.
kappa()
computes by default (an estimate of) the 2norm
condition number of a matrix or of the R
matrix of a QR
decomposition, perhaps of a linear fit. The 2norm condition number
can be shown to be the ratio of the largest to the smallest
nonzero singular value of the matrix.
rcond()
computes an approximation of the reciprocal
condition number, see the details.
Usage
kappa(z, ...)
## Default S3 method:
kappa(z, exact = FALSE,
norm = NULL, method = c("qr", "direct"),
inv_z = solve(z),
triangular = FALSE, uplo = "U", ...)
## S3 method for class 'lm'
kappa(z, ...)
## S3 method for class 'qr'
kappa(z, ...)
.kappa_tri(z, exact = FALSE, LINPACK = TRUE, norm = NULL, uplo = "U", ...)
rcond(x, norm = c("O","I","1"), triangular = FALSE, uplo = "U", ...)
Arguments
z , x 
a numeric or complex matrix or a result of

exact 
logical. Should the result be exact (up to small rounding error) as opposed to fast (but quite inaccurate)? 
norm 
character string, specifying the matrix norm with respect
to which the condition number is to be computed, see the function

method 
a partially matched character string specifying the method to be used;

inv_z 
for 
triangular 
logical. If true, the matrix used is just the upper or
lower triangular part of 
uplo 
character string, either 
LINPACK 
logical. If true and 
... 
further arguments passed to or from other methods;
for 
Details
For kappa()
, if exact = FALSE
(the default) the
condition number is estimated by a cheap approximation to the 1norm of
the triangular matrix R
of the qr(x)
decomposition
z = QR
. However, the exact 2norm calculation (via
svd
) is also likely to be quick enough.
Note that the approximate 1 and Infnorm condition numbers via
method = "direct"
are much faster to
calculate, and rcond()
computes these reciprocal
condition numbers, also for complex matrices, using standard LAPACK
routines.
Currently, also the kappa*()
functions compute these
approximations whenever exact
is false, i.e., by default.
kappa
and rcond
are different interfaces to
partly identical functionality.
.kappa_tri
is an internal function called by kappa.qr
and
kappa.default
; tri
is for triangular and its methods
only consider the upper or lower triangular part of the matrix, depending
on uplo = "U"
or "L"
, where "U"
was internally hard
wired before R 4.4.0.
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.
Value
The condition number, kappa
, or an approximation if
exact = FALSE
.
Author(s)
The design was inspired by (but differs considerably from) the S function of the same name described in Chambers (1992).
Source
The LAPACK routines DTRCON
and ZTRCON
and the LINPACK
routine DTRCO
.
LAPACK and LINPACK are from https://netlib.org/lapack/ and https://netlib.org/linpack/ and their guides are listed in the references.
References
Anderson. E. and ten others (1999)
LAPACK Users' Guide. Third Edition. SIAM.
Available online at
https://netlib.org/lapack/lug/lapack_lug.html.
Chambers, J. M. (1992) Linear models. Chapter 4 of Statistical Models in S eds J. M. Chambers and T. J. Hastie, Wadsworth & Brooks/Cole.
Dongarra, J. J., Bunch, J. R., Moler, C. B. and Stewart, G. W. (1978) LINPACK Users Guide. Philadelphia: SIAM Publications.
See Also
norm
;
svd
for the singular value decomposition and
qr
for the QR
one.
Examples
kappa(x1 < cbind(1, 1:10)) # 15.71
kappa(x1, exact = TRUE) # 13.68
kappa(x2 < cbind(x1, 2:11)) # high! [x2 is singular!]
hilbert < function(n) { i < 1:n; 1 / outer(i  1, i, `+`) }
sv9 < svd(h9 < hilbert(9))$ d
kappa(h9) # pretty high; by default {exact=FALSE, method="qr"} :
kappa(h9) == kappa(qr.R(qr(h9)), norm = "1")
all.equal(kappa(h9, exact = TRUE), # its definition:
max(sv9) / min(sv9),
tolerance = 1e12) ## the same (typically down to 2.22e16)
kappa(h9, exact = TRUE) / kappa(h9) # 0.677 (i.e., rel.error = 32%)
## Exact kappa for rectangular matrix
## panmagic.6npm1(7) :
pm7 < rbind(c( 1, 13, 18, 23, 35, 40, 45),
c(37, 49, 5, 10, 15, 27, 32),
c(24, 29, 41, 46, 2, 14, 19),
c(11, 16, 28, 33, 38, 43, 6),
c(47, 3, 8, 20, 25, 30, 42),
c(34, 39, 44, 7, 12, 17, 22),
c(21, 26, 31, 36, 48, 4, 9))
kappa(pm7, exact=TRUE, norm="1") # no problem for square matrix
m76 < pm7[,1:6]
(m79 < cbind(pm7, 50:56, 63:57))
## MoorePenrose inverse { ~= MASS::ginv(); differing tol (value & meaning)}:
## pinv := p(seudo) inv(erse)
pinv < function(X, s = svd(X), tol = 64*.Machine$double.eps) {
if (is.complex(X))
s$u < Conj(s$u)
dx < dim(X)
## X = U D V' ==> Result = V {1/D} U'
pI < function(u,d,v) tcrossprod(v, u / rep(d, each = dx[1L]))
pos < (d < s$d) > max(tol * max(dx) * d[1L], 0)
if (all(pos))
pI(s$u, d, s$v)
else if (!any(pos))
array(0, dX[2L:1L])
else { # some pos, some not:
i < which(pos)
pI(s$u[, i, drop = FALSE], d[i],
s$v[, i, drop = FALSE])
}
}
## rectangular
kappa(m76, norm="1")
try( kappa(m76, exact=TRUE, norm="1") )# error in solve().. must be square
## ==> use pseudoinverse instead of solve() for rectangular {and norm != "2"}:
iZ < pinv(m76)
kappa(m76, exact=TRUE, norm="1", inv_z = iZ)
kappa(m76, exact=TRUE, norm="M", inv_z = iZ)
kappa(m76, exact=TRUE, norm="I", inv_z = iZ)
iX < pinv(m79)
kappa(m79, exact=TRUE, norm="1", inv_z = iX)
kappa(m79, exact=TRUE, norm="M", inv_z = iX)
kappa(m79, exact=TRUE, norm="I", inv_z = iX)