7.6 模型估计和选择

估计ETS模型

除了通过最小化误差平方和来估计参数,另一种方法是通过 最大化“似然函数” 。似然函数是特定模型产生数据的概率。因此,一个好的模型通常有着较大的似然函数。对于一个加性误差模型,最大化似然函数(假定误差正态分布)与最小化误差平方和给出的结果相同。但是,乘性误差模型会得到不同的结果。在本节中,我们将通过最大化似然函数来估计平滑参数 \(\alpha\)\(\beta\)\(\gamma\)\(\phi\),以及初始状态 \(\ell_0\)\(b_0\)\(s_0,s_{-1},\dots,s_{-m+1}\)

平滑参数可以采用的可能值会受到限制。传统上,参数被限制在0和1之间,以便于将方程解释为加权平均值。也就是说, \(0< \alpha,\beta^*,\gamma^*,\phi<1\)。对于状态空间模型,我们设置了 \(\beta=\alpha\beta^*\)\(\gamma=(1-\alpha)\gamma^*\)。因此,传统的限制转换为 \(0< \alpha <1\)\(0 < \beta < \alpha\) 以及 \(0< \gamma < 1-\alpha\)。在实践中,阻尼参数 \(\phi\) 通常会受到进一步限制,以防止模型估计过程中的计算困难。在R中,它被限制为 \(0.8<\phi<0.98\)

查看参数的另一种方法是通过考虑状态空间模型的数学性质。参数受到限制是为了防止很久之前的观测值对当前的预测产生持续的影响。这导致了对参数的一些 可容性 限制 ,这种限制通常(但不总是)比一般领域 (Hyndman et al., 2008, p. Ch10) 的限制要小。例如,对于ETS(A,N,N)模型,通常的参数区域是 \(0< \alpha <1\),但可容性区域是 \(0< \alpha <2\)。对于ETS(A,A,N)模型,通常的参数区域为 \(0<\alpha<1\)\(0<\beta<\alpha\),但可容性区域是 \(0<\alpha<2\)\(0<\beta<4-2\alpha\)

模型选择

ETS统计框架的一大优点是信息准则可用于模型选择。可以使用 5.5 节中介绍的AIC,AIC\(_{\text{c}}\) 和BIC来确定哪个ETS模型最适合一个给定的时间序列。

对于ETS模型,赤池信息标准(AIC)被定义为 \[ \text{AIC} = -2\log(L) + 2k, \] 其中 \(L\) 是模型的似然函数, \(k\) 是已估计的参数个数和初始状态的总和(包括残差的方差)。

针对小样本偏差修正的AIC(AIC\(_\text{c}\))被定义为 \[ \text{AIC}_{\text{c}} = \text{AIC} + \frac{k(k+1)}{T-k-1}, \] 另外,贝叶斯信息准则(BIC)被定义为 \[ \text{BIC} = \text{AIC} + k[\log(T)-2]。 \]

(误差,趋势,季节性)的三种组合可能会导致数值计算上的困难。具体来说,可能导致这种不稳定性的模型是ETS(A,N,M),ETS(A,A,M),和ETS(A,A\(_d\),M),因为在状态方程中分母可能接近于零。在选择模型时,我们通常不会考虑这些特定组合。

当数据严格为正值时,带有乘性误差的模型 很有用,但当数据包含零或负值时,模型在数值上并不稳定。因此,如果时间序列不是严格为正,那么乘性误差模型将不被考虑。在这种情况下,只有六个完全加性模型将会被应用。

R中的 ets() 函数

在R中,可以使用forecast包里的 ets() 函数来估计这些模型。 与 ses()holt()hw() 函数不同, ets() 函数不会生成预测值,而是估计模型参数值并返回有关拟合模型的信息。默认情况下,它使用AICc来选择适当的模型,但也可以选择其他信息准则。

下面的R代码显示了此函数可以采用的最重要参数及其默认值。如果只指定了时间系列,并且所有其他参数都保留其默认值,则将自动选择适当的模型。我们在代码下面对参数值进行了解释。请参阅帮助文件以获取更完整的说明。

ets(y, model="ZZZ", damped=NULL, alpha=NULL, beta=NULL, gamma=NULL,
    phi=NULL, lambda=NULL, biasadj=FALSE, additive.only=FALSE,
    restrict=TRUE, allow.multiplicative.trend=FALSE)

y : 要预测的时间序列。

model : 一个由三个字母组成代码,表示使用ETS分类和符号表示的用来估计的模型。可能的输入值为“N”表示无,“A”表示加性,“M”表示乘性,或“Z”表示自动选择。如果任何一个输入为“Z”,则根据信息准则选择该成分。默认值ZZZ表示所有成分都使用信息准则进行选择。

damped : 如果 damped=TRUE,则将使用衰减趋势(A或M)。如果 damped=FALSE,那么将使用无衰减趋势。如果 ‘damped=NULL’(默认值),那么将选择衰减或无衰减趋势,具体取决于哪个模型具有最小的信息准则值。

alphabetagammaphi : 可以使用这些参数来指定平滑参数的值。如果它们设置为 NULL(每个参数的默认值),则对这些参数进行估计。

lambda : Box-Cox转换参数。如果 lambda=NULL(默认值),它将被忽略。否则,时间序列将在模型估计之前进行变换。当 lambda 不是 NULL 时, additive.only 被设置为 TRUE

biasadj : 如果取值为 TRUElambda 不是 NULL,那么反向变换的拟合值和预测值将被进行偏差调整。

additive.only : 如果 additive.only=TRUE,则只考虑含有加性成分的模型。否则,将考虑所有模型。

restrict : 如果 restrict=TRUE(默认值),那么在模型选择中不会考虑导致计算困难的模型。

allow.multiplicative.trend : 乘性趋势模型也使用,但本书并未涉及。将此参数设置为 TRUE 以允许考虑这些模型。

使用 ets 对象

ets() 函数将返回 ets 类的对象。有许多R函数可以方便地处理 ets 对象,将其中几个描述如下。

coef() : 返回所有拟合参数值。

accuracy() : 返回在训练数据上计算出的准确性度量。

summary() : 打印一些关于拟合模型的摘要信息。

autoplot()plot() : 生成各成分的时序图。

residuals() : 返回估计模型的残差。

fitted() : 返回训练数据的一步预测值。

simulate() : 将模拟来自拟合模型的未来样本路径。

forecast() : 计算将在下一节所讲的点预测值和预测区间。

示例:澳大利亚国际游客过夜的天数

我们现在使用ETS统计框架来预测2016–2019年期间国际旅客在澳大利亚过夜的天数。我们让 ets() 函数通过最小化AICc来选择模型。

aust <- window(austourists, start=2005)
fit <- ets(aust)
summary(fit)
#> ETS(M,A,M) 
#> 
#> Call:
#>  ets(y = aust) 
#> 
#>   Smoothing parameters:
#>     alpha = 0.1908 
#>     beta  = 0.0392 
#>     gamma = 2e-04 
#> 
#>   Initial states:
#>     l = 32.3679 
#>     b = 0.9281 
#>     s = 1.022 0.9628 0.7683 1.247
#> 
#>   sigma:  0.0383
#> 
#>   AIC  AICc   BIC 
#> 224.9 230.2 240.9 
#> 
#> Training set error measures:
#>                   ME  RMSE  MAE     MPE  MAPE   MASE
#> Training set 0.04837 1.671 1.25 -0.1846 2.693 0.4095
#>                ACF1
#> Training set 0.2006

所选模型是ETS(M,A,M): \[\begin{align*} y_{t} &= (\ell_{t-1} + b_{t-1})s_{t-m}(1 + \varepsilon_t)\\ \ell_t &= (\ell_{t-1} + b_{t-1})(1 + \alpha \varepsilon_t)\\ b_t &=b_{t-1} + \beta(\ell_{t-1} + b_{t_1})\varepsilon_t\\ s_t &= s_{t-m}(1+ \gamma \varepsilon_t). \end{align*}\] 参数估计值为 \(\hat\alpha=0.1908\)\(\hat\beta=0.03919\),和 \(\hat\gamma=0.0001917\)。输出中还返回初始状态 \(\ell_0\)\(b_0\)\(s_{0}\)\(s_{-1}\)\(s_{-2}\)\(s_{-3}\) 的估计值。将这些数值与表 7.5 中提供的相应的Holt-Winters乘性季节性方法得到的值进行比较。ETS(M,A,M)模型与乘性Holt-Winters方法给出的点预测不同,因为参数估计的方式不同。使用 ets() 函数时,默认的估计方法是最大化似然函数而不是最小化平方和。

7.9 显示了状态随时间的变化情况,而图 7.11 显示了模型生成的点预测和预测区间。 \(\beta\)\(\gamma\) 的较小取值表示斜率和季节性分量随时间的变化很小。较窄的预测区间表明该时间序列由于较强的趋势和季节性而相对容易预测。

autoplot(fit)
状态的估计随时间变化的图形表示。

图 7.9: 状态的估计随时间变化的图形表示。

由于这个模型具有乘性误差,所以残差并不等于一步训练误差。残差由 \(\hat{\varepsilon}_t\) 给出,而一步训练误差定义为 \(y_t - \hat{y}_{t|t-1}\)。我们可以使用 residuals() 函数来获得这它们的值。

cbind('Residuals' = residuals(fit),
      'Forecast errors' = residuals(fit, type='response')) %>%
  autoplot(facet=TRUE) + xlab("年份") + ylab("")+
  theme(text = element_text(family = "STHeiti"))+
  theme(plot.title = element_text(hjust = 0.5))
ETS(M,A,M)模型的残差和一步预测误差。

图 7.10: ETS(M,A,M)模型的残差和一步预测误差。

residuals() 函数中的 type 参数被用来区分残差和预测误差。默认值是 type ='innovation',它给出了常规的残差。

参考文献

Hyndman, R. J., Koehler, A. B., Ord, J. K., & Snyder, R. D. (2008). Forecasting with exponential smoothing: The state space approach. Berlin: Springer-Verlag. http://www.exponentialsmoothing.net