Sweep in matrix operations in R not working as intended

Question

I am currently trying to efficiently implement a new square matrix, M, as a function of two square matrices of the same dimension, W and V, as follows: F_ij = W_ii·V_ji. I.e., the i^th row of F is the i^th diagonal element of W times the i^th column of V.

There is a function in R, sweep, that allows one to sweep a statistic across a margin of an array. I.e., doing

mean.att <- apply(attitude, 2, mean)
sweep(data.matrix(attitude), 2, mean.att)

will de-mean the columns (2^nd dimension) of the attitude matrix (taken from the help). Also, one can supply a function, like FUN="*". In this case, the columns of attitude matrix will be multiplied by their respective medians.

So the intended code to produce the matrix F would be

Fij <- matrix(NA, ncol=N, nrow=N)
for (i in 1:N) {
  Fij[i, ] <- w[i, i] * v[, i]
}

Apparently, since the forte of R is vectorisation, this can be done with a sweep operation: I want to multiply the columns (FUN="*") of V by the diagonal of W, i.e.

Fij2 <- sweep(v, 2, diag(w), FUN="*")

However, whenever I check if Fij==Fij2, they are not!

MWE:

set.seed(1)
w <- matrix(rnorm(16), nrow=4)
v <- matrix(rnorm(16), nrow=4)
Fij <- matrix(NA, ncol=4, nrow=4)
for (i in 1:4) { # This loop can and should be vectorised
  Fij[i, ] <- w[i, i] * v[, i]
}
Fij2 <- sweep(v, 2, diag(w), FUN="*")
Fij
Fij2

The diagonal elements are equal, but the off-diagonal ones are not.

I should be very grateful if someone clarified which of the implementations of F_ij is false!

Ralf Stubner Ralf Stubner · Accepted Answer · 2018-02-20T18:17:47

The implementation using sweep gives the transposed result. You can see this easily if you use simpler input matrices:

w <- matrix(1:4, nrow = 2)
v <- matrix(5:8, nrow = 2)
Fij <- matrix(NA, ncol=2, nrow=2)
for (i in 1:2) {
  Fij[i, ] <- w[i, i] * v[, i]
}
Fij2 <- sweep(v, 2, diag(w), FUN="*")
Fij
#>      [,1] [,2]
#> [1,]    5    6
#> [2,]   28   32
Fij2
#>      [,1] [,2]
#> [1,]    5   28
#> [2,]    6   32
t(Fij2)
#>      [,1] [,2]
#> [1,]    5    6
#> [2,]   28   32

So Fij matches your description, as does t(Fij2).

Sweep in matrix operations in R not working as intended

2 Answers