library(pracma)
## Warning: package 'pracma' was built under R version 3.4.4
In this problem, we’ll verify using R that SVD and Eigenvalues are related as worked out in the weekly module. Given a 2x3 matrix \(A\)
\[\mathbf{A} = \left[\begin{array} {rrr} 1 & 2 & 3\\ -1 & 0 & 4 \end{array}\right]\]
Write code in R to compute \(X=AA^{T}\) and \(Y=A^{T}A\). Then, compute the eigenvalues and eigenvectors of \(X\) and \(Y\) using the built-in commands in R.
First, we create matrix \(A\)
cat("2x3 Matrix A","\n")
## 2x3 Matrix A
(A <- matrix(c(1,2,3,
-1,0,4),nrow=2, byrow=TRUE))
## [,1] [,2] [,3]
## [1,] 1 2 3
## [2,] -1 0 4
First we multiply matrix \(A\) to matrix \(A\) transpose and assign the values to matrix \(X\)
(X = A%*%t(A))
## [,1] [,2]
## [1,] 14 11
## [2,] 11 17
Using the eigen function in base R we will obtain the eigenvalues and eigenvectors of matrix \(X\).
eigen(X)
## eigen() decomposition
## $values
## [1] 26.601802 4.398198
##
## $vectors
## [,1] [,2]
## [1,] 0.6576043 -0.7533635
## [2,] 0.7533635 0.6576043
Next we compute the following equation \(Y=A^{T}A\) which assigns matrix \(A\) transpose multiplied by matrix \(A\) to matrix \(Y\).
(Y = t(A)%*%A)
## [,1] [,2] [,3]
## [1,] 2 2 -1
## [2,] 2 4 6
## [3,] -1 6 25
Following suit with matrix \(X\), we will now obtain the eigenvalues and eigenvectors of matrix \(Y\).
eigen(Y)
## eigen() decomposition
## $values
## [1] 2.660180e+01 4.398198e+00 -6.098625e-16
##
## $vectors
## [,1] [,2] [,3]
## [1,] -0.01856629 0.6727903 0.7396003
## [2,] 0.25499937 0.7184510 -0.6471502
## [3,] 0.96676296 -0.1765824 0.1849001
Then, compute the left-singular, singular values, and right-singular vectors of \(A\) using the svd command. 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 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\). Your code should compute all these vectors and scalars and store them in variables.
Using the svd function we find the following:
d a vector containing the singular values of matrix \(A\)
u a matrix whose columns contain the left singular vectors of matrix \(A\), as one can see the eigenvectors for matrix \(X\) matches this output.
v a matrix whose columns contain the right singular vectors of matrix \(A\), similar to the findings of matrix \(X\) the eigenvectors for matrix \(Y\) match those from this output.
svd(A)
## $d
## [1] 5.157693 2.097188
##
## $u
## [,1] [,2]
## [1,] -0.6576043 -0.7533635
## [2,] -0.7533635 0.6576043
##
## $v
## [,1] [,2]
## [1,] 0.01856629 -0.6727903
## [2,] -0.25499937 -0.7184510
## [3,] -0.96676296 0.1765824
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)\) where \(A\) is a matrix and \(B\) is its inverse and \(AxB=I\). The off-diagonal elements of \(I\) should be close to zero, if not zero. Likewise, the diagonal elements should be close to 1, if not 1. Small numerical precision errors are acceptable but the function myinverse should be correct and must use co-factors and determinant of \(A\) to compute the inverse.
First we create a full rank Matrix named matrix \(A\)
cat("Matrix A","\n")
## Matrix A
(A <- matrix(c(1,2,3,
-1,0,1,
0,1,-2),nrow=3, byrow=TRUE))
## [,1] [,2] [,3]
## [1,] 1 2 3
## [2,] -1 0 1
## [3,] 0 1 -2
cat("\n","The Rank of Matrix A","\n")
##
## The Rank of Matrix A
qr(A)$rank
## [1] 3
cat("\n","Inverse of Matrix A","\n")
##
## Inverse of Matrix A
inv(A)
## [,1] [,2] [,3]
## [1,] 0.125 -0.875 -0.25
## [2,] 0.250 0.250 0.50
## [3,] 0.125 0.125 -0.25
The formula to find the inverse of this matrix is \(A^{-1}=C^{T}/|A|\) which is the inverse of matrix \(A\) being equaled to the co-factor matrix transpose which is divided by the determinant of matrix \(A\).
According to the Inverse of Matrix - Singular Value Decomposition lecture, a co-factor is the determinant of the min matrix sub-matrix.
myinverse=function(A){
#remove the first row of matrix A
C1 <- A[-1,]
#remove the first column of matrix A
C <- C1[,-1]
#transpose the sub matrix
t(C)
#run the determinant of matrix A
det(A)
#divide matrix C transpose by the determinant of matrix A
t(C)/det(A)}
After running my function, it seems that some of the variables I arrived are are true. Yet the inv feature from the pracma provided output for the entire 3x3 matrix while following the lecture we were instructed to calculate the inverse based off using the co-factors which provided a 2x2 matrix.
A <- matrix(c(1,2,3,
-1,0,1,
0,1,-2),nrow=3, byrow=TRUE)
myinverse(A)
## [,1] [,2]
## [1,] 0.000 -0.125
## [2,] -0.125 0.250