Python Day56

Task：
1.假设检验基础知识
a.原假设与备择假设
b.P值、统计量、显著水平、置信区间
2.白噪声
a.白噪声的定义
b.自相关性检验：ACF检验和Ljung-Box 检验
c.偏自相关性检验：PACF检验
3.平稳性
a.平稳性的定义
b.单位根检验
4.季节性检验
a.ACF检验
b.序列分解：趋势+季节性+残差
记忆口诀：p越小，落在置信区间外，越拒绝原假设。

1. 假设检验基础知识

a. 原假设与备择假设

• 原假设 (Null Hypothesis, $H_0$): 通常是我们想要否定或反驳的陈述。在统计学中，它常常代表“没有效应”、“没有差异”或“没有关系”。
• 备择假设 (Alternative Hypothesis, $H_1$ or $H_a$): 与原假设相反的陈述。它代表我们希望找到证据支持的陈述，即“存在效应”、“存在差异”或“存在关系”。

示例：

• 原假设 ($H_0$): 该药物对降低血压没有影响。
• 备择假设 ($H_1$): 该药物可以降低血压。

b. P值、统计量、显著水平、置信区间

• 统计量 (Statistic): 从样本数据计算出来的量，用于估计总体的参数。它是检验统计量（test statistic）的基础，用于衡量样本数据与原假设之间的偏差。
• P值 (P-value): 在原假设为真的前提下，观察到当前样本数据或更极端数据的概率。
  • 如果 P值较小 (通常 < 显著水平): 说明当前样本数据不太可能在原假设为真的情况下出现，因此我们有理由拒绝原假设，接受备择假设。
  • 如果 P值较大 (通常 >= 显著水平): 说明当前样本数据在原假设为真的情况下出现的可能性较高，我们没有足够的证据拒绝原假设。
• 显著水平 (Significance Level, $\alpha$): 在做出拒绝原假设的决定之前，我们愿意承担的犯第一类错误（即错误地拒绝一个真实的 $H_0$）的最大概率。常见的显著水平有 0.05 (5%)，0.01 (1%) 和 0.10 (10%)。
• 置信区间 (Confidence Interval): 针对总体参数的一个区间估计。它表示我们有一定信心（例如 95%）认为真实的总体参数落在这个区间内。
  • 置信区间与 P 值和显著水平的关系: 如果显著水平为 $\alpha$，那么 $(1-\alpha )$ 的置信区间包含了我们想要检验的某个特定值（例如，在检验均值是否等于某个值时，检查该值是否在置信区间内），那么我们就不能拒绝原假设。反之，如果置信区间不包含该特定值，我们就可以拒绝原假设。

Python 代码示例 (使用 SciPy 进行 t 检验)：

假设我们想检验一组样本数据的均值是否显著不等于 0。

import numpy as np
from scipy import stats

# 假设的样本数据
sample_data = np.array([0.5, 0.2, 0.8, -0.1, 0.6, 0.3, 0.9, 0.0, 0.7, 0.4])

# 原假设 H0: 样本均值等于 0
# 备择假设 H1: 样本均值不等于 0

# 使用独立样本 t 检验 (one-sample t-test)
# 检验样本均值是否与指定的均值 (这里是 0) 有显著差异
t_statistic, p_value = stats.ttest_1samp(sample_data, 0)

# 设置显著水平
alpha = 0.05

print(f"T-statistic: {t_statistic:.4f}")
print(f"P-value: {p_value:.4f}")
print(f"Significance level (alpha): {alpha}")

# 决策规则
if p_value < alpha:
    print("由于 P 值 < 显著水平，拒绝原假设 H0。")
    print("样本均值与 0 存在显著差异。")
else:
    print("由于 P 值 >= 显著水平，未能拒绝原假设 H0。")
    print("没有足够的证据表明样本均值与 0 存在显著差异。")

# 计算 95% 置信区间
# 注意：ttest_1samp 函数直接返回 p_value，但为了演示置信区间，我们可以手动计算或使用 stats.t.interval
# mean = np.mean(sample_data)
# sem = stats.sem(sample_data) # Standard error of the mean
# df = len(sample_data) - 1
# confidence_interval = stats.t.interval(0.95, df, loc=mean, scale=sem)

# 更直接的方式是利用 ttest_1samp 的结果和均值来构建置信区间
mean_sample = np.mean(sample_data)
std_err_sample = stats.sem(sample_data)
# 对于双侧检验，置信区间是 mean ± t_critical * std_err
# t_critical 可以通过 ppf (percent point function) 找到
t_critical = stats.t.ppf(1 - alpha/2, len(sample_data) - 1)
confidence_interval = (mean_sample - t_critical * std_err_sample, mean_sample + t_critical * std_err_sample)

print(f"95% Confidence Interval for the mean: ({confidence_interval[0]:.4f}, {confidence_interval[1]:.4f})")

# 通过置信区间判断
if 0 < confidence_interval[0] or 0 > confidence_interval[1]:
    print("置信区间不包含 0，这与 P 值 < alpha 的结论一致。")
else:
    print("置信区间包含 0，这与 P 值 >= alpha 的结论一致。")

2. 白噪声

a. 白噪声的定义

白噪声是指所有频率分量在幅度上相等，且相位是随机的信号。在时间序列分析中，一个纯白噪声序列具有以下特征：

1. 均值为零: 序列的期望值为 0。
2. 方差恒定: 序列的方差不随时间变化。
3. 无自相关性: 序列中任何两个不同时间点的观测值之间没有任何线性关系。即在任何滞后阶数 $k>0$ 下，自协方差都为零，自相关系数也为零。

满足以上三个条件的序列是纯白噪声。如果一个序列仅满足后两个条件（方差恒定且无自相关性），则称为宽平稳白噪声。在实际应用中，我们通常关注的是宽平稳白噪声。

Python 代码示例 (生成和可视化白噪声)：

import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns

# 设置随机种子以确保结果可复现
np.random.seed(42)

# 生成一个纯白噪声序列
# 均值为 0，标准差为 1 的正态分布随机数是典型的白噪声实现方式
white_noise = np.random.normal(loc=0, scale=1, size=200)

# 打印序列的前几个值
print("White Noise Sequence (first 10 values):")
print(white_noise[:10])

# 可视化白噪声序列
plt.figure(figsize=(12, 6))
plt.plot(white_noise, label='White Noise')
plt.title('White Noise Sequence')
plt.xlabel('Time')
plt.ylabel('Value')
plt.legend()
plt.grid(True)
plt.show()

b. 自相关性检验：ACF 检验和 Ljung-Box 检验

自相关函数 (Autocorrelation Function, ACF) 衡量了时间序列在不同滞后阶数下的相关性。对于白噪声序列，其 ACF 应该在滞后 0 处等于 1 (自身相关)，在所有其他滞后阶数 $k>0$ 处都接近于零，并且统计上不显著。

• ACF 检验: 通过绘制 ACF 图，观察在各个滞后阶数下，ACF 值是否在置信区间（通常是围绕零的虚线）内。如果大部分滞后阶数的 ACF 值都在置信区间内，则认为序列可能没有显著的自相关性。
• Ljung-Box 检验: 这是一个统计检验，用于检验序列是否存在整体的自相关性。
  • 原假设 ($H_0$): 序列在指定的滞后阶数 $K$ 之前的自相关性都不显著。
  • 备择假设 ($H_1$): 序列在指定的滞后阶数 $K$ 之前存在显著的自相关性。
  • Ljung-Box 统计量是基于前 $K$ 个自相关系数计算得出的。P 值表示在原假设为真时，观察到当前统计量或更极端值的概率。如果 P 值小于显著水平，则拒绝原假设。

Python 代码示例 (ACF 检验和 Ljung-Box 检验)：

import numpy as np
import matplotlib.pyplot as plt
from statsmodels.graphics.tsaplots import plot_acf
from statsmodels.stats.diagnostic import acorr_ljungbox
import seaborn as sns

# 使用前面生成的白噪声序列
# white_noise = np.random.normal(loc=0, scale=1, size=200) # 如果未运行前一段代码

# 1. ACF 检验 (可视化)
plt.figure(figsize=(12, 6))
plot_acf(white_noise, lags=20, title='ACF of White Noise')
plt.show()

# 2. Ljung-Box 检验
# 设定检验的滞后阶数 K
# 通常选择一个与时间序列长度相关的阶数，例如 20 或 40
k_lags = 20

ljung_box_result = acorr_ljungbox(white_noise, lags=[k_lags], return_df=True)

print(f"\nLjung-Box Test for White Noise (up to lag {k_lags}):")
print(ljung_box_result)

# 对 Ljung-Box 检验结果进行解释
alpha = 0.05
p_value_ljungbox = ljung_box_result.iloc[0]['lb_pvalue'] # 获取指定滞后阶数的 p 值

if p_value_ljungbox < alpha:
    print(f"P值 ({p_value_ljungbox:.4f}) < 显著水平 ({alpha})，拒绝原假设。")
    print("序列存在显著的自相关性。")
else:
    print(f"P值 ({p_value_ljungbox:.4f}) >= 显著水平 ({alpha})，未能拒绝原假设。")
    print("序列不存在显著的自相关性 (符合白噪声特性)。")

# 再次验证 ACF 图，看滞后阶数 20 时的 ACF 值是否在置信区间内
# 如果 ACF 图显示在 20 滞后阶数时 ACF 值在置信区间内，并且 Ljung-Box 检验不拒绝原假设，
# 则进一步支持其为白噪声的结论。

c. 偏自相关性检验：PACF 检验

偏自相关函数 (Partial Autocorrelation Function, PACF) 衡量了在剔除了中间滞后观测值的影响后，当前时间点观测值与某个滞后时间点观测值之间的直接相关性。

对于白噪声序列，PACF 的行为与 ACF 类似：在滞后 0 处等于 1，在所有其他滞后阶数 $k>0$ 处都接近于零，并且统计上不显著。

Python 代码示例 (PACF 检验)：

import numpy as np
import matplotlib.pyplot as plt
from statsmodels.graphics.tsaplots import plot_pacf
import seaborn as sns

# 使用前面生成的白噪声序列
# white_noise = np.random.normal(loc=0, scale=1, size=200) # 如果未运行前一段代码

# PACF 检验 (可视化)
plt.figure(figsize=(12, 6))
plot_pacf(white_noise, lags=20, title='PACF of White Noise')
plt.show()

# PACF 图的解释：
# 同样，观察 PACF 值是否在置信区间（通常是围绕零的虚线）内。
# 对于白噪声，所有的 PACF 值（除滞后0外）都应该在置信区间内。

3. 平稳性

a. 平稳性的定义

时间序列的平稳性是时间序列分析中的一个核心概念。一个平稳的时间序列的统计性质不随时间变化。具体来说，主要有以下几种类型的平稳性：

1. 严平稳 (Strict Stationarity): 如果一个序列的联合概率分布不随时间移动而改变，则该序列是严平稳的。这意味着对于任意 $t_1,\dots ,t_k$ 和任意 $h$，以及任意函数 $f$，都有：
   $P(X_{t_1}\le x_1,\dots ,X_{t_k}\le x_k)=P(X_{t_1+h}\le x_1,\dots ,X_{t_k+h}\le x_k)$
   严平稳是比较强的条件，在实际中较难检验。

2. 弱平稳 (Weak Stationarity) 或协方差平稳 (Covariance Stationarity):

   • 均值恒定: $E[X_t]=\mu$ 对所有 $t$ 都成立，即均值不随时间变化。
   • 方差恒定: $Var(X_t)=E[(X_t-\mu {)}^2]={\sigma }^2$ 对所有 $t$ 都成立，即方差不随时间变化。
   • 滞后协方差仅依赖于滞后阶数: 对于任意 $t$ 和滞后阶数 $k$，协方差 $Cov(X_t,X_{t+k})=E[(X_t-\mu )(X_{t+k}-\mu )]$ 仅依赖于 $k$，而不依赖于 $t$。

在大多数时间序列分析中，我们关注的是弱平稳，因为它是模型（如 ARIMA 模型）能够有效工作的基本假设。如果一个序列不是弱平稳的，通常可以通过差分等方法将其转化为平稳序列。

b. 单位根检验

单位根 (Unit Root) 是指时间序列模型（如 AR(1) 模型）的特征方程的根等于 1。如果一个时间序列存在单位根，那么它通常不是平稳的，而是非平稳的，并且表现出随机游走 (random walk) 的特征，即序列会持续受到“惊喜”的影响。

单位根检验是用来检验时间序列是否包含单位根（即是否非平稳）的统计方法。常见的单位根检验方法包括：

1. 增广迪基-福勒检验 (Augmented Dickey-Fuller, ADF):

   • 原假设 ($H_0$): 序列包含单位根（即非平稳）。
   • 备择假设 ($H_1$): 序列不包含单位根（即平稳）。
   • ADF 检验通过对序列进行回归分析来估计单位根的存在，其统计量是基于一个特定的检验统计量。P 值用于判断是否拒绝原假设。

2. 菲利普斯-佩龙检验 (Phillips-Perron, PP):

   • 与 ADF 检验类似，但它通过非参数的方法来处理序列的异方差性和自相关性，对模型的误差项形式不作要求。
   • 原假设 ($H_0$): 序列包含单位根（即非平稳）。
   • 备择假设 ($H_1$): 序列不包含单位根（即平稳）。

Python 代码示例 (ADF 检验)：

import numpy as np
import matplotlib.pyplot as plt
from statsmodels.tsa.stattools import adfuller
import pandas as pd
import seaborn as sns

# 1. 生成一个非平稳序列 (例如，随机游走)
np.random.seed(42)
random_walk = np.cumsum(np.random.normal(loc=0, scale=1, size=200))

# 打印序列的前几个值
print("Random Walk Sequence (first 10 values):")
print(random_walk[:10])

# 可视化随机游走序列
plt.figure(figsize=(12, 6))
plt.plot(random_walk, label='Random Walk (Non-stationary)')
plt.title('Random Walk Sequence')
plt.xlabel('Time')
plt.ylabel('Value')
plt.legend()
plt.grid(True)
plt.show()

# 2. 进行 ADF 检验
# adfuller 函数返回 ADF 统计量，P 值，滞后阶数等
# 默认情况下，ADF 检验的原假设是 "序列存在单位根" (非平稳)
# 备择假设是 "序列不存在单位根" (平稳)

# 可以指定 maxlag 来确定用于 ADF 回归的滞后阶数
# 可以指定 autolag 来自动选择最大滞后阶数
# 可以指定 regression ('c'为常数项, 'ct'为常数项+趋势项, 'ctt'为常数项+线性趋势项+二次趋势项, 'nc'为无常数项和趋势项)
# 对于大多数时间序列，我们通常会包含常数项 'c'。如果序列明显有趋势，可以考虑 'ct'。
adf_result = adfuller(random_walk, regression='c')

print("\nAugmented Dickey-Fuller (ADF) Test for Random Walk:")
print(f'ADF Statistic: {adf_result[0]:.4f}')
print(f'P-value: {adf_result[1]:.4f}')
print('Critical Values:')
for key, value in adf_result[4].items():
    print(f'  {key}: {value:.4f}')

# 解释 ADF 检验结果
alpha = 0.05
p_value_adf = adf_result[1]

if p_value_adf < alpha:
    print(f"\nP值 ({p_value_adf:.4f}) < 显著水平 ({alpha})，拒绝原假设。")
    print("序列不存在单位根，是平稳的。")
else:
    print(f"\nP值 ({p_value_adf:.4f}) >= 显著水平 ({alpha})，未能拒绝原假设。")
    print("序列存在单位根，是非平稳的。")

# 可以与平稳序列进行对比
# 生成一个平稳序列 (白噪声)
# white_noise = np.random.normal(loc=0, scale=1, size=200) # 如果未运行前一段代码
adf_result_wn = adfuller(white_noise, regression='c')

print("\nAugmented Dickey-Fuller (ADF) Test for White Noise:")
print(f'ADF Statistic: {adf_result_wn[0]:.4f}')
print(f'P-value: {adf_result_wn[1]:.4f}')

if adf_result_wn[1] < alpha:
    print(f"P值 ({adf_result_wn[1]:.4f}) < 显著水平 ({alpha})，拒绝原假设。")
    print("白噪声序列不存在单位根，是平稳的。")
else:
    print(f"P值 ({adf_result_wn[1]:.4f}) >= 显著水平 ({alpha})，未能拒绝原假设。")
    print("白噪声序列存在单位根，是非平稳的。")

4. 季节性检验

季节性是指时间序列在固定的周期内重复出现的模式，例如每年的销售额在特定月份出现高峰。

a. ACF 检验 (用于季节性)

ACF 图是识别季节性的一个重要工具。对于一个具有明显季节性的序列，其 ACF 图会在季节性周期（滞后阶数）的倍数处显示出显著的、且周期性重复的峰值。

例如，如果一个序列具有月度季节性（周期为 12），那么在滞后 12, 24, 36 等阶数处，ACF 值会非常高且显著。

Python 代码示例 (ACF 检验用于季节性)：

假设我们有一个具有年季节性的数据（例如月度数据，周期为 12）。

import numpy as np
import matplotlib.pyplot as plt
from statsmodels.graphics.tsaplots import plot_acf
import pandas as pd
import seaborn as sns

# 生成一个带有季节性的合成序列
np.random.seed(42)
n_samples = 120 # 10 年的月度数据

# 基础趋势 (可选)
trend = np.linspace(0, 10, n_samples)

# 季节性成分 (周期为 12)
seasonality_amplitude = 5
seasonal_component = seasonality_amplitude * np.sin(np.linspace(0, 3 * np.pi, n_samples) * 2 * np.pi / 12) # 制造更平滑的季节性

# 随机噪声
noise = np.random.normal(loc=0, scale=2, size=n_samples)

# 合成序列
seasonal_series = trend + seasonal_component + noise

# 可视化合成序列
plt.figure(figsize=(12, 6))
plt.plot(seasonal_series, label='Seasonal Series')
plt.title('Seasonal Time Series Data')
plt.xlabel('Time (Months)')
plt.ylabel('Value')
plt.legend()
plt.grid(True)
plt.show()

# 绘制 ACF 图以识别季节性
# 对于月度数据 (周期 12)，我们通常会绘制到滞后 24 或 36
plt.figure(figsize=(12, 6))
plot_acf(seasonal_series, lags=36, title='ACF of Seasonal Series')
plt.show()

# 解释 ACF 图:
# 注意在滞后 12, 24, 36 等阶数处是否出现显著的峰值。
# 如果出现，表明存在季节性。

b. 序列分解：趋势 + 季节性 + 残差

序列分解是一种将时间序列分解为三个主要成分的方法：

• 趋势 (Trend): 序列的长期方向性运动。
• 季节性 (Seasonal): 周期性的重复模式。
• 残差 (Residual): 除去趋势和季节性后剩余的部分，通常认为是随机噪声。

序列分解有助于理解时间序列的组成部分，并可以用于移除季节性（例如，进行季节性调整）。分解模型可以是加法模型或乘法模型：

• 加法模型: $Y_t=T_t+S_t+R_t$ (当季节性幅度随时间变化不大时)
• 乘法模型: $Y_t=T_t\times S_t\times R_t$ (当季节性幅度随时间呈比例变化时)

Python 代码示例 (使用 Statsmodels 进行序列分解)：

import numpy as np
import matplotlib.pyplot as plt
from statsmodels.tsa.seasonal import seasonal_decompose
import pandas as pd
import seaborn as sns

# 使用前面生成的带有季节性的合成序列
# seasonal_series = ... # 如果未运行前一段代码

# 确定季节性周期
# 对于月度数据，周期是 12
seasonal_period = 12

# 进行序列分解 (这里使用加法模型，如果数据变化不明显，也可以尝试乘法模型)
# 需要 pandas Series 作为输入
decomposition = seasonal_decompose(pd.Series(seasonal_series), model='additive', period=seasonal_period)

# 提取各成分
trend_component = decomposition.trend
seasonal_component_decomp = decomposition.seasonal
residual_component = decomposition.resid

# 绘制分解结果
fig, (ax1, ax2, ax3, ax4) = plt.subplots(4, 1, figsize=(12, 10), sharex=True)

ax1.plot(seasonal_series, label='Original')
ax1.set_ylabel('Original')
ax1.legend()
ax1.grid(True)

ax2.plot(trend_component, label='Trend')
ax2.set_ylabel('Trend')
ax2.legend()
ax2.grid(True)

ax3.plot(seasonal_component_decomp, label='Seasonal')
ax3.set_ylabel('Seasonal')
ax3.legend()
ax3.grid(True)

ax4.plot(residual_component, label='Residual')
ax4.set_ylabel('Residual')
ax4.legend()
ax4.grid(True)

plt.xlabel('Time (Months)')
plt.tight_layout()
plt.show()

# 解释分解结果：
# 趋势图显示了数据的长期方向。
# 季节性图显示了周期性的模式，应该与我们之前合成的 sin wave 形状相似。
# 残差图应该看起来像随机噪声，并且我们也可以对其进行 ACF 和 Ljung-Box 检验来确认其是否为白噪声。

# 对残差进行白噪声检验 (可选)
print("\nTesting Residuals for White Noise:")
# ACF
plt.figure(figsize=(10, 4))
plot_acf(residual_component.dropna(), lags=20, title='ACF of Residuals') # dropna() 移除分解过程中产生的 NaN 值
plt.show()

# Ljung-Box Test
k_lags_resid = 20
ljung_box_result_resid = acorr_ljungbox(residual_component.dropna(), lags=[k_lags_resid], return_df=True)
print(f"Ljung-Box Test for Residuals (up to lag {k_lags_resid}):")
print(ljung_box_result_resid)

# 解释残差的白噪声检验结果：
# 如果 ACF 图显示残差在置信区间内，且 Ljung-Box 检验 P 值大于显著水平，
# 则说明分解是有效的，残差可以被认为是白噪声。