September 5th, 2019
One of the basic laws of matrixes is that AB != BA, matrix multiplication is not commutative — order matters. This is generally true with some exceptions (see question 2 in this problem set).
A <- matrix(c(2, 2, 1, 4, 2, 3, 1, 0, 1),nrow=3)
B <- t(A)
A %*% B == B %*% A
## [,1] [,2] [,3]
## [1,] FALSE FALSE FALSE
## [2,] FALSE FALSE FALSE
## [3,] FALSE FALSE FALSE
A <- matrix(c(2, 2, 1, 4, 2, 3, 1, 0, 1, 2, 3, 1),nrow=3)
B <- t(A)
#We get an error since the dimansions are incompatible
#A %*% B == B %*% A
#Error in B * A : non-conformable arrays
#Error in A * B : non-conformable arrays
t(A) x A and A x t(A) are equal when the matrix is symetric. A symmetric matrix is a square matrix that is equal to its transpose.
A <- matrix(c(2, 1, 1, 1, 2, 1, 1, 1, 2),nrow=3)
B <- t(A)
A
## [,1] [,2] [,3]
## [1,] 2 1 1
## [2,] 1 2 1
## [3,] 1 1 2
B
## [,1] [,2] [,3]
## [1,] 2 1 1
## [2,] 1 2 1
## [3,] 1 1 2
(A %*% B) == (B %*% A)
## [,1] [,2] [,3]
## [1,] TRUE TRUE TRUE
## [2,] TRUE TRUE TRUE
## [3,] TRUE TRUE TRUE
Matrix factorization is a very important problem. There are supercomputers built just to do matrix factorizations. Every second you are on an airplane, matrices are being factorized. Radars that track flights use a technique called Kalman filtering. At the heart of Kalman Filtering is a Matrix Factorization operation. Kalman Filters are solving linear systems of equations when they track your flight using radars.
Write an R function to factorize a square matrix A into LU or LDU, whichever you prefer. Please submit your response in an R Markdown document using our class naming convention, E.g. LFulton_Assignment2_PS2.png
You don’t have to worry about permuting rows of A and you can assume that A is less than 5x5, if you need to hard-code any variables in your code. If you doing the entire assignment in R, then please submit only one markdown document for both the problems.
factorization <- function(A) {
#Let's assign our matrix to 'upper' so that we can calculate U from that
upper <- A
#Let's create an identity matrix of the same dimension and call it 'lower' so that we can use that as a starting point for L calculation
lower <- diag(nrow = dim(A)[1])
#The following loops move through the matrix and calculate the Upper and Lower Triangular Matrixes:
for (i in 2:dim(A)[1]){ #
for (j in 1:(i-1)){
lower[i,j] <- upper[i,j] / upper[j,j]
upper[i, ] <- -lower[i,j] * upper[j, ] + upper[i, ]
}
}
print("Upper Triangular Martix:")
print(upper)
print("Lower Triangular Martix:")
print(lower)
print("Original Matrix A:")
print(A)
print("LU:")
print(lower%*%upper)
print("Comparing A and LU")
print(A==lower%*%upper)
#If we want to return the U and L rather than print - we would use the following code:
#return(list(lower,upper))
}
A <- matrix(c(2, 4, -4, 1, -4, 3, -6, -9, 5), nrow=3, ncol=3, byrow = TRUE)
factorization(A)
## [1] "Upper Triangular Martix:"
## [,1] [,2] [,3]
## [1,] 2 4 -4.0
## [2,] 0 -6 5.0
## [3,] 0 0 -4.5
## [1] "Lower Triangular Martix:"
## [,1] [,2] [,3]
## [1,] 1.0 0.0 0
## [2,] 0.5 1.0 0
## [3,] -3.0 -0.5 1
## [1] "Original Matrix A:"
## [,1] [,2] [,3]
## [1,] 2 4 -4
## [2,] 1 -4 3
## [3,] -6 -9 5
## [1] "LU:"
## [,1] [,2] [,3]
## [1,] 2 4 -4
## [2,] 1 -4 3
## [3,] -6 -9 5
## [1] "Comparing A and LU"
## [,1] [,2] [,3]
## [1,] TRUE TRUE TRUE
## [2,] TRUE TRUE TRUE
## [3,] TRUE TRUE TRUE
B <- matrix(c(-1, 3, 4, 2, 4, 2, -9, 3, -1), nrow=3, ncol=3, byrow = TRUE)
factorization(B)
## [1] "Upper Triangular Martix:"
## [,1] [,2] [,3]
## [1,] -1 3 4
## [2,] 0 10 10
## [3,] 0 0 -13
## [1] "Lower Triangular Martix:"
## [,1] [,2] [,3]
## [1,] 1 0.0 0
## [2,] -2 1.0 0
## [3,] 9 -2.4 1
## [1] "Original Matrix A:"
## [,1] [,2] [,3]
## [1,] -1 3 4
## [2,] 2 4 2
## [3,] -9 3 -1
## [1] "LU:"
## [,1] [,2] [,3]
## [1,] -1 3 4
## [2,] 2 4 2
## [3,] -9 3 -1
## [1] "Comparing A and LU"
## [,1] [,2] [,3]
## [1,] TRUE TRUE TRUE
## [2,] TRUE TRUE TRUE
## [3,] TRUE TRUE TRUE