时间序列分析R语言程序Word文档格式.docx
《时间序列分析R语言程序Word文档格式.docx》由会员分享,可在线阅读,更多相关《时间序列分析R语言程序Word文档格式.docx(17页珍藏版)》请在冰豆网上搜索。
168]#矩阵转置转向量
plot(tdat1.3,type='
l'
#例2.2续
acf(tdat1.3)#把字去掉
pacf(tdat1.3)
#例2.3绘制1949——1998年北京市每年最高气温序列时序图
Data1.4=read.csv("
\\Users\\Administrator\\Desktop\\附录1.4.csv"
header=T)
plot(Data1.4,type='
##不会定义坐标轴
#例2.3续
tdat1.4=Data1.4[,2]
a1.4=acf(tdat1.4)
Box.test(tdat1.4,type="
Ljung-Box"
lag=6)
lag=12)
#例2.4随机产生1000个服从标准正态分布的白噪声序列观察值,并绘制时序图
Data2.4=rnorm(1000,0,1)
Data2.4
plot(Data2.4,type='
#例2.4续
a2.4=acf(Data2.4)
Box.test(Data2.4,type="
#例2.5对1950——1998年北京市城乡居民定期储蓄所占比例序列的平稳性与纯随机性进行检验
Data1.5=read.csv("
\\Users\\Administrator\\Desktop\\附录1.5.csv"
plot(Data1.5,type='
xlim=c(1950,2010),ylim=c(60,100))
tdat1.5=Data1.5[,2]
a1.5=acf(tdat1.5)
#白噪声检验
Box.test(tdat1.5,type="
#例2.5续选择合适的ARMA模型拟合序列
acf(tdat1.5)
pacf(tdat1.5)
#根据自相关系数图和偏自相关系数图可以判断为AR
(1)模型
#例2.5续P81口径的求法在文档上
#P83
arima(tdat1.5,order=c(1,0,0),method="
ML"
)#极大似然估计
ar1=arima(tdat1.5,order=c(1,0,0),method="
summary(ar1)
ev=ar1$residuals
acf(ev)
pacf(ev)
#参数的显著性检验
t1=0.6914/0.0989
p1=pt(t1,df=48,lower.tail=F)*2
#ar1的显著性检验
t2=81.5509/1.7453
p2=pt(t2,df=48,lower.tail=F)*2
#残差白噪声检验
Box.test(ev,type="
lag=6,fitdf=1)
lag=12,fitdf=1)
#例2.5续P94预测及置信区间
predict(arima(tdat1.5,order=c(1,0,0)),n.ahead=5)
tdat1.5.fore=predict(arima(tdat1.5,order=c(1,0,0)),n.ahead=5)
U=tdat1.5.fore$pred+1.96*tdat1.5.fore$se
L=tdat1.5.fore$pred-1.96*tdat1.5.fore$se
plot(c(tdat1.5,tdat1.5.fore$pred),type="
l"
col=1:
2)
lines(U,col="
blue"
lty="
dashed"
lines(L,col="
#例3.1.1例3.5例3.5续
#方法一plot.ts(arima.sim(n=100,list(ar=0.8)))
#方法二
x0=runif
(1)
x=rep(0,1500)
x[1]=0.8*x0+rnorm
(1)
for(iin2:
length(x))
{x[i]=0.8*x[i-1]+rnorm
(1)}
plot(x[1:
100],type="
acf(x)
pacf(x)
##拟合图没有画出来
#例3.1.2
x[1]=-1.1*x0+rnorm
(1)
{x[i]=-1.1*x[i-1]+rnorm
(1)}
#例3.1.3
方法一
plot.ts(arima.sim(n=100,list(ar=c(1,-0.5))))
x1=runif
(1)
x[1]=x1
x[2]=x1-0.5*x0+rnorm
(1)
for(iin3:
{x[i]=x[i-1]-0.5*x[i-2]+rnorm
(1)}
#例3.1.4
x[2]=x1+0.5*x0+rnorm
(1)
{x[i]=x[i-1]+0.5*x[i-2]+rnorm
(1)}
又一个式子
x[2]=-x1-0.5*x0+rnorm
(1)
{x[i]=-x[i-1]-0.5*x[i-2]+rnorm
(1)}
#均值和方差
smu=mean(x)
svar=var(x)
#例3.2求平稳AR
(1)模型的方差例3.3
mu=0
mvar=1/(1-0.8^2)#书上51页
#总体均值方差
cat("
populationmeanandvarare"
c(mu,mvar),"
\n"
#样本均值方差
samplemeanandvarare"
#例题3.4
svar=(1+0.5)/((1-0.5)*(1-1-0.5)*(1+1-0.5))
#例题3.6MA模型自相关系数图截尾和偏自相关系数图拖尾
#3.6.1
法一:
x=arima.sim(n=1000,list(ma=-2))
plot.ts(x,type='
法二
x=rep(0:
1000)
for(iin1:
{x[i]=rnorm[i]-2*rnorm[i-1]}
plot(x,type='
#3.6.2
x=arima.sim(n=1000,list(ma=-0.5))
{x[i]=rnorm[i]-0.5*rnorm[i-1]}
##错误于rnorm[i]:
类别为'
closure'
的对象不可以取子集
#3.6.3
x=arima.sim(n=1000,list(ma=c(-4/5,16/25)))
法二:
{x[i]=rnorm[i]-4/5*rnorm[i-1]+16/25*rnorm[i-2]}
##错误于x[i]=rnorm[i]-4/5*rnorm[i-1]+16/25*rnorm[i-2]:
##更换参数长度为零
#例3.6续根据书上64页来判断
#例3.7拟合ARMA(1,1)模型,x(t)-0.5x(t-1)=u(t)-0.8*(u-1),并直观观察该模型自相关系数和偏自相关系数的拖尾性。
#法一:
x=rep(0,1000)
x[1]=0.5*x0+rnorm
(1)-0.8*rnorm
(1)
{x[i]=0.5*x[i-1]+rnorm
(1)-0.8*rnorm
(1)}
##图和书上不一样
#法二
x=arima.sim(n=1000,list(ar=0.5,ma=-0.8))
#图和书上一样
#例3.8选择合适的ARMA模型拟合加油站57天的OVERSHORT序列
Data1.6=read.csv("
\\Users\\Administrator\\Desktop\\附录1.6.csv"
tdat1.6=as.vector(t(as.matrix(Data1.6)))[1:
57]
plot(tdat1.6,type='
acf(tdat1.6)
pacf(tdat1.6)#把字去掉
arima(tdat1.6,order=c(0,0,1),method="
CSS"
)#最小二乘估计
ma1=arima(tdat1.6,order=c(0,0,1),method="
summary(ma1)
ev=ma1$residuals
##错误于arima(tdat1.6,order=c(0,0,1),method="
):
##'
x'
必需为数值
#例3.9选择合适的ARMA模型拟合1880——1985年全球气温改变差值差分序列
##没有数据
#例3.10例3.11例3.12##矩估计
#例3.13对等时间间隔的连续70