## Friday, September 7, 2012

### Matrix operations in R

# There are many ways of inputing matrices into R

# Cbind will bind vectors or other matrices together by adding on the columns of one to the other.
A = cbind(c(2,3),c(3,2),c(1,3))

# by the way, the c() function binds a set of elements seperated by commas into a vector.

A

# Rbind does the same as cbind but uses rows instead
B = cbind(c(-1, 1, 4), c(0, 1, 5), c(3, -1, 1))

B

# We can also create a matrix by using the matrix command.
# Data input in this manner is read from a vector into columns.

C = matrix(c(1,4,3,2),nrow = 2, ncol=2)

C

# An array can accompish the same effect as a matrix.

# The largest difference is that an array can take on more than two dimensions.

D = array(c(10,4,5,2), dim=c(2,2))

D

# One need not populate an matrix/array at creation when specifying matrix or array.

E = matrix( NA, nrow = 4, ncol = 4)

E
# We can see the matrix is filled with empty values

# We can replace individual elements by specifying their positions
E[1,1] = 3
E[2,1] = 0
E[3,1] = 4
E[4,1] = -1

# As well as entire columns or rows
E[,2] = c(0,2,3,2)

# Or subsection of the matrix with another matrix
E[,3:4] = cbind(c(2,3,2,1),c(-1,2,1,0))

E

# Sometimes you do not need to specify every element of the matrix if there are common elements
F = matrix(0, nrow=3, ncol=3)

F[1,1] = 4
F[2,2] = 3
F[3,3] = 6

F
# F is a diagnol matrix which is populated first with 0s then filled with the diagnol values.

# Sometimes we start with a data set and want to convert it to a matrix
G = data.frame(score1=c(5,10,-7),score2=c(-1,17,3),score3=c(0,5,2))

G

# Now let's convert it to a matrix
G = data.matrix(G)

G

# Now let's do some matrix operations with the matrices we have defined:

# a. C x D

C
D

C %*% D

# Which is different from element wise multiplication which is

C * D

# b. A x B

A
B

A %*%B

# c. B x A

# d. E x E' = E times the transpose of E

E
E %*% t(E)

# e. F^-1 or F inverse.  We know the inverse of F exists because it is a diagnol matrix with poisitive diagnol non-zero elements.

solve(F)

require(MASS)
ginv(F)

# Neither commands work for this particular application.

# However, it is easy to see that in a diagnol matrix the inverse is just.

# First we will make Finv the same as F
Finv = F

# This is an interesting bit of R functionality

# On the right the diag command is retrieving the diagnol of F and doing a element wise 1/x operation.

# The diag on the left is specifying target values to be replaced.
diag(Finv) = 1/diag(F)

# This might seem a little odd.  Let's try this kind of thing on Z
Z <- Finv

Finv

diag(Z) <- 23
Z
# In this case because 23 is a single number that can be duplicated throughout the diagnol of Z, there is no problem.

# We can check that this is really the inverse of F

Finv %*% F

# Or equally good

F %*% Finv

# Interestingly element by element multiplation yeilds the same result in this example

F * Finv

# f. Rank(D)

D
qr(D)$rank # F only has rank 1. This can be seen by dividing col 1 by col 2. D[,1]/D[,2] # Both numbers are two meaning item [1,1] = [1,2]*2 and [2,1]=[2,2]*2 # Interestingly this trick works for rows as well: D[1,]/D[2,] # g. B + G B G B + G # h. B - G B - G # Addition and subtraction is element by element in matrix notation # i. Rank(C) C qr(C)$rank

# C does not suffer from linearity problems
C[1,]/C[2,]

# j. -G
G
-G

# k. Trace(E)

# Trace is the sum of diagnol elements

E
sum(diag(E))

# l. Rank(F)

# Because F is a diagnol matrix (with no zero diagnol elements) the rank must be equal to the min of the dimensions, 3.
qr(F)\$rank

# kF, where k = -7

-7*F

# n. BF

B
F
B %*% F

# o. FB

F %*% B
# With matrices, BF != FB

# p. determinate of C or |C|
C
det(C)

# q. |D|
D
det(D)
# Which looks kind of funny but that is because R usings search algorithms to find the determinate and they are not exact.

# But calculating the determinate of a 2d matrix is easy:
D[1,1]*D[2,2]-D[1,2]*D[2,1]

# r. | CD |
C; D
H = C %*% D
det( H )
# Which also happens to be zero

H[1,1]*H[2,2]-H[1,2]*H[2,1]

# Sorry for the repetitive examples.  I had homework for a class and thought I might as well turn it into a post that someone might find useful.