P13
1.2 附录B中的B1，B2分别是北京地区1985~2000年的月平均气温和降水量数据，其中缺少1989年的数据，B2还缺少1995年1月数据。
（1）用简单方法补齐1989年的数据和B2中1995年1月的数据，给出季节项的周期； （2）对1990~2000年的两种数据个给出一种计算趋势项、季节项和随机项的公式； （3）利用（2）的公式对所述的数据进行时间序列的分解计算，用数据图列出结果； （4）用（2）中的结果补充1989年的数据 解：（1）对于B1：
使用分段趋势法有：（将趋势项定为年平均值）
趋势项：
B1减去趋势项剩下的部分为季节项和随机项。剩下的部分对每月取平均作为季节项，则有：
季节项：
余下的部分则为随机项：
对于B2：
同理：
（4）
对于B1：1989年的数据为当年的趋势项加季节项
对于B2：989年的数据为当月的趋势项加季节项
1.3 附录B中的B6是1973~1978年美国在意外事故中的死亡人数。 利用至少两种方法对该时间序列进行分解，要求如下： （1）画出数据图，给出数据的周期T；
数据图如下：
使用R中decompose函数进行时序分解如下：
可以看出趋势项大致为二次曲线，因此可以考虑采用二次曲线来拟合趋势项。
（2）给出趋势项、季节项和随机项的计算公式；
方法一：二次曲线
认为趋势项满足二次线性方程，由最小二乘公式计算出趋势项如图：
B6减去趋势项剩下的部分为季节项和随机项。剩下的部分对每月取平均作为季节项，剩下的部分即为随机项。
方法二：回归方法（多元）
认为趋势项满足多元线性方程，由最小二乘公式，计算出，再计算得出趋势项、季节项和随机项：
（3）画出趋势项、季节项和随机项的数据图； 方法一： 趋势项：
季节项：
随机项：
方法二：
趋势项、季节项：
随机项：
（4）对1979年的意外死亡人数做出预测。 方法一：
预测：1979年的数据为当月的趋势项加季节项：
预测图：
方法二：
预测：1979年的数据为回归方程结果：
预测图：
程序代码：
########### #B1
B1=matrix(c(-4.7,-3.7,-3.6,-2.9,-4.9,-2.3,-1.1,-3.7,-1.6,-0.7,-2.2,-3.8,-3.9,-1.6,-6.4,-1.9,-1.8,0.1,-1.4,-0.6,0.1,1.8,1.6,0.8,2.1,-0.4,1.3,2.4,2.1,-1.5,3.4,6.9,4.1,4.4,7.6,4.4,6.7,8.1,5.6,7.7,6.2,8.7,7.6,4.7,8.0,14.8,15.0,13.5,15.0,13.7,13.9,15.5,14.0,17.3,14.7,14.3,14.5,15.0,14.4,14.6,19.5,21.3,19.9,20.1,19.6,19.9,20.5,21.5,21.0,19.8,21.6,20.0,19.9,19.3,20.4,24.2,25.3,23.3,24.9,24.8,24.1,23.5,25.4,26.8,24.3,25.4,24.6,23.6,25.3,26.7,25.5,25.1,26.6,25.8,25.6,25.9,26.8,25.2,27.7,25.9,25.5,28.2,26.5,28.0,29.6,25.0,24.5,24.8,24.4,25.4,27.1,24.6,25.2,26.5,25.4,23.9,26.6,25.1,25.5,25.7,18.6,19.8,21.0,21.2,20.2,20.4,20.5,21.3,21.1,19.0,20.7,18.6,22.2,20.9,21.8,13.8,11.4,13.7,14.1,15.3,13.8,12.2,13.9,14.1,14.5,12.8,14.0,14.8,12.9,12.6,3.8,3.4,3.9,6.9,6.4,4.6,3.4,3.7,6.4,7.7,4.2,5.4,4.0,5.9,3.0,-3.6,-1.7,-0.3,-0.2,-0.8,-1.8,-0.3,-0.8,-1.4,-0.4,0.9,-1.5,0.1,-0.7,-0.6),ncol=12) B1.mean=apply(B1,2,mean)
B1=rbind(B1[1:4,],B1.mean,B1[5:15,],deparse.level =0) B1.trend=apply(B1,1,mean)#trend
plot(rep(c(t(B1.trend)),rep(12,16)),type=\"l\B1.temp=B1-B1.trend
B1.season=apply(B1.temp,2,mean)
B1.season=B1.season-sum(B1.season)/12#seasonal plot(rep(c(t(B1.season)),16),type=\"l\B1.rand=t(t(B1.temp)-B1.season)#random plot(c(t(B1.rand)),type=\"h\
# transform matrix into time series B1=c(t(B1))
B1.ts=ts(data=B1,start=1985,frequency=12) plot(B1.ts)
# use decompose function Decom.B1=decompose(B1.ts) plot(Decom.B1)
########### #B2
B2=c(1.5,7.5,7.8,13.6,24.5,32.0,289.5,297.7,38.4,3.8,4.6,0.1,0.0,6.0,17.0,1.0,5.0,203.0,163.0,143.0,114.0,4.0,6.0,4.0,4.3,2.4,13.0,41.8,64.6,91.2,130.9,246.5,46.2,4.1,35.4,3.5,0.9,1.3,8.9,8.2,37.4,61.8,278.7,204.0,48.8,22.8,0.0,0.5,4.7,21.6,40.5,59.7,119.6,4.0,223.0,157.0,63.1,0.3,3.6,0.2,0.3,0.8,25.1,17.1,214.6,236.3,198.0,124.7,72.0,12.2,1.0,4.7,0.7,0.0,3.4,10.5,52.8,69.4,153.9,141.4,54.5,38.1,16.7,0.1,3.7,1.5,0.3,16.9,8.6,39.2,206.4,158.5,18.3,9.9,43.4,0.0,0.0,5.0,0.0,1.9,66.0,23.6,459.2,214.2,15.2,10.3,12.7,5.1,NA,1.7,6.6,5.3,45.6,68.9,195.6,119.9,116.3,9.6,0.2,2.8,0.2,0.0,11.0,6.2,1.8,55.1,307.4,250.0,32.9,30.8,2.6,2.9,4.9,0.0,10.6,17.4,41.5,35.5,139.8,83.2,44.1,43.0,2.1,8.8,1.3,26.3,4.3,54.7,61.5,142.9,247.9,114.4,4.7,61.8,11.3,0.6,0.0,0.0,5.2,33.6,32.4,23.8,62.7,63.5,44.5,3.9,9.5,0.7,11.9,0.0,8.8,18.3,37.7,19.0,61.5,150.5,18.4,35.2,9.7,0.1) B2=matrix(B2,ncol=12,byrow=TRUE)
#fill the missing value
B2.mean=apply(B2,2,mean,na.rm=TRUE) B2=rbind(B2[1:4,],B2.mean,B2[5:15,]) B2=c(t(B2))
B2[is.na(B2)]=B2.mean[1]
# transform matrix into time series
B2.ts=ts(data=B2,start=1985,frequency=12) plot(B2.ts)
# use decompose function Decom.B2=decompose(B2.ts) plot(Decom.B2)
# least-squares fitting
Y=matrix(c(rep(1,192),1:192),ncol=192,byrow=TRUE) Result=solve(Y%*%t(Y))%*%Y%*%B2 B2.trend=c(t(Y)%*%Result)#trend plot(B2.trend,type=\"l\B2.temp=B2-B2.trend
B2.temp=matrix(B2.temp,ncol=12,byrow=TRUE) B2.season=apply(B2.temp,2,mean)
B2.season=B2.season-sum(B2.season)/12#seasonal plot(rep(c(t(B2.season)),16),type=\"l\B2.rand=t(t(B2.temp)-B2.season)#random plot(c(t(B2.rand)),type=\"h\
########### #(4) #B1
B1.1989=B1.trend[5]+B1.season #B2
B2.1989=B2.trend[(4*12+1):(4*12+12)]+B2.season
################# #Question 1.3
B6=c(9007,8106,8928,9137,10017,10826,11317,10744,9713,9938,9161,8927,7750,6981,8038,8422,8714,9512,10120,9823,8743,9192,8710,8680,8162,7306,8124,7870,9387,9556,10093,9620,8285,8433,8160,8034,7717,7461,7776,7925,8634,8945,10078,9179,8037,8488,7874,8647,7792,6957,7726,8106,8890,9299,10625,9302,8314,8850,8265,8796,7836,6892,7791,8129,9115,9434,10484,9827,9110,9070,8633,9240) B6=matrix(B6,ncol=12,byrow=TRUE) B6.mean=apply(B6,2,mean) B6=c(t(B6))
# transform matrix into time series
B6.ts=ts(data=B6,start=1973,frequency=12) plot(B6.ts)
# use decompose function Decom.B6=decompose(B6.ts) plot(Decom.B6)
###########
# least-squares fitting(quadratic)
Y=matrix(c(rep(1,72),1:72,1:72*1:72),ncol=72,byrow=TRUE) Result=solve(Y%*%t(Y))%*%Y%*%B6
B6.trend=c(t(Y)%*%Result)#trend plot(B6.trend,type=\"l\B6.temp=B6-B6.trend
B6.temp=matrix(B6.temp,ncol=12,byrow=TRUE) B6.season=apply(B6.temp,2,mean)
B6.season=B6.season-sum(B6.season)/12#seasonal plot(rep(c(t(B6.season)),6),type=\"l\B6.rand=t(t(B6.temp)-B6.season)#random plot(c(t(B6.rand)),type=\"h\
# predict 1979
Y.1979=matrix(c(rep(1,12),73:84,73:84*73:84),ncol=12,byrow=TRUE) B6.1979_1=c(t(Y.1979)%*%Result)+B6.season
plot(1:72,B6,type=\"l\lines(73:84,B6.1979_1,col=\"BLUE\
###########
# least-squares fitting(multivariate) Y2=matrix(rep(diag(12),6),ncol=72) Y2=rbind(t(1:72),Y2,deparse.level =0) Result2=solve(Y2%*%t(Y2))%*%Y2%*%B6
B6.trend2=c(t(Y2)%*%Result2)#trend+seasonal plot(B6.trend2,type=\"l\B6.rand2=B6-B6.trend2#random
B6.rand2=matrix(B6.rand2,ncol=12,byrow=TRUE) plot(c(t(B6.rand2)),type=\"h\
# predict 1979
Y2.1979=rbind(t(73:84),diag(12),deparse.level =0) B6.1979_2=c(t(Y2.1979)%*%Result2)
plot(1:72,B6,type=\"l\lines(73:84,B6.1979_2,col=\"BLUE\
因篇幅问题不能全部显示，请点此查看更多更全内容