Problem set 1

A <- matrix(data = c(1,2,3,-1,0,4), nrow = 2, ncol = 3, byrow = TRUE)
A
##      [,1] [,2] [,3]
## [1,]    1    2    3
## [2,]   -1    0    4
X <- A %*% t(A)
Y <- t(A) %*% A
X
##      [,1] [,2]
## [1,]   14   11
## [2,]   11   17
Y
##      [,1] [,2] [,3]
## [1,]    2    2   -1
## [2,]    2    4    6
## [3,]   -1    6   25

Compute Eigenvalues and Eigenvectors of X and Y

eigen_X <- eigen(X) 
eigen_values_X <- round(eigen_X$values, 3) 
eigen_values_X
## [1] 26.602  4.398
eigen_Y <- eigen(Y) 
eigen_values_Y <- round(eigen_Y$values, 3) 
eigen_values_Y
## [1] 26.602  4.398  0.000
eigen_vectors_X <- eigen_X$vectors
eigen_vectors_X
##           [,1]       [,2]
## [1,] 0.6576043 -0.7533635
## [2,] 0.7533635  0.6576043
eigen_vectors_Y <- eigen_Y$vectors
eigen_vectors_Y
##             [,1]       [,2]       [,3]
## [1,] -0.01856629 -0.6727903  0.7396003
## [2,]  0.25499937 -0.7184510 -0.6471502
## [3,]  0.96676296  0.1765824  0.1849001

Compute left-singular, singular values, and right-singular vectors of A

svd_A <- svd(A)
left_singular_A <- svd_A$u
left_singular_A
##            [,1]       [,2]
## [1,] -0.6576043 -0.7533635
## [2,] -0.7533635  0.6576043
singular_values_A <- round(svd_A$d, 3)
singular_values_A
## [1] 5.158 2.097
right_singular_A <- svd_A$v
right_singular_A 
##             [,1]       [,2]
## [1,]  0.01856629 -0.6727903
## [2,] -0.25499937 -0.7184510
## [3,] -0.96676296  0.1765824

Examine the two sets of singular vectors and show that they are indeed eigenvectors of X and Y. In addition, the two non-zero eigenvalues (the 3rd value will be very close to zero, if not zero) of both X and Y are the same and are squares of the non-zero singular values of A.

round(svd(X)$d, 3)
## [1] 26.602  4.398
round(eigen_values_X, 3)
## [1] 26.602  4.398
round(svd(Y)$d, 3)
## [1] 26.602  4.398  0.000
round(eigen_values_Y, 3)
## [1] 26.602  4.398  0.000
squares <- round(singular_values_A ^ 2, 3)
squares
## [1] 26.605  4.397

Problem set 2

Using the procedure outlined in section 1 of the weekly handout, write a function to compute the inverse of a well-conditioned full-rank square matrix using co-factors. In order to compute the co-factors, you may use built-in commands to compute the determinant. Your function should have the following signature:

B = myinverse(A)

myinverse <- function(A){
  dimension_A <- nrow(A) 
  matrix <- matrix(NA, nrow = dimension_A, ncol = dimension_A) 
  for (i in 1:dimension_A) {
    for (j in 1:dimension_A) {
      matrix[i,j] <- det(A[-i,-j]) * (-1) ^ (i + j) 
    }
  }
  tranpose_matrix <- t(matrix) 
  inverse_matrix <- tranpose_matrix / det(A) 
  return(inverse_matrix)
}
A <- matrix(c(1,2,3,6,5,6,7,8,9), nrow = 3)
B <- round(myinverse(A), 5)
B 
##       [,1] [,2]     [,3]
## [1,] -0.25   -1  1.08333
## [2,]  0.50   -1  0.50000
## [3,] -0.25    1 -0.58333
round(A%*%B, 3)
##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0    1    0
## [3,]    0    0    1