spam-solve.Rd
backsolve
and forwardsolve
solve a system
of linear equations where the coefficient matrix
is upper or lower triangular.
solve
solves a linear system or computes the inverse
of a matrix if the right-hand-side is missing.
# S4 method for spam
solve(a, b, Rstruct=NULL, ...)
# S4 method for spam
backsolve(r, x, ...)
# S4 method for spam
forwardsolve(l, x, ...)
# S4 method for spam
chol2inv(x, ...)
symmetric positive definite matrix of class spam
or a Cholesky factor
as the result of a chol
call.
object of class spam
or spam.chol.
method returned by the function
chol
.
vector or regular matrix of right-hand-side(s) of a system of linear equations.
the Cholesky structure of a
.
further arguments passed to or from other methods, see ‘Details’ below.
We can solve A %*% x = b
by first computing the Cholesky decomposition A =
t(R)%*%R)
, then solving t(R)%*%y = b
for y
, and
finally solving R%*%x = y
for x
.
solve
combines chol
, a Cholesky decomposition of a
symmetric positive definite sparse matrix, with forwardsolve
and
then backsolve
.
In case a
is from a chol
call, then solve
is an
efficient way to calculate backsolve(a, forwardsolve( t(a), b))
.
However, for a.spam
and a.mat
from a chol
call
with a sparse and ordinary matrix, note that
forwardsolve( a.mat, b, transpose=T, upper.tri=T)
is equivalent to
forwardsolve( t(a.mat), b)
and backsolve(a.spam, forwardsolve(a.spam, b, transpose=T, upper.tri=T))
yields the desired result. But
backsolve(a.spam,forwardsolve(t(a.spam), resid))
is
wrong because t(a.spam)
is a spam
and not a
spam.chol.NgPeyton
object.
forwardsolve
and backsolve
solve a system of linear
equations where the coefficient matrix is lower (forwardsolve
) or
upper (backsolve
) triangular. Usually, the triangular matrix is
result from a chol
call and it is not required to transpose it
for forwardsolve
. Note that arguments of the default
methods k
, upper.tri
and transpose
do not have any
effects here.
Notice that it is more efficient to solve successively the linear equations (both triangular solves) than to implement these in the Fortran code.
If the right-hand-side in solve
is missing it will compute
the inverse of a matrix. For details about the specific Cholsesky
decomposition, see chol
.
Recall that the Cholesky factors are from ordered matrices.
chol2inv(x)
is a faster way to solve(x)
.
There is intentionally no S3 distinction between the classes
spam
and spam.chol.
method.
See references in chol
.
# Generate multivariate form a covariance inverse:
# (usefull for GRMF)
set.seed(13)
n <- 25 # dimension
N <- 1000 # sample size
Sigmainv <- .25^abs(outer(1:n,1:n,"-"))
Sigmainv <- as.spam( Sigmainv, eps=1e-4)
Sigma <- solve( Sigmainv) # for verification
iidsample <- array(rnorm(N*n),c(n,N))
mvsample <- backsolve( chol(Sigmainv), iidsample)
norm( var(t(mvsample)) - Sigma)
#> [1] 1.049394
# compare with:
mvsample <- backsolve( chol(as.matrix( Sigmainv)), iidsample, n)
#### ,n as patch
norm( var(t(mvsample)) - Sigma)
#> [1] 1.040029
# 'solve' step by step:
b <- rnorm( n)
R <- chol(Sigmainv)
norm( backsolve( R, forwardsolve( R, b))-
solve( Sigmainv, b) )
#> [1] 0
norm( backsolve( R, forwardsolve( R, diag(n)))- Sigma )
#> [1] 0
# 'update':
R1 <- update( R, Sigmainv + diag.spam( n))