1. Problem set 1

(1) Show that A^T A # AA^T in general. (Proof and demonstration.)

# Sample matrix
A <- matrix(c(1,0,1,3,1,0,6,4,0), nrow = 3, byrow = T)
A
##      [,1] [,2] [,3]
## [1,]    1    0    1
## [2,]    3    1    0
## [3,]    6    4    0
# Transpose of matrix A
t(A)
##      [,1] [,2] [,3]
## [1,]    1    3    6
## [2,]    0    1    4
## [3,]    1    0    0
# Transposed matrix A multiplied by  matrix A 
TAA <- t(A) %*% A

TAA
##      [,1] [,2] [,3]
## [1,]   46   27    1
## [2,]   27   17    0
## [3,]    1    0    1
# Matrix A multiplied by transposed matrix A (A^T)
ATA <- A %*% t(A)
ATA
##      [,1] [,2] [,3]
## [1,]    2    3    6
## [2,]    3   10   22
## [3,]    6   22   52
# So, A^T.A != A.A^T.

TAA == ATA
##       [,1]  [,2]  [,3]
## [1,] FALSE FALSE FALSE
## [2,] FALSE FALSE FALSE
## [3,] FALSE FALSE FALSE

(2) For a special type of square matrix A, we get AT A = AAT . Under what conditions

could this be true? (Hint: The Identity matrix I is an example of such a matrix).

# Sample matrix
A <- matrix(c(5,2,4,2,0,2,4,2,5), nrow = 3, byrow = T)
A
##      [,1] [,2] [,3]
## [1,]    5    2    4
## [2,]    2    0    2
## [3,]    4    2    5
# Transpose of matrix A
t(A)
##      [,1] [,2] [,3]
## [1,]    5    2    4
## [2,]    2    0    2
## [3,]    4    2    5
# Transposed matrix A multiplied by  matrix A 
TAA <- t(A) %*% A

TAA
##      [,1] [,2] [,3]
## [1,]   45   18   44
## [2,]   18    8   18
## [3,]   44   18   45
# Matrix A multiplied by transposed matrix A (A^T)
ATA <- A %*% t(A)
ATA
##      [,1] [,2] [,3]
## [1,]   45   18   44
## [2,]   18    8   18
## [3,]   44   18   45
# So, A^T.A != A.A^T.

TAA == ATA
##      [,1] [,2] [,3]
## [1,] TRUE TRUE TRUE
## [2,] TRUE TRUE TRUE
## [3,] TRUE TRUE TRUE
#Since the example is an identity matrix, the condition is TRUE

2. Problem set 2

Write an R function to factorize a square matrix A into LU or LDU

#Declare the function
factorize <- function(A) {

#Check size of matrix - first 2x2, 3x3
if (nrow(A) == 2) {

  #2x2 matrix

  #Subract multipler R21
  R21 <- -A[2,1] / A[1,1]
  L <- matrix(c(1,0,R21,1),nrow=2, byrow = T)
  U <- L %*% A
  #U

  print("###A###")
  print(A)
  print("###L###")
  print(L)
  print("###U###")
  print(U)

  L <- solve(L)
  #Inverse of L
  print("###Inverse of L ###")
  print(L)
  (L %*% U == A)

} else if (nrow(A) == 3) {

  #3x3 matrix

  #Subract multipler R21
  R21 <- -A[2,1] / A[1,1]
  mR21 <- matrix(c(1,0,0,R21,1,0,0,0,1), nrow=3, byrow = T)
  A2 <- mR21 %*% A
  #A2

  #Subract multipler E31
  R31 <- -A2[3,1] / A2[1,1]
  mR31 <- matrix(c(1,0,0,0,1,0,R31,0,1), nrow=3, byrow = T)
  A3 <- mR31 %*% A2
  #A3

  #Subtract multipler E32
  R32 <- -A3[3,2] / A3[2,2]
  mR32 <- matrix(c(1,0,0,0,1,0,0,R32,1), nrow=3, byrow = T)
  U <- mR32 %*% A3
  
  L <- solve(mR21) %*% solve(mR31) %*% solve(mR32)
  A <- L %*% U
  

  print("###A###")
  print(A)
  print("###L###")
  print(L)
  print("###U###")
  print(U)
  
  
  (L %*% U == A)

}}

#Check with 2x2 matrix

A <- matrix(c(2,3,6,8),nrow=2, byrow = T)
factorize(A)
## [1] "###A###"
##      [,1] [,2]
## [1,]    2    3
## [2,]    6    8
## [1] "###L###"
##      [,1] [,2]
## [1,]    1    0
## [2,]   -3    1
## [1] "###U###"
##      [,1] [,2]
## [1,]    2    3
## [2,]    0   -1
## [1] "###Inverse of L ###"
##      [,1] [,2]
## [1,]    1    0
## [2,]    3    1
##      [,1] [,2]
## [1,] TRUE TRUE
## [2,] TRUE TRUE
#Check with 3x3 matrix

A <- matrix(c(2,3,6,8,2,1,0,3,1),nrow=3, byrow = T)
factorize(A)
## [1] "###A###"
##      [,1] [,2] [,3]
## [1,]    2    3    6
## [2,]    8    2    1
## [3,]    0    3    1
## [1] "###L###"
##      [,1] [,2] [,3]
## [1,]    1  0.0    0
## [2,]    4  1.0    0
## [3,]    0 -0.3    1
## [1] "###U###"
##      [,1] [,2]  [,3]
## [1,]    2    3   6.0
## [2,]    0  -10 -23.0
## [3,]    0    0  -5.9
##      [,1] [,2] [,3]
## [1,] TRUE TRUE TRUE
## [2,] TRUE TRUE TRUE
## [3,] TRUE TRUE TRUE