Jak wydajniej obliczać toczną kowariancję

Question

Próbuję obliczyć toczną kowariancję między zestawem danych (każda kolumna mojej zmiennej x) i jedną inną (zmienną y) w R. Myślałem, że mogę użyć jednego z wnioskującego funkcje, ale nie może znaleźć w tym samym czasie dwóch zestawów wejść. Oto, co próbowałem:

set.seed(1) 
x<-matrix(rnorm(500),nrow=100,ncol=5) 
y<-rnorm(100) 
rollapply(x,width=5,FUN= function(x) {cov(x,y)}) 
z<-cbind(x,y) 
rollapply(z,width=5, FUN=function(x){cov(z,z[,6])}) 

Ale nikt nie robi tego, co bym chciał. Jednym rozwiązaniem znalazłem jest użycie dla pętli, ale zastanawiam się, czy mogę być bardziej wydajny w R niż:

dResult<-matrix(nrow=96,ncol=5) 
for(iLine in 1:96){ 
    for(iCol in 1:5){ 
     dResult[iLine,iCol]=cov(x[iLine:(iLine+4),iCol],y[iLine:(iLine+4)]) 
    } 
} 

co daje mi oczekiwany wynik:

head(dResult) 


      [,1]  [,2]  [,3]  [,4]  [,5] 
[1,] 0.32056460 0.05281386 -1.13283586 -0.01741274 -0.01464430 
[2,] -0.03246014 0.78631603 -0.34309778 0.29919297 -0.22243572 
[3,] -0.16239479 0.56372428 -0.27476604 0.39007645 0.05461355 
[4,] -0.56764687 0.09847672 0.11204244 0.78044096 -0.01980684 
[5,] -0.43081539 0.01904417 0.01282632 0.35550327 0.31062580 
[6,] -0.28890607 0.03967327 0.58307743 0.15055881 0.60704533 

Answer 1

set.seed(1) 
x<-as.data.frame(matrix(rnorm(500),nrow=100,ncol=5)) 
y<-rnorm(100) 

library(zoo) 

covResult = sapply(x,function(alpha) { 

cov_value = rollapply(cbind(alpha,y),width=5,FUN = function(beta) cov(beta[,1],beta[,2]),by.column=FALSE,align="right") 

return(cov_value) 

}) 

head(covResult) 
#    V1   V2   V3   V4   V5 
#[1,] 0.32056460 0.05281386 -1.13283586 -0.01741274 -0.01464430 
#[2,] -0.03246014 0.78631603 -0.34309778 0.29919297 -0.22243572 
#[3,] -0.16239479 0.56372428 -0.27476604 0.39007645 0.05461355 
#[4,] -0.56764687 0.09847672 0.11204244 0.78044096 -0.01980684 
#[5,] -0.43081539 0.01904417 0.01282632 0.35550327 0.31062580 
#[6,] -0.28890607 0.03967327 0.58307743 0.15055881 0.60704533 

Sprawdź również:

library(PerformanceAnalytics) 
?chart.rollingCorrelation 

Answer 2

Teraz ja” m działa kilka długich symulacji, więc nie można użyć R, ale należy to liczyć to powinno działać. Warstwa zewnętrzna stosowana przez kolumny pobiera kolumnę, przekazuje ją do rollapply, gdzie zostanie użyta do wykonania kowariancji okna przewijania z y. Miejmy nadzieję: D

apply(x,2,function(x) rollapply(x,width=5,function(z) cov(x,y))) 

Answer 3

Jest to rozwiązanie z rollapply() i sapply():

sapply(1:5, function(j) rollapply(1:100, 5, function(i) cov(x[i, j], y[i]))) 

myślę, że to jest bardziej czytelna i bardziej R-owski niż rozwiązania z for-pętle, ale sprawdziłem z microbenchmark i wydaje się wolniejsze.

Answer 4

Jeśli potrzebujesz czegoś szybciej i nie trzeba żadnej z innych niż domyślne argumenty cov, można użyć TTR::runCov. Zwróć uwagę, że domyślnie jest on wyposażony w wiodącą NA.

Różnica prędkości będzie miała większe znaczenie przy większych danych. Oto przykład, w jaki sposób z niego korzystać:

cov_joshua <- function() { 
    apply(x, 2, function(x, y) TTR::runCov(x, y, 5), y = y) 
} 

A tu porównanie z obecnie przyjętym odpowiedź za pomocą małego zestawu danych OP stanowił:

cov_osssan <- function() { 
    f <- function(b) cov(b[,1], b[,2]) 
    apply(x, 2, function(a) { 
    rollapplyr(cbind(a,y), width=5, FUN = f, by.column=FALSE) 
    }) 
} 
require(zoo) # for cov_osssan 
require(microbenchmark) 
set.seed(1) 
nr <- 100 
nc <- 5 
x <- matrix(rnorm(nc*nr),nrow=nr,ncol=nc) 
y <- rnorm(nr) 
microbenchmark(cov_osssan(), cov_joshua()) 
# Unit: milliseconds 
#   expr  min  lq median  uq  max neval 
# cov_osssan() 22.881253 24.569906 25.625623 27.44348 32.81344 100 
# cov_joshua() 5.841422 6.170189 6.706466 7.47609 31.24717 100 
all.equal(cov_osssan(), cov_joshua()[-(1:4),]) # rm leading NA 
# [1] TRUE 

Teraz, wykorzystując większy zestaw danych:

system.time(cov_joshua()) 
# user system elapsed 
# 2.117 0.032 2.158 
system.time(cov_osssan()) 
# ^C 
# Timing stopped at: 144.957 0.36 145.491 

Zmęczony czekaniem (po ~ 2,5 minuty) na ukończenie cov_osssan.