R - Wektoryzacja realizacji funkcji trójskładnikowej

Question

Mam trzy wektory X, Y i Z o równej długości n. Potrzebuję utworzyć tablicę n x n x n z funkcją f(X[i],Y[j],Z[k]). Prostym sposobem na to jest sekwencyjne przechodzenie przez każdy element każdego z 3 wektorów. Jednak czas wymagany do obliczenia tablicy rośnie wykładniczo z n. Czy istnieje sposób wdrożenia tego za pomocą operacji wektoryzacji?

EDYCJA: Jak wspomniano w komentarzach, dodałem prosty przykład tego, co jest potrzebne.

set.seed(1) 
X = rnorm(10) 
Y = seq(11,20) 
Z = seq(21,30) 

F = array(0, dim=c(length(X),length(Y),length(Z))) 
for (i in 1:length(X)) 
    for (j in 1:length(Y)) 
    for (k in 1:length(Z)) 
     F[i,j,k] = X[i] * (Y[j] + Z[k]) 

Dzięki.

Answer 1

Można używać zagnieżdżonych outer:

set.seed(1) 
X = rnorm(10) 
Y = seq(11,20) 
Z = seq(21,30) 

F = array(0, dim = c(length(X),length(Y),length(Z))) 
for (i in 1:length(X)) 
    for (j in 1:length(Y)) 
    for (k in 1:length(Z)) 
     F[i,j,k] = X[i] * (Y[j] + Z[k]) 

F2 <- outer(X, outer(Y, Z, "+"), "*") 

> identical(F, F2) 
[1] TRUE 

A microbenchmark tym expand.grid rozwiązania zaproponowanego przez Nicka K:

X = rnorm(100) 
Y = seq(1:100) 
Z = seq(101:200) 

forLoop <- function(X, Y, Z) { 
    F = array(0, dim = c(length(X),length(Y),length(Z))) 
    for (i in 1:length(X)) 
    for (j in 1:length(Y)) 
     for (k in 1:length(Z)) 
     F[i,j,k] = X[i] * (Y[j] + Z[k]) 
    return(F) 
} 

nestedOuter <- function(X, Y, Z) { 
    outer(X, outer(Y, Z, "+"), "*") 
} 

expandGrid <- function(X, Y, Z) { 
    df <- expand.grid(X = X, Y = Y, Z = Z) 
    G <- df$X * (df$Y + df$Z) 
    dim(G) <- c(length(X), length(Y), length(Z)) 
    return(G) 
} 

library(microbenchmark) 
mbm <- microbenchmark(
    forLoop = F1 <- forLoop(X, Y, Z), 
    nestedOuter = F2 <- nestedOuter(X, Y, Z), 
    expandGrid = F3 <- expandGrid(X, Y, Z), 
    times = 50L) 

> mbm 
Unit: milliseconds 
expr   min   lq  mean  median   uq  max neval 
forLoop 3261.872552 3339.37383 3458.812265 3388.721159 3524.651971 4074.40422 50 
nestedOuter 3.293461 3.36810 9.874336 3.541637 5.126789 54.24087 50 
expandGrid 53.907789 57.15647 85.612048 88.286431 103.516819 235.45443 50 

Answer 2

można użyć expand.grid następująco:

df <- expand.grid(X = X, Y = Y, Z = Z) 
G <- df$X * (df$Y + df$Z) 
dim(G) <- c(length(X), length(Y), length(Z)) 
all.equal(F, G) 

Jeśli miałeś funkcję wektorowania, to by działało równie dobrze. Jeśli nie, możesz użyć plyr :: daply.

Answer 3

Oto dodatkowa opcja, możliwa implementacja Rcpp (w przypadku, gdy lubisz swoje pętle). I nie był w stanie wyprzedzić @Juliens rozwiązanie chociaż (a może ktoś może), ale są one bardziej lub mniej mają taką samą taktowanie

library(Rcpp) 
cppFunction('NumericVector RCPP(NumericVector X, NumericVector Y, NumericVector Z){ 

      int nrow = X.size(), ncol = 3, indx = 0; 
      double temp(1) ; 
      NumericVector out(pow(nrow, ncol)) ; 
      IntegerVector dim(ncol) ; 

      for (int l = 0; l < ncol; l++){ 
       dim[l] = nrow; 
      }    

      for (int j = 0; j < nrow; j++) { 
       for (int k = 0; k < nrow; k++) { 
        temp = Y[j] + Z[k] ; 
        for (int i = 0; i < nrow; i++) { 
         out[indx] = X[i] * temp ; 
         indx += 1 ; 
        } 
       } 
      } 

      out.attr("dim") = dim; 
      return out; 
}') 

Validating

identical(RCPP(X, Y, Z), F) 
## [1] TRUE 

Szybka odniesienia

set.seed(123) 
X = rnorm(100) 
Y = 1:100 
Z = 101:200 

nestedOuter <- function(X, Y, Z) outer(X, outer(Y, Z, "+"), "*") 

library(microbenchmark) 
microbenchmark( 
    nestedOuter = nestedOuter(X, Y, Z), 
    RCPP = RCPP(X, Y, Z), 
    unit = "relative", 
    times = 1e4) 

# Unit: relative 
#  expr  min  lq  mean median  uq  max neval 
# nestedOuter 1.000000 1.000000 1.000000 1.000000 1.000000 1.0000000 10000 
#  RCPP 1.164254 1.141713 1.081235 1.100596 1.080133 0.7092394 10000