15.1.6 诊断:是白噪声序列吗?
时间序列分析主要是寻找一种方法来拟和这个序列(或者就像我们引言所讲,由白噪声序列来生成一个近似的序列)。我们可以以如下方向来建立模型:考察数据序列并将其转换成类似白噪声序列。为了检验我们的分析是否正确,需要检验模型残差是否为白噪声序列。至此,我们可以考察自相关函数图(平均而言,仅能有5%的相关系数线超过虚线,如果有更多,那么我们的分析或者说结果是有疑问的)。
<br />
z <- rnorm(200)<br />
op <- par(mfrow=c(2,1), mar=c(5,4,2,2)+.1)<br />
plot(ts(z))<br />
acf(z, main = "")<br />
par(op)<br />
图略
<br />
x <- diff(co2)<br />
y <- diff(x,lag=12)<br />
op <- par(mfrow=c(2,1), mar=c(5,4,2,2)+.1)<br />
plot(ts(y))<br />
acf(y, main="")<br />
par(op)<br />
图略
为了得到数值结果(p值),我们可以进行Box-Pierce或者Ljung-Box 检验(这些检验也称为“portmanteau 检验”):其思想是考察(加权)平均和一阶自相关系数,这些求和结果(渐近)服从卡方分布(Ljung-Box 检验较Box-Pierce检验更为稳健,小样本性质更好近似的卡方分布)。
> Box.test(z) # Box-Pierce
Box-Pierce test
X-squared = 0.014, df = 1, p-value = 0.9059
> Box.test(z, type="Ljung-Box")
Box-Ljung test
X-squared = 0.0142, df = 1, p-value = 0.9051
> Box.test(y)
Box-Pierce test
X-squared = 41.5007, df = 1, p-value = 1.178e-10
> Box.test(y, type=”Ljung”)
Box-Ljung test
X-squared = 41.7749, df = 1, p-value = 1.024e-10
<br />
op <- par(mfrow=c(2,1))<br />
plot.box.ljung <- function (z,k = 15,main = "p-value of the Ljung-Box test",ylab = "p-value") {<br />
p <- rep(NA, k)<br />
for (i in 1:k) { <br />
p[i] <- Box.test(z, i,<br />
type = "Ljung-Box")$p.value<br />
}<br />
plot(p,type = "h",ylim = c(0,1),lwd = 3,main = main,ylab = ylab)<br />
abline(h = c(0,.05),lty = 3)<br />
}<br />
plot.box.ljung(z, main="Random data")<br />
plot.box.ljung(y, main="diff(diff(co2),lag=12)")<br />
par(op)<br />
图略
还有其它的检验方法:McLeod-Li, Turning-point, difference-sign, rank 检验等等。我们也可以用DW检验,在回归问题中我们已经应用了DW检验。
TODO: check that I actually mention it.
library(car)
?durbin.watson
library(lmtest)
?dwtest
这里应用的同样的检验方法,但并非以同样的方式检验,结果也许不同:
> dwtest(LakeHuron ~ 1)
Durbin-Watson test
data: LakeHuron ~ 1
DW = 0.3195, p-value = < 2.2e-16
alternative hypothesis: true autocorrelation is greater than 0
> durbin.watson(lm(LakeHuron ~ 1))
lag Autocorrelation D-W Statistic p-value
1 0.8319112 0.3195269 0
Alternative hypothesis: rho != 0
<br />
op <- par(mfrow=c(2,1))<br />
library(lmtest)<br />
plot(LakeHuron,main = "Lake Huron")<br />
acf(LakeHuron,main = paste("Durbin-Watson: p =",signif( dwtest( LakeHuron ~ 1 ) $ p.value, 3 )))<br />
par(op)<br />
图略
<br />
n <- 200<br />
x <- rnorm(n)<br />
op <- par(mfrow=c(2,1))<br />
x <- ts(x)<br />
plot(x, main="White noise", ylab="")<br />
acf(x,main = paste("Durbin-Watson: p =",signif( dwtest( x ~ 1 ) $ p.value, 3)))<br />
par(op)<br />
图略
<br />
n <- 200<br />
x <- rnorm(n)<br />
op <- par(mfrow=c(2,1))<br />
y <- filter(x,.8,method="recursive")<br />
plot(y, main="AR(1)", ylab="")<br />
acf(y,main = paste("p =",signif( dwtest( y ~ 1 ) $ p.value, 3 )))<br />
par(op)<br />
图略
但是注意:默认的检验是检验自相关是零还是正的,如果自相关显著为负,结果会令人费解。
<br />
set.seed(1)<br />
n <- 200<br />
x <- rnorm(n)<br />
y <- filter(x, c(0,1), method="recursive")<br />
op <- par(mfrow=c(3,1), mar=c(2,4,2,2)+.1)<br />
plot(y,main = paste("one-sided DW test: p =",signif( dwtest ( y ~ 1 ) $ p.value, 3 )))<br />
acf( y, main="")<br />
pacf(y, main="")<br />
par(op)<br />
图略
<br />
op <- par(mfrow=c(3,1), mar=c(2,4,2,2)+.1)<br />
res <- dwtest( y ~ 1, alternative="two.sided")<br />
plot(y,main = paste("two-sided p =",signif( res$p.value, 3 )))<br />
acf(y, main="")<br />
pacf(y, main="")<br />
par(op)<br />
图略
还有其它检验,例如runs检验(它看起来像数据在run,一个run是对同一个符号的连续观测;这类检验主要用来进行时间序列的定性检验)。
library(tseries)
?runs.test
或者Cowles-Jones 检验
TODO
A sequence is a pair of consecutive returns of the same sign; a
reversal is a pair of consecutive returns of opposite signs.
Their ratio, the Cowles-Jones ratio,
number of sequences
CJ = ---------------------
number of reversals
should be around one IF the drift is zero.