Kronecker product
Definition
For any two matrices \(\textbf{A} = \{a_{ij}\}\) of dimensions \(m\times n\) and \(\textbf{B} = \{b_{ij}\}\) of dimensions \(p\times q\), the direct Kronecker product between them is the matrix of dimensions \(mp\times nq\) defined as the block matrix
\[ \textbf{A}\otimes\textbf{B} = \{a_{ij}\textbf{B}\} \]
This can be computed using the kronecker() function from
the ‘base’ R-package. Alternatively, the Kronecker()
function from the ‘tensorEVD’ R-package can be used.
Examples
Example 1
Simple scalar multiplication. Let \(a\) be a scalar and \(\textbf{B}\) any matrix. Then, computing their Kronecker product is the same as multiplying \(\textbf{B}\) by the scalar:
a <- 10
( B <- matrix(1:6, ncol=2) )## [,1] [,2]
## [1,] 1 4
## [2,] 2 5
## [3,] 3 6Kronecker(a, B)## [,1] [,2]
## [1,] 10 40
## [2,] 20 50
## [3,] 30 60# In this case, should be equal to Kronecker(B, a)
Kronecker(B, a)## [,1] [,2]
## [1,] 10 40
## [2,] 20 50
## [3,] 30 60Example 2
Block diagonal matrix. A block diagonal matrix can be formed using the Kronecker product with a diagonal matrix
D <- diag(1, 2)
Kronecker(D, B)## [,1] [,2] [,3] [,4]
## [1,] 1 4 0 0
## [2,] 2 5 0 0
## [3,] 3 6 0 0
## [4,] 0 0 1 4
## [5,] 0 0 2 5
## [6,] 0 0 3 6# Is not equal to
Kronecker(B, D) # this a 'striped' matrix## [,1] [,2] [,3] [,4]
## [1,] 1 0 4 0
## [2,] 0 1 0 4
## [3,] 2 0 5 0
## [4,] 0 2 0 5
## [5,] 3 0 6 0
## [6,] 0 3 0 6Example 3
Outer product. Let \(\textbf{a}\) and \(\textbf{b}\) be any two vectors. Then, the Kronecker product between \(\textbf{a}\) and the transpose of \(\textbf{b}\) form an outer product:
a <- c(1,2,3)
b <- c(4,5)
Kronecker(a, t(b))## [,1] [,2]
## [1,] 4 5
## [2,] 8 10
## [3,] 12 15# Should be equal to the product a b'
tcrossprod(a, b)## [,1] [,2]
## [1,] 4 5
## [2,] 8 10
## [3,] 12 15Performance
Here we compare tensorEVD’s Kronecker() function with
the kronecker() function from the ‘base’ R-package in the
computation of a Kronecker product of two general matrices \(\textbf{A}\) and \(\textbf{B}\).
# Simulating matrices A and B
m = 20; n = 20
p = 40; q = 30
A <- matrix(rnorm(m*n), ncol=n)
B <- matrix(rnorm(p*q), ncol=q)
# Making the Kronecker product
K1 <- kronecker(A, B)
K2 <- Kronecker(A, B)
# Should be equal
all.equal(K1,K2)## [1] TRUEBenchmark
Here we compare these methods in terms of computational speed in
scenarios with small and large matrices. We included in the benchmark
the kronecker.prod() function from the ‘fastmatrix’
R-package. The following benchmark was performed using the code provided
in this script
run on a Linux environment based on the following system settings:
- Machine: Intel(R) Xeon(R) Gold 6148 CPU @ 2.40GHz
- Memory: 64 GB in RAM
- R version 4.1.1 (2021-08-10)
Subsetting a Kronecker
Selecting specific rows and columns from the Kronecker can be done by pre- and post- multiplication with incidence matrices, for instance
\[ \textbf{R}(\textbf{A}\otimes\textbf{B})\textbf{C}' \]
where \(\textbf{R}\) is an incidence matrix for rows and \(\textbf{C}\) is an incidence matrix for columns.
This sub-matrix can be obtained by indexing rows or columns using,
for instance, integer vectors rows and cols as
Kronecker(A, B)[rows,cols]. However, this approach computes
first the Kronecker product then makes the sub-setting.
This approach can be inefficient if a relatively small number of rows
and columns are to be selected. The Kronecker() function
can derive this sub-matrix directly from \(\textbf{A}\) and \(\textbf{B}\) on the fly without forming the
whole Kronecker product using rows and cols as
arguments. For example,
dm <- c(nrow(A)*nrow(B), ncol(A)*ncol(B)) # dimension of the Kronecker
# Subsetting a matrix with 30% of rows/columns
rows <- sample(seq(dm[1]), 0.3*dm[1])
cols <- sample(seq(dm[2]), 0.3*dm[2])
K1 <- Kronecker(A, B)[rows,cols]
K2 <- Kronecker(A, B, rows=rows, cols=cols)
dim(K1) # small size## [1] 240 180all.equal(K1, K2)## [1] TRUEBenchmark
Here we show some benchmark results on the time performance of
subsetting a Kronecker matrix using
Kronecker(A, B, rows, cols) and
Kronecker(A, B)[rows,cols], for small (30% of the
rows/columns) and large (300% of the rows/columns) sub-matrices
scenarios, from a Kronecker product matrix of dimension \(5000\times 5000\). The code to perform this
benchmark is provided in this script.
Extras
Inplace calculation
If \(\textbf{A}\) is a matrix and \(\textbf{B}\) is a scalar, the resulting Kronecker will be of the same dimensions as \(\textbf{A}\); therefore, we could overwrite the result on the memory occupied by \(\textbf{A}\).
Usually, assigning an output to an object will occupy a different memory address than inputs:
B <- rnorm(1) # B is a scalar
K <- Kronecker(A, B)
c(K=pryr::address(K), A=pryr::address(A))## K A
## "0x7fd4a1e12200" "0x7fd4a1277c00"The parameter inplace can be used to overwrite the
output at the same address as the input:
A <- Kronecker(A, B, inplace=TRUE)
c(A=pryr::address(A)) # output address remain unchanged## A
## "0x7fd4a1277c00"all.equal(K, A) # contains the desired result## [1] TRUEMaking dimension names
Row and column names for the Kronecker product can be retrieved using
the make.dimnames argument. Attribute dimnames
of the Kronecker will be produced by crossing rownames and
colnames of input \(\textbf{A}\) with those of input \(\textbf{B}\). For instance,
A <- matrix(1:9, ncol=3)
B <- matrix(10*(1:4), ncol=2)
dimnames(A) <- list(paste("week",1:3), month.abb[1:3])
dimnames(B) <- list(c("office","home"), LETTERS[1:2])
Kronecker(A, B, make.dimnames=TRUE)## Jan:A Jan:B Feb:A Feb:B Mar:A Mar:B
## week 1:office 10 30 40 120 70 210
## week 1:home 20 40 80 160 140 280
## week 2:office 20 60 50 150 80 240
## week 2:home 40 80 100 200 160 320
## week 3:office 30 90 60 180 90 270
## week 3:home 60 120 120 240 180 360