Python实现时间序列数据的特征工程：利用小波变换与经验模态分解（EMD）

大家好，今天我们来聊聊时间序列数据特征工程中两个非常有用的技术：小波变换（Wavelet Transform）和经验模态分解（Empirical Mode Decomposition, EMD）。时间序列数据分析和预测在各个领域都至关重要，而有效的特征工程能够显著提升模型的性能。小波变换和EMD都是信号分解技术，它们可以将原始时间序列分解成多个具有不同频率和时间尺度的子序列，从而提取出隐藏在数据中的有用信息。

一、时间序列特征工程的重要性

在深入探讨小波变换和EMD之前，我们先来简单回顾一下时间序列特征工程的重要性。机器学习模型，特别是深度学习模型，通常需要丰富的特征才能达到良好的预测效果。对于时间序列数据，直接使用原始数据往往不够，需要进行适当的特征提取。常见的特征包括：

统计特征: 均值、方差、标准差、最大值、最小值、中位数、偏度、峰度等。
时域特征: 自相关系数（ACF）、偏自相关系数（PACF）。
频域特征: 功率谱密度（PSD）、频谱峰值。
时频域特征: 短时傅里叶变换（STFT）、小波变换。

今天我们重点介绍后两类中的小波变换以及EMD。

二、小波变换（Wavelet Transform）

小波变换是一种时频分析方法，它将信号分解成一系列小波函数，这些小波函数在时间和频率上都具有局部性。与傅里叶变换相比，小波变换能够更好地捕捉信号的非平稳特性。

2.1 基本原理

小波变换的核心思想是使用一组称为小波的基函数，通过伸缩和平移来分析信号。小波函数 $psi(t)$ 必须满足以下条件：

有限能量: $int_{-infty}^{infty} |psi(t)|^2 dt < infty$
零均值: $int_{-infty}^{infty} psi(t) dt = 0$

常见的小波函数包括Daubechies小波、Haar小波、Morlet小波等。连续小波变换（CWT）定义如下：

$W(a, b) = frac{1}{sqrt{a}} int_{-infty}^{infty} x(t) psi^*(frac{t-b}{a}) dt$

其中：

$x(t)$ 是输入信号。
$psi(t)$ 是母小波函数。
$a$ 是尺度参数（scale），控制小波的伸缩。
$b$ 是平移参数（translation），控制小波的平移。
$psi^*(t)$ 是 $psi(t)$ 的复共轭。

离散小波变换（DWT）是CWT的离散化版本，它使用离散的尺度和平移参数。DWT通常采用Mallat算法实现，该算法通过滤波器组将信号分解成近似系数（Approximation coefficients）和细节系数（Detail coefficients）。

2.2 Python实现

Python中可以使用pywt库来进行小波变换。首先，安装pywt库：

pip install pywavelets

下面是一个使用Daubechies 4小波对时间序列数据进行两层分解的例子：

import pywt
import numpy as np
import matplotlib.pyplot as plt

# 生成示例时间序列数据
np.random.seed(42)
time = np.arange(0, 10, 0.01)
signal = np.sin(2 * np.pi * time) + 0.5 * np.cos(6 * np.pi * time) + 0.1 * np.random.randn(len(time))

# 小波分解
wavelet = 'db4'  # Daubechies 4小波
level = 2       # 分解层数
coeffs = pywt.wavedec(signal, wavelet, level=level)

# 提取近似系数和细节系数
cA2, cD2, cD1 = coeffs  # 近似系数cA2，细节系数cD2和cD1

# 可视化结果
plt.figure(figsize=(12, 8))

plt.subplot(4, 1, 1)
plt.plot(time, signal)
plt.title('Original Signal')

plt.subplot(4, 1, 2)
plt.plot(cA2)
plt.title('Approximation Coefficients (cA2)')

plt.subplot(4, 1, 3)
plt.plot(cD2)
plt.title('Detail Coefficients (cD2)')

plt.subplot(4, 1, 4)
plt.plot(cD1)
plt.title('Detail Coefficients (cD1)')

plt.tight_layout()
plt.show()

# 特征提取：计算各层系数的统计特征
features = []
features.extend([np.mean(cA2), np.std(cA2), np.max(cA2), np.min(cA2)]) # 近似系数的统计特征
features.extend([np.mean(cD2), np.std(cD2), np.max(cD2), np.min(cD2)]) # 细节系数的统计特征
features.extend([np.mean(cD1), np.std(cD1), np.max(cD1), np.min(cD1)]) # 细节系数的统计特征

print("Extracted Features:", features)

在这个例子中，我们首先生成一个包含正弦波和噪声的示例时间序列数据。然后，使用pywt.wavedec函数对信号进行小波分解，得到近似系数和细节系数。最后，我们可视化了原始信号以及各层分解的系数，并计算了各层系数的均值、标准差、最大值和最小值作为特征。

2.3 小波变换特征工程的应用

小波变换在时间序列特征工程中有很多应用，例如：

去噪: 细节系数通常包含噪声信息，可以通过阈值处理去除噪声。
特征提取: 近似系数可以捕捉信号的整体趋势，细节系数可以捕捉信号的局部波动。
异常检测: 分析细节系数的异常波动可以检测时间序列中的异常点。
信号压缩: 通过保留重要的系数，可以实现信号的压缩。

三、经验模态分解（Empirical Mode Decomposition, EMD）

经验模态分解（EMD）是一种自适应的信号分解方法，它可以将非线性、非平稳信号分解成一系列本征模态函数（Intrinsic Mode Functions, IMFs）。

3.1 基本原理

EMD的核心思想是迭代地从信号中提取出具有不同频率尺度的IMF。 IMF必须满足以下两个条件：

在整个数据长度上，极值点的数量和过零点的数量必须相等或最多相差一个。
在任意时刻，由局部极大值点定义的上包络线和由局部极小值点定义的下包络线的均值必须为零。

EMD的分解步骤如下：

找到信号的所有局部极大值点和局部极小值点。
通过三次样条插值分别拟合上包络线和下包络线。
计算上下包络线的均值 $m(t)$。
从原始信号中减去均值包络线，得到残差 $h(t) = x(t) – m(t)$。
判断残差 $h(t)$ 是否满足IMF的条件。如果满足，则将 $h(t)$ 作为一个IMF，并从原始信号中减去它。如果不满足，则将 $h(t)$ 作为新的信号，重复步骤1-4。
重复上述步骤，直到残差信号满足停止条件（例如，残差信号的能量低于某个阈值，或者残差信号是单调函数）。

3.2 Python实现

Python中可以使用PyEMD库来进行EMD分解。首先，安装PyEMD库：

pip install EMD-signal

下面是一个使用EMD对时间序列数据进行分解的例子：

from PyEMD import EMD
import numpy as np
import matplotlib.pyplot as plt

# 生成示例时间序列数据
np.random.seed(42)
time = np.arange(0, 10, 0.01)
signal = np.sin(2 * np.pi * time) + 0.5 * np.cos(6 * np.pi * time) + 0.1 * np.random.randn(len(time))

# EMD分解
emd = EMD()
imfs = emd(signal)

# 可视化结果
plt.figure(figsize=(12, 8))

plt.subplot(len(imfs) + 1, 1, 1)
plt.plot(time, signal)
plt.title('Original Signal')

for i, imf in enumerate(imfs):
    plt.subplot(len(imfs) + 1, 1, i + 2)
    plt.plot(time, imf)
    plt.title(f'IMF {i+1}')

plt.tight_layout()
plt.show()

# 特征提取：计算各IMF的统计特征
features = []
for i, imf in enumerate(imfs):
    features.extend([np.mean(imf), np.std(imf), np.max(imf), np.min(imf)])
print("Extracted Features:", features)

在这个例子中，我们首先生成一个包含正弦波和噪声的示例时间序列数据。然后，使用PyEMD.EMD()类对信号进行EMD分解，得到一系列IMF。最后，我们可视化了原始信号以及分解得到的IMF，并计算了各IMF的均值、标准差、最大值和最小值作为特征。

3.3 EMD特征工程的应用

EMD在时间序列特征工程中有很多应用，例如：

去噪: 高频IMF通常包含噪声信息，可以通过去除高频IMF来去噪。
特征提取: 不同的IMF反映了信号在不同频率尺度上的波动，可以提取各IMF的统计特征、能量特征等作为特征。
趋势分析: 低频IMF可以反映信号的长期趋势。
故障诊断: 分析IMF的异常波动可以诊断设备的故障。

四、小波变换 vs EMD

小波变换和EMD都是强大的信号分解工具，但它们也存在一些差异：

特性	小波变换	EMD
基函数	预定义的，例如Daubechies、Haar、Morlet等	自适应的，根据信号本身确定
分解方式	基于尺度和平移	基于局部均值包络线
适用性	适用于分析平稳信号和非平稳信号，但需要选择合适的小波基函数	适用于分析非线性、非平稳信号
计算复杂度	较低	较高
参数选择	需要选择小波基函数和分解层数	无需选择基函数，但可能需要设置停止条件

总的来说，小波变换是一种线性变换，它使用预定义的基函数来分解信号。EMD是一种非线性、自适应的分解方法，它根据信号本身来确定基函数。在实际应用中，可以根据信号的特性选择合适的方法。

五、实际案例：股票价格预测

让我们来看一个实际案例：使用小波变换和EMD对股票价格进行特征工程，并用于股票价格预测。

import yfinance as yf
import pywt
from PyEMD import EMD
import numpy as np
import pandas as pd
from sklearn.model_selection import train_test_split
from sklearn.ensemble import RandomForestRegressor
from sklearn.metrics import mean_squared_error

# 1. 数据获取
ticker = "AAPL"  # 苹果公司股票代码
data = yf.download(ticker, start="2020-01-01", end="2023-01-01")
close_prices = data['Close'].values

# 2. 特征工程 (小波变换和EMD)
def wavelet_features(data, wavelet='db4', level=3):
    coeffs = pywt.wavedec(data, wavelet, level=level)
    features = []
    for coeff in coeffs:
        features.extend([np.mean(coeff), np.std(coeff), np.max(coeff), np.min(coeff)])
    return features

def emd_features(data):
    emd = EMD()
    imfs = emd(data)
    features = []
    for imf in imfs:
        features.extend([np.mean(imf), np.std(imf), np.max(imf), np.min(imf)])
    return features

# 3. 创建数据集
sequence_length = 10  # 使用过去10天的数据预测未来一天
X, y = [], []
for i in range(len(close_prices) - sequence_length):
    # 原始数据特征
    seq_x = close_prices[i:i + sequence_length]
    # 小波变换特征
    wavelet_x = wavelet_features(seq_x)
    # EMD特征
    emd_x = emd_features(seq_x)
    # 将所有特征合并
    combined_x = np.concatenate([seq_x, wavelet_x, emd_x])
    X.append(combined_x)
    y.append(close_prices[i + sequence_length])

X = np.array(X)
y = np.array(y)

# 4. 划分训练集和测试集
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# 5. 模型训练 (随机森林)
model = RandomForestRegressor(n_estimators=100, random_state=42)  # 使用100棵树的随机森林
model.fit(X_train, y_train)

# 6. 模型评估
y_pred = model.predict(X_test)
rmse = np.sqrt(mean_squared_error(y_test, y_pred))
print(f"Root Mean Squared Error: {rmse}")

# 7. 可视化预测结果
plt.figure(figsize=(12, 6))
plt.plot(y_test, label="Actual Prices")
plt.plot(y_pred, label="Predicted Prices")
plt.xlabel("Time")
plt.ylabel("Stock Price")
plt.title("Stock Price Prediction using Wavelet and EMD Features")
plt.legend()
plt.show()

这个案例展示了如何使用yfinance库获取股票数据，然后使用小波变换和EMD提取特征，最后使用随机森林模型进行预测。请注意，这只是一个简单的示例，实际的股票价格预测需要更复杂的模型和特征工程。

六、总结与展望

今天我们学习了时间序列特征工程中两个重要的技术：小波变换和经验模态分解。小波变换通过预定义的基函数将信号分解成不同尺度上的系数，而EMD则是一种自适应的分解方法，可以根据信号本身来确定基函数。这两种方法都可以有效地提取时间序列数据的特征，并应用于各种实际问题中，例如去噪、特征提取、异常检测和预测等。

在未来，我们可以进一步研究以下几个方向：

结合深度学习: 将小波变换和EMD与深度学习模型结合，例如使用卷积神经网络（CNN）或循环神经网络（RNN）来学习分解后的系数的特征。
自适应参数选择: 研究如何自动选择小波基函数和EMD的停止条件，以提高特征提取的效率和准确性。
多模态数据融合: 将时间序列数据与其他类型的数据（例如文本数据、图像数据）融合，以提高预测的准确性。

七、核心思想回顾

小波变换和EMD为时间序列特征工程提供了强有力的工具，通过分解信号，提取多尺度特征，有效提升模型的预测能力。实践中需要根据数据特点选择合适的技术，并结合具体应用场景进行优化。

更多IT精英技术系列讲座，到智猿学院

Python实现时间序列数据的特征工程：利用小波变换与经验模态分解（EMD）

发表回复 取消回复

发表回复取消回复