Problem 1

x <- matrix(c(9,2,6,5,8,12,8,6,4,10,3,4,0,2,1),5,3)
colnames(x) <- c("x1", "x2", "x3")

1.A.

x_bar

x_bar <- apply(x,2,mean)
x_bar
## x1 x2 x3 
##  6  8  2

1.b.

Sn, the symmetric matrix The diagonal we have already found [6 8 2] S12 = S21 S13 = S31 S23 = S32

s12 <- (1/5) * ((9-6)*(12-8) + (2-6)*(8-8) + (6-6)*(6-8) + (5-6)*(4-8) + (8-6)*(10-8))
s13 <- (1/5) * ((9-6)*(3-2) + (2-6)*(4-2) + (6-6)*(0-2) + (5-6)*(2-2) + (8-6)*(1-2))
s23 <- (1/5) * ((12-8)*(3-2) + (8-8)*(4-2) + (6-8)*(0-2) + (4-8)*(2-2) + (10-8)*(1-2))
Sn <- matrix(c(6,s12,s13,s12,8,s23,s13,s23,2),3,3)
Sn
##      [,1] [,2] [,3]
## [1,]  6.0  4.0 -1.4
## [2,]  4.0  8.0  1.2
## [3,] -1.4  1.2  2.0

1.c.

Rn r12 = r21 r13 = r31 r23 = r32

r12 <- s12/(sqrt(6)*sqrt(8)) 
r13 <- s13/(sqrt(6)*sqrt(2)) 
r23 <- s23/(sqrt(8)*sqrt(2)) 
R <- matrix(c(1,r12,r13,r12,1,r23,r13,r23,1),3,3)
R
##            [,1]      [,2]       [,3]
## [1,]  1.0000000 0.5773503 -0.4041452
## [2,]  0.5773503 1.0000000  0.3000000
## [3,] -0.4041452 0.3000000  1.0000000

Problem 2

A <- matrix(c(-1,4,3,2),2,2)

B <- matrix(c(4,1,-2,-3,2,0),3,2)

2.a.

5A

5*A
##      [,1] [,2]
## [1,]   -5   15
## [2,]   20   10

2.b.

BA

B %*% A
##      [,1] [,2]
## [1,]  -16    6
## [2,]    7    7
## [3,]    2   -6

2.c.

A’B’ Yes, A’B’ = (BA)’

t(A) %*% t(B)
##      [,1] [,2] [,3]
## [1,]  -16    7    2
## [2,]    6    7   -6

2.d.

AB is undefined because A = 2x2 and B = 2x3. The outter values must match in order to perform matrix multiplication.

Problem 3

A <- matrix(c(5,2,3,4,-3,7,4,1,2),3,3)
B <- matrix(c(1,0,1,0,1,2,1,0,3),3,3)

3.a.

tr(A) and tr(B) - flip along a mirror down the diagonal

sum(diag(A))
## [1] 4
sum(diag(B))
## [1] 5

3.b.

tr(A + B) Yes, tr(A + B) DOES = tr(A) + tr(B)

sum(diag(A) + diag(B))
## [1] 9
sum(diag(A)) + sum(diag(B))
## [1] 9

3.c.

|A| and |B|

det(A)
## [1] 23
det(B)
## [1] 2

3.d.

|AB| Yes, |AB| DOES = |A|*|B|

det(A%*%B)
## [1] 46
det(A)%*%det(B)
##      [,1]
## [1,]   46

3.e.

rank of A is 3. rank of B is 3.

x <-qr(A)
x$rank
## [1] 3
y <-qr(B)
y$rank
## [1] 3

Problem 4

A <- matrix(c(1,2,1,-1,1,-1,1,0,-1),3,3)

4.a.

Yes, since the dot product between all the columns is 0, the matrix can be said to be mutually orthogonal.

A1 <- c(1,2,1)
A2 <- c(-1,1,-1)
A3 <- c(1,0,-1)

A1%*%A2
##      [,1]
## [1,]    0
A2%*%A1
##      [,1]
## [1,]    0
A1%*%A3
##      [,1]
## [1,]    0
A3%*%A1
##      [,1]
## [1,]    0
A3%*%A2
##      [,1]
## [1,]    0
A2%*%A3
##      [,1]
## [1,]    0

4.b.

aa <- A1 / sqrt(sum(A1^2)) 
ab <- A2 / sqrt(sum(A2^2))
ac <- A3 / sqrt(sum(A3^2))

C <- as.matrix(data.frame(aa, ab, ac))

C
##             aa         ab         ac
## [1,] 0.4082483 -0.5773503  0.7071068
## [2,] 0.8164966  0.5773503  0.0000000
## [3,] 0.4082483 -0.5773503 -0.7071068

4.c.

Yes, CC’ = C’C = I so C is an orthogonal matrix.

C%*%t(C) #rounding issue here, but this is essentially = I 
##              [,1] [,2]         [,3]
## [1,] 1.000000e+00    0 2.220446e-16
## [2,] 0.000000e+00    1 0.000000e+00
## [3,] 2.220446e-16    0 1.000000e+00
t(C)%*%C
##    aa ab ac
## aa  1  0  0
## ab  0  1  0
## ac  0  0  1

Problem 5

A <- matrix(c(1,2,2,-2),2,2)
A
##      [,1] [,2]
## [1,]    1    2
## [2,]    2   -2

5.a.

Eigen values and vectors shown below

x <-eigen(A)
x$values
## [1]  2 -3
x$vectors
##            [,1]       [,2]
## [1,] -0.8944272 -0.4472136
## [2,] -0.4472136  0.8944272

5.b.

A^-1

solve(A)
##           [,1]       [,2]
## [1,] 0.3333333  0.3333333
## [2,] 0.3333333 -0.1666667

5.c.

tr(A) = lambda1 + lambda2 = 2 + -3 = -1 this is true.

sum(diag(A))
## [1] -1

5.d.

|A| = lambda1lambda2 = 2-3 = -6

det(A)
## [1] -6

5.e.

A is not positive definite. All eigen values must be positive.

library(matrixcalc)
is.positive.definite(A, tol=1e-8)
## [1] FALSE

Problem 6

A <- matrix(c(9,-2,-2,6),2,2)
A
##      [,1] [,2]
## [1,]    9   -2
## [2,]   -2    6

6.a.

eigen(A)$values
## [1] 10  5

6.b.

library(expm)
## Loading required package: Matrix
## 
## Attaching package: 'expm'
## The following object is masked from 'package:Matrix':
## 
##     expm
a_12 <- sqrtm(A)
a_12
##            [,1]       [,2]
## [1,]  2.9770357 -0.3704839
## [2,] -0.3704839  2.4213099

6.c.

sqrtm(A)%*%sqrtm(A)
##      [,1] [,2]
## [1,]    9   -2
## [2,]   -2    6

6.d.

"%^%" <- function(S, power) 
   with(eigen(S), vectors %*% (values^power * t(vectors))) 
a_neg12 <- A%^%(-0.5)
a_neg12
##            [,1]       [,2]
## [1,] 0.34242493 0.05239433
## [2,] 0.05239433 0.42101643

6.e.

As it is shown below, \(A^{1/2}A^{-1/2} = A^{-1/2}A^{1/2} = I\)

a_12 %*% a_neg12 
##               [,1]          [,2]
## [1,]  1.000000e+00 -2.775558e-17
## [2,] -1.387779e-16  1.000000e+00
a_neg12 %*% a_12
##               [,1]          [,2]
## [1,]  1.000000e+00 -5.551115e-17
## [2,] -1.110223e-16  1.000000e+00

7.a.

A symmetric matrix is equal to its transpose. So, we take the transpose of A’A.

(A’A)‘= A’(A’)’ = A’A

7.b.

y = Ax y’y = x’A’Ax

First, we know A’A is symmetric.

Second, take \(x \in \mathbb{R}^{nx1}\). We wish to prove that \(x'A'Ax \geq 0\).

(x’A’Ax)‘= x’A’(A’)‘(x’)‘= x’A’Ax = (Ax)’(Ax)

Now we notice that (Ax) \(\in \mathbb{R}^{nx1}\).

The dot product of \((Ax)'(Ax) = (Ax) \dot (Ax) = ||Ax|| \geq 0\)

Thus, A’A is nonnegative definite.