r语言|统计学--基于R(第3版)(基于R应用的统计学丛书)作者(贾俊平 习题答案 第十一章)

11.1

#11.1 load("C:/exercise/ch11/exercise11_1.RData") exercise11_1 #采用指数平滑法预测2016年的PPI,并对实际值和预测值绘图进行比较 #PPI的简单指数平滑预测 exercise11_1<-ts(exercise11_1,start=2006) cpiforecast<-HoltWinters(exercise11_1[,2],beta=FALSE,gamma=FALSE) cpiforecast #历史数据的拟合值 cpiforecast$fitted #绘制观测值和拟合值图 plot(exercise11_1[,2],type='o',xlab="时间",ylab="PPI") lines(exercise11_1[,1][-1],cpiforecast$fitted[,1],type="o",lty=2,col="blue") legend(x="topleft",legend=c("观测值","拟合值"),lty=1:2,col=c(1,4)) #获得样本外的预测值(2016年) library(forecast) cpiforecast1<-forecast(cpiforecast,h=1) cpiforecast1 #实际值和预测值(2016年)图 plot(cpiforecast1,type="o",xlab="时间",ylab="PPI",main="")

11.2
#11.2 load("C:/exercise/ch11/exercise11_2.RData") exercise11_2 #分别采用一元线性回归模型和Holt指数平滑模型预测2011年的原油产量,并绘制残差图对两种方法的预测效果进行比较 model1<-lm(原油产量~年份,data=https://www.it610.com/article/exercise11_2) summary(model1) #y^=-571251.54+293.90x1 h<--571251.54+293.90*2011 h#即,一元回归模型预测出2011年原油产量大概为19781.36 #确定模拟参数α和β以及模型系数a和b exercise11_2<-ts(exercise11_2,start=1990) grainforecast<-HoltWinters(exercise11_2[,2],gamma=FALSE) grainforecast #绘制拟合图 plot(exercise11_2[,2],type="o",xlab="年份",ylab="原油产量") lines(exercise11_2[,1][c(-1,-2)],grainforecast$fitted[,1],type="o",lty=2,col="blue") legend(x="topleft",legend=c("观测值","拟合值"),lty=1:2) #2011年原油产量 library(forecast) grainforecast1<-forecast(grainforecast,h=1) grainforecast1 #实际值和预测值图 plot(grainforecast1,type="o",xlab="年份",ylab="原油产量",main="") #预测值为20309.77 #残差比较 residual1<-grainforecast1$residuals residual2<-model1$residuals plot(as.numeric(residual1),ylim=c(-2000,2000),xlab="",ylab="残差",pch=1) points(residual2,pch=2,col="red") abline(h=0) legend(x="topright",legend=c("Holt模型预测的残差","一元线性回归预测的残差"),pch=1:2,col=c("black","red"),cex=0.8)

11.3
#11.3 load("C:/exercise/ch11/exercise11_3.RData") exercise11_3 #分别采用Holt指数平滑模型和指数模型预测2009年的财政收入,并对实际值和预测值绘图进行比较 #财政收入的Holt指数平滑预测 exercise11_3<-ts(exercise11_3,start=1991) grainforecast<-HoltWinters(exercise11_3[,2],gamma=FALSE) grainforecast #绘制拟合图 plot(exercise11_3[,2],type="o",xlab="年份",ylab="财政收入") lines(exercise11_3[,1][c(-1,-2)],grainforecast$fitted[,1],type="o",lty=2,col="blue") legend(x="topleft",legend=c("观测值","拟合值"),lty=1:2) #2009年财政收入预测 library(forecast) grainforecast1<-forecast(grainforecast,h=1) grainforecast1 #71943.02 #实际值和预测值图 plot(grainforecast1,type="o",xlab="年份",ylab="财政收入",main="") #指数模型预测 #指数曲线拟合 y<-log(exercise11_3[,2]) x<-1:18 fit<-lm(y~x) fit exp(7.8486) #历史数据及2009年财政收入的预测 predata<-exp(predict(fit,data.frame(x=1:19))) predata #65345.638 #各年的预测残差 predata<-exp(predict(fit,data.frame(x=1:18))) predata<-ts(predata,start=1991) residuals<-exercise11_3[,2]-predata residuals #实际值和预测值的比较图 predata<-exp(predict(fit,data.frame(x=1:19))) plot(1991:2009,predata,type="o",lty=2,col="blue",xlab="年份",ylab="财政收入") points(exercise11_3[,2],type="o",pch=19) legend(x="topleft",legend=c("观测值","预测值"),lty=1:2) abline(v=2008,lty=6,col="gray") ##残差图 plot(1991:2008,residuals,type="o",lty=2,xlab="年份",ylab="residuals") abline(h=0)

11.4
#10.4 load("C:/exercise/ch11/exercise11_4.RData") exercise11_4 #分别拟合回归直线,二次曲线,三次曲线,并绘制图形对结果进行比较 #回归直线 fit<-lm(盘价格~时间,data=https://www.it610.com/article/exercise11_4) summary(fit) #Y1^=374.16134-0.61373t #二次曲线 y<-exercise11_4[,2] t<-1:35 fit1<-lm(y~t+I(t^2)) fit1 #Y2^=381.64416-1.82715t+0.03371t^2 #三次曲线 fit2<-lm(y~t+I(t^2)+I(t^3)) fit2 #Y3^=372.561669+1.002958t-0.160088t^2+0.003589t^3 #图 predata<-predict(fit,data.frame(时间=1:36))#各年的预测值 predata res<-residuals(fit)#各年的预测残差 res predata1<-predict(fit1,data.frame(t=1:36))#二次曲线预测值 predata1 residual1<-fit1$residuals#二次曲线预测的残差 residual1 predata2<-predict(fit2,data.frame(t=1:36))#三次曲线预测值 predata2 residual2<-fit2$residuals#三次曲线预测的残差 residual2 #实际值和预测值曲线 predata<-predict(fit,data.frame(时间=1:36)) predata1<-predict(fit1,data.frame(t=1:36)) predata2<-predict(fit2,data.frame(t=1:36)) plot(1:36,predata2,type="o",lty=2,col="red",xlab="时间",ylim=c(340,400),ylab="盘价格") lines(1:36,predata1,type="o",lty=3,col="green",xlab="时间",ylim=c(340,400),ylab="盘价格") lines(1:36,predata,type="o",lty=4,col="blue") points(exercise11_4[,1],exercise11_4[,2],type="o",pch=19) legend(x="bottom",legend=c("观测值","回归直线","二次曲线","三次曲线"),lty=1:4,col=c("black","blue","green","red")) abline(v=35,lty=6)

11.5
#11.5 load("C:/exercise/ch11/exercise11_5.RData") exercise11_5 #分别使用Winter指数平滑模型和分解法预测2014年各月份的社会消费品零售总额,并绘制实际值和预测值图,分析预测效果 #Winter指数平滑模型 #把数据转化成以月份为周期的时间序列格式 exercise11_5<-ts(exercise11_5[,3],start=2009,frequency=12) exercise11_5 #确定模型参数的α、β和γ以及模型的a,b和s saleforecast<-HoltWinters(exercise11_5) saleforecast #Winter模型的拟合图 plot(exercise11_5,type="o",xlab="时间",ylab="零售总额") lines(saleforecast$fitted[,1],type="o",lty=2,col="blue") legend(x="topleft",legend=c("观测值","拟合值"),lty=1:2,col=c(1,4)) #Winter模型2014年的预测 library(forecast) saleforecast1<-forecast(saleforecast,h=12) saleforecast1 #Winter模型的预测图 plot(saleforecast1,type="o",xlab="时间",ylab="零售总额",main="") abline(v=2014,lty=6,col="gray") #分解预测 #计算月份指数 salecompose<-decompose(exercise11_5,type="multiplicative") names(salecompose) salecompose$seasonal #月份调整后的序列图 seasonaladjust<-exercise11_5/salecompose$seasonal plot(exercise11_5,xlab="时间",ylab="零售总额",type="o",pch=19) lines(seasonaladjust,lty=2,type="o",col="blue") legend(x="topleft",legend=c("零售总额","零售总额的月份调整"),lty=1:2) #月份调整后序列的线性模型 x<-1:60 fit<-lm(seasonaladjust~x) fit #最终预测值 predata<-predict(fit,data.frame(x=1:72))*rep(salecompose$seasonal[1:12],6) predata<-ts(predata,start=2009,frequency=12) predata #预测的残差 residuals1<-exercise11_5-predict(fit,data.frame(x=1:60))*salecompose$seasonal residuals1<-ts(residuals1,start=2009,frequency=12) round(residuals1,4) #实际值和预测值的比较图 predata<-predict(fit,data.frame(x=1:72))*rep(salecompose$seasonal[1:12],6) predata<-ts(predata,start=2009,frequency=12) plot(predata,type="o",lty=2,col="blue",xlab="时间",ylab="零售总额") lines(exercise11_5) legend(x="topleft",legend=c("实际零售总额","预测零售总额"),lty=1:2,col=c("black","blue"),cex=0.7) abline(v=2014,lty=6,col="grey") #分解预测与Winter模型预测残差的比较 residuals1<-exercise11_5-predict(fit,data.frame(x=1:60))*salecompose$seasonal saleforecast<-HoltWinters(exercise11_5) residuals2<-exercise11_5[-(1:12)]-saleforecast$fitted[,1] par(cex=0.7,mai=c(0.4,0.7,0.1,0.1)) plot(as.numeric(residuals1),xlab="",ylab="残差",ylim=c(-800,800),pch=1) abline(h=0) points(as.numeric(residuals2),pch=2,col="red") abline(h=0) legend(x="topright",legend=c("分解预测的残差","Winter模型预测的残差"),pch=1:2,col=1:2,cex=0.8)

【r语言|统计学--基于R(第3版)(基于R应用的统计学丛书)作者(贾俊平 习题答案 第十一章)】11.6
#11.6 #利用11.4的数据,分别采用m=5和m=10对收盘价格进行平滑,并绘制出实际值与平滑值的图形进行比较 load("C:/exercise/ch11/exercise11_4.RData") exercise11_4 exercise11_4<-ts(exercise11_4,start=1) library(DescTools) ma5<-MoveAvg(exercise11_4[,2],order=5,align="center",endrule="keep") ma10<-MoveAvg(exercise11_4[,2],order=10,align="center",endrule="NA") y1<-exercise11_4[,1] y2<-exercise11_4[,2] data.frame("时间"=y1,"盘价格"=y2,ma5,ma10) #绘制实际值和平滑值的比较图形 par(mai=c(0.7,0.7,0.1,0.1),cex=0.8) plot(exercise11_4[,2],type="o",xlab="时间",ylab="盘价格") lines(ma5,type="o",lty=2,col="red") lines(ma10,type="o",lty=2,col="blue") legend(x="topleft",legend=c("盘价格","ma5","ma10"),lty=c(1,2,6),col=c("black","red","blue"),cex=0.8)


    推荐阅读