2. Evaluation of determinants
Michael Friendly
2024-10-05
Source:vignettes/a2-det-ex2.Rmd
a2-det-ex2.Rmd
This example shows two classical ways to find the determinant, of a square matrix. They each work by reducing the problem to a series of smaller ones which can be more easily calculated.
1. Calculate det()
by cofactor expansion
Set up a matrix, and find its determinant (so we know what the answer should be).
## [1] 50
Find cofactors of row 1 elements
The cofactor of element is the signed determinant of what is left when row i, column j of the matrix are deleted. NB: In R, negative subscripts delete rows or columns.
## 18 == 18
## -8 == -8
## -6 == -6
det() = product of row with cofactors
In symbols:
rowCofactors()
is a convenience function, that
calculates these all together
rowCofactors(A, 1)
## [1] 18 -8 -6
Voila: Multiply row 1 times the cofactors of its elements. NB: In R, this multiplication gives a matrix.
A[1,] %*% rowCofactors(A, 1)
## [,1]
## [1,] 50
all.equal( det(A), c(A[1,] %*% rowCofactors(A, 1)) )
## [1] TRUE
2. Finding det()
by Gaussian elimination
(pivoting)
This example follows Green and Carroll, Table 2.2. Start with a 4 x 4
matrix,
,
and save det(M)
.
M <- matrix(c(2, 3, 1, 2,
4, 2, 3, 4,
1, 4, 2, 2,
3, 1, 0, 1), nrow=4, ncol=4, byrow=TRUE)
(dsave <- det(M))
## [1] 15
# ### 'pivot' on the leading diagonal element, M[1,1]:
det()
will be the product of the ‘pivots’, the leading
diagonal elements. This step reduces row 1 and column 1 to 0, so it may
be discarded. NB: In R, dropping a row/column can change a matrix to a
vector, so we use drop = FALSE
inside the subscript.
(d <- M[1,1])
## [1] 2
#-- Reduce row 1, col 1 to 0
(M[1,] <- M[1,, drop=FALSE] / M[1, 1])
## [,1] [,2] [,3] [,4]
## [1,] 1 1.5 0.5 1
(M <- M - M[,1] %*% M[1,, drop=FALSE])
## [,1] [,2] [,3] [,4]
## [1,] 0 0.0 0.0 0
## [2,] 0 -4.0 1.0 0
## [3,] 0 2.5 1.5 1
## [4,] 0 -3.5 -1.5 -2
#-- Drop first row and column
M <- M[-1, -1]
#-- Accumulate the product of pivots
d <- d * M[1, 1]
Repeat, reducing new row, col 1 to 0
(M[1,] <- M[1,, drop=FALSE] / M[1,1])
## [,1] [,2] [,3]
## [1,] 1 -0.25 0
(M <- M - M[,1] %*% M[1,, drop=FALSE])
## [,1] [,2] [,3]
## [1,] 0 0.000 0
## [2,] 0 2.125 1
## [3,] 0 -2.375 -2
M <- M[-1, -1]
d = d * M[1, 1]