霍尔特-温特斯模型


# 此单元格将不会被渲染,但用于隐藏警告并限制显示的行数。

import warnings
warnings.filterwarnings("ignore")

import logging
logging.getLogger('statsforecast').setLevel(logging.ERROR)

import pandas as pd
pd.set_option('display.max_rows', 6)

使用 Holt Winters ModelStatsforecast 的逐步指南。

在这一过程当中,我们将熟悉主要的 StatsForecast 类和一些相关的方法,如 StatsForecast.plotStatsForecast.forecastStatsForecast.cross_validation 等。

目录

引言

霍尔特-温特模型,也称为三重指数平滑法,是一种广泛应用于时间序列分析的预测技术。它由查尔斯·霍尔特和彼得·温特于1960年开发,作为对霍尔特的双重指数平滑方法的改进。

霍尔特-温特模型用于预测显示趋势和季节性的时间序列的未来值。该模型使用三个平滑参数,一个用于估计趋势,另一个用于估计时间序列的水平或基准水平,第三个用于估计季节性。这些参数分别称为α、β和γ。

霍尔特-温特模型是霍尔特的双重指数平滑方法的扩展,后者仅使用两个平滑参数来估计时间序列的趋势和基准水平。霍尔特-温特模型通过增加一个用于季节性的第三个平滑参数来提高预测的准确性。

霍尔特-温特模型的主要优点之一是易于实施,并且不需要大量历史数据即可生成准确的预测。此外,该模型高度适应,可以定制以适应各种具有季节性的时间序列。

然而,霍尔特-温特模型也有一些局限性。例如,该模型假设时间序列是平稳的且季节性是恒定的。如果时间序列不是平稳的或季节性不恒定,霍尔特-温特模型可能不是最合适的选择。

总的来说,霍尔特-温特模型是时间序列分析中一种有用且广泛使用的技术,尤其是在序列预期呈现恒定趋势和季节性时。

霍尔特-温特斯法

霍尔特-温特斯季节性方法包括预测方程和三个平滑方程 — 一个用于水平 \(\ell_{t}\) ,一个用于趋势 \(b_t\) ,一个用于季节性成分 \(s_t\) ,以及相应的平滑参数 \(\alpha\)\(\beta^*\)\(\gamma\) 。我们用 \(m\) 来表示季节性的周期,即一年中的季节数量。例如,对于季度数据 \(m=4\) ,对于月度数据 \(m=12\)

该方法有两种变体,主要在于季节性成分的性质。当季节性变化在整个系列中大致保持恒定时,优先使用加法方法;当季节性变化与系列水平成正比时,优先使用乘法方法。在加法方法中,季节性成分以观察系列的绝对值表示,而在水平方程中,系列通过减去季节性成分进行季节调整。在每一年内,季节性成分的总和将近似为零。在乘法方法中,季节性成分以相对值(百分比)表示,系列通过除以季节性成分进行季节调整。在每一年内,季节性成分的总和将近似为 \(m\)

霍尔特-冬季加法法

霍尔特-冬季加法法是一种时间序列预测技术,通过纳入加法季节性成分来扩展霍尔特-冬季法。该方法适用于显示随时间变化的季节性模式的时间序列数据。

霍尔特-冬季加法法使用三个平滑参数 - α(alpha),β(beta)和γ(gamma) - 来估计时间序列的水平、趋势和季节成分。α参数控制水平成分的平滑,β参数控制趋势成分的平滑,γ参数控制加法季节成分的平滑。

预测过程包括三个步骤:首先,利用平滑参数和历史数据估计水平、趋势和季节成分;其次,使用这些成分预测时间序列的未来值;最后,使用加法因子调整预测值以考虑季节成分。

霍尔特-冬季加法法的一个优点是能够处理具有加法季节性成分的时间序列数据,这在许多实际应用中是常见的。该方法也易于实施,并且可以扩展以处理具有变化季节模式的时间序列数据。

然而,该方法也存在一些局限性。它假设季节性模式是加法的,这可能并非所有时间序列的情况。此外,该方法需要足够的历史数据,以准确估计平滑参数和季节成分。

总体而言,霍尔特-冬季加法法是一种强大且广泛使用的预测技术,可以用于生成对具有加法季节性成分的时间序列数据的准确预测。该方法易于实施,并可以扩展以处理具有变化季节模式的时间序列数据。

加法法的成分形式为:

\[\begin{align*} \hat{y}_{t+h|t} &= \ell_{t} + hb_{t} + s_{t+h-m(k+1)} \\ \ell_{t} &= \alpha(y_{t} - s_{t-m}) + (1 - \alpha)(\ell_{t-1} + b_{t-1})\\ b_{t} &= \beta^*(\ell_{t} - \ell_{t-1}) + (1 - \beta^*)b_{t-1}\\ s_{t} &= \gamma (y_{t}-\ell_{t-1}-b_{t-1}) + (1-\gamma)s_{t-m}, \end{align*}\]

其中 \(k\)\((h-1)/m\) 的整数部分,这确保用于预测的季节性指数估计来自样本的最后一年。水平方程显示了季节调整观察值 \((y_{t} - s_{t-m})\) 和非季节性预测 \((\ell_{t-1}+b_{t-1})\) 在时间 \(t\) 的加权平均。趋势方程与霍尔特的线性方法一致。季节方程显示当前季节指数 \((y_{t}-\ell_{t-1}-b_{t-1})\) 和上一年同一季节的季节指数(即 \(m\) 个时间段之前)的加权平均。

季节成分的方程通常表示为

\[s_{t} = \gamma^* (y_{t}-\ell_{t}) + (1-\gamma^*)s_{t-m}.\]

如果我们从上述成分形式的平滑方程中替换 \(\ell_{t}\),我们得到

\[s_{t} = \gamma^*(1-\alpha) (y_{t}-\ell_{t-1}-b_{t-1}) + [1-\gamma^*(1-\alpha)]s_{t-m},\]

这与我们在此指定的季节成分的平滑方程相同,其中 \(\gamma=\gamma^*(1-\alpha)\)。通常的参数限制是 \(0\le\gamma^*\le1\),这转化为 \(0\le\gamma\le 1-\alpha\)

Holt-Winters的乘法方法

Holt-Winters的乘法方法使用三个平滑参数 - alpha (α), beta (β), 和 gamma (γ) - 来估计时间序列的水平、趋势和季节性成分。alpha参数控制水平成分的平滑,beta参数控制趋势成分的平滑,而gamma参数控制乘法季节性成分的平滑。

预测过程包括三个步骤:首先,使用平滑参数和历史数据估计水平、趋势和季节性成分;其次,这些成分用于预测时间序列的未来值;最后,使用乘法因子调整预测值以考虑季节性成分。

Holt-Winters的乘法方法的一个优点是它能够处理具有乘法季节性成分的时间序列数据,这在许多现实世界的应用中是很常见的。该方法实现简单,并且可以扩展以处理具有变化季节模式的时间序列数据。

然而,该方法也存在一些局限性。它假设季节性模式是乘法的,这可能并不适用于所有时间序列。此外,该方法需要足够的历史数据,以准确估计平滑参数和季节性成分。

总体来说,Holt-Winters的乘法方法是一种强大且广泛使用的预测技术,可以用于生成对具有乘法季节性成分的时间序列数据的准确预测。该方法易于实现,并且可以扩展以处理具有变化季节模式的时间序列数据。

在乘法版本中,季节性平均值为1。如果季节变化随着水平的增加而增加,则使用乘法方法。

\[\begin{align*} \hat{y}_{t+h|t} &= (\ell_{t} + hb_{t})s_{t+h-m(k+1)} \\ \ell_{t} &= \alpha \frac{y_{t}}{s_{t-m}} + (1 - \alpha)(\ell_{t-1} + b_{t-1})\\ b_{t} &= \beta^*(\ell_{t}-\ell_{t-1}) + (1 - \beta^*)b_{t-1} \\ s_{t} &= \gamma \frac{y_{t}}{(\ell_{t-1} + b_{t-1})} + (1 - \gamma)s_{t-m}. \end{align*}\]

ETS 分类中的数学模型

我希望读者能够更清晰地理解ETS框架是如何建立在时间序列分解的理念之上的。通过引入不同的组件、定义其类型,并添加它们的更新方程,我们可以构建更能捕捉时间序列关键特征的模型。但我们也应该考虑组件随时间变化的潜在可能性。“转移”或“状态”方程应该反映这种变化;它们解释了水平、趋势或季节性组件是如何演变的。

正如第2.2节所讨论的,根据不同类型的组件及其交互关系,我们最终得到了分类中的30个模型。表1和表2从数学上总结了在图1和图2中图示的30个ETS模型,呈现了测量和转移方程的公式。

表1:加法误差ETS模型

非季节性 加法 乘法
无趋势 \(\begin{aligned} &y_{t} = l_{t-1} + \epsilon_t \\ &l_t = l_{t-1} + \alpha \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = l_{t-1} + s_{t-m} + \epsilon_t \\ &l_t = l_{t-1} + \alpha \epsilon_t \\ &s_t = s_{t-m} + \gamma \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = l_{t-1} s_{t-m} + \epsilon_t \\ &l_t = l_{t-1} + \alpha \frac{\epsilon_t}{s_{t-m}} \\ &s_t = s_{t-m} + \gamma \frac{\epsilon_t}{l_{t-1}} \end{aligned}\)
加法 \(\begin{aligned} &y_{t} = l_{t-1} + b_{t-1} + \epsilon_t \\ &l_t = l_{t-1} + b_{t-1} + \alpha \epsilon_t \\ &b_t = b_{t-1} + \beta \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = l_{t-1} + b_{t-1} + s_{t-m} + \epsilon_t \\ &l_t = l_{t-1} + b_{t-1} + \alpha \epsilon_t \\ &b_t = b_{t-1} + \beta \epsilon_t \\ &s_t = s_{t-m} + \gamma \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = (l_{t-1} + b_{t-1}) s_{t-m} + \epsilon_t \\ &l_t = l_{t-1} + b_{t-1} + \alpha \frac{\epsilon_t}{s_{t-m}} \\ &b_t = b_{t-1} + \beta \frac{\epsilon_t}{s_{t-m}} \\ &s_t = s_{t-m} + \gamma \frac{\epsilon_t}{l_{t-1} + b_{t-1}} \end{aligned}\)
加法阻尼 \(\begin{aligned} &y_{t} = l_{t-1} + \phi b_{t-1} + \epsilon_t \\ &l_t = l_{t-1} + \phi b_{t-1} + \alpha \epsilon_t \\ &b_t = \phi b_{t-1} + \beta \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = l_{t-1} + \phi b_{t-1} + s_{t-m} + \epsilon_t \\ &l_t = l_{t-1} + \phi b_{t-1} + \alpha \epsilon_t \\ &b_t = \phi b_{t-1} + \beta \epsilon_t \\ &s_t = s_{t-m} + \gamma \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = (l_{t-1} + \phi b_{t-1}) s_{t-m} + \epsilon_t \\ &l_t = l_{t-1} + \phi b_{t-1} + \alpha \frac{\epsilon_t}{s_{t-m}} \\ &b_t = \phi b_{t-1} + \beta \frac{\epsilon_t}{s_{t-m}} \\ &s_t = s_{t-m} + \gamma \frac{\epsilon_t}{l_{t-1} + \phi b_{t-1}} \end{aligned}\)
乘法 \(\begin{aligned} &y_{t} = l_{t-1} b_{t-1} + \epsilon_t \\ &l_t = l_{t-1} b_{t-1} + \alpha \epsilon_t \\ &b_t = b_{t-1} + \beta \frac{\epsilon_t}{l_{t-1}} \end{aligned}\) \(\begin{aligned} &y_{t} = l_{t-1} b_{t-1} + s_{t-m} + \epsilon_t \\ &l_t = l_{t-1} b_{t-1} + \alpha \epsilon_t \\ &b_t = b_{t-1} + \beta \frac{\epsilon_t}{l_{t-1}} \\ &s_t = s_{t-m} + \gamma \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = l_{t-1} b_{t-1} s_{t-m} + \epsilon_t \\ &l_t = l_{t-1} b_{t-1} + \alpha \frac{\epsilon_t}{s_{t-m}} \\ &b_t = b_{t-1} + \beta \frac{\epsilon_t}{l_{t-1}s_{t-m}} \\ &s_t = s_{t-m} + \gamma \frac{\epsilon_t}{l_{t-1} b_{t-1}} \end{aligned}\)
乘法阻尼 \(\begin{aligned} &y_{t} = l_{t-1} b_{t-1}^\phi + \epsilon_t \\ &l_t = l_{t-1} b_{t-1}^\phi + \alpha \epsilon_t \\ &b_t = b_{t-1}^\phi + \beta \frac{\epsilon_t}{l_{t-1}} \end{aligned}\) \(\begin{aligned} &y_{t} = l_{t-1} b_{t-1}^\phi + s_{t-m} + \epsilon_t \\ &l_t = l_{t-1} b_{t-1}^\phi + \alpha \epsilon_t \\ &b_t = b_{t-1}^\phi + \beta \frac{\epsilon_t}{l_{t-1}} \\ &s_t = s_{t-m} + \gamma \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = l_{t-1} b_{t-1}^\phi s_{t-m} + \epsilon_t \\ &l_t = l_{t-1} b_{t-1}^\phi + \alpha \frac{\epsilon_t}{s_{t-m}} \\ &b_t = b_{t-1}^\phi + \beta \frac{\epsilon_t}{l_{t-1}s_{t-m}} \\ &s_t = s_{t-m} + \gamma \frac{\epsilon_t}{l_{t-1} b_{t-1}} \end{aligned}\)

表2:乘法误差ETS模型

非季节性 加法 乘法
无趋势 \(\begin{aligned} &y_{t} = l_{t-1}(1 + \epsilon_t) \\ &l_t = l_{t-1}(1 + \alpha \epsilon_t) \end{aligned}\) \(\begin{aligned} &y_{t} = (l_{t-1} + s_{t-m})(1 + \epsilon_t) \\ &l_t = l_{t-1} + \alpha \mu_{y,t} \epsilon_t \\ &s_t = s_{t-m} + \gamma \mu_{y,t} \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = l_{t-1} s_{t-m}(1 + \epsilon_t) \\ &l_t = l_{t-1}(1 + \alpha \epsilon_t) \\ &s_t = s_{t-m}(1 + \gamma \epsilon_t) \end{aligned}\)
加法 \(\begin{aligned} &y_{t} = (l_{t-1} + b_{t-1})(1 + \epsilon_t) \\ &l_t = (l_{t-1} + b_{t-1})(1 + \alpha \epsilon_t) \\ &b_t = b_{t-1} + \beta \mu_{y,t} \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = (l_{t-1} + b_{t-1} + s_{t-m})(1 + \epsilon_t) \\ &l_t = l_{t-1} + b_{t-1} + \alpha \mu_{y,t} \epsilon_t \\ &b_t = b_{t-1} + \beta \mu_{y,t} \epsilon_t \\ &s_t = s_{t-m} + \gamma \mu_{y,t} \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = (l_{t-1} + b_{t-1}) s_{t-m}(1 + \epsilon_t) \\ &l_t = (l_{t-1} + b_{t-1})(1 + \alpha \epsilon_t) \\ &b_t = b_{t-1} + \beta (l_{t-1} + b_{t-1}) \epsilon_t \\ &s_t = s_{t-m} (1 + \gamma \epsilon_t) \end{aligned}\)
加法阻尼 \(\begin{aligned} &y_{t} = (l_{t-1} + \phi b_{t-1})(1 + \epsilon_t) \\ &l_t = (l_{t-1} + \phi b_{t-1})(1 + \alpha \epsilon_t) \\ &b_t = \phi b_{t-1} + \beta \mu_{y,t} \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = (l_{t-1} + \phi b_{t-1} + s_{t-m})(1 + \epsilon_t) \\ &l_t = l_{t-1} + \phi b_{t-1} + \alpha \mu_{y,t} \epsilon_t \\ &b_t = \phi b_{t-1} + \beta \mu_{y,t} \epsilon_t \\ &s_t = s_{t-m} + \gamma \mu_{y,t} \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = (l_{t-1} + \phi b_{t-1}) s_{t-m}(1 + \epsilon_t) \\ &l_t = l_{t-1} + \phi b_{t-1} (1 + \alpha \epsilon_t) \\ &b_t = \phi b_{t-1} + \beta (l_{t-1} + \phi b_{t-1}) \epsilon_t \\ &s_t = s_{t-m}(1 + \gamma \epsilon_t) \end{aligned}\)
乘法 \(\begin{aligned} &y_{t} = l_{t-1} b_{t-1} (1 + \epsilon_t) \\ &l_t = l_{t-1} b_{t-1} (1 + \alpha \epsilon_t) \\ &b_t = b_{t-1} (1 + \beta \epsilon_t) \end{aligned}\) \(\begin{aligned} &y_{t} = (l_{t-1} b_{t-1} + s_{t-m})(1 + \epsilon_t) \\ &l_t = l_{t-1} b_{t-1} + \alpha \mu_{y,t} \epsilon_t \\ &b_t = b_{t-1} + \beta \frac{\mu_{y,t}}{l_{t-1}} \epsilon_t \\ &s_t = s_{t-m} + \gamma \mu_{y,t} \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = l_{t-1} b_{t-1} s_{t-m} (1 + \epsilon_t) \\ &l_t = l_{t-1} b_{t-1} (1 + \alpha \epsilon_t) \\ &b_t = b_{t-1} (1 + \beta \epsilon_t) \\ &s_t = s_{t-m} (1 + \gamma \epsilon_t) \end{aligned}\)
乘法阻尼 \(\begin{aligned} &y_{t} = l_{t-1} b_{t-1}^\phi (1 + \epsilon_t) \\ &l_t = l_{t-1} b_{t-1}^\phi (1 + \alpha \epsilon_t) \\ &b_t = b_{t-1}^\phi (1 + \beta \epsilon_t) \end{aligned}\) \(\begin{aligned} &y_{t} = (l_{t-1} b_{t-1}^\phi + s_{t-m})(1 + \epsilon_t) \\ &l_t = l_{t-1} b_{t-1}^\phi + \alpha \mu_{y,t} \epsilon_t \\ &b_t = b_{t-1}^\phi + \beta \frac{\mu_{y,t}}{l_{t-1}} \epsilon_t \\ &s_t = s_{t-m} + \gamma \mu_{y,t} \epsilon_t \end{aligned}\) \(\begin{aligned} &y_{t} = l_{t-1} b_{t-1}^\phi s_{t-m} (1 + \epsilon_t) \\ &l_t = l_{t-1} b_{t-1}^\phi \left(1 + \alpha \epsilon_t\right) \\ &b_t = b_{t-1}^\phi \left(1 + \beta \epsilon_t\right) \\ &s_t = s_{t-m} \left(1 + \gamma \epsilon_t\right) \end{aligned}\)

从统计学的角度来看,表1和表2中的公式对应于“真实模型”,它们解释了潜在数据所依据的模型,但在它们的构建和估计中,\(\epsilon_t\)被估计值\(e_t\)替代(根据错误类型不同而有所不同),时间序列组件和光滑参数也用它们的估计值替代(例如,用\(\hat \alpha\)代替\(\alpha\))。然而,如果这些模型参数的值是已知的,就可以从这些模型中预测点估计,或得出条件h步前的预期。

模型选择

Holt Winters统计框架的一个巨大优势是可以使用信息准则进行模型选择。这里可以使用 AIC, AIC_cBIC 来确定哪一种 Holt Winters 模型最适合给定的时间序列。

对于 Holt Winters 模型,赤池信息量准则(AIC)定义为

\[\text{AIC} = -2\log(L) + 2k,\]

其中 \(L\) 是模型的似然值,\(k\) 是已估计的参数和初始状态的总数(包括残差方差)。

修正小样本偏差的 AICAIC_c)定义为

\[AIC_c = AIC + \frac{2k(k+1)}{T-k-1}\]

而贝叶斯信息量准则(BIC)为

\[\text{BIC} = \text{AIC} + k[\log(T)-2]\]

在(误差、趋势、季节性)三种组合中,某些组合可能导致数值上的困难。具体而言,可能导致此类不稳定性的模型包括 ETS(A,N,M), ETS(A,A,M)ETS(A,Ad,M),这是由于状态方程中可能接近于零的值导致的除法运算。我们通常在选择模型时不会考虑这些特定组合。

当数据严格为正时,乘法误差模型是有用的,但当数据包含零或负值时,它们在数值上并不稳定。因此,如果时间序列不是严格为正,则不考虑乘法误差模型。在这种情况下,仅应用六种完全加性模型。

加载库和数据

Tip

需要使用Statsforecast。有关安装的说明,请参见说明

接下来,我们导入绘图库并配置绘图样式。

import matplotlib.pyplot as plt
import seaborn as sns
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf
plt.style.use('grayscale') # 五三八 灰度 经典
plt.rcParams['lines.linewidth'] = 1.5
dark_style = {
    'figure.facecolor': '#008080',  # #212946
    'axes.facecolor': '#008080',
    'savefig.facecolor': '#008080',
    'axes.grid': True,
    'axes.grid.which': 'both',
    'axes.spines.left': False,
    'axes.spines.right': False,
    'axes.spines.top': False,
    'axes.spines.bottom': False,
    'grid.color': '#000000',  #2A3459
    'grid.linewidth': '1',
    'text.color': '0.9',
    'axes.labelcolor': '0.9',
    'xtick.color': '0.9',
    'ytick.color': '0.9',
    'font.size': 12 }
plt.rcParams.update(dark_style)


from pylab import rcParams
rcParams['figure.figsize'] = (18,7)

读取数据

df=pd.read_csv("https://raw.githubusercontent.com/Naren8520/Serie-de-tiempo-con-Machine-Learning/main/Data/ads.csv")
df.head()
Time Ads
0 2017-09-13T00:00:00 80115
1 2017-09-13T01:00:00 79885
2 2017-09-13T02:00:00 89325
3 2017-09-13T03:00:00 101930
4 2017-09-13T04:00:00 121630

StatsForecast 的输入始终是一个长格式的数据框,包含三列:unique_id、ds 和 y:

  • unique_id(字符串、整数或类别)表示系列的标识符。

  • ds(日期戳)列应为 Pandas 期望的格式,理想情况下是 YYYY-MM-DD 格式的日期或 YYYY-MM-DD HH:MM:SS 格式的时间戳。

  • y(数字)表示我们希望进行预测的测量值。

df["unique_id"]="1"
df.columns=["ds", "y", "unique_id"]
df.head()
ds y unique_id
0 2017-09-13T00:00:00 80115 1
1 2017-09-13T01:00:00 79885 1
2 2017-09-13T02:00:00 89325 1
3 2017-09-13T03:00:00 101930 1
4 2017-09-13T04:00:00 121630 1
print(df.dtypes)
ds           object
y             int64
unique_id    object
dtype: object

我们可以看到我们的时间变量 (ds) 是以对象格式表示的,我们需要将其转换为日期格式。

df["ds"] = pd.to_datetime(df["ds"])

使用 plot 方法探索数据

使用 StatsForecast 类的 plot 方法绘制一些系列。此方法从数据集中打印一个随机系列,并且对于基本的探索性数据分析(EDA)非常有用。

from statsforecast import StatsForecast

StatsForecast.plot(df)

增强型迪基-福勒检验

增强型迪基-福勒(ADF)检验是一种统计检验,旨在确定时间序列数据中是否存在单位根。单位根可能导致时间序列分析结果不可预测。在单位根检验中形成一个原假设,以确定时间序列数据受到趋势影响的强度。接受原假设意味着我们接受时间序列数据不平稳的证据。拒绝原假设或接受替代假设,则表明我们接受时间序列数据是由平稳过程生成的证据。这个过程也被称为平稳趋势。ADF检验统计量的值为负。较低的ADF值表示对原假设的拒绝更为强烈。

增强型迪基-福勒检验是一种常见的统计检验,用于测试给定时间序列是否平稳。我们可以通过定义原假设和替代假设来实现这一目标。

  • 原假设:时间序列是非平稳的。它给出一个与时间有关的趋势。

  • 替代假设:时间序列是平稳的。换句话说,序列不依赖于时间。

  • ADF或t统计量 < 临界值:拒绝原假设,时间序列是平稳的。

  • ADF或t统计量 > 临界值:未能拒绝原假设,时间序列是非平稳的。

from statsmodels.tsa.stattools import adfuller

def Augmented_Dickey_Fuller_Test_func(series , column_name):
    print (f'Dickey-Fuller test results for columns: {column_name}')
    dftest = adfuller(series, autolag='AIC')
    dfoutput = pd.Series(dftest[0:4], index=['Test Statistic','p-value','No Lags Used','Number of observations used'])
    for key,value in dftest[4].items():
        dfoutput['Critical Value (%s)'%key] = value
    print (dfoutput)
    if dftest[1] <= 0.05:
        print("Conclusion:====>")
        print("Reject the null hypothesis")
        print("The data is stationary")
    else:
        print("Conclusion:====>")
        print("The null hypothesis cannot be rejected")
        print("The data is not stationary")
Augmented_Dickey_Fuller_Test_func(df["y"],'Ads')
Dickey-Fuller test results for columns: Ads
Test Statistic         -7.089634e+00
p-value                 4.444804e-10
No Lags Used            9.000000e+00
                            ...     
Critical Value (1%)    -3.462499e+00
Critical Value (5%)    -2.875675e+00
Critical Value (10%)   -2.574304e+00
Length: 7, dtype: float64
Conclusion:====>
Reject the null hypothesis
The data is stationary

自相关图

自相关函数

定义1.\(\{x_t;1 ≤ t ≤ n\}\) 是来自 \(\{X_t\}\) 的大小为 n 的时间序列样本。 1. \(\bar x = \sum_{t=1}^n \frac{x_t}{n}\) 称为 \(\{X_t\}\) 的样本均值。 2. \(c_k =\sum_{t=1}^{n−k} (x_{t+k}- \bar x)(x_t−\bar x)/n\) 被称为 \(\{X_t\}\) 的样本自协方差函数。 3. \(r_k = c_k /c_0\) 被认为是 \(\{X_t\}\) 的样本自相关函数。

关于这个定义,注意以下几点:

  • 和大多数文献一样,本指南使用 ACF 来表示样本自相关函数和自相关函数。ACF 的含义可以在上下文中轻易辨识。

  • 显然,c0 是 \(\{X_t\}\) 的样本方差。此外,\(r_0 = c_0/c_0 = 1\) 并且对于任何整数 \(k, |r_k| ≤ 1\)

  • 当我们计算任何固定长度 \(n\) 的样本序列的 ACF 时,不能对大 k 的 \(r_k\) 值给予过多信任,因为随着 \(k\) 增大,可用于计算 \(r_k\)\((x_{t +k }, x_t)\) 对的数量减少。一条经验法则是对于 \(k > n/3\) 不要估计 \(r_k\),另一个是 \(n ≥ 50, k ≤ n/4\)。无论如何,保持谨慎始终是个好主意。

  • 我们还可以通过定义1计算非平稳时间序列样本的 ACF。在这种情况下,然而,ACF 或 \(r_k\) 随着 \(k\) 增加非常缓慢或几乎不减。

  • 将 ACF \((r_k)\) 绘制与滞后 \(k\) 的关系非常简单,但在分析时间序列样本时非常有帮助。这样的 ACF 图称为相关图。

  • 如果 \(\{X_t\}\) 是平稳的且 \(E(X_t)=0\) 并且所有 \(k \neq 0\)\(\rho_k =0\),即它是一个白噪声序列,则 \(r_k\) 的抽样分布在渐近上是正态分布的,均值为 0,方差为 \(1/n\)。因此,\(r_k\) 落在区间 \([−1.96/\sqrt{n}, 1.96/\sqrt{n}]\) 的概率约为 95%。

现在我们可以总结一下: (1) 如果时间序列图明显显示出趋势和/或季节性,那么它肯定是非平稳的; (2) 如果 ACF \(r_k\) 随着滞后 \(k\) 的增加非常缓慢或几乎不减,那么该时间序列也应该是非平稳的。

fig, axs = plt.subplots(nrows=1, ncols=2)

plot_acf(df["y"],  lags=30, ax=axs[0],color="fuchsia")
axs[0].set_title("Autocorrelation");

plot_pacf(df["y"],  lags=30, ax=axs[1],color="lime")
axs[1].set_title('Partial Autocorrelation')

plt.show();

时间序列的分解

如何分解时间序列以及原因是什么?

在时间序列分析中,为了预测新值,了解过去的数据非常重要。更正式地说,了解值随时间变化的模式是非常重要的。有许多原因可能导致我们的预测值出现偏差。基本上,时间序列由四个组成部分构成。这些组成部分的变化导致时间序列模式的改变。这些组成部分是:

  • 水平(Level): 这是随时间变化的主要平均值。
  • 趋势(Trend): 趋势是导致时间序列中出现上升或下降模式的值。
  • 季节性(Seasonality): 这是在时间序列中短时间内发生的周期性事件,导致时间序列中出现短期的上升或下降模式。
  • 残差/噪声(Residual/Noise): 这些是时间序列中的随机变化。

这些组成部分随时间的结合导致了时间序列的形成。大多数时间序列包含水平和噪声/残差,而趋势或季节性则是可选值。

如果季节性和趋势是时间序列的一部分,那么将会对预测值产生影响。因为预测的时间序列的模式可能与之前的时间序列不同。

时间序列中组成部分的结合可以分为两种类型: * 加法型 * 乘法型

加法型时间序列

如果时间序列的组成部分是通过相加形成的,那么该时间序列被称为加法型时间序列。通过可视化,我们可以说,如果时间序列的上升或下降模式在整个序列中是相似的,则该时间序列是加法型的。任何加法型时间序列的数学函数可以表示为: \[y(t) = level + Trend + seasonality + noise\]

乘法型时间序列

如果时间序列的组成部分是通过相乘形成的,那么该时间序列被称为乘法型时间序列。通过可视化,如果时间序列随着时间的推移呈现指数增长或衰退,则该时间序列可以视为乘法型时间序列。乘法型时间序列的数学函数可以表示为: \[y(t) = Level * Trend * seasonality * Noise\]

加性

from statsmodels.tsa.seasonal import seasonal_decompose 
a = seasonal_decompose(df["y"], model = "additive", period=24)
a.plot();

乘法性

from statsmodels.tsa.seasonal import seasonal_decompose 
a = seasonal_decompose(df["y"], model = "Multiplicative", period=24)
a.plot();

将数据分为训练集和测试集

让我们将数据分成几个集合

  1. 用于训练我们的 Holt Winters Model 的数据。
  2. 用于测试我们模型的数据。

对于测试数据,我们将使用最后30个小时的数据来测试和评估我们模型的性能。

train = df[df.ds<='2017-09-20 17:00:00'] 
test = df[df.ds>'2017-09-20 17:00:00'] 
train.shape, test.shape
((186, 3), (30, 3))

现在让我们绘制训练数据和测试数据。

sns.lineplot(train,x="ds", y="y", label="Train", linestyle="--",linewidth=2)
sns.lineplot(test, x="ds", y="y", label="Test", linewidth=2, color="yellow")
plt.title("Ads watched (hourly data)");
plt.show()

Holt-Winters 方法的实现与 StatsForecast

要更进一步了解 Holt-Winters 模型 函数的参数,下面列出了相关信息。如需更多信息,请访问 文档

season_length : int
    每单位时间的观察次数。例如:12 个月数据。
error_type : str
    ETS 模型的误差类型。可以是加法(A)或乘法(M)。
alias : str
    模型的自定义名称。
prediction_intervals : Optional[ConformalIntervals]
    计算符合预测区间的信息。
    默认情况下,该模型将计算原生的预测
    区间。

加载库

from statsforecast import StatsForecast
from statsforecast.models import HoltWinters

实例化模型

导入并实例化模型。设置参数有时可能会比较棘手。Rob Hyndmann 的这篇关于 季节周期 的文章对于 season_length 会很有帮助。

在这种情况下,我们将测试模型的两个替代方案,一个是加法模型,一个是乘法模型。

season_length = 24 # 每小时数据 
horizon = len(test) # 预测数量

models = [HoltWinters(season_length=season_length, error_type="A", alias="Add"),
          HoltWinters(season_length=season_length, error_type="M", alias="Multi")]

我们通过实例化一个新的StatsForecast对象并使用以下参数来拟合模型:

models:一个模型列表。选择您想要的模型并导入它们。

  • freq: 一个字符串,指示数据的频率。 (查看pandas可用的频率。)

  • n_jobs: n_jobs: int,用于并行处理的作业数量,使用-1表示所有核心。

  • fallback_model: 如果模型失败要使用的模型。

任何设置都传递给构造函数。然后您调用它的fit方法并传入历史数据框。

sf = StatsForecast(df=df,
                   models=models,
                   freq='H', 
                   n_jobs=-1)

拟合模型

sf.fit()
StatsForecast(models=[Add,Multi])

让我们看看我们霍尔特-温特斯模型的结果。我们可以通过以下指令观察它:

result=sf.fitted_[0,0].model_
print(result.keys())
print(result['fit'])
dict_keys(['loglik', 'aic', 'bic', 'aicc', 'mse', 'amse', 'fit', 'residuals', 'components', 'm', 'nstate', 'fitted', 'states', 'par', 'sigma2', 'n_params', 'method', 'actual_residuals'])
results(x=array([ 2.44373452e-02,  1.38129692e-03,  3.10806480e-02,  8.90103969e-01,
        1.24049460e+05, -8.68954810e+01, -3.89677433e+04, -2.62776804e+04,
       -1.53306670e+04,  1.04556014e+04,  3.40799540e+04,  3.61115645e+04,
        2.95135877e+04,  2.70191012e+04,  2.61439737e+04,  2.73344984e+04,
        2.70239387e+04,  2.94934401e+04,  2.12740663e+04,  4.02181550e+03,
       -7.08633837e+03, -1.29309441e+04, -1.23216723e+04, -9.32780505e+03,
       -3.92528658e+03, -1.18781713e+03, -2.43122964e+04, -3.56913114e+04,
       -4.32670579e+04]), fn=5078.452315148393, nit=1000, simplex=array([[ 2.40487996e-02,  1.40159178e-03,  2.98715521e-02,
         8.89367428e-01,  1.23345215e+05, -8.56751447e+01,
        -3.94103267e+04, -2.58132744e+04, -1.46476758e+04,
         1.04203903e+04,  3.50465498e+04,  3.59654863e+04,
         2.93863399e+04,  2.77527548e+04,  2.63918127e+04,
         2.74059591e+04,  2.59885977e+04,  2.94583381e+04,
         2.19996508e+04,  4.22020655e+03, -7.12092088e+03,
        -1.29452569e+04, -1.21734962e+04, -9.28051739e+03,
        -4.09832793e+03, -1.22159944e+03, -2.40047026e+04,
        -3.59264081e+04, -4.37141281e+04],
       [ 2.30584659e-02,  1.29870085e-03,  3.23295869e-02,
         8.95699018e-01,  1.23226757e+05, -8.51134029e+01,
        -3.94188347e+04, -2.66463805e+04, -1.46040706e+04,
         1.00508287e+04,  3.48121192e+04,  3.66574880e+04,
         2.95156329e+04,  2.75214239e+04,  2.61766232e+04,
         2.70638403e+04,  2.59821671e+04,  2.92220715e+04,
         2.15353216e+04,  4.11207780e+03, -7.26888692e+03,
        -1.31025981e+04, -1.24323871e+04, -9.46450899e+03,
        -3.87523173e+03, -1.29619824e+03, -2.48262400e+04,
        -3.65297643e+04, -4.26697834e+04],
       [ 2.38566990e-02,  1.38393698e-03,  2.86485845e-02,
         8.91732054e-01,  1.22667885e+05, -8.63064081e+01,
        -3.95610795e+04, -2.64459058e+04, -1.50401318e+04,
         1.03748608e+04,  3.46772905e+04,  3.61258281e+04,
         3.03028728e+04,  2.71476029e+04,  2.59766785e+04,
         2.83781299e+04,  2.62171593e+04,  2.87088506e+04,
         2.11399641e+04,  4.26251638e+03, -6.86374732e+03,
        -1.30572982e+04, -1.25961133e+04, -9.28515605e+03,
        -3.92065232e+03, -1.29918039e+03, -2.41322902e+04,
        -3.62733255e+04, -4.35070867e+04],
       [ 2.27798959e-02,  1.34604573e-03,  3.14497882e-02,
         8.97672321e-01,  1.23457441e+05, -8.39505753e+01,
        -3.84752143e+04, -2.70271535e+04, -1.51565030e+04,
         1.01558921e+04,  3.56374371e+04,  3.57659522e+04,
         2.92214649e+04,  2.64450530e+04,  2.59612798e+04,
         2.79127971e+04,  2.65311803e+04,  2.95053263e+04,
         2.20459941e+04,  4.27175853e+03, -7.11671275e+03,
        -1.28483284e+04, -1.20024908e+04, -9.14008142e+03,
        -3.94099609e+03, -1.27913818e+03, -2.43917933e+04,
        -3.59036544e+04, -4.37145499e+04],
       [ 2.29779847e-02,  1.35988450e-03,  3.20546039e-02,
         8.97392076e-01,  1.23999353e+05, -8.57416139e+01,
        -3.86796082e+04, -2.61112159e+04, -1.46410300e+04,
         1.00447411e+04,  3.46746838e+04,  3.56015924e+04,
         2.92653819e+04,  2.72512111e+04,  2.64957723e+04,
         2.82119577e+04,  2.61112530e+04,  2.90487409e+04,
         2.15387214e+04,  4.05550604e+03, -6.88742058e+03,
        -1.28509450e+04, -1.21550165e+04, -9.26733160e+03,
        -4.09280472e+03, -1.31561243e+03, -2.42237989e+04,
        -3.58787629e+04, -4.36947890e+04],
       [ 2.34931397e-02,  1.36066473e-03,  2.92309535e-02,
         8.94631252e-01,  1.23220366e+05, -7.99821563e+01,
        -3.97471338e+04, -2.62212988e+04, -1.54426040e+04,
         1.01804753e+04,  3.45628830e+04,  3.61783359e+04,
         2.83321227e+04,  2.75443313e+04,  2.65554748e+04,
         2.85866656e+04,  2.65994443e+04,  2.99352483e+04,
         2.12904274e+04,  4.19579939e+03, -7.04756866e+03,
        -1.35378259e+04, -1.21517302e+04, -9.38313255e+03,
        -3.93098135e+03, -1.23125948e+03, -2.43161758e+04,
        -3.56808651e+04, -4.30934906e+04],
       [ 2.34226303e-02,  1.36919062e-03,  3.05044515e-02,
         8.93651006e-01,  1.22119669e+05, -8.45830782e+01,
        -3.83859519e+04, -2.66545546e+04, -1.46663235e+04,
         1.03044558e+04,  3.37863355e+04,  3.58787384e+04,
         2.98614652e+04,  2.79521053e+04,  2.65483349e+04,
         2.80130243e+04,  2.67030580e+04,  2.95716350e+04,
         2.13867299e+04,  4.20532353e+03, -7.19231268e+03,
        -1.29869745e+04, -1.22214265e+04, -9.35200831e+03,
        -3.99601887e+03, -1.25191151e+03, -2.45210288e+04,
        -3.55315031e+04, -4.33561528e+04],
       [ 2.33686092e-02,  1.38777500e-03,  3.05842986e-02,
         8.95765424e-01,  1.24435621e+05, -8.44258777e+01,
        -3.95388846e+04, -2.70585507e+04, -1.42279741e+04,
         1.04108267e+04,  3.41935238e+04,  3.61927880e+04,
         2.96134192e+04,  2.72521230e+04,  2.63894944e+04,
         2.73577833e+04,  2.63681631e+04,  2.89873337e+04,
         2.13122122e+04,  4.27459659e+03, -7.10258507e+03,
        -1.26644787e+04, -1.22537168e+04, -9.49233627e+03,
        -4.01010478e+03, -1.23803020e+03, -2.39978907e+04,
        -3.60888141e+04, -4.29553133e+04],
       [ 2.34816233e-02,  1.37991019e-03,  3.06089532e-02,
         8.94874645e-01,  1.23617905e+05, -8.40833189e+01,
        -3.85116704e+04, -2.63905726e+04, -1.43692871e+04,
         1.05551899e+04,  3.51276038e+04,  3.58524187e+04,
         2.94603806e+04,  2.74703366e+04,  2.55213418e+04,
         2.83361549e+04,  2.60687975e+04,  2.86600785e+04,
         2.14422052e+04,  4.10932930e+03, -7.16315898e+03,
        -1.34079834e+04, -1.23877734e+04, -9.20876927e+03,
        -3.95489699e+03, -1.23349778e+03, -2.47612959e+04,
        -3.48066942e+04, -4.33211841e+04],
       [ 2.41575092e-02,  1.39907614e-03,  3.02357059e-02,
         8.91102368e-01,  1.23794366e+05, -8.48796907e+01,
        -3.90608848e+04, -2.64189926e+04, -1.45758290e+04,
         9.86006031e+03,  3.44455393e+04,  3.62279168e+04,
         2.93705791e+04,  2.75920864e+04,  2.50923293e+04,
         2.84262381e+04,  2.67984741e+04,  2.89640135e+04,
         2.16286886e+04,  4.29540433e+03, -7.25171917e+03,
        -1.29078390e+04, -1.19948534e+04, -9.29963339e+03,
        -3.99568567e+03, -1.24285574e+03, -2.43650156e+04,
        -3.65231527e+04, -4.28366363e+04],
       [ 2.36961761e-02,  1.38273987e-03,  3.10961396e-02,
         8.93478329e-01,  1.22419456e+05, -8.54817124e+01,
        -3.83003011e+04, -2.62799796e+04, -1.47684198e+04,
         1.03088809e+04,  3.52182251e+04,  3.52679730e+04,
         2.92136442e+04,  2.74002342e+04,  2.61909267e+04,
         2.78569566e+04,  2.68610704e+04,  2.89820187e+04,
         2.15006485e+04,  4.28929006e+03, -6.97610950e+03,
        -1.28970410e+04, -1.24884156e+04, -9.67725756e+03,
        -3.89959668e+03, -1.20865626e+03, -2.39696630e+04,
        -3.60623458e+04, -4.25044368e+04],
       [ 2.33721259e-02,  1.37576033e-03,  3.07287778e-02,
         8.95381214e-01,  1.23361917e+05, -8.59292902e+01,
        -3.97494538e+04, -2.59881823e+04, -1.42733538e+04,
         1.01850859e+04,  3.48288372e+04,  3.50586273e+04,
         2.97507897e+04,  2.71064540e+04,  2.67061201e+04,
         2.78353166e+04,  2.65264771e+04,  2.98277569e+04,
         2.14489843e+04,  4.23889638e+03, -7.08556474e+03,
        -1.31498142e+04, -1.22079975e+04, -9.28165029e+03,
        -3.89871487e+03, -1.25831398e+03, -2.48016766e+04,
        -3.59037067e+04, -4.20714772e+04],
       [ 2.26415137e-02,  1.34430820e-03,  3.23341885e-02,
         8.99285053e-01,  1.22878310e+05, -8.53335700e+01,
        -3.92562952e+04, -2.59721752e+04, -1.48816794e+04,
         1.00898476e+04,  3.50587088e+04,  3.62440964e+04,
         2.98438930e+04,  2.79605178e+04,  2.61637574e+04,
         2.76187001e+04,  2.60357145e+04,  2.90427564e+04,
         2.07810067e+04,  4.42036847e+03, -7.03406014e+03,
        -1.26056030e+04, -1.19217451e+04, -8.98726114e+03,
        -3.90645105e+03, -1.24678256e+03, -2.42136819e+04,
        -3.58509789e+04, -4.43751915e+04],
       [ 2.31801372e-02,  1.37749747e-03,  2.93541139e-02,
         8.95770558e-01,  1.22858350e+05, -8.80052223e+01,
        -3.86617013e+04, -2.59278016e+04, -1.47751488e+04,
         1.03042927e+04,  3.41708552e+04,  3.56988813e+04,
         2.91682837e+04,  2.67894217e+04,  2.58456196e+04,
         2.70119460e+04,  2.68253248e+04,  2.93127849e+04,
         2.12248716e+04,  4.22651962e+03, -7.12436510e+03,
        -1.29474406e+04, -1.21520227e+04, -9.53356225e+03,
        -3.96352288e+03, -1.27666223e+03, -2.44328890e+04,
        -3.56137046e+04, -4.35671367e+04],
       [ 2.30363008e-02,  1.34895700e-03,  3.31052656e-02,
         8.96994359e-01,  1.23002355e+05, -8.50280976e+01,
        -3.91939881e+04, -2.63829796e+04, -1.46326214e+04,
         1.04634298e+04,  3.45295261e+04,  3.54928661e+04,
         2.95432569e+04,  2.75345263e+04,  2.50916467e+04,
         2.68369024e+04,  2.68635715e+04,  3.03912706e+04,
         2.11714956e+04,  4.28329750e+03, -6.90261152e+03,
        -1.30645780e+04, -1.22326849e+04, -9.08641469e+03,
        -4.04851644e+03, -1.32494294e+03, -2.38092257e+04,
        -3.57905204e+04, -4.35116477e+04],
       [ 2.44373452e-02,  1.38129692e-03,  3.10806480e-02,
         8.90103969e-01,  1.24049460e+05, -8.68954810e+01,
        -3.89677433e+04, -2.62776804e+04, -1.53306670e+04,
         1.04556014e+04,  3.40799540e+04,  3.61115645e+04,
         2.95135877e+04,  2.70191012e+04,  2.61439737e+04,
         2.73344984e+04,  2.70239387e+04,  2.94934401e+04,
         2.12740663e+04,  4.02181550e+03, -7.08633837e+03,
        -1.29309441e+04, -1.23216723e+04, -9.32780505e+03,
        -3.92528658e+03, -1.18781713e+03, -2.43122964e+04,
        -3.56913114e+04, -4.32670579e+04],
       [ 2.34225937e-02,  1.37504026e-03,  3.20762187e-02,
         8.94715930e-01,  1.22649715e+05, -8.51815678e+01,
        -3.92349107e+04, -2.60372086e+04, -1.46713893e+04,
         1.03940388e+04,  3.42200364e+04,  3.54151977e+04,
         2.86168968e+04,  2.72219207e+04,  2.65021110e+04,
         2.77555263e+04,  2.61144978e+04,  2.88679128e+04,
         2.14504642e+04,  4.26955607e+03, -7.33178601e+03,
        -1.28672736e+04, -1.23366734e+04, -9.27407389e+03,
        -3.98684306e+03, -1.26689130e+03, -2.42014964e+04,
        -3.60431164e+04, -4.35735281e+04],
       [ 2.30245102e-02,  1.33633606e-03,  3.17433016e-02,
         8.95865886e-01,  1.23029482e+05, -8.29602756e+01,
        -3.96755895e+04, -2.63099175e+04, -1.48152395e+04,
         1.05185870e+04,  3.52063359e+04,  3.52312467e+04,
         2.99077028e+04,  2.71533181e+04,  2.66237164e+04,
         2.78798999e+04,  2.74977313e+04,  2.80935898e+04,
         2.09427481e+04,  4.09547938e+03, -7.09751151e+03,
        -1.29766286e+04, -1.18250113e+04, -9.54750290e+03,
        -3.97635927e+03, -1.32445291e+03, -2.43437156e+04,
        -3.60393575e+04, -4.30058401e+04],
       [ 2.35957957e-02,  1.31866210e-03,  3.04466963e-02,
         8.94704385e-01,  1.23506905e+05, -8.81033940e+01,
        -3.92910011e+04, -2.67861535e+04, -1.43723619e+04,
         1.02596397e+04,  3.42680222e+04,  3.49198652e+04,
         2.96607616e+04,  2.71335556e+04,  2.62749412e+04,
         2.78456008e+04,  2.67384571e+04,  2.90090597e+04,
         2.21391413e+04,  4.25656010e+03, -6.97489434e+03,
        -1.31391614e+04, -1.21060497e+04, -9.02606877e+03,
        -3.96765488e+03, -1.25601904e+03, -2.41311166e+04,
        -3.64834934e+04, -4.39441979e+04],
       [ 2.36751650e-02,  1.36543428e-03,  2.91691871e-02,
         8.94846111e-01,  1.22793505e+05, -8.69820381e+01,
        -3.90199341e+04, -2.60601993e+04, -1.47352545e+04,
         1.01516410e+04,  3.51794657e+04,  3.57287886e+04,
         2.98594202e+04,  2.80104436e+04,  2.60402555e+04,
         2.73922257e+04,  2.62342050e+04,  2.91820141e+04,
         2.17442081e+04,  4.06104487e+03, -7.27342325e+03,
        -1.29440707e+04, -1.20623302e+04, -9.39257599e+03,
        -4.11688432e+03, -1.22491978e+03, -2.42825476e+04,
        -3.61681582e+04, -4.31967736e+04],
       [ 2.36356407e-02,  1.40249143e-03,  3.11055579e-02,
         8.94345588e-01,  1.22524471e+05, -8.18301894e+01,
        -3.97067939e+04, -2.64168124e+04, -1.48209981e+04,
         1.02443748e+04,  3.46517106e+04,  3.60741202e+04,
         2.94602626e+04,  2.75042300e+04,  2.60176469e+04,
         2.71780405e+04,  2.62591460e+04,  2.90272777e+04,
         2.21503102e+04,  4.14004508e+03, -6.53769777e+03,
        -1.28541437e+04, -1.20646212e+04, -9.61019431e+03,
        -4.01591048e+03, -1.26480647e+03, -2.51893997e+04,
        -3.57918334e+04, -4.39872769e+04],
       [ 2.31857120e-02,  1.35188810e-03,  3.04128998e-02,
         8.96437230e-01,  1.23431237e+05, -8.47006111e+01,
        -3.99591747e+04, -2.66301164e+04, -1.49368574e+04,
         1.04096179e+04,  3.42782255e+04,  3.52664081e+04,
         2.93964291e+04,  2.71483335e+04,  2.55844710e+04,
         2.85473780e+04,  2.62439803e+04,  2.95137214e+04,
         2.16999982e+04,  4.21616087e+03, -7.13754628e+03,
        -1.20070878e+04, -1.24988683e+04, -9.28901448e+03,
        -3.95465776e+03, -1.31514929e+03, -2.43986957e+04,
        -3.54573986e+04, -4.29463528e+04],
       [ 2.37859054e-02,  1.39168374e-03,  3.09656129e-02,
         8.93949631e-01,  1.22688446e+05, -8.27726967e+01,
        -3.85019011e+04, -2.61638638e+04, -1.47580029e+04,
         1.02622252e+04,  3.44676069e+04,  3.60705301e+04,
         2.96185616e+04,  2.73176309e+04,  2.65596714e+04,
         2.75277449e+04,  2.62402989e+04,  2.92453590e+04,
         2.17648536e+04,  4.23176214e+03, -7.23546281e+03,
        -1.29448992e+04, -1.21012800e+04, -9.30023489e+03,
        -3.83220485e+03, -1.35293557e+03, -2.39553815e+04,
        -3.50787839e+04, -4.32983892e+04],
       [ 2.31755551e-02,  1.36972669e-03,  3.10073566e-02,
         8.95714770e-01,  1.23159605e+05, -8.46149210e+01,
        -3.94983330e+04, -2.62406849e+04, -1.47819273e+04,
         1.01373745e+04,  3.46650362e+04,  3.63296076e+04,
         2.87720185e+04,  2.77427498e+04,  2.67639279e+04,
         2.72905992e+04,  2.74381779e+04,  2.87693629e+04,
         2.18008833e+04,  4.23717627e+03, -7.10957534e+03,
        -1.29620482e+04, -1.25000167e+04, -9.03037734e+03,
        -3.96084115e+03, -1.27278527e+03, -2.43366281e+04,
        -3.52451140e+04, -4.32780284e+04],
       [ 2.36133672e-02,  1.35985554e-03,  3.11513958e-02,
         8.94275972e-01,  1.23447227e+05, -7.88721079e+01,
        -3.84923510e+04, -2.59996897e+04, -1.50010097e+04,
         1.02508651e+04,  3.47205430e+04,  3.55858360e+04,
         3.04894948e+04,  2.62117116e+04,  2.68062560e+04,
         2.76499265e+04,  2.64629471e+04,  2.93271473e+04,
         2.09505773e+04,  4.35292061e+03, -7.25190305e+03,
        -1.28526873e+04, -1.24771698e+04, -8.93110886e+03,
        -4.29429884e+03, -1.25765558e+03, -2.46098978e+04,
        -3.55435309e+04, -4.37079615e+04],
       [ 2.36439940e-02,  1.33889715e-03,  3.01748889e-02,
         8.94229932e-01,  1.23290474e+05, -8.80150842e+01,
        -3.88097928e+04, -2.63137034e+04, -1.44499856e+04,
         1.03428185e+04,  3.49849551e+04,  3.63993785e+04,
         2.87588452e+04,  2.74069044e+04,  2.61451007e+04,
         2.81270537e+04,  2.64494936e+04,  2.92261510e+04,
         2.06734305e+04,  4.38913406e+03, -6.90340334e+03,
        -1.27939382e+04, -1.19989815e+04, -9.20233589e+03,
        -4.06279150e+03, -1.35076768e+03, -2.50712735e+04,
        -3.61164361e+04, -4.28106292e+04],
       [ 2.38312121e-02,  1.35348442e-03,  3.11372755e-02,
         8.92879959e-01,  1.23435424e+05, -8.75236266e+01,
        -3.98877150e+04, -2.63840774e+04, -1.41741810e+04,
         1.00040238e+04,  3.50671919e+04,  3.62991187e+04,
         2.94645434e+04,  2.66796707e+04,  2.60886751e+04,
         2.78872141e+04,  2.65673983e+04,  2.94202012e+04,
         2.15059513e+04,  4.29327789e+03, -6.99260131e+03,
        -1.29350082e+04, -1.21356464e+04, -9.63679687e+03,
        -3.99062189e+03, -1.21073235e+03, -2.41469724e+04,
        -3.49995291e+04, -4.46667304e+04],
       [ 2.28820531e-02,  1.34709017e-03,  3.12930885e-02,
         8.98042023e-01,  1.23363231e+05, -8.40520054e+01,
        -3.88532804e+04, -2.57281498e+04, -1.43468660e+04,
         1.04202224e+04,  3.44137444e+04,  3.67956433e+04,
         2.94663774e+04,  2.69333507e+04,  2.59652714e+04,
         2.82039633e+04,  2.67333645e+04,  2.91868638e+04,
         2.16875935e+04,  4.18880730e+03, -7.00016622e+03,
        -1.28357717e+04, -1.21595618e+04, -9.33735990e+03,
        -3.90969568e+03, -1.27052486e+03, -2.42933318e+04,
        -3.69259218e+04, -4.41040622e+04],
       [ 2.38049951e-02,  1.39030825e-03,  3.06405290e-02,
         8.92994409e-01,  1.24289098e+05, -8.40140819e+01,
        -3.85337082e+04, -2.61684441e+04, -1.46993627e+04,
         1.01634568e+04,  3.47937273e+04,  3.51012725e+04,
         2.94226408e+04,  2.77672823e+04,  2.62946187e+04,
         2.75564784e+04,  2.68482111e+04,  2.92367476e+04,
         2.13296806e+04,  4.24587345e+03, -7.14882503e+03,
        -1.30628882e+04, -1.22409211e+04, -9.52582674e+03,
        -3.83384250e+03, -1.23631839e+03, -2.45485642e+04,
        -3.59539920e+04, -4.45739752e+04],
       [ 2.37778011e-02,  1.37778734e-03,  3.02599768e-02,
         8.92655976e-01,  1.22304944e+05, -8.16288959e+01,
        -3.89852188e+04, -2.65123655e+04, -1.43437856e+04,
         1.02586834e+04,  3.55133375e+04,  3.57277510e+04,
         2.93369637e+04,  2.69961210e+04,  2.64681386e+04,
         2.77244225e+04,  2.68049260e+04,  2.95848184e+04,
         2.10648643e+04,  4.03875969e+03, -7.09865270e+03,
        -1.28342462e+04, -1.23781572e+04, -9.07575255e+03,
        -4.05235581e+03, -1.28589801e+03, -2.43611569e+04,
        -3.65053195e+04, -4.43511137e+04]]))

现在让我们可视化我们模型的拟合值。

如我们所见,上述结果以字典形式输出,为了从字典中提取每个元素,我们将使用 .get() 函数来提取元素,然后将其保存到 pd.DataFrame() 中。

residual=pd.DataFrame(result.get("residuals"), columns=["residual Model"])
residual
residual Model
0 -2012.193015
1 -699.563899
2 1246.127401
... ...
213 1355.434583
214 4926.937284
215 2467.814647

216 rows × 1 columns

import scipy.stats as stats

fig, axs = plt.subplots(nrows=2, ncols=2)

residual.plot(ax=axs[0,0])
axs[0,0].set_title("Residuals");

sns.distplot(residual, ax=axs[0,1]);
axs[0,1].set_title("Density plot - Residual");

stats.probplot(residual["residual Model"], dist="norm", plot=axs[1,0])
axs[1,0].set_title('Plot Q-Q')

plot_acf(residual,  lags=35, ax=axs[1,1],color="fuchsia")
axs[1,1].set_title("Autocorrelation");

plt.show();

预测方法

如果您想在生产环境中提高速度,特别是在处理多个序列或模型时,我们建议使用 StatsForecast.forecast 方法,而不是 .fit.predict

主要区别在于 .forecast 不存储拟合值,并且在分布式环境中具有很好的可扩展性。

预测方法接受两个参数:下一个 h(时间跨度)和 level

  • h (int): 表示预测未来 h 步。在这种情况下,为 30 小时。

  • level (list of floats): 此可选参数用于概率预测。设置预测区间的水平(或置信百分位)。例如,level=[90] 意味着模型期望真实值在该区间内的概率为 90%。

此处的预测对象是一个新的数据框,其中包含模型名称和 y hat 值的列,以及不确定性区间的列。根据您的计算机,此步骤大约需要 1 分钟。

Y_hat = sf.forecast(horizon, fitted=True)
Y_hat
ds Add Multi
unique_id
1 2017-09-22 00:00:00 75020.843750 75016.414062
1 2017-09-22 01:00:00 73714.437500 74739.921875
1 2017-09-22 02:00:00 81196.570312 80907.437500
... ... ... ...
1 2017-09-23 03:00:00 92962.234375 92432.453125
1 2017-09-23 04:00:00 115738.320312 112928.710938
1 2017-09-23 05:00:00 113103.835938 111366.179688

30 rows × 3 columns

通过预测方法,我们还可以从模型中提取拟合值并进行图形可视化,使用以下指令可以实现这一点。

values=sf.forecast_fitted_values()
values.head()
ds y Add Multi
unique_id
1 2017-09-13 00:00:00 80115.0 82127.195312 77955.531250
1 2017-09-13 01:00:00 79885.0 80584.562500 77934.929688
1 2017-09-13 02:00:00 89325.0 88078.875000 84123.835938
1 2017-09-13 03:00:00 101930.0 99432.601562 96152.687500
1 2017-09-13 04:00:00 121630.0 122571.562500 117567.875000
StatsForecast.plot(values)

添加95%置信区间与预测方法

sf.forecast(h=horizon, level=[95])
ds Add Add-lo-95 Add-hi-95 Multi Multi-lo-95 Multi-hi-95
unique_id
1 2017-09-22 00:00:00 75020.843750 56622.953125 93418.734375 75016.414062 62760.343750 87272.476562
1 2017-09-22 01:00:00 73714.437500 55310.503906 92118.367188 74739.921875 62483.855469 86995.992188
1 2017-09-22 02:00:00 81196.570312 62786.074219 99607.062500 80907.437500 68651.367188 93163.500000
... ... ... ... ... ... ... ...
1 2017-09-23 03:00:00 92962.234375 74304.234375 111620.242188 92432.453125 79563.328125 105301.585938
1 2017-09-23 04:00:00 115738.320312 97069.117188 134407.531250 112928.710938 100059.578125 125797.835938
1 2017-09-23 05:00:00 113103.835938 94423.414062 131784.265625 111366.179688 98497.054688 124235.304688

30 rows × 7 columns

Y_hat=Y_hat.reset_index()
Y_hat
unique_id ds Add Multi
0 1 2017-09-22 00:00:00 75020.843750 75016.414062
1 1 2017-09-22 01:00:00 73714.437500 74739.921875
2 1 2017-09-22 02:00:00 81196.570312 80907.437500
... ... ... ... ...
27 1 2017-09-23 03:00:00 92962.234375 92432.453125
28 1 2017-09-23 04:00:00 115738.320312 112928.710938
29 1 2017-09-23 05:00:00 113103.835938 111366.179688

30 rows × 4 columns

Y_hat1 = pd.concat([df,Y_hat],  keys=['unique_id', 'ds'])
Y_hat1
ds y unique_id Add Multi
unique_id 0 2017-09-13 00:00:00 80115.0 1 NaN NaN
1 2017-09-13 01:00:00 79885.0 1 NaN NaN
2 2017-09-13 02:00:00 89325.0 1 NaN NaN
... ... ... ... ... ... ...
ds 27 2017-09-23 03:00:00 NaN 1 92962.234375 92432.453125
28 2017-09-23 04:00:00 NaN 1 115738.320312 112928.710938
29 2017-09-23 05:00:00 NaN 1 113103.835938 111366.179688

246 rows × 5 columns

fig, ax = plt.subplots(1, 1)
plot_df = pd.concat([train, Y_hat1]).set_index('ds').tail(300)
plot_df['y'].plot(ax=ax, linewidth=2)
plot_df[ "Add"].plot(ax=ax, linewidth=2, color="yellow")
plot_df[ "Multi"].plot(ax=ax, linewidth=2, color="red")
ax.set_title(' Forecast', fontsize=22)
ax.set_ylabel('Ads watched (hourly data)', fontsize=20)
ax.set_xlabel('Monthly [t]', fontsize=20)
ax.legend(prop={'size': 15})
ax.grid(True)

带置信区间的预测方法

要生成预测,请使用预测方法。

预测方法接受两个参数:预测下一个 h(表示时间范围)和 level

  • h (int): 表示预测向未来 h 步的值。在这种情况下,表示 30 小时。

  • level (list of floats): 此可选参数用于概率预测。设置您预测区间的水平(或置信百分位数)。例如,level=[95] 意味着模型期望真实值在该区间内出现的概率为 95%。

此处的预测对象是一个新的数据框,其中包含模型名称及 y hat 值的列,以及不确定性区间的列。

此步骤应少于 1 秒。

sf.predict(h=horizon)
ds Add Multi
unique_id
1 2017-09-22 00:00:00 75020.843750 75016.414062
1 2017-09-22 01:00:00 73714.437500 74739.921875
1 2017-09-22 02:00:00 81196.570312 80907.437500
... ... ... ...
1 2017-09-23 03:00:00 92962.234375 92432.453125
1 2017-09-23 04:00:00 115738.320312 112928.710938
1 2017-09-23 05:00:00 113103.835938 111366.179688

30 rows × 3 columns

forecast_df = sf.predict(h=horizon, level=[80,95]) 
forecast_df
ds Add Add-lo-95 Add-lo-80 Add-hi-80 Add-hi-95 Multi Multi-lo-95 Multi-lo-80 Multi-hi-80 Multi-hi-95
unique_id
1 2017-09-22 00:00:00 75020.843750 56622.953125 62991.109375 87050.578125 93418.734375 75016.414062 62760.343750 67002.601562 83030.218750 87272.476562
1 2017-09-22 01:00:00 73714.437500 55310.503906 61680.750000 85748.117188 92118.367188 74739.921875 62483.855469 66726.109375 82753.734375 86995.992188
1 2017-09-22 02:00:00 81196.570312 62786.074219 69158.593750 93234.546875 99607.062500 80907.437500 68651.367188 72893.625000 88921.242188 93163.500000
... ... ... ... ... ... ... ... ... ... ... ...
1 2017-09-23 03:00:00 92962.234375 74304.234375 80762.421875 105162.054688 111620.242188 92432.453125 79563.328125 84017.789062 100847.125000 105301.585938
1 2017-09-23 04:00:00 115738.320312 97069.117188 103531.187500 127945.460938 134407.531250 112928.710938 100059.578125 104514.039062 121343.375000 125797.835938
1 2017-09-23 05:00:00 113103.835938 94423.414062 100889.359375 125318.312500 131784.265625 111366.179688 98497.054688 102951.507812 119780.851562 124235.304688

30 rows × 11 columns

我们可以使用pandas函数pd.concat()将预测结果与历史数据连接起来,然后能够使用该结果进行绘图。

pd.concat([df, forecast_df]).set_index('ds')
y unique_id Add Add-lo-95 Add-lo-80 Add-hi-80 Add-hi-95 Multi Multi-lo-95 Multi-lo-80 Multi-hi-80 Multi-hi-95
ds
2017-09-13 00:00:00 80115.0 1 NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
2017-09-13 01:00:00 79885.0 1 NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
2017-09-13 02:00:00 89325.0 1 NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
... ... ... ... ... ... ... ... ... ... ... ... ...
2017-09-23 03:00:00 NaN NaN 92962.234375 74304.234375 80762.421875 105162.054688 111620.242188 92432.453125 79563.328125 84017.789062 100847.125000 105301.585938
2017-09-23 04:00:00 NaN NaN 115738.320312 97069.117188 103531.187500 127945.460938 134407.531250 112928.710938 100059.578125 104514.039062 121343.375000 125797.835938
2017-09-23 05:00:00 NaN NaN 113103.835938 94423.414062 100889.359375 125318.312500 131784.265625 111366.179688 98497.054688 102951.507812 119780.851562 124235.304688

246 rows × 12 columns

df_plot= pd.concat([df, forecast_df]).set_index('ds')
df_plot
y unique_id Add Add-lo-95 Add-lo-80 Add-hi-80 Add-hi-95 Multi Multi-lo-95 Multi-lo-80 Multi-hi-80 Multi-hi-95
ds
2017-09-13 00:00:00 80115.0 1 NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
2017-09-13 01:00:00 79885.0 1 NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
2017-09-13 02:00:00 89325.0 1 NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN
... ... ... ... ... ... ... ... ... ... ... ... ...
2017-09-23 03:00:00 NaN NaN 92962.234375 74304.234375 80762.421875 105162.054688 111620.242188 92432.453125 79563.328125 84017.789062 100847.125000 105301.585938
2017-09-23 04:00:00 NaN NaN 115738.320312 97069.117188 103531.187500 127945.460938 134407.531250 112928.710938 100059.578125 104514.039062 121343.375000 125797.835938
2017-09-23 05:00:00 NaN NaN 113103.835938 94423.414062 100889.359375 125318.312500 131784.265625 111366.179688 98497.054688 102951.507812 119780.851562 124235.304688

246 rows × 12 columns

None现在让我们可视化我们的预测结果和时间序列的历史数据,同时绘制我们在进行95%置信度预测时获得的置信区间

def plot_forecasts(y_hist, y_pred, models):
    _, ax = plt.subplots(1, 1, figsize = (20, 7))
    df_plot = pd.concat([y_hist, y_pred]).set_index('ds').tail(12*10)
    df_plot[['y'] + models].plot(ax=ax, linewidth=3 , )
    colors = ['black', "lime"]
    for model, color in zip(models, colors):
        ax.fill_between(df_plot.index, 
                        df_plot[f'{model}-lo-95'], 
                        df_plot[f'{model}-hi-95'],
                        alpha=.35,
                        color=color,
                        label=f'HoltWinter{model}-level-95')
    for model, color in zip(models, colors):
        ax.fill_between(df_plot.index, 
                        df_plot[f'{model}-lo-80'], 
                        df_plot[f'{model}-hi-80'],
                        alpha=.35,
                        color=color,
                        label=f'HoltWinter{model}-level-80')
    ax.set_title('', fontsize=22)
    ax.set_ylabel("Ads watched (hourly data)", fontsize=20)
    ax.set_xlabel('Hourly', fontsize=20)
    ax.legend(prop={'size': 20})
    ax.grid(True)
    plt.show()
plot_forecasts(df, forecast_df, models=["Add"])

plot_forecasts(df, forecast_df, models=["Multi"])

让我们使用 Statsforecast 中的 plot 函数绘制相同的图,如下所示。

sf.plot(df, forecast_df)

交叉验证

在之前的步骤中,我们使用历史数据来预测未来。然而,为了评估其准确性,我们还希望了解模型在过去的表现。为了评估模型在数据上的准确性和稳健性,执行交叉验证。

对于时间序列数据,交叉验证是通过在历史数据上定义一个滑动窗口并预测随后的时期来进行的。这种形式的交叉验证使我们能够更好地估计模型在更广泛时间实例上的预测能力,同时保持训练集中的数据连续,这也是我们模型所要求的。

以下图表描述了这种交叉验证策略:

执行时间序列交叉验证

时间序列模型的交叉验证被认为是最佳实践,但大多数实现都非常慢。statsforecast库将交叉验证实现为分布式操作,使得执行过程不那么耗时。如果您有大型数据集,您还可以使用Ray、Dask或Spark在分布式集群中执行交叉验证。

在这种情况下,我们希望评估每个模型在最后5个月的表现(n_windows=),每两个月进行一次预测(step_size=12)。根据您的计算机,这一步骤应该大约需要1分钟。

StatsForecast类的cross_validation方法接受以下参数。

  • df: 训练数据框

  • h (int): 表示预测的未来h步。在这种情况下,是12个月。

  • step_size (int): 每个窗口之间的步长。换句话说:您希望多频繁地运行预测过程。

  • n_windows(int): 用于交叉验证的窗口数量。换句话说:您想评估过去多少次的预测过程。

crossvalidation_df = sf.cross_validation(df=df,
                                         h=horizon,
                                         step_size=30,
                                         n_windows=3)

crossvaldation_df对象是一个新的数据框,包括以下列:

  • unique_id: 索引。如果您不喜欢使用索引,只需运行 crossvalidation_df.resetindex()
  • ds: 日期戳或时间索引
  • cutoff: n_windows 的最后一个日期戳或时间索引。
  • y: 真实值
  • model: 包含模型名称和拟合值的列。
crossvalidation_df
ds cutoff y Add Multi
unique_id
1 2017-09-18 06:00:00 2017-09-18 05:00:00 99440.0 134578.328125 133820.109375
1 2017-09-18 07:00:00 2017-09-18 05:00:00 97655.0 133548.781250 133734.000000
1 2017-09-18 08:00:00 2017-09-18 05:00:00 97655.0 134798.656250 135216.046875
... ... ... ... ... ...
1 2017-09-21 21:00:00 2017-09-20 17:00:00 103080.0 103021.726562 103086.851562
1 2017-09-21 22:00:00 2017-09-20 17:00:00 95155.0 89544.054688 90028.406250
1 2017-09-21 23:00:00 2017-09-20 17:00:00 80285.0 78090.210938 78823.953125

90 rows × 5 columns

模型评估

我们现在可以使用适当的准确性度量来计算预测的准确性。在这里我们将使用均方根误差(RMSE)。为此,我们首先需要 install datasetsforecast,这是一个由 Nixtla 开发的Python库,包含计算RMSE的函数。

%%capture
!pip install datasetsforecast

计算RMSE的函数需要两个参数:

  1. 实际值。
  2. 预测值,在本例中为 霍尔特-温特斯模型
from datasetsforecast.losses import mse, mae, rmse

def evaluate_cross_validation(df, metric):
    models = df.drop(columns=['ds', 'cutoff', 'y']).columns.tolist()
    evals = []
    for model in models:
        eval_ = df.groupby(['unique_id', 'cutoff']).apply(lambda x: metric(x['y'].values, x[model].values)).to_frame() # 计算每个唯一标识符、模型和截止时间的损失。
        eval_.columns = [model]
        evals.append(eval_)
    evals = pd.concat(evals, axis=1)
    evals = evals.groupby(['unique_id']).mean(numeric_only=True) # 对每种模型和unique_id组合的所有截止值的误差指标取平均值。
    evals['best_model'] = evals.idxmin(axis=1)
    return evals
evaluation_df = evaluate_cross_validation(crossvalidation_df, rmse)

evaluation_df
Add Multi best_model
unique_id
1 12742.019531 13012.641602 Add

致谢

我们要感谢 Naren Castellon 为写这篇教程所做的贡献。

参考文献

  1. 黄昌权 • 阿拉·佩图基娜. 斯普林格系列 (2022). 使用Python进行应用时间序列分析与预测.
  2. 伊万·斯韦顿科夫. 增强动态自适应模型 (ADAM) 的预测与分析
  3. 詹姆斯·D·汉密尔顿. 时间序列分析 普林斯顿大学出版社, 新泽西州普林斯顿, 第1版, 1994.
  4. Nixtla 参数.
  5. Pandas 可用频率.
  6. 罗布·J·汉德曼和乔治·阿萨纳索普洛斯 (2018). “预测原则与实践,时间序列交叉验证”..
  7. 季节性周期-罗布·J·汉德曼.

Give us a ⭐ on Github