This file can be seen fully rendered on RPubs here.

Problem Set 1:

1. Show that \(A^TA\neq AA^T\)
# First create matrix A
A <- matrix(c(8,9,10,0,1,2,11,12,13), nrow = 3)
A
##      [,1] [,2] [,3]
## [1,]    8    0   11
## [2,]    9    1   12
## [3,]   10    2   13
# Then transpose matrix A using t()
t_A <- t(A)
t_A
##      [,1] [,2] [,3]
## [1,]    8    9   10
## [2,]    0    1    2
## [3,]   11   12   13
# Multiply the transposed matrix A by matrix A
ATA <- t_A %*% A
ATA
##      [,1] [,2] [,3]
## [1,]  245   29  326
## [2,]   29    5   38
## [3,]  326   38  434
# Multiply A by transposed matrix A (AA^T)
AAT <- A %*% t_A
AAT
##      [,1] [,2] [,3]
## [1,]  185  204  223
## [2,]  204  226  248
## [3,]  223  248  273

This shows that generally \(A^TA\neq AA^T\). You can check this using logical equality test:

# Check using logical comparison
ATA == AAT
##       [,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 \(A^TA = AA^T\). Under what conditions could this be true?

For the identity matrix, its multiples, and symetric matrices \(A^TA = AA^T\). This is true because a symmetric matrix and its transpose are equal, meaning you are just multiplying a matrix by itself.

# Show this using 2x2 identiy matrix
id <- matrix(c(1,0,0,1), nrow = 2)
id
##      [,1] [,2]
## [1,]    1    0
## [2,]    0    1
# Transpose the id matrix
t_id <- t(id)
# try both
id1 <- t_id %*% id
id2 <- id %*% t_id
# check using logic
id1 == id2
##      [,1] [,2]
## [1,] TRUE TRUE
## [2,] TRUE TRUE
3. Write an R function to factorize a square matrix \(A\) into LU or LDU

This function will factorize a 2x2 and a 3x3 matrix into LU

f_LU <- function(A) {

#Check dimensions of the matrix - 2x2 or 3x3

# 2x2
if (nrow(A) == 2) {

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

  L <- solve(L)
  
  print ("------ A ------")
  print (A)
  print ("------ L ------")
  print(L)
  print ("------ U ------")
  print(U)

  (L %*% U == A)

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

  #Subract multipler
  e_21 <- -A[2,1] / A[1,1]
  me_21 <- matrix(c(1,0,0,e_21,1,0,0,0,1), nrow=3)
  A2 <- me_21 %*% A

  #Subract multipler
  e_31 <- -A2[3,1] / A2[1,1]
  me_31 <- matrix(c(1,0,0,0,1,0,e_31,0,1), nrow=3)
  A3 <- me_31 %*% A2

  #Subtract multipler
  e_32 <- -A3[3,2] / A3[2,2]
  me_32 <- matrix(c(1,0,0,0,1,0,0,e_32,1), nrow=3)
  U <- me_32 %*% A3

  L <- solve(me_21) %*% solve(me_31) %*% solve(me_32)

  A <- L %*% U
  
  print ("------ A ------")
  print(A)
  print ("------ L ------")
  print(L)
  print ("------ U ------")
  print(U)

  (L %*% U == A)
 
  }
}

Check with a 2x2 and a 3x3 matrix

# 2x2
A2x2 <- matrix(c(2,3,4,5), nrow = 2)
f_LU(A2x2)
## [1] "------ A ------"
##      [,1] [,2]
## [1,]    2    4
## [2,]    3    5
## [1] "------ L ------"
##      [,1] [,2]
## [1,]    1  1.5
## [2,]    0  1.0
## [1] "------ U ------"
##      [,1] [,2]
## [1,] -2.5 -3.5
## [2,]  3.0  5.0
##      [,1] [,2]
## [1,] TRUE TRUE
## [2,] TRUE TRUE
# 3x3
A3x3 <- matrix(c(8,9,10,0,1,2,11,12,13), nrow = 3)
f_LU(A3x3)
## [1] "------ A ------"
##      [,1] [,2] [,3]
## [1,]    8    0   11
## [2,]    9    1   12
## [3,]   10    2   13
## [1] "------ L ------"
##      [,1]  [,2]      [,3]
## [1,]    1 1.125 -2.455882
## [2,]    0 1.000  2.000000
## [3,]    0 0.000  1.000000
## [1] "------ U ------"
##           [,1]      [,2]      [,3]
## [1,]  44.93382  8.286765  58.67647
## [2,] -11.00000 -3.000000 -14.00000
## [3,]  10.00000  2.000000  13.00000
##      [,1] [,2] [,3]
## [1,] TRUE TRUE TRUE
## [2,] TRUE TRUE TRUE
## [3,] TRUE TRUE TRUE