Site de Jean-Michel RICHER

Maître de Conférences en Informatique à l'Université d'Angers

Ce site est en cours de reconstruction certains liens peuvent ne pas fonctionner ou certaines images peuvent ne pas s'afficher.


8.1. Operations on Matrices

8.1.1. Define a matrix

You can use:

  • diag() to define an identity matrix
  • or the matrix() function for which you can specify the number or rows and columns
> id <- diag(4)  # create identity matrix of size 4x4
> id
     [,1] [,2] [,3] [,4]
[1,]    1    0    0    0
[2,]    0    1    0    0
[3,]    0    0    1    0
[4,]    0    0    0    1
 
> m1 <- matrix(1:10, nrow=2, ncol=5)
> m1
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    3    5    7    9
[2,]    2    4    6    8   10
> m2 <- matrix(c(1,-1,2,-2,3,-3), nrow=2, ncol=3)
> m2
     [,1] [,2] [,3]
[1,]    1    2    3
[2,]   -1   -2   -3
 
> print(readLines("matrix.Rdata"))
[1] "1  2  3  4  5  6  7  8  9 10 11 12"
> m3 <- matrix(scan("matrix.Rdata"), nrow=4, ncol=3)
Read 12 items
> m3
     [,1] [,2] [,3]
[1,]    1    5    9
[2,]    2    6   10
[3,]    3    7   11
[4,]    4    8   12

8.1.2. cbind, rbind

Combine matrices horizontally (cbind) or vertically (rbind)

> x1 <- c(1,2,3)
> x2 <- c(4,5,6)
> cbind(x1,x2)
     x1 x2
[1,]  1  4
[2,]  2  5
[3,]  3  6
> rbind(x1,x2)
   [,1] [,2] [,3]
x1    1    2    3
x2    4    5    6
 
> m <- matrix(10:15, nrow=3, ncol=2)
> m
     [,1] [,2]
[1,]   10   13
[2,]   11   14
[3,]   12   15
> cbind(m,x1)
           x1
[1,] 10 13  1
[2,] 11 14  2
[3,] 12 15  3
 

8.1.3. Sum, difference, product, division

Use +, -, %*% and %/%.

> m1 <- matrix(1:16, nrow=4, ncol=4)
> m2 <- matrix(16:1, nrow=4, ncol=4)
> m1 + m2
     [,1] [,2] [,3] [,4]
[1,]   17   17   17   17
[2,]   17   17   17   17
[3,]   17   17   17   17
[4,]   17   17   17   17
> m1 - m2
     [,1] [,2] [,3] [,4]
[1,]  -15   -7    1    9
[2,]  -13   -5    3   11
[3,]  -11   -3    5   13
[4,]   -9   -1    7   15
> m1 %*% m2
     [,1] [,2] [,3] [,4]
[1,]  386  274  162   50
[2,]  444  316  188   60
[3,]  502  358  214   70
[4,]  560  400  240   80
> m1 %/% m2
     [,1] [,2] [,3] [,4]
[1,]    0    0    1    3
[2,]    0    0    1    4
[3,]    0    0    1    7
[4,]    0    0    2   16

8.1.4. Transpose

Use the t() function.

> m1 <- matrix(1:15, nrow=3, ncol=5)
> m1
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    4    7   10   13
[2,]    2    5    8   11   14
[3,]    3    6    9   12   15
> t(m1)
     [,1] [,2] [,3]
[1,]    1    2    3
[2,]    4    5    6
[3,]    7    8    9
[4,]   10   11   12
[5,]   13   14   15

8.1.5. rowSums, rowMeans, colSums, colMeans

> m1 <- matrix(1:15, nrow=3, ncol=5)
> m1
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    4    7   10   13
[2,]    2    5    8   11   14
[3,]    3    6    9   12   15
> rowSums(m1)
[1] 35 40 45
> colSums(m1)
[1]  6 15 24 33 42
> rowMeans(m1)
[1] 7 8 9

8.1.6. Equations

We want to solve the following equation $y = m * x$ where $x$ and $y$ are vectors and $m$ is a matrix. We only know $m$ and $y$ and we want to find $x$. For this we can use the solve() function:

$[\table -8;13;31;] = [\table 1, 2, -1; 3, 4, 7; -5, 1, 6 ] × [\table x_1; x_2; x_3]$

> m <- matrix(c(1, 2, -1, 3, 4, 7, -5, 1, 6), nrow=3, ncol=3)
> y <- c(-8,13,31)
> x <- solve(m,y)
     [,1]
[1,]    1
[2,]    2
[3,]    3
 
> x
[1] 1 2 3

or we could compute $x = m^{-1} × y$, where $m^{-1}$ is the inverse of $m$. To do this you will need the ginv function from the MASS package.

# do this once !
> install.packages('MASS')
...
> inv_m <- MASS::ginv(m)  # compute inverse of matrix m
> inv_m
           [,1]         [,2]        [,3]
[1,] -0.1517857  0.473214286 -0.20535714
[2,]  0.1160714 -0.008928571  0.09821429
[3,] -0.1607143  0.089285714  0.01785714
>  x <- inv_m %*% y
> x
     [,1]
[1,]    1
[2,]    2
[3,]    3