применяя перекрестное произведение (кронекер) много раз

Предположим, у меня есть список list0 длины D, где каждый элемент представляет собой матрицу N x T.

Я пытаюсь создать продукт Кронекера, ряд за рядом, который делает следующее.

для (я в 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

r list function apply cross-product
person wolfsatthedoor    schedule 28.04.2018    source источник
comment
Редактировать: Мое простое решение заключалось в том, чтобы хранить каждый крон в отдельном списке, затем повторять его после D›2, а затем вставлять в окончательную матрицу. Хотя, возможно, есть более разумное решение.   -  person wolfsatthedoor    schedule 28.04.2018

Ответы (3)


arrow_upward
2
arrow_downward

lapply дает список i-й строки матриц, а Reduce объединяет их вместе. Затем sapply собирает это в окончательную матрицу.

N <- nrow(list0[[1]])
sapply(1:N, function(i) Reduce("%x%", init = 1, lapply(rev(list0), "[", i, TRUE)))
person G. Grothendieck    schedule 28.04.2018

arrow_upward
1
arrow_downward

Похоже, массив может помочь вам здесь:

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
person ms609    schedule 28.04.2018

arrow_upward
1
arrow_downward

Вот решение с использованием 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

Честно говоря, я удивлен, что цикл настолько медленнее.

person John Fox    schedule 28.04.2018
comment
Масштабируется ли это, если D › 4? - person ms609; 29.04.2018