I would like to create a block diagonal matrix with the diagonal blocks being repeated a certain number of times, with the off-diagonal blocks being all zero matrices. For example, suppose we start with a matrix with the following:
> diag.matrix
[,1] [,2] [,3] [,4] [,5]
[1,] 1.0 0.5 0.5 0.5 0.5
[2,] 0.5 1.0 0.5 0.5 0.5
[3,] 0.5 0.5 1.0 0.5 0.5
[4,] 0.5 0.5 0.5 1.0 0.5
[5,] 0.5 0.5 0.5 0.5 1.0
I would like this matrix to be the diagonal block matrix so that in the end I have something like:
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
[1,] 1.0 0.5 0.5 0.5 0.5 0.0 0.0 0.0 0.0 0.0
[2,] 0.5 1.0 0.5 0.5 0.5 0.0 0.0 0.0 0.0 0.0
[3,] 0.5 0.5 1.0 0.5 0.5 0.0 0.0 0.0 0.0 0.0
[4,] 0.5 0.5 0.5 1.0 0.5 0.0 0.0 0.0 0.0 0.0
[5,] 0.5 0.5 0.5 0.5 1.0 0.0 0.0 0.0 0.0 0.0
[6,] 0.0 0.0 0.0 0.0 0.0 1.0 0.5 0.5 0.5 0.5
[7,] 0.0 0.0 0.0 0.0 0.0 0.5 1.0 0.5 0.5 0.5
[8,] 0.0 0.0 0.0 0.0 0.0 0.5 0.5 1.0 0.5 0.5
[9,] 0.0 0.0 0.0 0.0 0.0 0.5 0.5 0.5 1.0 0.5
[10,] 0.0 0.0 0.0 0.0 0.0 0.5 0.5 0.5 0.5 1.0
Here, we have the same block matrix above repeated twice in the block diagonals. If I wanted to do this efficiently an arbitrary number of times, is there a way to do it? thanks!
1) kronecker If M
is your matrix and k
is the number of times you want it repeated then:
kronecker(diag(k), M)
For example,
M <- matrix(0.5, 5, 5) + diag(0.5, 5)
k <- 2
kronecker(diag(k), M)
giving:
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
[1,] 1.0 0.5 0.5 0.5 0.5 0.0 0.0 0.0 0.0 0.0
[2,] 0.5 1.0 0.5 0.5 0.5 0.0 0.0 0.0 0.0 0.0
[3,] 0.5 0.5 1.0 0.5 0.5 0.0 0.0 0.0 0.0 0.0
[4,] 0.5 0.5 0.5 1.0 0.5 0.0 0.0 0.0 0.0 0.0
[5,] 0.5 0.5 0.5 0.5 1.0 0.0 0.0 0.0 0.0 0.0
[6,] 0.0 0.0 0.0 0.0 0.0 1.0 0.5 0.5 0.5 0.5
[7,] 0.0 0.0 0.0 0.0 0.0 0.5 1.0 0.5 0.5 0.5
[8,] 0.0 0.0 0.0 0.0 0.0 0.5 0.5 1.0 0.5 0.5
[9,] 0.0 0.0 0.0 0.0 0.0 0.5 0.5 0.5 1.0 0.5
[10,] 0.0 0.0 0.0 0.0 0.0 0.5 0.5 0.5 0.5 1.0
1a) %x% The last line of code can alternately be written as:
diag(k) %x% M
2) Matrix::bdiag Another possibility if you want to save space is to create a sparse matrix of class "dgMCatrix"
. It does not store the zero values. See ?bdiag
:
library(Matrix)
bdiag(replicate(k, M, simplify = FALSE))
giving:
10 x 10 sparse Matrix of class "dgCMatrix"
[1,] 1.0 0.5 0.5 0.5 0.5 . . . . .
[2,] 0.5 1.0 0.5 0.5 0.5 . . . . .
[3,] 0.5 0.5 1.0 0.5 0.5 . . . . .
[4,] 0.5 0.5 0.5 1.0 0.5 . . . . .
[5,] 0.5 0.5 0.5 0.5 1.0 . . . . .
[6,] . . . . . 1.0 0.5 0.5 0.5 0.5
[7,] . . . . . 0.5 1.0 0.5 0.5 0.5
[8,] . . . . . 0.5 0.5 1.0 0.5 0.5
[9,] . . . . . 0.5 0.5 0.5 1.0 0.5
[10,] . . . . . 0.5 0.5 0.5 0.5 1.0
2b) Diagonal or to create a sparse matrix of class "dgTMatrix"
:
Diagonal(k) %x% M
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With