有没有办法在 ARIMA 模型上使用 R 中的 strucchange 包?我一直找不到。非常感谢。
ARIMA 模型上的结构包
该包strucchange
需要将线性模型的公式作为输入传递给lm
. 我认为没有一种直接的方法可以将包与 function 一起使用arima
。我不知道任何其他实现此功能的 R 包,但我可以提供一些可能对您的目的有所帮助的基本指南。
您可以基于残差的累积平方和 (CUMSUM) 以及基于不同子样本中模型参数的 F 检验来执行一些诊断。
让我们以下面的模拟 AR 过程为例x
。前 50 个观测值是从 AR(1) 模型生成的,接下来的 100 个观测值是从 AR(2) 模型生成的:
set.seed(135)
x1 <- arima.sim(model = list(order = c(1,0,0), ar = -0.2), n = 50)
x2 <- arima.sim(model = list(order = c(2,0,0), ar = c(0.3, 0.5)), n = 100)
x <- ts(c(x1, x2))
CUMSUM 方法:一旦将 AR 模型拟合到整个系列,CUMSUM 过程可以得到如下:
fit <- arima(x, order = c(2,0,0), include.mean = FALSE)
e <- residuals(fit)
sigma <- sqrt(fit$sigma2)
n <- length(x)
cs <- cumsum(e) / sigma
作为参考,可以在strucchange
基于 OLS 的 CUSUM 测试的包中获得置信限。为此,我们可以创建一个类对象efp
并绘制它:
require(strucchange)
retval <- list()
retval$coefficients <- coef(fit)
retval$sigma <- sigma
retval$process <- cs
retval$type.name <- "OLS-based CUSUM test"
retval$lim.process <- "Brownian bridge"
retval$datatsp <- tsp(x)
class(retval) <- c("efp")
plot(retval)
置信限仅供参考,我不确定它们是否是在这种情况下进行正式测试的正确值。无论如何,序列中的突然变化或转变cs
可以解释为在那个时间点附近发生了一些事情的迹象,可能是结构变化。在图中,我们观察到在观察 50 左右,我们在数据生成过程中引入了变化。
F 检验:另一种方法基于 F 检验统计量,计算如下:
rss <- sum(residuals(fit)^2)
sigma2 <- fit$sigma2
stats <- rep(NA, n)
for (i in seq.int(20, n-20))
{
fit1 <- arima(x[seq(1,i)], order = c(2,0,0), include.mean = FALSE)
fit2 <- arima(x[seq(i+1,n)], order = c(2,0,0), include.mean = FALSE)
ess <- sum(c(residuals(fit1), residuals(fit2))^2)
stats[i] <- (rss - ess)/sigma2
}
与 CUMSUM 图类似,F 统计量图可能揭示结构变化的存在。根据卡方分布可以获得 95% 的置信限。
plot(stats)
abline(h = qchisq(0.05, df = length(coef(fit)), lower.tail = FALSE), lty = 2, col = "red")
如果与每个统计量相关的最小 p 值低于显着性水平,例如 0.05,那么我们可以怀疑此时存在结构变化。在观察 50 发生的这个模拟系列中,当 AR 系数在数据生成过程中发生变化时:
which.min(1 - pchisq(stats, df = 2))
#[1] 50
您可以在您可能已经知道的包装插图和其中的参考资料中找到更多详细信息。strucchange
我已经写了一篇关于使用 R 中的 strucchange 包检测结构中断的博客。这非常简单 - 这是大纲:
# assuming you have a 'ts' object in R
# 1. install package 'strucchange'
# 2. Then write down this code:
library(strucchange)
# store the breakdates
bp_ts <- breakpoints(ts)
# this will give you the break dates and their confidence intervals
summary(bp_ts)
# store the confidence intervals
ci_ts <- confint(bp_ts)
## to plot the breakpoints with confidence intervals
plot(ts)
lines(bp_ts)
lines(ci_ts)
我博客中使用的时间序列数据恰好是一个 ARIMA(0,1,1) 过程。如果您想验证这一点,请查看我的Github 存储库。
如果您想检查结构中断的存在,我建议您:
使用包进行异常值分析
tsoutliers
。通过使用 函数tso
,您可以检查您的模型是否有孤立的尖峰(加性异常值)、平均水平的突然变化(水平偏移)、需要几个周期才能消失的尖峰(瞬态变化)或模型的创新(干预异常值)。通过 QLR 检验分析参数的稳定性(参见 Andrews,2006 中的临界值)。该检验最适合以下情况:i) 中断日期未知,ii) 滞后变量是解释变量 e/或 iii) 误差是异方差/自相关的。
要运行此测试,您必须在模型中添加假人并制作循环,其中会计算多次回归,并且假人的值会随着时间的推移而变化,以捕获可能的中断。
#Define the window to be tested. In general, the first and last 15% observations are excluded
window_breaks <- seq(inic_anomes, last_anomes, 1/12)
# Loop to run the regressions and compute the test statistic
for(i in 1:length(window_breaks)) {
# Set up dummy variable
D <- time(y) > window_breaks[i]
# Estimate model with dummy
model <- lm(y ~ x + D + D*x)
# Compute and save the F-statistic
Fstats[i] <- linearHypothesis(model,
c("D", "x:D"),
vcov = kernHAC)$F[2]
}
QLR <- max(Fstats) #if QLR < critical value, there is no break