1. Problem set 1
(1) Show that ATA=/ AAT in general. (Proof and demonstration.)
A<- matrix(c(1, 2, 3, 4, 5, 6, 7, 8, 9), nrow=3, ncol=3)
AT<- matrix(c(1, 4, 7, 2, 5, 8, 3, 6, 9), nrow=3, ncol=3)
A
## [,1] [,2] [,3]
## [1,] 1 4 7
## [2,] 2 5 8
## [3,] 3 6 9
AT
## [,1] [,2] [,3]
## [1,] 1 2 3
## [2,] 4 5 6
## [3,] 7 8 9
AAT =/ATA
A%*%AT
## [,1] [,2] [,3]
## [1,] 66 78 90
## [2,] 78 93 108
## [3,] 90 108 126
AT%*%A
## [,1] [,2] [,3]
## [1,] 14 32 50
## [2,] 32 77 122
## [3,] 50 122 194
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).
A square matrix in which all the main diagonal elements are 1’s and all the remaining elements are 0’s is called an Identity Matrix.
In this case AAT=ATA
A <- matrix(c(1, 2, 3, 4, 5, 6, 7, 8, 9), nrow=3, ncol=3)
AT <- matrix(c(1, 2, 3, 4, 5, 6, 7, 8, 9), nrow=3, ncol=3)
A
## [,1] [,2] [,3]
## [1,] 1 4 7
## [2,] 2 5 8
## [3,] 3 6 9
AT
## [,1] [,2] [,3]
## [1,] 1 4 7
## [2,] 2 5 8
## [3,] 3 6 9
A%*%AT
## [,1] [,2] [,3]
## [1,] 30 66 102
## [2,] 36 81 126
## [3,] 42 96 150
AT%*%A
## [,1] [,2] [,3]
## [1,] 30 66 102
## [2,] 36 81 126
## [3,] 42 96 150
2. Problem set 2 Write an R function to factorize a square matrix A into LU or LDU, whichever you prefer.You don’t have to worry about permuting rows of A and you can assume that A is less than 5x5, if you need to hard-code any variables in your code.
lud_function <- function(a,dia){
#initiate identity matrix
l = diag(dia)
u = a
#check if a[1,1] equal zero move it following row
for (i in 1:nrow(u)){
if (u[1,1] == 0){
if (i+1 <= nrow(u)){
swap = u[i+1,]
u[i+1,] = u[1,]
u[1,] = swap
if (u[1,1] != 0){
break
}
}
}
}
k = 0
#read columns
for (j in 1:ncol(a)){
k = k + 1
#read rows
for (i in 1:nrow(u)){
#if row > column and value for the row is not equal to zero
if (i > j){
if (u[i,j] != 0){
#initiate the sign of the factor
s = -1
if (u[i,j] < 0){
if (u[k,j] < 0){
#when both numbers are negative keep sign positive
s = 1
}
}
f = s * u[i,j]/u[k,j]
swap1 = s * (u[k,] * u[i,j])/u[k,j]
#get current row
swap2 = u[i,]
u[i,] = swap1 + swap2
l[i,j] = -1*f
}
}
}
}
return(list(l=l,u=u))
}
a = matrix(c(2,6,6,4,2, 3,1,3,2,1, 4,3,1,4,2, 1,1,2,7,4, 3,2,5,8,2), nrow=5, ncol=5)
m = lud_function(a,5)
#Lower Matrix
m$l
## [,1] [,2] [,3] [,4] [,5]
## [1,] 1 0.00 0.0000000 0.000000 0
## [2,] 3 1.00 0.0000000 0.000000 0
## [3,] 3 -0.75 1.0000000 0.000000 0
## [4,] 2 -0.50 -0.4788732 1.000000 0
## [5,] 1 -0.25 -0.2394366 0.678392 1
#Upper Matrix
m$u
## [,1] [,2] [,3] [,4] [,5]
## [1,] 2 3.000000 4.000000 1.000000e+00 3.0000000
## [2,] 0 -8.000000 -9.000000 -2.000000e+00 -7.0000000
## [3,] 0 -12.000000 -17.750000 -2.500000e+00 -9.2500000
## [4,] 0 -13.746479 -17.000000 2.802817e+00 -5.9295775
## [5,] 0 2.452261 3.032663 -2.220446e-16 -0.9422111
#L*U
m$l%*%m$u
## [,1] [,2] [,3] [,4] [,5]
## [1,] 2 3 4 1 3
## [2,] 6 1 3 1 2
## [3,] 6 3 1 2 5
## [4,] 4 2 4 7 8
## [5,] 2 1 2 4 2
a #The original matrix matches the lower * the upper
## [,1] [,2] [,3] [,4] [,5]
## [1,] 2 3 4 1 3
## [2,] 6 1 3 1 2
## [3,] 6 3 1 2 5
## [4,] 4 2 4 7 8
## [5,] 2 1 2 4 2