Wzór na autokorelację w R vs. Excel

Próbuję dowiedzieć się, w jaki sposób R oblicza autokorelację lag-k (najwyraźniej jest to ta sama formuła stosowana przez Minitab i SAS), dzięki czemu mogę ją porównać do korzystania z funkcji CORREL programu Excel zastosowanej do serii i jej wersji z opóźnieniem k. R i Excel (przy użyciu CORREL) dają nieco inne wartości autokorelacji.

Chciałbym również dowiedzieć się, czy jedno obliczenie jest bardziej poprawne niż drugie.

Dokładne równanie podano w: Venables, WN i Ripley, BD (2002) Modern Applied Statistics with S. Fourth Edition. Springer-Verlag. Dam ci przykład:

### simulate some data with AR(1) where rho = .75
xi <- 1:50
yi <- arima.sim(model=list(ar=.75), n=50)

### get residuals
res <- resid(lm(yi ~ xi))

### acf for lags 1 and 2
cor(res[1:49], res[2:50])      ### not quite how this is calculated by R
cor(res[1:48], res[3:50])      ### not quite how this is calculated by R

### how R calculates these
acf(res, lag.max=2, plot=F)

### how this is calculated by R
### note: mean(res) = 0 for this example, so technically not needed here
c0 <- 1/50 * sum( (res[1:50] - mean(res)) * (res[1:50] - mean(res)) ) 
c1 <- 1/50 * sum( (res[1:49] - mean(res)) * (res[2:50] - mean(res)) ) 
c2 <- 1/50 * sum( (res[1:48] - mean(res)) * (res[3:50] - mean(res)) ) 
c1/c0
c2/c0

I tak dalej (np. res[1:47]I res[4:50]dla opóźnienia 3).

Naiwnym sposobem obliczania automatycznej korelacji (i ewentualnie tego, co wykorzystuje Excel) jest utworzenie 2 kopii wektora, a następnie usunięcie pierwszego n elementów z pierwszej kopii i ostatnich n elementów z drugiej kopii (gdzie n jest opóźnieniem obliczam z). Następnie przekaż te 2 wektory do funkcji, aby obliczyć korelację. Ta metoda jest OK i da rozsądną odpowiedź, ale ignoruje fakt, że 2 porównywane wektory są tak naprawdę miarami tego samego.

Udoskonalona wersja (jak pokazuje Wolfgang) jest podobną funkcją do regularnej korelacji, z tym wyjątkiem, że wykorzystuje cały wektor do obliczenia średniej i wariancji.