R语言时间序列作业.docx
《R语言时间序列作业.docx》由会员分享,可在线阅读,更多相关《R语言时间序列作业.docx(17页珍藏版)》请在冰豆网上搜索。
![R语言时间序列作业.docx](https://file1.bdocx.com/fileroot1/2023-1/5/c0d38f62-2975-4af4-b4c5-87f948f43536/c0d38f62-2975-4af4-b4c5-87f948f435361.gif)
R语言时间序列作业
2016年第二学期时间序列分析及应用R语言课后作业
第三章趋势
3.4(a)data(hours);plot(hours,ylab='MonthlyHours',type='o')
画出时间序列图
(b)data(hours);plot(hours,ylab='MonthlyHours',type='l')
type='o'表示每个数据点都叠加在曲线上;type='b'表示在曲线上叠加数据点,但是该数据点附近是断开的;type='l'表示只显示各数据点之间的连接线段;type='p'只想显示数据点。
points(y=hours,x=time(hours),pch=as.vector(season(hours)))
3.10(a)data(hours);hours.lm=lm(hours~time(hours)+I(time(hours)^2));summary(hours.lm)
用最小二乘法拟合二次趋势,结果显示如下:
Call:
lm(formula=hours~time(hours)+I(time(hours)^2))
Residuals:
Min1QMedian3QMax
-1.00603-0.25431-0.022670.228840.98358
Coefficients:
EstimateStd.ErrortvaluePr(>|t|)
(Intercept)-5.122e+051.155e+05-4.4334.28e-05***
time(hours)5.159e+021.164e+024.4314.31e-05***
I(time(hours)^2)-1.299e-012.933e-02-4.4284.35e-05***
---
Signif.codes:
0‘***’0.001‘**’0.01‘*’0.05‘.’0.1‘’1
Residualstandarderror:
0.423on57degreesoffreedom
MultipleR-squared:
0.5921,AdjustedR-squared:
0.5778
F-statistic:
41.37on2and57DF,p-value:
7.97e-12
(b)plot(y=rstudent(hours.lm),x=as.vector(time(hours)),type='l',ylab='StandardizedResiduals')
points(y=rstudent(hours.lm),x=as.vector(time(hours)),pch=as.vector(season(hours)))
标准残差的时间序列,应用月度绘图标志。
(为了更容易识别季节性)
带季节性图标的的残差-时间图
(c)runs(rstudent(hours.lm))
对标准差进行游程检验
$pvalue
[1]0.00012
$observed.runs
[1]16
$expected.runs
[1]30.96667
$n1
[1]31
$n2
[1]29
$k
[1]0
结果解释:
P值为0.00012,表明非随机性是合理的。
(d)acf(rstudent(hours.lm))
标准残差的样本自相关函数
季节均值模型残差的样本自相关系数
(e)qqnorm(rstudent(hours.lm));qqline(rstudent(hours.lm))
(QQ图)
正态性可以通过正态得分或者分位数-分位数(QQ)图来检验。
此处的直线型图形支持了该模型中随机项是正态分布的假设。
hist(rstudent(hours.lm),xlab='StandardizedResiduals')
标准残差的直方图(季节均值模型的标准残差直方图)
shapiro.test(rstudent(hours.lm))
正态性检验(Shapiro-Wilk检验)本质是:
计算残差与相应的正态分位数之间的相关系数。
相关性越小,就越有理由否定正态性。
Shapiro-Wilknormalitytest
data:
rstudent(hours.lm)
W=0.99385,p-value=0.9909
根据上面的检验结果,我们不能拒绝模型的随机项是正态分布的假设。
第四章平稳时间序列模型
4.4
第五章非平稳时间序列模型
5.1(a)ARMA(2,1)p=2,q=1,参数值φ和θ
φ1=1φ2=-0.25Θ1=0.1
(b)IMA(2,0)p=2,d=1,q=0,参数值φ和θ
(c)ARMA(2,2)p=2,q=2,参数值φ和θ
φ1=0.5φ2=-0.5Θ1=0.5Θ2=-0.25
5.7(a)A:
AR
(2)φ1=0.9φ2=0.09
B:
IMA(1,1)Θ1=0.1
(b)一个是固定的一个是不固定的。
5.11(a)data(winnebago);win.graph(width=6.5,height=3,pointsize=8)
plot(winnebago,type='o',ylab='WinnebagoMonthlySales')
时间序列图
表明公司的休闲车的销量在逐渐增加。
(b)plot(log(winnebago),type='o',ylab='Log(MonthlySales)')
取对数之后的时间序列图
仍然呈现增加的趋势,但是比没有取对数之前增加的缓慢一些。
(c)percentage=na.omit((winnebago-zlag(winnebago))/zlag(winnebago))
win.graph(width=3,height=3,pointsize=8)
plot(x=diff(log(winnebago))[-1],y=percentage[-1],ylab='PercentageChange',xlab='DifferenceofLogs')
cor(diff(log(winnebago))[-1],percentage[-1])
[1]0.9646886
结果显示:
0.96认为一致。
第六章模型识别
6.29(a)set.seed(762534);series=arima.sim(n=60,list(ar=0.4,ma=0.6))
phi=0.4;theta=0.6;ACF=ARMAacf(ar=phi,ma=-theta,lag.max=10)
plot(y=ACF[-1],x=1:
10,xlab='Lag',ylab='ACF',type='h',ylim=c(-.2,.2));abline(h=0)
(b)acf(series)
第八章模型诊断
8.7(a)data(hare);model=arima(sqrt(hare),order=c(3,0,0))
win.graph(width=6.5,height=3,pointsize=8);acf(rstandard(model))
残差的样本自相关图
(b)LB.test(model,lag=9)
Box-Ljungtest
data:
residualsfrommodel
X-squared=6.2475,df=6,p-value=0.396
JB统计量结果表明不拒绝误差项的独立性。
(c)对残差进行检验。
runs(rstandard(model))
$pvalue
[1]0.602
$observed.runs
[1]18
$expected.runs
[1]16.09677
$n1
[1]13
$n2
[1]18
$k
[1]0
P值为0.602,不拒绝误差项的独立性。
(d)win.graph(width=3,height=3,pointsize=8)
qqnorm(residuals(model))
残差的正态QQ图
QQ图看出有一点小的曲率,但是这种现象可能是俩个极端值造成的。
(e)对残差的正态性进行shapiro-wilk检验
shapiro.test(residuals(model))
Shapiro-Wilknormalitytest
data:
residuals(model)
W=0.93509,p-value=0.06043
结果表明,我们不会拒绝通常意义水平的正态性。
第九章预测
9.2(a)
所以,
(b)从上式中看出
,
(c)
即2008年预测的95%预测极限为7.67-13.33.
(d)因为有,
所以,
第十章季节模型
10.12(a)data(boardings);series=boardings[,1]
plot(series,type='l',ylab='LightRail&BusBoardings')
points(series,x=time(series),pch=as.vector(season(series)))
(b)acf(as.vector(series),ci.type='ma')
滞后期为1,5,6,12,时候,存在显著的自相关。
(c)model=arima(series,order=c(0,0,3),seasonal=list(order=c(1,0,0),period=12));model
Call:
arima(x=series,order=c(0,0,3),seasonal=list(order=c(1,0,0),period=12))
Coefficients:
ma1ma2ma3sar1intercept
0.72900.61160.29500.877612.5455
s.e.0.11860.11720.11180.05070.0354
sigma^2estimatedas0.0006542:
loglikelihood=143.54,aic=-277.09
所有的变量都是显著的。
(d)model2=arima(series,order=c(0,0,4),seasonal=list(order=c(1,0,0),period=12));model2
Call:
arima(x=series,order=c(0,0,4),seasonal=list(order=c(1,0,0),period=12))
Coefficients:
ma1ma2ma3ma4sar1intercept
0.72770.66860.42440.14140.891812.5459
s.e.0.12120.13270.16810.12280.04450.0419
sigma^2estimatedas0.0006279:
loglikelihood=144.22,aic=-276.45
模型2中AIC=-276.45,模型1中的AIC=-277.09,AIC越小越好,所以模型是过度拟合的。
第12章异方差时间序列模型
12.1library(TSA)
data(CREF)
r.cref=diff(log(CREF))*100
win.graph(width=4.875,height=2.5,pointsize=8)
plot(abs(r.cref))
win.graph(width=4.875,height=2.5,pointsize=8)
plot(r.cref^2)
12.9(a)>data(google)
>plot(google)
收益率数据的时间序列图
>acf(google)
>pacf(google)
根据ACF和PACF可以得知,无自相关。
(b)计算google日收益率均值。
>t.test(google,alternative='greater')
OneSamplet-test
data:
google
t=2.5689,df=520,p-value=0.00524
alternativehypothesis:
truemeanisgreaterthan0
95percentconfidenceinterval:
0.000962967Inf
sampleestimates:
meanofx
0.002685589
x的均值为0.002685589,备择假设为:
均值异于0,根据P值显示,0.00524接受备择假设。
(c)McLeod-Li检验ARCH效应。
>win.graph(width=4.875,height=3,pointsize=8)
>McLeod.Li.test(y=google)
根据图显示,所有滞后值在5%的水平上均显著。
说明数据具有ARCH特征。
(d)识别GARCH模型,估计识别的模型并对拟合的模型进行模型诊断检验。
>eacf(google^2)
取值平方的样本EACF
AR/MA
012345678910111213
0xxoooooooxooox
1xooooooooxooox
2xooooooooxooox
3xxxooooooxooox
4xxxooooooooooo
5xxxooooooooooo
6xxxxoooooooooo
7oxxooxoooooooo
>eacf(abs(google))
绝对值的样本EACF
AR/MA
012345678910111213
0xxxoooxooxooxx
1xoooooooooooox
2xxooooooooooox
3xxxoooooooooox
4xoxooooooooooo
5xoxoxooooooooo
6oxxxxxoooooooo
7xoxxxoxooooooo
根据上图,得知设定GARCH(1,1)模型。
Google日收益率的平方值相应的样本EACF也得知模型GARCH(1,1)符合。
检验模型:
>m1=garch(x=google,order=c(1,1))
>summary(m1)
Call:
garch(x=google,order=c(1,1))
Model:
GARCH(1,1)
Residuals:
Min1QMedian3QMax
-3.64597-0.464860.082320.653795.73937
Coefficient(s):
EstimateStd.ErrortvaluePr(>|t|)
a05.058e-051.232e-054.1064.03e-05***
a11.264e-012.136e-025.9203.21e-09***
b17.865e-013.578e-0221.980<2e-16***
---
Signif.codes:
0‘***’0.001‘**’0.01‘*’0.05‘.’0.1‘’1
DiagnosticTests:
JarqueBeraTest
data:
Residuals
X-squared=223.86,df=2,p-value<2.2e-16
Box-Ljungtest
data:
Squared.Residuals
X-squared=0.00066246,df=1,p-value=0.9795
模型诊断:
对google日收益率拟合的GARCH(1,1)模型的标准残差
>plot(residuals(m1),type='h',ylab='StandardizedResiduals')
标准残差的QQ正态得分图
>win.graph(width=2.5,height=2.5,pointsize=8)
>qqnorm(residuals(m1));qqline(residuals(m1))
如果模型识别正确,那么残差应该是近似独立同分布的。
从上qq图中显示并不明确。
所以我们进行广义混合检验。
绘制google日收益率的GARCH(1,1)模型的标准残差平方的样本ACF
acf(residuals(m1)^2,na.action=na.omit)
从图形中得出总体印象是残差平方序列不相关。
进行广义混合检验得到的p值
gBox(m1,method='squared')
结果显示错误!
!
!
不清楚为什么!
!
!
在这里继续使用绝对标准残差重新对模型进行检验。
绝对标准残差的样本ACF
acf(abs(residuals(m1)),na.action=na.omit)
gBox(m1,method='absolute')
(e)绘制并评论估计的条件方差的时间序列图
图中显示出来有几个时期的波动率较高。
(f)对拟合模型的标准残差绘制QQ图
qqnorm(residuals(m1));qqline(residuals(m1))
残差看起来并不正常。
(g)构造b1的95%的置信区间。
(0.7865-1.96*0.03578,0.7865+1.96*0.03578)=(0.7164,0.8566)
欢迎您的下载,
资料仅供参考!
致力为企业和个人提供合同协议,策划案计划书,学习资料等等
打造全网一站式需求