在 R 中从头开始计算自相关函数

Question

我想在 R 中了解如何从头开始计算自相关函数。当 y 是一个 ts 对象时，如何使用 cor(x=y, y=lag(x=y, k=2)) 得到 ACF？

我已经尝试了 use 参数 [use = "complete.obs" # c("everything", "all.obs", "complete.obs", "na.or.complete", "pairwise.complete.obs")] 的所有选项。

set.seed(1)
y <-ts(data = rnorm(20), start = c(2010, 1), frequency = 4)
y

# Qtr1        Qtr2        Qtr3        Qtr4
# 2010  0.91897737  0.78213630  0.07456498 -1.98935170
# 2011  0.61982575 -0.05612874 -0.15579551 -1.47075238
# 2012 -0.47815006  0.41794156  1.35867955 -0.10278773
# 2013  0.38767161 -0.05380504 -1.37705956 -0.41499456
# 2014 -0.39428995 -0.05931340  1.10002537  0.76317575

acf(x=y, plot=FALSE)
# Autocorrelations of series ‘y’, by lag
# 
# 0.00   0.25   0.50   0.75   1.00   1.25   1.50   1.75   2.00   2.25   2.50   2.75   3.00   3.25 
# 1.000 -0.122 -0.185 -0.049  0.147 -0.283 -0.255  0.212  0.097 -0.120 -0.181  0.286 -0.063  0.094 

cor(
      x   = y
    , y   = lag(x=y, k=2)
    , use = "complete.obs" # c("everything", "all.obs", "complete.obs", "na.or.complete", "pairwise.complete.obs")
    )

# [1] 1

Answer 1

根据@khashaa 提供的 link 得到了正确答案。谢谢

acf(x=y, plot=FALSE)
# Autocorrelations of series ‘y’, by lag
# 
# 0.00   0.25   0.50   0.75   1.00   1.25   1.50   1.75   2.00   2.25   2.50   2.75   3.00   3.25 
# 1.000 -0.122 -0.185 -0.049  0.147 -0.283 -0.255  0.212  0.097 -0.120 -0.181  0.286 -0.063  0.094 

sum((y-mean(y))*(lag(x=y, k=0)-mean(y)))/sum((y-mean(y))^2)
[1] 1
sum((y-mean(y))*(lag(x=y, k=1)-mean(y)))/sum((y-mean(y))^2)
[1] -0.1222859
sum((y-mean(y))*(lag(x=y, k=2)-mean(y)))/sum((y-mean(y))^2)
[1] -0.1852114
sum((y-mean(y))*(lag(x=y, k=3)-mean(y)))/sum((y-mean(y))^2)
[1] -0.04940401

Answer 2

看来你有两个问题。

第一个问题：cor(x, lag(x, k = k))将永远是1，因为这两个系列在时间上没有对齐。在使用 cor

之前，您需要使用 ts.union 或 cbind.ts

X <- ts.union(yt = y, yt2 = lag(x = y, k =  2))
head(X)
##            yt      yt2
## [1,]       NA -0.62645
## [2,]       NA  0.18364
## [3,] -0.62645 -0.83563
## [4,]  0.18364  1.59528
## [5,] -0.83563  0.32951
## [6,]  1.59528 -0.82047
tail(X)
##              yt      yt2
## [17,]  1.124931 -0.01619
## [18,] -0.044934  0.94384
## [19,] -0.016190  0.82122
## [20,]  0.943836  0.59390
## [21,]  0.821221       NA
## [22,]  0.593901       NA

问题是没有时间索引，原始数据是相同的（有时间偏移）。你可以自己测试一下

x1 <- as.vector(y)
x2 <- as.vector(lag(y, k = 2))
all.equal(x1, x2)
## [1] TRUE

因此，如果要计算时间序列与其滞后之间的相关系数，可以使用 X（使用 ts.union 创建）

cor(X[, 1], X[, 2], use = "complete.obs")
## [1] -0.19018

结果仍然不同于 acf(y, plot = FALSE)$acf[3]

acf(y, plot = FALSE)$acf[3]
## [1] -0.18521

这给我们带来了第二个原因为什么你不能使用cor来计算acf：

Acf 的数学定义假设至少有一个 second-order 平稳性（每个滞后的均值和方差相等）：

但如果您使用 cor 的标准实现，对于每个序列，将计算不同的均值和方差值（滞后值将与原始序列不同）。

c0 <- var(y)
m <- mean(y)
n <- length(y)
ct <- sum((X[, 1] - m) * (X[, 2] - m), na.rm = TRUE) / (n - 1)
(rt <- ct / c0)
## [1] -0.18521

all.equal(rt, acf(y, plot = FALSE)$acf[3])
## [1] TRUE

在 R 中从头开始计算自相关函数

Calculating Autocorrelation Function from scratch in R

statistics

r