假设我有一个长度 D
的列表 list0
,其中每个元素都是矩阵N x T
。
我正在尝试创建一个kronecker产品,逐行进行以下操作。
for(i in 1:n({
dummy[,i] <- list0[[D]][i,] %x% ...( (list0[[2]][i,] %x% list0[[1]][i,]))
}
有人知道应用此功能的最聪明方法吗?以下是我手动输入它的示例,但我想要任意d。
set.seed(1)
N = 2
T = 3
D = 4
dummy = matrix(0,(T)^D,N)
list0 = list()
for(d in 1:D) {
list0[[d]] <- matrix(rnorm(N*T,0,1),N,T)
}
for(i in 1:N){
dummy[,i] <- list0[[4]][i,] %x% (list0[[3]][i,] %x% (list0[[2]][i,] %x% list0[[1]][i,]))
}
head(dummy)
[,1] [,2]
[1,] 0.15578313 -0.1783412
[2,] 0.20779959 -1.5492222
[3,] -0.08194020 0.7967800
[4,] 0.18402067 0.0737661
[5,] 0.24546573 0.6407946
[6,] -0.09679284 -0.3295669
lapply
给出了矩阵的ITH行和 Reduce
kroneckers的列表。然后将sapply
组装到最终矩阵中。
N <- nrow(list0[[1]])
sapply(1:N, function(i) Reduce("%x%", init = 1, lapply(rev(list0), "[", i, TRUE)))
看起来数组可能在这里为您提供帮助:
set.seed(1)
arr0 <- array(rnorm(N*T*D, 0, 1), c(N, T, D))
result <- apply(arr0, 1, function (slice) {
xx <- slice[, 1]
for (i in seq_len(D)[-1]) xx <- xx %x% slice[, i]
xx
})
输出:
> head(result)
[,1] [,2]
[1,] 0.15578313 -0.178341215
[2,] 0.17432718 -0.234865849
[3,] 0.01414475 0.597377689
[4,] -0.28208921 -0.003618330
[5,] -0.31566842 -0.004765147
[6,] -0.02561305 0.012120078
这是使用sapply()
:
library(microbenchmark)
microbenchmark(
loop={
dummy <- matrix(0, (T)^D, N)
for(i in 1:N){
dummy[, i] <- list0[[4]][i, ] %x% (list0[[3]][i, ] %x%
(list0[[2]][i, ] %x% list0[[1]][i, ]))
}
},
sapply={
dummy2 <- sapply(1:N, function(i) list0[[4]][i,] %x% (list0[[3]][i,] %x%
(list0[[2]][i,] %x% list0[[1]][i,])))
}
)
Unit: microseconds
expr min lq mean median uq max neval cld
loop 5014.211 5190.6955 5578.7469 5320.988 5505.053 9268.179 100 b
sapply 199.230 212.0995 278.1589 229.025 244.364 4927.115 100 a
all.equal(dummy, dummy2)
[1] TRUE
我坦率地说,循环速度较慢。