"Diagonalize" each row of a matrix

Question

I have an n x p matrix that looks like this:

n = 100
p = 10    
x <- matrix(sample(c(0,1), size = p*n, replace = TRUE), n, p)

I want to create an n x p x p array A whose kth item along the 1st dimension is a p x p diagonal matrix containing the elements of x[k,]. What is the most efficient way to do this in R? I'm looking for a way that uses outer (or some other vectorized approach) rather than one of the apply functions.

Solution using lapply:

A <- aperm(simplify2array(lapply(1:nrow(x), function(i) diag(x[i,]))), c(3,2,1))

I'm looking for something more efficient than this.

Thanks.

lmo lmo · Accepted Answer · 2016-07-28T20:38:01

As a starting point, here is a humble for loop method with pre-allocation of the matrix.

# pre-allocate matrix of desired size
myArray <- array(0, dim=c(ncol(x), ncol(x), nrow(x)))
# fill in array
for(i in seq_len(nrow(x))) myArray[,,i] <- diag(x[i,])

It should run relatively fast. On my machine, for a 1000 X 100 matrix, the lapply method took 0.87 seconds, while the for loop (including the array pre-allocation) took 0.25 seconds to transform the matrix into to your desired array. So the for loop was about 3.5 times faster.

transpose your original matrix

Note also that row operations on R matrices tend to be slower than column operations. This is because matrices are stored in memory by column. If you transpose your matrix, and perform the operation this way, the time to complete the operation on 100X1000 matrix drops to 0.14, half that of the first for loop, and 7 times faster than the lapply method.

"Diagonalize" each row of a matrix

1 Answers