Come posso vectorize questo codice in R? Magari con la applicano function ()?

Question

Sono veramente fatica a replicare l'output della funzione dist () nel codice R senza l'utilizzo di 1 o 2 per loop. (Se vi state chiedendo perché sto facendo questo, è così che io possa giocare con il calcolo della distanza, e anche per migliorare le mie competenze R - quindi per favore solo soluzioni che coinvolgono R)

Sommario: matrice viene passato a dist (), che calcola la distanza euclidea row-wise ed emette una matrice piena distanza della distanza tra ciascuna riga (ad esempio la distanza tra le righe 1 e 50 sarà in distancematrix [1, 50] e distancematrix [50, 1]). Gli sguardi codice veloce come questo:

distancematrix <- as.matrix(dist(myMatrix, method="euclidean", diag = T, upper = T))

ho prodotto con successo la stessa uscita in R utilizzando il seguente codice:

for (i in 1:nrow(myMatrix)) {
  for (j in 1:nrow(myMatrix)) {
    distancematrix[i, j] <- sum(abs(myMatrix[i,] - myMatrix[j,])) 
  }
}

Tuttavia, utilizzando due nidificato per cicli è molto più lento rispetto all'utilizzo di dist (). Ho letto un sacco su come utilizzare apply () per ottimizzare il più lento per i cicli, ma non sono stato in grado di ottenere la mia testa intorno ad esso finora. Credo che almeno uno dei cicli for è sicuramente evitabile semplicemente l'output di un vettore, e trattare con esso alla fine. Tuttavia, non posso per la vita di me capire come rimuovere sia per i cicli.

Qualcuno ha qualche idea?

Answer 1

Prima di tutto va notato che il codice che avete inviato in realtà non replicare l'uscita della funzione dist, perché la linea:

distancematrix[i, j] <- sum(abs(myMatrix[i,] - myMatrix[j,]))

non calcola la distanza euclidea; dovrebbe essere:

distancematrix[i, j] <- sqrt(sum((myMatrix[i,] - myMatrix[j,]) ^ 2))

Qui ci sono due soluzioni che si basano su apply. Essi sono semplificate, in particolare non sfruttano la simmetria della matrice distanza (che, se ritenuto, porterebbe ad un aumento di velocità 2 volte). In primo luogo, generare alcuni dati di test:

# Number of data points
N <- 2000
# Dimensionality
d <- 10
# Generate data
myMatrix = matrix(rnorm(N * d), nrow = N)

Per convenienza, definire:

# Wrapper for the distance function
d_fun <- function(x_1, x_2) sqrt(sum((x_1 - x_2) ^ 2))

Il primo approccio è una combinazione di apply e sapply:

system.time(
    D_1 <- 
        apply(myMatrix, 1, function(x_i) 
            sapply(1:nrow(myMatrix), function(j) d_fun(x_i, myMatrix[j, ]))
    )
)

   user  system elapsed 
 14.041   0.100  14.001 

mentre il secondo utilizza solo apply (ma andare oltre gli indici, che sono abbinati expand.grid):

system.time(
    D_2 <- 
        matrix(apply(expand.grid(i = 1:nrow(myMatrix), j = 1:nrow(myMatrix)), 1, function(I) 
            d_fun(myMatrix[I[["i"]], ], myMatrix[I[["j"]], ])
        )
    )
)

   user  system elapsed 
 39.313   0.498  39.561 

Tuttavia, come previsto entrambi sono molto più lento rispetto dist:

system.time(
    distancematrix <- as.matrix(
        dist(myMatrix, method = "euclidean", diag = T, upper = T)
    )
)

   user  system elapsed 
  0.337   0.054   0.388