《统计学实验》课件和数据第15章时间序列分析R

2026/2/2《统计学实验》第15章时间序列分析15-1统计学实验

—SPSS和R软件应用与实例

主编：费宇2026/2/2《统计学实验》第15章时间序列分析15-2第15章时间序列分析2026/2/2《统计学实验》第15章时间序列分析15-3一、实验目的1.运用R中plot()函数绘制时间序列图，acf()和pacf()函数作时间序列数据的自相关和偏自相关图检验序列的平稳性，HoltWinters()滤波函数实现时间序列数据的指数平滑；一、实验目的运用arima()函数进行ARIMA模型的识别、估计和预测，用dpose()函数作时间序列的季节分解；掌握运用R软件中的函数对时间序列数据进行分析的基本操作过程，并能读懂R输出的结果。2026/2/2《统计学实验》第15章时间序列分析15-42026/2/2《统计学实验》第15章时间序列分析15-5二、实验环境系统软件Windows2000或WindowsXP或Windows7；统计软件R2.13.2或更高版本。2026/2/2《统计学实验》第15章时间序列分析15-6三、实验内容检验时间序列的平稳性根据自相关系数和偏自相关系数来初步判断ARMA模型的阶数p和q对ARIMA模型进行估计、预测和诊断识别时间序列数据的季节变动，能看出其季节波动趋势，学会剔除季节因素的方法2026/2/2《统计学实验》第15章时间序列分析15-7第15章时间序列分析15.1时间序列的图形化观察15.2指数平滑15.3时间序列的Box-Jenkin

模型15.4时间序列的季节分解2026/2/2《统计学实验》第15章时间序列分析15-815.1时间序列的图形化观察【例15.1】(数据文件为li15.1.txt)选取我国从1978年至2008年国内生产总值(GDP，亿元)数值，总31个观测值，数据来源“国泰安宏观经济数据库”。试对该数据进行分析，判断其是否具有平稳性。2026/2/2《统计学实验》第15章时间序列分析15-9【统计理论】时间序列分析中一个重要的概念是平稳性，它是时间序列建模的一个重要基础。如果一个随机过程的均值和方差在时间过程上都是常数，并且在任何两时期的协方差仅依赖于这两个时期间的间隔，而不依赖于计算这个协方差的实际时间，就称它是宽平稳的。2026/2/2《统计学实验》第15章时间序列分析15-10【统计理论】时间序列平稳性检验的常用方法有：

(1)时序图依据均值和方差的统计意义，一个平稳序列时间序列的实现大致是由在某一水平线附近等幅波动的点构成的。平稳序列的时序图应该显示出序列始终围绕一个常数值波动，且波动的范围不大。2026/2/2《统计学实验》第15章时间序列分析15-11【统计理论】实践中真正具有平稳性的时间序列并不多见，通常都会表现出非平稳性。非平稳性的表现形式多种多多样，主要特征有：趋势性、异方差性、等。在间序列数据中，差分可以使时间序列达到平稳性目的。2026/2/2《统计学实验》第15章时间序列分析15-12假设

为一时间序列，则该序列的一阶差分定义为

，其中

表示一阶差分算子。为了达到平稳性的目的，一阶差分是最常用的方法，有时我们甚至需要对其进行多次一阶差分。【统计理论】2026/2/2《统计学实验》第15章时间序列分析15-13【统计理论】

(2)自相关和偏自相关系数所谓自相关是指序列与其自身经过某些阶数滞后形成的序列之间存在某种程度的相关性。随着

的增加，自相关函数下降且趋于0。从下降的速度看，平稳序列要比非平稳序列快得多。2026/2/2《统计学实验》第15章时间序列分析15-14【软件操作】首先，用plot()函数绘制时间序列图setwd("D:/R-Statistics/data/chap-15")#设定工作路径x=read.table("li15.1.txt",header=T)#从li15.1.txt中读入样本数据xdx=diff(x[,2])#对DGP作一阶差分变换【软件操作】plot(x[,1],x[,2],type=“1",lty=1,col=1,xlab="时间",ylab="")lines(x[2:31,1],dx,lty=2,col=1)legend(1980,250000,c("GDP","GDP_1"),col=c(1,1),text.col="black",lty=c(1,2))2026/2/2《统计学实验》第15章时间序列分析15-152026/2/2《统计学实验》第15章时间序列分析15-16【运行结果】图15.1例15.1中GDP及GDP_1的时序图【软件操作】然后，采用acf()和pacf()函数作相关图x=read.table("li15.1.txt",header=T)#从li15.1.txt中读入样本数据xGDP=x[,2]#从数据文件中读取GDP数据acf(GDP)#GDP序列的自相关函数条形图2026/2/2《统计学实验》第15章时间序列分析15-17【运行结果】2026/2/2《统计学实验》第15章时间序列分析15-18图15.2例15.1中GDP的自相关函数的条形图2026/2/2《统计学实验》第15章时间序列分析15-1915.2指数平滑【例15.2】(数据文件为li15.2.txt)选取1990年至2004年我国月度社会消费品零售总额数据(亿元)，利用指数平滑方法试对该数据进行分析。

2026/2/2《统计学实验》第15章时间序列分析15-20【统计理论】当利用过去数据的加权平均预测将来的数值时(这一过程成为平滑)，对序列中较近的数据给予较大的权重，远期的数据给于较小的权重。令

，则

的预测值为：即利用过去的观测值进行加权平均预测时间序列在

时刻的数值。【软件操作】使用HoltWinters()滤波函数可以作时间序列数据的指数平滑x=read.table("li15.2.txt",header=T)#从li15.2.txt中读入样本数据xtx=ts(x,frequency=12,start=c(1990,1),end=c(2004,8))#添加合适的时间格式2026/2/2《统计学实验》第15章时间序列分析15-21【软件操作】m=HoltWinters(tx,gamma=FALSE,beta=FALSE,seasonal="mult")#参数gamma=FALSE,滤波函数将作指数平滑，参数beta为Holt-Winters滤波函数的参数，参数seasonal为滤波函数中季节趋势，seasonal="mult"表示季节效应为乘法模型。2026/2/2《统计学实验》第15章时间序列分析15-22【软件操作】ts.plot(tx)#作时间序列tx的点图lines(fitted(m)[,1],lty=2)#fitted(m)为拟合的序列legend(1990,4200,c("原序列","指数平滑后拟合值"),text.col="black",lty=c(1,2))2026/2/2《统计学实验》第15章时间序列分析15-232026/2/2《统计学实验》第15章时间序列分析15-24【运行结果】图15.3原序列和通过指数平滑后拟合值和观测值的时间序列图2026/2/2《统计学实验》第15章时间序列分析15-2515.3Box-Jenkin

模型

15.3.1ARMA模型的识别和估计【例15.3】(数据文件为li15.3.txt)该数据是由一个AR(1)模型模拟得到的数据，其中模拟参数，试对该数据进行模型识别和参数估计。

(数据来源：费宇等,《统计学》第7章,高等教育出版社,2010)2026/2/2《统计学实验》第15章时间序列分析15-26【统计理论】(1)自回归(AR)模型①模型形式假设

为观测时间序列，

；对

的统计建模，一个自然的想法是用其过去某个时段上的数据对其进行回归，写成模型形式，即其中

为模型参数，

是相互独立的均值为零，方差为

的随机误差。由于它是对自身变量的回归，模型(15.4)称为p阶自回归模型，记为。2026/2/2《统计学实验》第15章时间序列分析15-27【统计理论】②识别条件当使用ARMA模型拟合时间序列数据时，关键的问题是如何确定其阶数，即即AR(p)中的p，MA(q)中的q，以及ARMA(p,q)中的p，q。通过某种方法确定这些阶数的过称成为ARMA模型的模型识别。识别模型阶数的两个重要统计量是自相关函数(autocorrelationfunction，简记为ACF)和偏自相关函数(partialautocorrelationfunction，简记为PACF)。2026/2/2《统计学实验》第15章时间序列分析15-28【统计理论】对于AR(p)模型，自相关函数

随着p的增加呈现指数衰减或震荡性衰减；而偏相关函数

是

阶截尾的，即当时，

等于零。以AR(1)为例，此时

，由于，当

时，

呈现指数衰减，而当

时，

呈或震荡性衰减。偏相关函数

，因此是一阶截尾的。2026/2/2《统计学实验》第15章时间序列分析15-29【统计理论】假设

为后移(滞后)算子，即，

，令则(15.4)可以写为：

2026/2/2《统计学实验》第15章时间序列分析15-30【统计理论】有几种信息准则也可用来决定AR过程的阶p，它们都基于似然函数。例如，著名的AIC准则定义为：其中

是残差平方和的估计。在实际应用中，事先给定一个正整数

，对

计算

，然后选择阶，使

达最小值。2026/2/2《统计学实验》第15章时间序列分析15-31【统计理论】③平稳条件假设

表示

的一个

阶多项式，则AR(p)模型平稳的充分必要条件是：的根全部落在单位圆之外。最简单的AR模型是AR(1)，模型为：

。此时，平稳模型的充要条件是

。2026/2/2《统计学实验》第15章时间序列分析15-32【统计理论】④模型意义仅通过时间序列变量的自身历史观测值来反映有关因素对预测目标的影响和作用，影响序列变化的主要因素是时间序列在不同时期的取值。2026/2/2《统计学实验》第15章时间序列分析15-33【统计理论】(2)某著名企业平均(MA)模型①模型形式对时间序列

，一个滑动平均模型假设

可以由该时刻和该时刻之前的随机误差线性表示，即

我们称(15.6)为

阶滑动平均模型，记为MA(q)。2026/2/2《统计学实验》第15章时间序列分析15-34【统计理论】滑动平均模型

(15.6)也可以表为，其中

。由于

只有有限项，因此式(15.6)总满足平稳性条件，是一个平稳过程。2026/2/2《统计学实验》第15章时间序列分析15-35【统计理论】②识别条件对于MA(p)模型，自相关函数

是

阶截尾的，即当

时，

接近于零。而偏相关系数

随着的增加呈现指数衰减或震荡性衰减，并逐渐趋于零。以MA(1)为例，此时

，

，即一阶截尾；而偏相关系数

，因为

，偏相关系数将出现指数衰减。2026/2/2《统计学实验》第15章时间序列分析15-36【统计理论】③可逆条件MA模型中的一个重要内容是研究它的可逆性。令

表示z的q阶多项式，则MA(q)模型满足可逆性的充分必要条件是：

的根全部落在单位圆之外。一阶滑动平均模型MA(1)的可逆条件为：

。当满足可逆条件时，MA(q)模型往往可以转化为AR(p)模型。2026/2/2《统计学实验》第15章时间序列分析15-37【统计理论】④模型意义用过去各个时期的随机干扰或预测误差的线性组合来表达当前预测值。2026/2/2《统计学实验》第15章时间序列分析15-38【统计理论】(3)自回归滑动平均(ARMA)模型①模型形式当模型既有自回归项，又有滑动平均项时，就得到一个混合模型，称为自回归滑动平均模型：

ARMA(p,q)模型的平稳条件为：的根全部落在单位圆之外。2026/2/2《统计学实验》第15章时间序列分析15-39【统计理论】②识别条件偏相关函数

和自相关函数

均不截尾，但较快收敛到0，则该时间序列可能是ARMA（p,q）。2026/2/2《统计学实验》第15章时间序列分析15-40【统计理论】

③模型阶数ARMA模型阶数的确定很难从ACF和PACF获得充分的信息，TsayandTiao(1984)提出了一种称为EACF(extendedautocorrelationfunction)方法来识别ARMA(p,q)模型阶数。2026/2/2《统计学实验》第15章时间序列分析15-41【统计理论】在模型识别过程中，ACF及PACF能提供ARMA模型阶数的重要信息，但在建模过程中，我们还应该结合模型拟合统计量，如

、AIC及BIC，进行综合分析。其中一种选择模型阶数的方法是计算对应各种不同阶数的AIC或BIC统计量，其中对应AIC或BIC数值较小的模型即可选为较好的模型。2026/2/2《统计学实验》第15章时间序列分析15-42【软件操作】首先，通过acf()或pacf()函数确定模型的阶数x=read.table("li15.3.txt",header=T)#从li15.3.txt中读入样本数据xacf(x)#作自相关图，拖尾pacf(x)#作偏自相关图，1阶截尾2026/2/2《统计学实验》第15章时间序列分析15-43【运行结果】图15.4（a）

例15.3中数据的ACF图【运行结果】2026/2/2《统计学实验》第15章时间序列分析15-44图15.4（b）

例15.3中数据的PACF图【软件操作】然后，用arima()函数对模型进行参数估计arima(x,order=c(1,0,0),method="ML",include.mean=F)#用AR(1)模型拟合参数order=控制模型的阶数，按顺序分别为AR的阶数、差分的步数、MA的阶数，参数method控制模型参数估计方法，method="ML"表示估计方法为最大似然估计，参数include.mean=F表示模型不含截距项2026/2/2《统计学实验》第15章时间序列分析15-45【软件操作】a=arima(x,order=c(1,0,0),method="ML",include.mean=F)#把估计的结果赋于变量a，以便于通过输出对象a进行残差分析r=a$residuals#用r来保存残差Box.test(r,type="Ljung-Box")#检验模型拟合残差是否存在序列相关,参数type控制检验的方法，还可以选择"Box-Pierce"检验统计量2026/2/2《统计学实验》第15章时间序列分析15-46【运行结果】Call:arima(x=x,order=c(1,0,0),include.mean=F,method="ML")Coefficients:ar1

0.6783s.e.0.0327sigma^2estimatedas0.9788:loglikelihood=-704.41,aic=1412.822026/2/2《统计学实验》第15章时间序列分析15-47【运行结果】残差检验的Ljung-Box检验结果为：

Box-Ljungtestdata:rX-squared=0.119,df=1,p-value=0.73012026/2/2《统计学实验》第15章时间序列分析15-482026/2/2《统计学实验》第15章时间序列分析15-4915.3.2ARIMA模型的识别和估计【例7.4】(数据文件为li7.4.txt)使用例7.1中1978年至2005年国内生产总值(GDP，亿元)数据，2006-2008年用于做预测。使用ARIMA模型进行分析，给出2006-2008年的预测值。

2026/2/2《统计学实验》第15章时间序列分析15-50【统计理论】在际问题中，时间序列可能是非平稳的，此时不能简单的用ARMA模型进行拟合。一个想法是将序列进行差分，使其达到平稳性要求，再用ARMA模型进行分析。2026/2/2《统计学实验》第15章时间序列分析15-51【统计理论】假设零均值时间序列

经过

次差分之后可以转化为一个平稳时间序列

，则拟合ARMA(p,q)模型后，有ARIMA(p,d,q)中

的平稳性条件和可逆性条件与ARMA模型是一致的。当预先确定好差分阶数

，ARIMA(p,d,q)模型的参数估计与ARMA(p,q)模型是一样的。2026/2/2《统计学实验》第15章时间序列分析15-52【软件操作】对于ARIMA模型，应先计算该序列的0(即原始数据)，1，2阶差分后的ACF和PACF，以确定序列的差分阶数x=read.table("li15.4.txt",header=T)tx=ts(x,frequency=1,start=c(1978),end=c(2005))x1=diff(tx,lag=1)x2=diff(tx,lag=2)【软件操作】layout(matrix(1:3,1,3))#把原序列的0，1，2阶差分后的ACF条形图平铺在同一acf(tx,main="GDP")#作原始数据的自相关图acf(x1,main="GDP_1")#1阶差分后的自相关图acf(x2,main="GDP_2")#2阶差分后的自相关图2026/2/2《统计学实验》第15章时间序列分析15-53【软件操作】layout(matrix(1:3,1,3))pacf(tx,main="GDP")#作原始数据的偏自相关图pacf(x1,main="GDP_1")#1阶差分后的偏自相关图pacf(x2,main="GDP_2")#2阶差分后的偏自相关图2026/2/2《统计学实验》第15章时间序列分析15-542026/2/2《统计学实验》第15章时间序列分析15-55【运行结果】图15.5GDP数据0,1,2阶差分后的ACF图【运行结果】2026/2/2《统计学实验》第15章时间序列分析15-56图15.6GDP数据0,1,2阶差分后的PACF图【软件操作】对于未差分数据，从其ACF和PACF的结构来看似乎符合一个AR(1)模型。但原序列不平稳，因此可以考虑MA(1)和ARMA(1,1)。如考虑一阶差分后的数据，其PACF在一阶截尾，因此我们选用ARIMA(1,1,0)、ARIMA(1,1,1)和ARIMA(0,1,2)进行拟合。接着我们对模型预设各种不同的阶数，并进行参数估计和Ljun-Box检验，通过输出结果中的统计量最终确定模型阶数2026/2/2《统计学实验》第15章时间序列分析15-57【软件操作】model1=arima(tx,order=c(1,0,0),method="ML",include.mean=F)model2=arima(tx,order=c(0,0,1),method="ML",include.mean=F)model3=arima(tx,order=c(1,1,0),method="ML",include.mean=F)model4=arima(tx,order=c(1,1,1),method="ML",include.mean=F)2026/2/2《统计学实验》第15章时间序列分析15-58【软件操作】model5=arima(tx,order=c(0,1,2),method="ML",include.mean=F)model6=arima(tx,order=c(1,2,0),method="ML",include.mean=F)AIC(model1,model2,model3,model4,model5,model6)Box.test(model1$residuals,type="Ljung-Box")2026/2/2《统计学实验》第15章时间序列分析15-59【软件操作】Box.test(model2$residuals,type="Ljung-Box")Box.test(model3$residuals,type="Ljung-Box")Box.test(model4$residuals,type="Ljung-Box")Box.test(model5$residuals,type="Ljung-Box")Box.test(model6$residuals,type="Ljung-Box")2026/2/2《统计学实验》第15章时间序列分析15-60【运行结果】各个预设模型的拟合残差的Ljung-Box检验结果汇总如下：dfAICp-valueofBox-Ljungtestmodel12600.67270.000model22678.57040.000model32504.96670.0503model43501.58560.8861model53518.94160.041model62481.32570.73272026/2/2《统计学实验》第15章时间序列分析15-61【软件操作】计算ARIMA(1,1,1)的参数估计arima(tx,order=c(1,1,1),method="ML",include.mean=F)2026/2/2《统计学实验》第15章时间序列分析15-62【运行结果】Call:arima(x=tx,order=c(1,1,1),include.mean=F,method="ML")Coefficients:ar1ma10.95430.6364s.e.0.05160.1394sigma^2estimatedas4748365:loglikelihood=-247.79,aic=501.592026/2/2《统计学实验》第15章时间序列分析15-63【软件操作】依据残差的平稳性特征再次验证模型阶数的合理性：layout(matrix(1:2,1,2))acf(model4$residuals,main="reisdual")pacf(model4$residuals,main="reisdual")2026/2/2《统计学实验》第15章时间序列分析15-64【运行结果】2026/2/2《统计学实验》第15章时间序列分析15-65图15.7拟合ARIMA(1,1,1)模型残差的ACF和PACF图【软件操作】最后，对2006-2008年的GDP进行预测：tx