高级平稳性理论与检验方法

概述

平稳性（Stationarity）是时间序列分析中最重要的概念之一，它描述了时间序列统计性质的稳定性。许多经典的时间序列模型（如 ARMA）都假设序列是平稳的。本章将深入讨论平稳性的定义、检验方法以及非平稳序列的处理策略。¹

---

一、平稳性的定义

1.1 严平稳性（Strict Stationarity）

严平稳性要求序列的所有统计性质在时间平移下完全不变：

定义：对于任意整数 $k$ 和任意正整数 $m$，如果

$$F_{X_1,X_2,\dots ,X_m}(x_1,x_2,\dots ,x_m)=F_{X_{1+k},X_{2+k},\dots ,X_{m+k}}(x_1,x_2,\dots ,x_m)$$

其中 $F$ 是联合分布函数，则称时间序列 $\{X_t\}$ 是严平稳（Strictly Stationary）的。

含义：无论选择哪一段时间，序列的分布都是相同的。这是一种非常强的条件，实际中几乎无法验证。

1.2 弱平稳性（Weak/Covariance Stationarity）

弱平稳性只要求序列的前二阶矩（均值、方差、协方差）不随时间变化：

定义：如果时间序列 $\{X_t\}$ 满足：

1. 均值平稳：$\mathbb{E}[X_t]=\mu$，对所有 $t$ 成立
2. 方差有限：$\mathbb{E}[X_t^2]<\infty$
3. 协方差平稳：$\text{Cov}(X_t,X_{t+h})=\gamma (h)$，只依赖于滞后 $h$，不依赖于 $t$

则称 $\{X_t\}$ 是弱平稳（Weakly Stationary）或宽平稳（Wide-sense Stationary）的。

1.3 平稳性条件对比

特性	严平稳	弱平稳
条件强度	非常强	中等
分布要求	所有阶联合分布	仅前二阶矩
实际应用	理论分析	实际建模
高斯过程	等价	等价

1.4 自协方差函数

弱平稳序列的自协方差函数 $\gamma (h)$ 满足以下性质：

$$\gamma (0)=\text{Var}(X_t)\ge 0$$ $$\gamma (-h)=\gamma (h)\text{（对称性）}$$ $$|\gamma (h)|\le \gamma (0)\text{（非负定性）}$$

自相关函数（ACF）：

$$\rho (h)=\frac{\gamma (h)}{\gamma (0)}$$

import numpy as np
import matplotlib.pyplot as plt
 
def plot_acf_properties():
    """演示 ACF 的性质"""
    np.random.seed(42)
    n = 200
    
    # 生成平稳 AR(1) 序列
    phi = 0.7
    x = np.zeros(n)
    x[0] = np.random.randn()
    for t in range(1, n):
        x[t] = phi * x[t-1] + np.random.randn()
    
    # 计算 ACF
    def acf(x, max_lag=50):
        n = len(x)
        x_centered = x - np.mean(x)
        acf_values = np.correlate(x_centered, x_centered, mode='full')
        acf_values = acf_values[n-1:n+max_lag]
        return acf_values / acf_values[0]
    
    lags, acf_values = range(51), acf(x, 50)
    
    plt.figure(figsize=(12, 4))
    plt.stem(lags, acf_values, use_line_collection=True)
    plt.axhline(y=0, color='black', linestyle='-', linewidth=0.5)
    plt.axhline(y=1.96/np.sqrt(n), color='red', linestyle='--', label='95% CI')
    plt.axhline(y=-1.96/np.sqrt(n), color='red', linestyle='--')
    plt.xlabel('Lag')
    plt.ylabel('ACF')
    plt.title('Autocorrelation Function of Stationary AR(1)')
    plt.grid(True)
    plt.show()
 
plot_acf_properties()

---

二、常见非平稳类型

2.1 确定性趋势

$$X_t=\mu +\beta t+Y_t$$

其中 $Y_t$ 是平稳序列，$\beta \ne 0$ 表示存在确定性趋势。

示例：

def deterministic_trend():
    """确定性趋势序列"""
    t = np.arange(200)
    trend = 0.5 * t  # 线性趋势
    noise = np.random.randn(200) * 0.5
    x = trend + noise
    
    plt.figure(figsize=(12, 4))
    plt.plot(x)
    plt.plot(trend, 'r--', label='Trend')
    plt.title('Series with Deterministic Trend')
    plt.xlabel('Time')
    plt.ylabel('Value')
    plt.legend()
    plt.show()
 
deterministic_trend()

2.2 单位根过程

单位根过程是另一类常见的非平稳序列：

$$X_t=\rho X_{t-1}+{\epsilon }_t,|\rho |=1$$

当 $\rho =1$ 时，称为随机游走：

$$X_t=X_{t-1}+{\epsilon }_t=X_0+\sum _{i=1}^t{\epsilon }_i$$

随机游走的特性：

$$\text{Var}(X_t)=t\cdot {\sigma }^2$$

方差随时间增长，这解释了为什么随机游走难以预测。

def random_walk():
    """随机游走序列"""
    np.random.seed(42)
    n = 200
    eps = np.random.randn(n)
    x = np.zeros(n)
    x[0] = 0
    for t in range(1, n):
        x[t] = x[t-1] + eps[t]
    
    plt.figure(figsize=(12, 4))
    plt.plot(x)
    plt.title('Random Walk')
    plt.xlabel('Time')
    plt.ylabel('Value')
    plt.grid(True)
    plt.show()
    
    # 验证方差增长
    variances = [np.var(x[:i]) for i in range(10, n, 10)]
    times = list(range(10, n, 10))
    print(f"方差随时间增长: {variances}")
 
random_walk()

2.3 季节性非平稳

$$X_t=S_t+Y_t$$

其中 $S_t=S_{t+s}$ 是周期为 $s$ 的季节成分。

---

三、单位根检验

3.1 ADF 检验

Augmented Dickey-Fuller（ADF）检验是最常用的单位根检验：

原假设：序列存在单位根（非平稳）

备择假设：序列是平稳的

检验形式：

$$\Delta X_t=\alpha +\beta t+\gamma X_{t-1}+\sum _{i=1}^p{\delta }_i\Delta X_{t-i}+{\epsilon }_t$$

其中：

• $\alpha$：截距（漂移项）
• $\beta t$：时间趋势
• $\gamma$：单位根系数
• 滞后项 $\Delta X_{t-i}$ 用于消除自相关

from statsmodels.tsa.stattools import adfuller, kpss
 
def adf_test(x, regression='c', autolag='AIC'):
    """
    ADF 单位根检验
    
    Args:
        x: 时间序列
        regression: 检验类型
            'c' - 仅截距
            'ct' - 截距和时间趋势
            'n' - 无截距无趋势
        autolag: 自动选择滞后阶数的方法
    Returns:
        adf_stat, p_value, used_lag, critical_values
    """
    result = adfuller(x, regression=regression, autolag=autolag)
    
    print("=" * 50)
    print("Augmented Dickey-Fuller Test")
    print("=" * 50)
    print(f"Test Statistic: {result[0]:.4f}")
    print(f"P-Value: {result[1]:.4f}")
    print(f"Lags Used: {result[2]}")
    print(f"Observations Used: {result[3]}")
    print("\nCritical Values:")
    for key, value in result[4].items():
        print(f"  {key}: {value:.4f}")
    
    # 解读
    print("\n" + "=" * 50)
    if result[1] < 0.05:
        print("结论: 拒绝原假设，序列是平稳的")
    else:
        print("结论: 无法拒绝原假设，序列是非平稳的")
    print("=" * 50)
    
    return result
 
 
# 示例
np.random.seed(42)
 
# 平稳序列
stationary = np.random.randn(200)
 
# 非平稳序列（随机游走）
non_stationary = np.cumsum(np.random.randn(200))
 
print("平稳序列 ADF 检验:")
adf_test(stationary)
 
print("\n非平稳序列 ADF 检验:")
adf_test(non_stationary)

3.2 KPSS 检验

KPSS 检验与 ADF 互补：

原假设：序列是平稳的（或趋势平稳）

备择假设：序列存在单位根

def kpss_test(x, regression='c', nlags='auto'):
    """
    KPSS 检验
    
    Args:
        regression: 'c' - 水平平稳, 'ct' - 趋势平稳
    """
    result = kpss(x, regression=regression, nlags=nlags)
    
    print("=" * 50)
    print("KPSS Test")
    print("=" * 50)
    print(f"Test Statistic: {result[0]:.4f}")
    print(f"P-Value: {result[1]:.4f}")
    print(f"Lags Used: {result[2]}")
    print("\nCritical Values:")
    for key, value in result[3].items():
        print(f"  {key}: {value:.4f}")
    
    print("\n" + "=" * 50)
    if result[1] < 0.05:
        print("结论: 拒绝原假设，序列是非平稳的")
    else:
        print("结论: 无法拒绝原假设，序列是平稳的")
    print("=" * 50)
    
    return result
 
 
# 综合检验
print("平稳序列:")
adf_test(stationary)
kpss_test(stationary)
 
print("\n非平稳序列:")
adf_test(non_stationary)
kpss_test(non_stationary)

3.3 PP 检验

Phillips-Perron 检验是 ADF 的非参数版本：

from statsmodels.tsa.stattools import phillips_perron
 
def pp_test(x, regression='c'):
    """Phillips-Perron 检验"""
    result = phillips_perron(x, regression=regression)
    
    print("=" * 50)
    print("Phillips-Perron Test")
    print("=" * 50)
    print(f"Test Statistic: {result[0]:.4f}")
    print(f"P-Value: {result[1]:.4f}")
    print(f"Lags Used: {result[2]}")
    print("\nCritical Values:")
    for key, value in result[3].items():
        print(f"  {key}: {value:.4f}")
    
    return result

3.4 综合检验策略

def comprehensive_stationarity_test(x, significance_level=0.05):
    """
    综合平稳性检验
    
    结合 ADF、KPSS、PP 三种检验结果
    """
    results = {}
    
    # ADF 检验
    adf_result = adfuller(x)
    results['ADF'] = {
        'statistic': adf_result[0],
        'p_value': adf_result[1],
        'stationary': adf_result[1] < significance_level
    }
    
    # KPSS 检验
    kpss_result = kpss(x, regression='c', nlags='auto')
    results['KPSS'] = {
        'statistic': kpss_result[0],
        'p_value': kpss_result[1],
        'stationary': kpss_result[1] >= significance_level  # KPSS原假设相反
    }
    
    # PP 检验
    pp_result = phillips_perron(x)
    results['PP'] = {
        'statistic': pp_result[0],
        'p_value': pp_result[1],
        'stationary': pp_result[1] < significance_level
    }
    
    # 综合判断
    n_stationary = sum(r['stationary'] for r in results.values())
    
    print("\n" + "=" * 50)
    print("综合检验结果")
    print("=" * 50)
    for test_name, result in results.items():
        status = "平稳" if result['stationary'] else "非平稳"
        print(f"{test_name}: {status} (p={result['p_value']:.4f})")
    
    print(f"\n投票结果: {n_stationary}/3 支持平稳")
    print("=" * 50)
    
    return results
 
# 综合检验示例
comprehensive_stationarity_test(stationary)
comprehensive_stationarity_test(non_stationary)

---

四、非平稳序列的处理

4.1 差分

一阶差分：

$$\Delta X_t=X_t-X_{t-1}$$

ADF 检验的差分形式：

如果原序列是单位根过程，一阶差分后应该平稳：

def difference_test():
    """检验差分后序列的平稳性"""
    np.random.seed(42)
    
    # 原始非平稳序列
    y = np.cumsum(np.random.randn(200))
    
    print("原始序列:")
    adf_test(y)
    
    # 一阶差分
    dy = np.diff(y)
    
    print("\n一阶差分后:")
    adf_test(dy)
    
    # 二阶差分（如需要）
    ddy = np.diff(dy)
    
    print("\n二阶差分后:")
    adf_test(ddy)
 
difference_test()

4.2 去趋势

from scipy import signal
 
def detrend_series(x, method='linear'):
    """
    去趋势
    
    Args:
        x: 时间序列
        method: 'linear' 或 'polynomial'
    """
    t = np.arange(len(x))
    
    if method == 'linear':
        # 线性去趋势
        coeffs = np.polyfit(t, x, 1)
        trend = np.polyval(coeffs, t)
    else:
        # 多项式去趋势
        coeffs = np.polyfit(t, x, 2)
        trend = np.polyval(coeffs, t)
    
    detrended = x - trend
    
    return detrended, trend
 
 
def high_pass_filter(x, cutoff=0.1, fs=1.0):
    """
    高通滤波去趋势
    """
    b, a = signal.butter(3, cutoff, btype='high', fs=fs)
    filtered = signal.filtfilt(b, a, x)
    return filtered
 
 
# 示例
t = np.arange(200)
trend = 0.5 * t + 0.01 * t**2
x = trend + np.random.randn(200) * 0.5
 
detrended, estimated_trend = detrend_series(x, method='linear')
 
plt.figure(figsize=(12, 6))
plt.subplot(2, 1, 1)
plt.plot(x, label='Original')
plt.plot(estimated_trend, 'r--', label='Estimated Trend')
plt.title('Original Series with Trend')
plt.legend()
 
plt.subplot(2, 1, 2)
plt.plot(detrended, label='Detrended')
plt.title('Detrended Series')
plt.legend()
 
plt.tight_layout()
plt.show()

4.3 季节差分

对于有季节性的序列：

$${\Delta }_s{X}_t=X_t-X_{t-s}$$

def seasonal_difference(x, s):
    """季节差分"""
    return x[s:] - x[:-s]
 
 
# 示例
np.random.seed(42)
t = np.arange(200)
seasonal = 10 * np.sin(2 * np.pi * t / 52)  # 年度季节
trend = 0.5 * t
noise = np.random.randn(200) * 0.5
x = trend + seasonal + noise
 
# 季节差分（假设年度周期）
x_diff = seasonal_difference(x, 52)
 
print("原始序列:")
adf_test(x)
 
print("\n季节差分后:")
adf_test(x_diff)

---

五、协整理论

5.1 协整的定义

Engle-Granger 协整：如果两个或多个非平稳序列的线性组合是平稳的，则称这些序列是协整的。

定义：如果序列 $\{X_t\}$ 和 $\{Y_t\}$ 都是 $I(1)$（一阶单整）的，但存在某个线性组合 $Z_t=X_t-\beta Y_t$ 是 $I(0)$（平稳的），则称 $\{X_t\}$ 和 $\{Y_t\}$ 是协整的，记为 $\text{CI}(1,1)$。

5.2 协整检验

from statsmodels.tsa.stattools import coint
 
def cointegration_test(x, y):
    """
    Engle-Granger 协整检验
    
    原假设：序列之间不存在协整关系
    """
    result = coint(x, y)
    
    print("=" * 50)
    print("Engle-Granger Cointegration Test")
    print("=" * 50)
    print(f"Test Statistic: {result[0]:.4f}")
    print(f"P-Value: {result[1]:.4f}")
    print("\nCritical Values:")
    for key, value in result[2].items():
        print(f"  {key}: {value:.4f}")
    
    print("\n" + "=" * 50)
    if result[1] < 0.05:
        print("结论: 拒绝原假设，序列之间存在协整关系")
    else:
        print("结论: 无法拒绝原假设，序列之间不存在协整关系")
    print("=" * 50)
    
    return result
 
 
# 示例：检验两个相关序列
np.random.seed(42)
x = np.cumsum(np.random.randn(200))
y = 0.8 * x + np.random.randn(200) * 0.5  # 与x相关
 
cointegration_test(x, y)

5.3 误差修正模型（ECM）

协整序列可以用 ECM 建模：

$$\Delta X_t=\alpha (X_{t-1}-\beta Y_{t-1})+\sum _{i=1}^p{\gamma }_i\Delta X_{t-i}+{\epsilon }_t$$

其中 $(X_{t-1}-\beta Y_{t-1})$ 是误差修正项。

from statsmodels.tsa.api import VARDynamicSUR
 
def error_correction_model(x, y, maxlag=5):
    """
    简单的误差修正模型估计
    
    使用OLS分两步估计：
    1. 估计协整回归
    2. 估计 ECM
    """
    import statsmodels.api as sm
    
    # 第一步：协整回归
    X_with_const = sm.add_constant(y)
    model_coin = sm.OLS(x, X_with_const).fit()
    beta = model_coin.params[1]
    residuals = model_coin.resid
    
    # 第二步：ECM
    dx = np.diff(x)
    dy = np.diff(y)
    ec_term = residuals[:-1]  # 滞后误差修正项
    
    # 构建 ECM 数据
    X_ecm = np.column_stack([
        ec_term,
        dx[:-1],
        dy[:-1]
    ])
    X_ecm = sm.add_constant(X_ecm)
    
    model_ecm = sm.OLS(dx[1:], X_ecm).fit()
    
    print("ECM 回归结果:")
    print(model_ecm.summary())
    
    return model_ecm, beta
 
 
# 示例
np.random.seed(42)
x = np.cumsum(np.random.randn(200))
y = 0.8 * x + np.random.randn(200) * 0.3
 
model, beta = error_correction_model(x, y)

---

六、参考

---

相关阅读

• 时间序列分解理论
• 时序频域分析
• 循环神经网络

Footnotes

1. Hamilton, J. D. (1994). Time Series Analysis. Princeton University Press. ↩