季节时间序列模型
在某些时间序列中，存在明显的周期性变化。这种周期是由于季节性变化（包括季度、月度、周度等变化）或其他一些固有因素引起的。这类序列称为季节性序列。比如一个地区的气温值序列（每隔一小时取一个观测值）中除了含有以天为周期的变化，还含有以年为周期的变化。在经济领域中，季节性序列更是随处可见。如季度时间序列、月度时间序列、周度时间序列等。处理季节性时间序列只用以上介绍的方法是不够的。描述这类序列的模型之一是季节时间序列模型(seasonal ARIMA model)，用SARIMA表示。较早文献也称其为乘积季节模型（multiplicative seasonal model）。
设季节性序列（月度、季度、周度等序列都包括其中）的变化周期为s，即时间间隔为s的观测值有相似之处。首先用季节差分的方法消除周期性变化。季节差分算子定义为，
(s = 1- Ls
若季节性时间序列用yt表示，则一次季节差分表示为
(s yt = (1- Ls) yt = yt - yt - s
对于非平稳季节性时间序列，有时需要进行D次季节差分之后才能转换为平稳的序列。在此基础上可以建立关于周期为s的P阶自回归Q阶移动平均季节时间序列模型（注意P、Q等于2时，滞后算子应为(Ls)2 = L2s。
(P (Ls) (sDyt = (Q (Ls) ut ()
对于上述模型，相当于假定ut是平稳的、非自相关的。
当ut非平稳且存在ARMA成分时，则可以把ut描述为
(p (L) (dut = (q (L) vt ()
其中vt为白噪声过程，p, q分别表示非季节自回归、移动平均算子的最大阶数，d表示ut的一阶（非季节）差分次数。由上式得
ut = (p-1(L) (-d (q (L) vt ()
把 () 式代入 () 式，于是得到季节时间序列模型的一般表达式。
(p(L) (P(Ls) ((d(sDyt) = (q(L) (Q(Ls) vt ()
其中下标P, Q, p, q分别表示季节与非季节自回归、移动平均算子的最大滞后阶数，d, D分别表示非季节和季节性差分次数。上式称作 (p, d, q) ( (P, D, Q)s 阶季节时间序列模型或乘积季节模型。
保证（(d(sDyt）具有平稳性的条件是(p(L)(P(Ls) = 0的根在单位圆外；保证（(d(sDyt）具有可逆性的条件是(q (L)(Q (Ls) = 0的根在单位圆外。
当P = D = Q = 0时，SARIMA模型退化为ARIMA模型；从这个意义上说，ARIMA模型是SARIMA模型的特例。当P = D = Q = p = q = d = 0时，SARIMA模型退化为白噪声模型。
(1, 1, 1) ( (1, 1, 1)12 阶月度SARIMA模型表达为
(1- (1 L) (1- (1 L12) ( (12 yt = (1+(1 L) (1+(1 L12) vt
( (12 yt具有平稳性的条件是 ( (1 ( < 1，( (1 ( < 1，( (12 yt具有可逆性的条件是 ( (1 ( < 1，( (1 ( < 1。
设log(Yt) = yt，变量( (12 yt在EViews中用DLOG(Y,1,12)表示（这样表示的好处是EViews可以直接预测到Y），上式的EViews估计命令是
DLOG(Y,1,12) AR(1) SAR(12) MA(1) SMA(12)
(0, 1, 1) ( (0, 1, 1)12 阶月度SARIMA模型表达为
( (12 yt = (1+ (1 L) (1+ (1 L12) vt ()
() 式的EViews估计命令是
DLOG(Y,1,12) MA(1) SMA(12)
由() 式得
((12 yt = (1+(1 L) (1+(1 L12 ) vt = vt +(1 L vt +(1 L12vt + (1 (1 L13vt
= vt +(1 vt –1 +(1 vt – 12 + (1 (1 vt – 13
上式对应的EViews估计命令是
DLOG(Y,1,12) MA(1) MA(12) MA(13)
模型表达式是
((12 yt = vt +(1 vt –1 +(12 vt – 12 + (13 vt – 13
这是一个非季节模型表达式。以上两个EViews估计命令是等价的，都是估计MA(13)模型。
注意：唯一不同点是上式对vt – 13的系数没有约束，而对季节模型来说，相当于增加了一个约束条件，(13 =(1 (1。
进一步化简
( (yt – yt - 12) = vt +(1 vt –1 +(1 vt – 12 + (1 (1 vt – 13
( yt – ( yt - 12 = vt +(1 vt –1 +(1 vt – 12 + (1 (1 vt – 13
用于预测的模型型式是
yt = yt -1 + yt - 12 – yt – 13 + vt +(1 vt –1 +(1 vt – 12 + (1 (1 vt – 13 ()
从上式可以看出SARIMA模型可以展开为ARIMA模型。
对乘积季节模型的季节阶数，即周期长度s的识别可以通过对实际问题的分析、时间序列图以及时间序列的相关图和偏相关图分析得到。
以相关图和偏相关图为例，如果相关图和偏相关图不是呈线性衰减趋势，而是在变化周期的整倍数时点上出现绝对值相当大的峰值并呈振荡式变化，就可以认为该时间序列可以用SARIMA模型描述。
建立SARIMA模型，
（1）首先要确定d, D。通过差分和季节差分把原序列变换为一个平稳的序列。令
xt = (d(sD yt
（2）然后用xt 建立 (p (L) (P (Ls) xt = (q (L) (Q (Ls) vt模型。
注意：
（1）用对数的季节时间序列数据建模时通常D不会大于1，P和Q不会大于3。
（2）乘积季节模型参数的估计、检验与前面介绍的估计、检验方法相同。利用乘积季节模型预测也与上面介绍的预测方法类似。
季节时间序列建模案例
案例1：（文件名：b2c3，5b2c3）北京市1978:1~1989:12社会商品零售额月度数据（yt，单位：亿元人民币）曲线见图，数据见表。yt与时间呈指数关系且存在递增型异方差。对数的社会商品零售额月度数据（Ln yt）曲线见图。Lnyt与时间近似呈线性关系（异方差问题也得到抑制）。
图 yt 图 Lnyt
通过Lnyt的相关图和偏相关图（见图）可以看到Lnyt是一个非平稳序列（相关图衰减很慢）且Lnyt与其12倍数的滞后期存在自回归关系。
图 Lnyt的相关图（下）和偏相关图（上）
对Lnyt进行一阶差分，得(Lnyt（图）。图是对Lnyt进行2次一阶差分的结果，序列(2Ln yt是过度差分序列。从 (Lnyt的相关图和偏相关图（图）可以看到，通过差分 (Lnyt的平稳性得到很大改进，但与其12倍数的滞后期存在显著的自相关关系。
图 (Ln yt 图 (2Ln yt
图 (Lnyt的相关图（下）和偏相关图（上）
对Lnyt进行一次季节性差分（或12阶差分），得 (12 Lnyt（图）。从 (12 Lnyt的相关图和偏相关图（图）可以看到 (12 Lnyt仍然是非平稳的。
图 (12 Lnyt，（EViews：DLOG(Y,0,12)）
图 (12 Lnyt的相关图（下）和偏相关图（上）
对Lnyt进行一阶差分和一阶季节性差分，得((12 Lnyt（见图）。从xt 的相关图和偏相关图（见图）可以看到((12 Lnyt近似为一个平稳过程。
图 ( (12 Lnyt = xt，（EViews：DLOG(Y,1,12)）
图 ((12 Lnyt的相关图（下）和偏相关图（上）
用1978:1~1989:11期间数据，估计yt 的 (1, 1, 1) ( (1, 1, 0)12阶季节时间序列模型，得结果如下：
(1+ L) (1 + L12) ((12Lnyt = (1+ L) vt ()
() () ()
R2 = , . = , Q36 = , ((36-2-1) = 44
EViews估计命令是
DLOG(Y,1,12) AR(1) SAR(12) MA(1)
EViews输出结果见图。
注意：
（1）仔细对照（）式和图输出结果，不要把自回归系数估计值的符号写错。通过自回归特征根倒数可知，把表达式中的算子写作(1+ L)是正确的。通过移动平均特征根倒数可知，把表达式中的算子写作(1+ L) 是正确的。
（2）表达式中，季节和非季节因子（特征多项式）之间是相乘关系。
（3）在EViews估计命令中把变量写作DLOG(Y,1,12)的好处是可以直接对yt和((12 Lnyt预测。
模型残差序列的相关与偏相关图如图。
图 EViews估计结果
图模型残差序列的相关与偏相关图
对于((12 Lnyt来，模型参数全部有显著性，Q36 = < ((36-2-1) = 44。两种检验通过。见输出结果（），对于((12 Lnyt，模型共有14个特征根。
图 D12DLnyt的实际与预测序列 图 yt的实际与预测序列
对1989年第12月份yt进行样本外1期预测，结果如图。
图 EViews预测结果
预测误差是
( ==
用1978:1~1989:12期间数据得EViews (0, 1, 1) ( (0, 1, 1)12 模型估计结果如下
( (12 Lnyt = (1- L) (1 - L12) vt ()
(- ) (- )
R2 = , DW = , F = , . = , Q36 = , ( (36-2) = 44
模型参数全部有显著性，Q36 = < ( (36-2) = 44。两种检验通过。上式变换为，
(Ln yt – ( Lnyt - 12 = vt - vt –1 - vt – 12 + vt – 13
Lnyt = Lnyt -1 +Ln yt - 12 –Lnyt – 13 + vt - vt –1 - vt – 12 + vt – 13 ()
（）式也是一个可以选用的模型。
表 北京市社会商品零售额（yt）月度数据（单位：亿元人民币，1978:1~1989:12）
年：月
yt
年：月
yt
年：月
yt
年：月
yt
年：月
yt
1978:01
1980:06
1982:11
1985:04
1987:09
1978:02
1980:07
1982:12
1985:05
1987:10
1978:03
1980:08
1983:01
1985:06
1987:11
1978:04
1980:09
1983:02
1985:07
1987:12
1978:05
1980:10
1983:03
1985:08
1988:01
1978:06
1980:11
1983:04
1985:09
1988:02
1978:07
1980:12
1983:05
1985:10
1988:03
1978:08
1981:01
1983:06
1985:11
1988:04
1978:09
1981:02
1983:07
1985:12
1988:05
1978:10
1981:03
1983:08
1986:01
1988:06
1978:11
1981:04
1983:09
1986:02
1988:07
1978:12
1981:05
1983:10
1986:03
1988:08
1979:01
1981:06
1983:11
1986:04
1988:09
1979:02
1981:07
1983:12
1986:05
1988:10
1979:03
1981:08
1984:01
1986:06
1988:11
1979:04
1981:09
1984:02
1986:07
1988:12
1979:05
1981:10
1984:03
1986:08
1989:01
1979:06
1981:11
1984:04
1986:09
1989:02
1979:07
1981:12
1984:05
1986:10
1989:03
1979:08
1982:01
1984:06
1986:11
1989:04
1979:09
1982:02
1984:07
1986:12
1989:05
1979:10
1982:03
1984:08
1987:01
1989:06
1979:11
1982:04
1984:09
1987:02
1989:07
1979:12
1982:05
1984:10
1987:03
1989:08
1980:01
1982:06
1984:11
1987:04
1989:09
1980:02
1982:07
1984:12
1987:05
1989:10
1980:03
1982:08
1985:01
1987:06
1989:11
1980:04
1982:09
1985:02
1987:07
1989:12
1980:05
1982:10
1985:03
1987:08
案例2 香港季节GDPt数据的拟合（季节时间序列模型，file：HongKong）
1980:1~2002:4年香港季度GDPt序列曲线见图（数据见表，单位：港元）。1980~1997年GDPt随时间呈指数增长。1997年由于遭受东南亚金融危机的影响，经济发展处于停滞状态，1998~2002年底GDPt总量几乎没有增长。另一个特征是GDPt随时间呈递增型异方差。所以，用对数的季度GDPt数据（LnGDPt，曲线见图）建立季节时间序列模型。
图 GDPt 图 LnGDPt
通过LnGDPt的相关图和偏相关图（图）可以看到LnGDPt是一个非平稳序列（相关图衰减得很慢）。
图 LnGDPt的相关图和偏相关图
对LnGDPt进行一阶差分，得 DLnGDPt（见图）。DLnGDPt的平稳性得到很大改进，但其季节因素影响还很大。从 DLnGDPt的相关图和偏相关图（图）也可以明显地看到这个特征。若对LnGDPt直接进行一次季节差分（四阶差分），得D4LnGDPt见图。其波动性也很大。相关图和偏相关图见图。D2LnGDPt显然是过度差分序列（图）。
图 DLnGDPt ,（. = ） 图 DLnGDPt的相关图和偏相关图
图 D4LnGDPt,（. = ） 图 D4LnGDPt的相关图和偏相关图
图 D2LnGDPt ,（. = ）
在DLnGDPt的基础上进行一阶季节差分，或在D4LnGDPt基础上进行一阶非季节差分，得 D4DLnGDPt（图）。其相关图和偏相关图见图。D4DLnGDPt中已经基本消除了季节变化因素。在D4DLnGDPt的基础上建立时间序列模型。
图 D4DLnGDPt,（. = ） 图 D4DLnGDPt的相关和偏相关图
通过对D4DLnGDPt的相关和偏相关图分析，应该建立(2, 1, 2) ( (1, 1, 1)4 模型。EViews估计命令是
DLOG(GDP,1,4) C AR(1) AR(2) SAR(4) MA(1) MA(2) SMA(4)
用1980:1~2002:3的数据得估计结果如下（EViews输出结果如图）：
D4DLnGDPt = - + ut （1980:1~2002:3）
()
( L+ L2) (1 - L4) ut = (1 - L+ L2) (1 - L4) vt ()
() () () () () ()
R2 = , DW = , F = , Q36 = , ( (36-3-3-1) =
图 EViews估计结果
图模型()误差项的相关和偏相关图
注意：
（1）不要把自回归系数估计值的符号写错。不要把均值（- ）项表达错。EViews仍然是对(D4DLnGDPt + )建立(2, 1, 2) ( (1, 1, 1)4 阶季节时间序列模型，而不是对D4DLnGDPt建立季节时间序列模型。
（2）季节和非季节因子之间是相乘关系。
（3）在EViews估计命令中把变量写作DLOG(GDP,1,4)，好处是预测时可直接预测GDPt，也可以预测D4DLnGDPt。
模型参数全部有显著性，Q36 = <( (36-3-3-1) = 。两种检验通过。依据输出结果，对于D4DLnGDPt，模型共有12个特征根。4个实根，8个复根。
图 D4DLnGDPt的实际与预测序列 图 GDPt的实际与预测序列
对2002年第4季度GDPt进行样本外1期预测，结果如下：
预测误差是
( ==
表 香港季度GDPt数据（1980:1~2002:4，单位：港元）
年：月
GDPt /1011
年：月
GDPt /1011
年：月
GDPt /1011
年：月
GDPt /1011
1980:1
1985:4
1991:3
1997:2
1980:2
1986:1
1991:4
1997:3
1980:3
1986:2
1992:1
1997:4
1980:4
1986:3
1992:2
1998:1
1981:1
1986:4
1992:3
1998:2
1981:2
1987:1
1992:4
1998:3
1981:3
1987:2
1993:1
1998:4
1981:4
1987:3
1993:2
1999:1
1982:1
1987:4
1993:3
1999:2
1982:2
1988:1
1993:4
1999:3
1982:3
1988:2
1994:1
1999:4
1982:4
1988:3
1994:2
2000:1
1983:1
1988:4
1994:3
2000:2
1983:2
1989:1
1994:4
2000:3
1983:3
1989:2
1995:1
2000:4
1983:4
1989:3
1995:2
2001:1
1984:1
1989:4
1995:3
2001:2
1984:2
1990:1
1995:4
2001:3
1984:3
1990:2
1996:1
2001:4
1984:4
1990:3
1996:2
2002:1
1985:1
1990:4
1996:3
2002:2
1985:2
1991:1
1996:4
2002:3
1985:3
1991:2
1997:1
2002:4
PAGE 2
PAGE 12
