例2.5中,我们用下面的语句对拟合arima模型之后的残差进行了LB检验: a=arima(prop, order = c(1,0,0),method=\r=a$residuals
a=arima(prop, order = c(1,0,0),method=\r=a$residuals
#用r来保存残差
Box.test(r,type=\#对残差进行纯随机性检验
最后一句不完整,需要加上参数fitdf=1,修改为 Box.test(r,type=\,fitdf=1)
fitdf表示p+q,number of degrees of freedom to be subtracted if x is a series of residuals,当检验的序列是残差到时候,需要加上命令fitdf,表示减去的自由度。 运行Box.test(r,type=\后,显示的结果: Box.test(r,type=\ Box-Ljung test data: r
X-squared = 5.8661, df = 5, p-value = 0.3195
“df = 5”表示自由度为5,由于参数lag=6,所以是滞后6期的检验。
第四讲
# example4_1 拟合线性模型
x1=c(12.79,14.02,12.92,18.27,21.22,18.81,25.73,26.27,26.75,28.73,31.71,33.95) a=as.ts(x1) is.ts(a) ts.plot(a)
t=1:12 t
lm1=lm(a~t)
summary(lm1) # 返回拟合参数的统计量 coef(lm1) #返回被估计的系数 fitted(lm1) #返回模拟值 residuals(lm1) #返回残差值 fit1=as.ts(fitted(lm1))
ts.plot(a);lines(fit1,col=\拟合图 #eg1
cs=ts(scan(\cs
ts.plot(cs) t=1:40
lm2=lm(cs~t)
summary(lm2) # 返回拟合参数的统计量 coef(lm2) #返回被估计的系数
fit2=as.ts(fitted(lm2)) #返回模拟值 residuals(lm2) #返回残差值
ts.plot(cs);lines(fit2,col=\拟合图
#example4_2 拟合非线性模型 t=1:14
x2=c(1.85,7.48,14.29,23.02,37.42,74.27,140.72,265.81,528.23,1040.27,2064.25,4113.73,8212.21,16405.95) x2
plot(t,x2)
m1=nls(x2~a*t+b^t,start=list(a=0.1,b=1.1),trace=T) summary(m1) # 返回拟合参数的统计量 coef(m1) #返回被估计的系数 fitted(m1) #返回模拟值 residuals(m1) #返回残差值
plot(t,x2);lines(t,fitted(m1)) #拟合图
#读取excel中读取文件,逗号分隔符 a=read.csv(\t=a$t x=a$x x
ts.plot(x)
m2=nls(x~a*t+b^t,start=list(a=0.1,b=1.1),trace=T) summary(m2) # 返回拟合参数的统计量 coef(m2) #返回被估计的系数 fitted(m2) #返回模拟值 residuals(m2) #返回残差值
plot(t,x);lines(t,fitted(m2)) #拟合图 #eg2
I<-scan(\I
x=ts(data=I,start=c(1991,1),f=12) #化为时间序列 x
plot.ts(x) t=1:130 t2=t^2
m3=lm(x~t+t2)
coef(m3) #返回被估计的系数
summary(m3) # 返回拟合参数的统计量
#去不显著的自变量 ,再次模拟 m4=lm(x~t2)
coef(m4) #返回被估计的系数
summary(m4) # 返回拟合参数的统计量 m2=fitted(m4) #返回模拟值 y=ts(data=m2,start=c(1991,1),f=12) y
ts.plot(x);lines(y)
#平滑法
#简单移动平均法 x=c(5,5.4,5.8,6.2) x
y=filter(x,rep(1/4,4),sides=1) y
#指数平滑 for(i in 1:3) {
x[1]=x[1]
x[i+1]=0.25*x[i+1]+0.75*x[i] }
#HoltWinters Filter
a=ts(read.csv(\a
m=HoltWinters(a,alpha=0.15,beta=0.1,gamma=FALSE,l.start=51259,b.start=4325) m
fitted(m) plot(m)
plot(fitted(m))
#综合
cs=ts(read.csv(\ #读取数据 cs
ts.plot(cs) #绘制时序图 cs.sea1=rep(0,12) cs.sea1
for(i in 1:12){
for(j in 1:8){
cs.sea1[i]=cs.sea1[i]+cs[i+12*(j-1)] } }
cs.sea=(cs.sea1/8)/(mean(cs)) cs.sea
cs.sea2=rep(cs.sea,8)
cs.sea2
x=cs/cs.sea2 x
plot(x) t=1:96
m1=lm(x~t) coef(m1)
summary(m1)
m=ts(fitted(m1),start=c(1993,1),f=12) ts.plot(x,type=\r=residuals(m1)
Box.test(r) #白噪声检验
第五讲
########################
#回顾 #例5.1
sha=ts(scan(\ts.plot(sha) diff(sha)
par(mfrow=c(2,1)) ts.plot(diff(sha)) acf(diff(sha))
#例5.2
car=ts(read.csv(\car
par(mfrow=c(3,1)) ts.plot(car) ts.plot(diff(car))
ts.plot(diff(car,differences=2))
#例5.3
milk=ts(scan(\milk
par(mfrow=c(3,1)) ts.plot(milk) ts.plot(diff(milk))
dm1=diff(diff(milk),lag=12) ts.plot(dm1) acf(dm1) #例5.5
x=ts(cumsum(rnorm(1000,0,100)))
ts.plot(x)
########################### #拟合ARIMA模型 #5.8.1
a=ts(scan(\par(mfrow=c(2,2)) ts.plot(a) da=diff(a) ts.plot(da) acf(da,20) pacf(da,20) Box.test(da,6)
fit1=arima(a,c(1,1,0),method=\predict(fit1,5)
#############################
incom=ts(read.csv(\incom
ts.plot(incom)
dincom=diff(incom) ts.plot(dincom)
acf(dincom,lag=18) #自相关图
Box.test(dincom,type=\白噪声检验 Box.test(dincom,type=\Box.test(dincom,type=\pacf(dincom,lag=18)
fit1=arima(dincom,order=c(0,0,1),method=\
fit2=arima(incom,order=c(0,1,1),xreg=1:length(incom),method=\ # 见http://www.stat.pitt.edu/stoffer/tsa2/Rissues.htm Box.test(fit2$resid,lag=6,type=\
fore=predict(fit2,10,newxreg=(length(incom)+1):(length(incom)+10)) #疏系数模型 #例5.8
w=ts(read.csv(\w=w[,1]
par(mfrow=c(2,2)) ts.plot(w) ts.plot(diff(w)) acf(diff(w),lag=18) pacf(diff(w),lag=18) dw=diff(w)
fit3=arima(dw,order=c(4,0,0),fixed=c(NA,0,0,NA,0),method=\Box.test(fit3$resid,lag=6,type=\Box.test(fit3$resid,lag=12,type=\