## Daily Archives: August 26, 2014

### Matrices in R

The natural next step after looking at R vectors is to examine matrices. Although a matrix, being a two dimensional grid of values, may seem the natural choice for representing tables of data, it is actually better to use R’s data frame for that purpose as it contains many more functions for manipulating the data. It’s best to use a matrix primarily for those operations you would normally perform on a mathematical matrix, such as matrix multiplication, inversion and so on. Here’s a simple matrix

```> x = matrix(1:10, 2, 5)
> x
[,1] [,2] [,3] [,4] [,5]
[1,]    1    3    5    7    9
[2,]    2    4    6    8   10
```

The matrix x is composed of the integers from 1 to 10, arranged in 2 rows and 5 columns. By default, R fills a matrix column-wise. We can refer to a particular element in a matrix by giving its row and column index (again, remember that R indexes start at 1, not 0), so x[2,3] is 6. If we specify just one of these indexes, we get a vector containing either a single row or single column:

```> x[,2]    # Column 2
 3 4
> x[1,]    # Row 1
 1 3 5 7 9
```

It’s important to note that a vector is not a one dimensional matrix in R; they are quite different beasts. If you do want a one dimensional matrix you can get it by including the option ‘drop = FALSE’:

```> x[,2, drop=FALSE]
[,1]
[1,]    3
[2,]    4
> x[1,,drop=FALSE]
[,1] [,2] [,3] [,4] [,5]
[1,]    1    3    5    7    9
```

The first expression gives a 2×1 matrix, while the second gives a 1×5 matrix. The word ‘drop’ means to drop the dimensions not used in the selection. The default is ‘drop=TRUE’ meaning that when you select a single row or column, the matrix reduces from 2 to 1 dimension, becoming a vector. We can even produce a 1×1 matrix by isolating a single element as in x[2, 4, drop=FALSE].

Sometimes it’s more convenient if we can refer to a row or column by a specific name, rather than by its index. Names can be any string:

```> colnames(x) = c("A","B","C","D","E")
> rownames(x) = c("Alpha", "Beta")
> x
A B C D  E
Alpha 1 3 5 7  9
Beta  2 4 6 8 10
```

With these names, we can refer to x[“Alpha”, “B”] to get the value 3.

The t() function returns the transpose of a matrix:

```  Alpha Beta
A     1    2
B     3    4
C     5    6
D     7    8
E     9   10
```

Addition, subtraction and component-wise multiplication are done with the +, – and * operators as usual. To do ‘real’ matrix multiplication, use the %*% operator, which works only if the number of columns in its left operand is equal to the number of rows in the right operand. We can multiply a matrix by its transpose (or a function of its transpose) in either order:

```> sqrt(t(x)) %*% x
A         B        C        D        E
A 3.828427  8.656854 13.48528 18.31371 23.14214
B 5.732051 13.196152 20.66025 28.12436 35.58846
C 7.135047 16.506163 25.87728 35.24839 44.61951
D 8.302606 19.250962 30.19932 41.14768 52.09603
E 9.324555 21.649111 33.97367 46.29822 58.62278
> x %*% t(x)
Alpha Beta
Alpha   165  190
Beta    190  220
```

A matrix inverse is found with the solve() function. The inverse of a matrix is defined only if the matrix is square and its rows and columns are linearly independent. If either of these conditions is violated, R will give an error.

```> xtx = x %*% t(x)
> xtxinv = solve(xtx)
> xtxinv
Alpha   Beta
Alpha  1.10 -0.950
Beta  -0.95  0.825
> xtxinv %*% xtx
Alpha Beta
Alpha     1    0
Beta      0    1
```

We’ve done a check that the matrix multiplied by its inverse gives the identity matrix.