一、前言

时间序列的平稳性本质上是指序列的统计特性在时间平移方面的不变性。时间序列研究的目的是为了通过对历史数据进行分析，预测未来的趋势。如果时间序列不是平稳的，那么由历史数据得到的统计特性对未来的预测就毫无意义。因此，确保时间序列平稳非常重要。

二、平稳性检验

• ADF：ADF检验（Augmented Dickey-Fuller Testing）是最常用的单位根检验方法之一。该方法假设时间序列服从有单位根的p阶自回归过程，即AR(p)过程，通过检验序列是否存在单位根来判断序列是否是平稳的。需要注意的是，该方法针对小样本的检验结果可能不够准确。
• PP：1988 年Phillips和Perron提出了一种非参数检验方法，它假设了序列的随机干扰服从无穷阶的MA过程。这种方法主要是为了解决残差项中潜在的序列相关和异方差问题。然而，这种方法在小样本情况下的检验结果与ADF类似，也不够准确。
• DF-GLS：DF-GLS检验是Eilliot（1996）针对具有确定趋势的AR模型提出的单位根检验方法，该方法使用广义最小二乘法检验数据进行一次“准差分”，然后利用准差分的数据对原序列进行退势处理，再利用ADF检验的模型形式对退势数据进行单位根检验。
• KPSS：KPSS检验是Kwiatkowski, Phillips, and Shin在1992年提出的趋势平稳检验方法，与之前介绍的三种检验方法都不同。KPSS检验的原假设是时间序列平稳或趋势平稳，备择假设是存在单位根。KPSS检验采用非参数修正方式解决趋势平稳零假设下的序列相关问题，其过程类似于PP检验。相对于以上三种检验方法，KPSS检验更加适合检验趋势平稳性。
• Zivot-Andrews：也叫Z-A检验，Zivot和Andrews于1992年提出了一种新型自回归平稳性检验，该检验通过计算基于移动窗口的ADF和PP检验得到的两个自相关函数之间的差异来评估序列平稳性。与其他针对自回归分析的平稳性检验不同，Zivot-Andrews检验考虑了季节性趋势以及观测值之间的相关程度，这是该方法的独特优势。
• Variance Ratio：即V-R检验，VR检验是一种基于Zivot-Andrews检验的时间序列平稳性检验方法。与其他检验方法相比，它更加有效且简单，使用了自回归残差和趋势分量来进行分析。检验的显著性水平为α=0.05，首先通过分解将原始序列分成自回归残差和线性趋势序列，然后使用一个统计量比较自回归方差比率和趋势方差比率，从而得出平稳性检验结果。
• RUR：RUR是2006年提出的一种时间序列检验方法，它从序列的运行范围进行构造，并推导出其极限分布。该检验方法具有很多理想特性，如对序列的单调变换不变性以及对重要参数变化的鲁棒性等，属于非参数统计方法。

这些检验方法对应的Python计算函数如下表所示：

以常用的ADF检验举例说明：

import numpy as np

# 生成随机数据
np.random.seed(123)
data = np.random.standard_normal(size=100)

# 使用arch的函数
from arch.unitroot import ADF
adf = ADF(data)
print("------Results of arch package Test------")
print(adf.summary().as_text())

# 使用statsmodels的函数 
from statsmodels.tsa.stattools import adfuller
def adf_test(timeseries):
    print("------Results of statsmodels package Test------")
    dftest = adfuller(timeseries, autolag="AIC")
    dfoutput = pd.Series(
        dftest[0:4],
        index=[
            "Test Statistic",
            "p-value",
            "#Lags Used",
            "Number of Observations Used",
        ],
    )
    for key, value in dftest[4].items():
        dfoutput["Critical Value (%s)" % key] = value
    print(dfoutput)

adf_test(data)

计算的结果如下，可以看到arch包会给出一个判断结果，statsmodels的函数则需要根据计算的结果如p-value < 0.05进行判断。

------Results of arch package Test------
   Augmented Dickey-Fuller Results   
=====================================
Test Statistic                 -9.822
P-value                         0.000
Lags                                0
-------------------------------------

Trend: Constant
Critical Values: -3.50 (1%), -2.89 (5%), -2.58 (10%)
Null Hypothesis: The process contains a unit root.
Alternative Hypothesis: The process is weakly stationary.
------Results of statsmodels package Test------
Test Statistic                -9.82
p-value                        0.00
#Lags Used                     0.00
Number of Observations Used   99.00
Critical Value (1%)           -3.50
Critical Value (5%)           -2.89
Critical Value (10%)          -2.58
dtype: float64

三、 平稳化处理

差分

• 一阶差分

通过去除时间序列中的一些变化特征来平稳化它的均值。这样做不仅消除了时间序列的趋势和季节性，同时还可能减小噪声对数据的影响。

• 二阶差分

有时候即使进行了一阶差分处理后的数据仍然不平稳，这时则需要再一次对数据进行差分来得到一个更加平稳的序列。

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

# 非平稳序列
np.random.seed(123)
data = np.cumsum(np.random.randint(-100, 100, size=100)) + 200

# 差分变换
df = pd.DataFrame(data)
df_diff1 = df.diff(1) # 1阶差分，步长为1
df_diff2 = df_diff1.diff(1) # 2阶差分，步长为1
axs.plot(data, label='data')
axs.plot(df_diff1, label='first-order')
axs.plot(df_diff2, label='second-order')
plt.legend()
plt.show()

平滑

平滑用于剔除时间序列中的随机波动，从而得出序列的基本轨迹，常用的平滑方法有：

• 简单移动平均

简单移动平均法（Simple Moving Average，SMA）是最基本的一种时间序列平滑处理方法，其原理是对于给定的时间窗口内的数据进行平均，将得到一个平滑后的数据点。

但是，简单移动平均法认为每个时期的数据同等重要，没有考虑数据的权重。而在实际应用中，往往需要根据数据的不同属性和意义来赋予不同的权重。因此，加权移动平均法就应运而生了。

• 加权移动平均

对序列的每个数据分别乘以一个权重系数，对近期数据给予较大的权重，对远期的数据给予较小的权重。它的平滑公式为：

用加权移动平均法求预测值，对近期的趋势反应较敏感，但如果一组数据有明显的季节性影响时，用加权移动平均法所得到的预测值可能会出现偏差。

• 一次指数平滑

一次指数平滑法是根据过去一段时间内的历史数据，对当前数据进行加权平均计算。它主要考虑到该序列的趋势发展，每个时刻的值是由历史时刻的观测值加权得出的，各观测值的权重随着时间的推移以幂函数的形式递减。该方法不仅弥补了简单平均法的不能体现各时期重要性的缺点，又弥补了加权平均法只能关注最近时期的缺点。它的平滑公式为：

平滑系数取值范围0~1，其值越接近于1，远期实际值对本期平滑值影响程度下降得越快；越接近于0，远期实际值对本期平滑值的影响程度下降得越慢。所以，当时间数列相对平稳时，可取较大值；当时间数列波动较大时，应取较小的值。

• 二次指数平滑

当时间序列的变动出现直线趋势时，用一次指数平滑法进行预测会存在明显的滞后偏差，此时可以对其进行二次指数平滑进行改进。平滑公式为：

预测第T期的公式为：

平滑方法的使用示例如下所示：

import numpy as np
import pandas as pd
from statsmodels.tsa.api import ExponentialSmoothing
import matplotlib.pyplot as plt

# 非平稳序列
np.random.seed(123)
data = np.cumsum(np.random.randint(-100, 100, size=100)) + 200
df = pd.DataFrame(data)

# 平滑处理
w = 5
SMA = df.rolling(window=w).mean()
WMA = df.rolling(window=w).apply(lambda x: x[::-1].cumsum().sum() * 2 / w / (w + 1))
EWA = df.ewm(halflife=w/2).mean()
EWA2 = ExponentialSmoothing(df, trend='additive', seasonal='additive', seasonal_periods=w).fit().fittedvalues
fig, ax = plt.subplots()
ax.plot(data, label='data')
ax.plot(SMA, label='SMA')
ax.plot(WMA, label='WMA')
ax.plot(EWA, label='EWA')
ax.plot(EWA2, label='EWA2')
plt.legend()
plt.show()

变换

• 对数变换

对数变换可以将指数增长的时间序列转化为线性增长，同时通过减小波动性去除了数据方差随时间变化的问题。在进行对数变换后，还可以结合差分法、平滑法或分解法进行进一步的处理，以获得更加平稳的时间序列。

• Box-Cov变换

Box-cox变换是一个包含对数变换和幂变换的参数化变换族，它旨在将数据分布映射到接近高斯分布，以实现方差稳定和最小化偏斜。需要注意的是，Box-cox变换只适用于严格正数据。其定义如下：

• Yeo-Johonson变换

Yeo-Johonson变换是一种广义幂变换方法，是Box-Cox变换在实数域的推广，可以明显改善数据的正态性和减小异方差性，其定义如下：

Python实现变换的示例代码如下：

import numpy as np
from scipy.stats import yeojohnson, boxcox

# 非平稳序列
np.random.seed(123)
data = np.cumsum(np.random.randint(1, 100, size=100))

# 变换
y_log = np.log(data)
y_bc, lam_bc = boxcox(data)
y_yj, lam_yj = yeojohnson(data)

# 数据分布对比
fig, ax = plt.subplots(1, 4)
ax[0].hist(data)
ax[0].set_title('data')
ax[1].hist(y_log)
ax[1].set_title('log')
ax[2].hist(y_bc)
ax[2].set_title('boxcox')
ax[3].hist(y_yj)
ax[3].set_title('yoejohnson')
plt.show()

分解

对于非季节性时间序列数据的平稳化，可以使用前面提到的几种方法，而对于季节性时间序列数据的平稳化，则需要使用季节性的时间序列分解法。具体来说，分解法将时间序列分解成三个部分：趋势部分、季节性部分和随机波动部分。其中，趋势部分表示时间序列的长期变化趋势，季节性部分表示时间序列中的周期性变化，而随机波动部分则代表噪声和非系统性因素。

常用的时间序列分解方法有STL分解，它的核心在于使用自适应的局部回归技术（LOESS）来估计非线性的季节和趋势，随后对原始序列与估计出的季节和趋势作差得到残差项。这个过程不断迭代直到残差项变得平稳。

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from statsmodels.tsa.seasonal import STL

# 生成季节性的非平稳数据
np.random.seed(2021)
trend = 3 * np.sin(np.linspace(0, 2*np.pi, 12*10))
seasonal = 5 * np.cos(np.linspace(0, 4*np.pi, 12*10))
residual = np.random.randn(12*10)
y = trend + seasonal + residual
data = pd.Series(y)

# STL分解
decom = STL(data, period=12).fit()
trend = decom.trend
seasonal = decom.seasonal
residual = decom.resid

fig, ax = plt.subplots()
ax.plot(data, label='data')
ax.plot(trend, label='Trend')
ax.plot(seasonal, label='Seasonal')
ax.plot(residual, label='Residual')
plt.legend()
plt.show()

四、总结

本文介绍了时间序列数据的平稳性检验的方法，总结了常用的平稳化处理手段并给出相应的Python计算实例，希望对读者有所帮助。