Problem Set 1

  1. Show that \(A^T\)A \(\neq\) A\(A^T\) in general. (Proof and demonstration)

Create a 3x3 matrix for demonstration

Mtx_3 <- matrix( c(1,3,5,2,4,6,1,5,7), nrow=3, byrow = TRUE)
Mtx_3
##      [,1] [,2] [,3]
## [1,]    1    3    5
## [2,]    2    4    6
## [3,]    1    5    7

Transpose the above matrix

Mtx_3T <- t(Mtx_3)
Mtx_3T
##      [,1] [,2] [,3]
## [1,]    1    2    1
## [2,]    3    4    5
## [3,]    5    6    7

Display the values of \(A^T\)A and A\(A^T\)

MTM <- Mtx_3T %*% Mtx_3
MTM
##      [,1] [,2] [,3]
## [1,]    6   16   24
## [2,]   16   50   74
## [3,]   24   74  110
MMT <- Mtx_3 %*% Mtx_3T  
MMT
##      [,1] [,2] [,3]
## [1,]   35   44   51
## [2,]   44   56   64
## [3,]   51   64   75

Comparing MTM and MMT

MTM == MMT
##       [,1]  [,2]  [,3]
## [1,] FALSE FALSE FALSE
## [2,] FALSE FALSE FALSE
## [3,] FALSE FALSE FALSE
  1. For a special type of square matrix A, we get \(A^T\)A = A\(A^T\) . Under what conditions could this be true?

Let’s verify it based on the hint provided using Identity matrix

Idm_3 <- diag(3)
Idm_3
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0    1    0
## [3,]    0    0    1

Transpose of Identity matrix and then multiply with the original Identity matrix

Idm_3T <- t(Idm_3)
Idm_3T
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0    1    0
## [3,]    0    0    1
MTM_I <- Idm_3T %*% Idm_3
MTM_I
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0    1    0
## [3,]    0    0    1
MMT_I <- Idm_3 %*% Idm_3T
MMT_I
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0    1    0
## [3,]    0    0    1

Comparing the above matrices

MTM_I == MMT_I
##      [,1] [,2] [,3]
## [1,] TRUE TRUE TRUE
## [2,] TRUE TRUE TRUE
## [3,] TRUE TRUE TRUE

Yes, if it is an Identity matrix, then the equation holds TRUE.

Now let’s try Diagonal Matrix

Dgm_3 <- matrix( c(1,0,0,0,5,0,0,0,7), nrow=3, byrow = TRUE)
Dgm_3
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0    5    0
## [3,]    0    0    7

Transpose the Diagonal Matrix

Dgm_3T <- t(Dgm_3)
Dgm_3T
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0    5    0
## [3,]    0    0    7

Repeating the calculations on the Diagomal matrix

DGMT_DGM <- Dgm_3T %*% Dgm_3
DGMT_DGM
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0   25    0
## [3,]    0    0   49
DGM_DGMT <- Dgm_3 %*% Dgm_3T
DGM_DGMT
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0   25    0
## [3,]    0    0   49
DGMT_DGM == DGM_DGMT
##      [,1] [,2] [,3]
## [1,] TRUE TRUE TRUE
## [2,] TRUE TRUE TRUE
## [3,] TRUE TRUE TRUE

The equation also holds TRUE for a Diagonal matrix. That means it will also hold true for a Scalar matrix.

Scalar_3 <- matrix( c(3,0,0,0,3,0,0,0,3), nrow=3, byrow = TRUE)
Scalar_3
##      [,1] [,2] [,3]
## [1,]    3    0    0
## [2,]    0    3    0
## [3,]    0    0    3

Problem Set 2

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.

Factorize a 3x3 (square) matrix into LU using R

library(matrixcalc)


function_LU <- function(M){
  #based on specs, only return if square matrix and don't worry about more than 5 x 5
  if(nrow(M)>= 5) {
    return("The matrix size should be less than 5.")
  }
  
  if(is.square.matrix(M) == FALSE) {
    return("Not a square matrix!")
  }
  
  #get number of rows from matrix
  n <- dim(M)[1]
  
  #create identity matrix L that's same size as M for elimination purposes
  L <- diag(n)
  
  #initiatlize variable
  f <- 0
  
  while (f <= dim(M)[1]){
     f = f+1
     i = f
     j = f
     while (i <= dim(M)[1]-1){
        c = M[i+1,j]/M[j,j]
        M[i+1,] = M[i+1,]-c*M[j,]
        L[i+1,j] = c
        i = 1+i
     }
  }
  return(L)
}

A <-  matrix( c(1,3,5,2,4,6,1,5,7), nrow=3, byrow = TRUE)
function_LU(A)
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    2    1    0
## [3,]    1   -1    1
luA <- lu.decomposition( A )
L <- luA$L
U <- luA$U
print( L )
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    2    1    0
## [3,]    1   -1    1
print( U )
##      [,1] [,2] [,3]
## [1,]    1    3    5
## [2,]    0   -2   -4
## [3,]    0    0   -2
print( L %*% U )
##      [,1] [,2] [,3]
## [1,]    1    3    5
## [2,]    2    4    6
## [3,]    1    5    7
print( A )
##      [,1] [,2] [,3]
## [1,]    1    3    5
## [2,]    2    4    6
## [3,]    1    5    7
B <- matrix( c( 2, -1, 3, -4, 6, 3, 4, -2, 8 ), nrow=3, byrow=TRUE )
function_LU(B)
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]   -2    1    0
## [3,]    2    0    1
luB <- lu.decomposition( B )
L <- luB$L
U <- luB$U
print( L )
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]   -2    1    0
## [3,]    2    0    1
print( U )
##      [,1] [,2] [,3]
## [1,]    2   -1    3
## [2,]    0    4    9
## [3,]    0    0    2
print( L %*% U )
##      [,1] [,2] [,3]
## [1,]    2   -1    3
## [2,]   -4    6    3
## [3,]    4   -2    8
print( B )
##      [,1] [,2] [,3]
## [1,]    2   -1    3
## [2,]   -4    6    3
## [3,]    4   -2    8