user34790
user34790

Reputation: 2036

Create a block circulant matrix in R

I am trying to create a block circulant matrix in R. The structure of a block circulant matrix is given below.

C0 C1 ... Cn-1
Cn-1 C0 C1 ... Cn-2
Cn-2 Cn-1 .... Cn-3

and so on

I have the blocks

C0 .... Cn-1

What is the easiest way to create the matrix. Is there a function already available?

Upvotes: 5

Views: 683

Answers (4)

mnel
mnel

Reputation: 115435

I think what you are looking for is circulant.matrix from the lgcp package.

If x is a matrix whose columns are the bases of the sub-blocks of a block circulant matrix, then this function returns the block circulant matrix of interest.

eg

x <- matrix(1:8,ncol=4)
 circulant(x)
#      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
# [1,]    1    2    3    4    5    6    7    8
# [2,]    2    1    4    3    6    5    8    7
# [3,]    7    8    1    2    3    4    5    6
# [4,]    8    7    2    1    4    3    6    5
# [5,]    5    6    7    8    1    2    3    4
# [6,]    6    5    8    7    2    1    4    3
# [7,]    3    4    5    6    7    8    1    2
# [8,]    4    3    6    5    8    7    2    1

Alternative approach

Here is a highly inefficient approach using kronecker and Reduce

bcirc <- function(list.blocks){
  P <- lapply(seq_along(list.blocks), function(x,y) x ==y, x = circulant(seq_along(list.blocks)))
  Reduce('+',Map(P = P, A=list.blocks, f = function(P,A) kronecker(P,A)))
  }

benchmarking with @flodel and @Ben Bolker

lbirary(microbenchmark)
microbenchmark(bcm(C), bcirc(C), bcMat(C))
Unit: microseconds
     expr       min         lq     median         uq       max neval
   bcm(C) 10836.719 10925.7845 10992.8450 11141.1240 21622.927   100
 bcirc(C)   444.983   455.7275   479.5790   487.0370   569.105   100
 bcMat(C)   288.558   296.4350   309.8945   348.4215  2190.231   100

Upvotes: 3

flodel
flodel

Reputation: 89097

Thanks for a challenging question! Here is a solution summing kronecker products of your matrices with sub- and super-diagonals.

Sample data, a list of matrices:

C <- lapply(1:3, matrix, nrow = 2, ncol = 2)

My solution:

bcm <- function(C) {
   require(Matrix)
   n <- length(C)
   Reduce(`+`, lapply((-n+1):(n-1),
                      function(i) kronecker(as.matrix(bandSparse(n, n, -i)),
                                            C[[1 + (i %% n)]])))
}
bcm(C)

#      [,1] [,2] [,3] [,4] [,5] [,6]
# [1,]    1    1    3    3    2    2
# [2,]    1    1    3    3    2    2
# [3,]    2    2    1    1    3    3
# [4,]    2    2    1    1    3    3
# [5,]    3    3    2    2    1    1
# [6,]    3    3    2    2    1    1

Upvotes: 6

Ben Bolker
Ben Bolker

Reputation: 226597

I don't know if this is particularly efficient, but as I interpret your question it does what you want.

rotList <- function(L,n) {
    if (n==0) return(L)
    c(tail(L,n),head(L,-n))
}
rowFun <- function(n,matList) do.call(rbind,rotList(matList,n))
bcMat <- function(matList) {
    n <- length(matList)
    do.call(cbind,lapply(0:(n-1),rowFun,matList))
}

Example:

bcMat(list(diag(3),matrix(1:9,nrow=3),matrix(4,nrow=3,ncol=3)))

Upvotes: 4

thelatemail
thelatemail

Reputation: 93938

Is something like this what you are looking for?

> vec <- 1:4
> sapply(rev(seq_along(vec)),function(x) c(tail(vec,x),head(vec,-x)) )

     [,1] [,2] [,3] [,4]
[1,]    1    2    3    4
[2,]    2    3    4    1
[3,]    3    4    1    2
[4,]    4    1    2    3

Upvotes: 0

Related Questions