Problem Set 1

  1. Show that \(A^{T} A ≠ A A^{T}\) in general. (Proof and demonstration).
Demonstration
library(magrittr)
## Warning: package 'magrittr' was built under R version 4.3.2
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.3.2
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
A <- matrix(seq(from=1,to=6), nrow=2, byrow=T)
A %*% t(A)  # 2 x 2 matrix
##      [,1] [,2]
## [1,]   14   32
## [2,]   32   77
t(A) %*% A  # 3 x 3 matrix
##      [,1] [,2] [,3]
## [1,]   17   22   27
## [2,]   22   29   36
## [3,]   27   36   45
  1. For a special type of square matrix A, we get \(A^{T} A = A A^{T}\) . Under what conditions could this be true? (Hint: The Identity matrix I is an example of such a matrix).
  • The only way this could be true is if the dimensions of matrix \(A\) are the same- meaning it has dimensions n * n. Then \(A^{T} A\) will have dimensions n * n, as will \(A A^{T}\)
  • Additionally, the transpose of \(A\) has to have the elements symmetrically reflected over the diagonal. That is, \(A\)i,j = \(A\)j,i. See demonstration below;
B <- matrix(seq(from=1,to=4), nrow=2, byrow=T)

B %*% t(B)  
##      [,1] [,2]
## [1,]    5   11
## [2,]   11   25
t(B) %*% B 
##      [,1] [,2]
## [1,]   10   14
## [2,]   14   20
# The dimensions are the same, but the products are different. Let's try with a different matrix

#--------------------------------------------------------------------------------
C <- matrix(c(2,5,5,2), nrow=2, byrow=T)
C %*% t(C) 
##      [,1] [,2]
## [1,]   29   20
## [2,]   20   29
t(C) %*% C 
##      [,1] [,2]
## [1,]   29   20
## [2,]   20   29
t(C) %*% C == C %*% t(C) # Gives true
##      [,1] [,2]
## [1,] TRUE TRUE
## [2,] TRUE TRUE

Problem Set 2

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.

# Let's perform LU decomposition
library(matrixcalc)

LU_matrix <- function(matrix){
  
  # Check if input matrix is actually a matrix
  if(is.matrix(matrix) != TRUE){
    message('Input must be a  matrix.')
    return(NA) 
  }
  
  # Check if the input matrix is square dimensional
  if(is.square.matrix(matrix) == FALSE){
  return('Matrix is not square')
  }

  # Derive L and U 
  n <- nrow(matrix)   # Number of rows in the matrix 
  U <- matrix 
  # Build lower triangular matrix with 0's below the diagonal
  L <- diag(n)   

  for (j in c(1:n)){
    for(i in c(2:n)){
      if(i > j){
        row <- U[j,]
        multiplier <- U[i, j] / row[j] # Gaussian multiplier from video
        U[i,] <- U[i,] - (multiplier * row)  
        L[i,j] <- multiplier
      }
    }
  }
  
  return(list(L = L, U = U))
}



# Test function on sample matrix
D <- matrix(c(5,3,8,2,4,6,1,9,2), nrow=3, byrow=TRUE)
D
##      [,1] [,2] [,3]
## [1,]    5    3    8
## [2,]    2    4    6
## [3,]    1    9    2
LU_matrix(D)
## $L
##      [,1] [,2] [,3]
## [1,]  1.0    0    0
## [2,]  0.4    1    0
## [3,]  0.2    3    1
## 
## $U
##      [,1] [,2] [,3]
## [1,]    5  3.0  8.0
## [2,]    0  2.8  2.8
## [3,]    0  0.0 -8.0
library(matrixcalc)

D <- matrix(c(5,3,8,2,4,6,1,9,2), nrow=3, byrow=TRUE)
lu.decomposition(D)
## $L
##      [,1] [,2] [,3]
## [1,]  1.0    0    0
## [2,]  0.4    1    0
## [3,]  0.2    3    1
## 
## $U
##      [,1] [,2] [,3]
## [1,]    5  3.0  8.0
## [2,]    0  2.8  2.8
## [3,]    0  0.0 -8.0