2016-10-01 5 views
-1

Ich habe eine n x n Symmetrix Toeplitz-Matrix T, ein Vektor v der Länge n, und ich möchte das Matrix-Vektor-Produkt T%*%v schnell berechnen. Gibt es ein Paket in R, das die Methode der schnellen Fourier-Transformation zur Berechnung T%*%v (oder eine andere Methode, falls eine existiert) verwenden kann? Zum Beispiel hat Matlab das Toeplitzmult-Paket.Toeplitz Matrix-Vektor Multiplikation in R

Antwort

1

Die folgende Funktion funktioniert. Beachten Sie, dass die Funktion ifft() die Bibliothek pracma benötigt.

toepmult <- function(A,v){ 
    n <- nrow(A) 
    x <- as.matrix(c(A[1,],0,A[1,][n:2])) 
    p <- c(v,rep(0,n)) 
    h <- as.vector(fft(p)*fft(x)) 
    out <- Re(pracma::ifft(h)[1:n]) 
    return(matrix(out,n)) 
} 

Für einen Vektor/Matrix der Größe 1000 die toepmult Funktion nimmt etwa 18% der Zeit A%*%v dauert.

A <- toeplitz(runif(1000)) 
v <- runif(1000) 
microb(A%*%v,toepmult(A,v),times=1000) 
#Unit: microseconds 
#   expr  min  lq  mean median  uq  max neval 
#  A %*% v 1515.858 1597.345 1809.3517 1693.533 1957.4350 3868.788 1000 
# toepmult(A, v) 185.901 215.395 331.2928 298.435 347.7335 4611.285 1000 
#[[1]] 
#   [,1]  [,2] 
#median 1693.533 298.435 
#ratio  1.000  0.176 
#diff  0.000 -1395.098 

Für einen Vektor/Matrix der Größe 10.000, die toepmult Funktion nimmt ungefähr 2,5% der Zeit A%*%v dauert.

A <- toeplitz(runif(10000)) 
v <- runif(10000) 
microb(A%*%v,toepmult(A,v),times=1000) 
#Unit: milliseconds 
#   expr  min   lq  mean  median   uq  max neval 
#  A %*% v 145.834304 160.395663 181.842779 170.396014 186.221449 495.2003 1000 
# toepmult(A, v) 2.802058 4.077408 4.990894 4.322707 4.911103 180.4926 1000 
#[[1]] 
#   [,1]  [,2] 
#median 170.396 4.323 
#ratio 1.000 0.025 
#diff  0.000 -166.073