前言

本文隶属于专栏《机器学习数学通关指南》，该专栏为笔者原创，引用请注明来源，不足和错误之处请在评论区帮忙指出，谢谢！

本专栏目录结构和参考文献请见《机器学习数学通关指南》

---

ima 知识库

知识库广场搜索：

知识库	创建人
机器学习	@Shockang
机器学习数学基础	@Shockang
深度学习	@Shockang

---

正文

自相关函数与功率谱密度：时域与频域的桥梁

📚 谱分析方法是研究信号频率成分的重要手段，是连接时域特性和频域分析的桥梁，在机器学习和数据科学中有着广泛应用。

🔄 一、自相关函数（Autocorrelation Function, ACF）

1.1 数学定义 📐

对于平稳随机过程 $X(t)$，自相关函数定义为：

$R_X(\tau )=E[X(t)X(t+\tau )]$

物理意义：描述信号在时间间隔 $\tau$ 上的自相似性。当 $\tau =0$ 时，$R_X(0)$ 表示信号的均方值（信号总功率）。

1.2 关键性质 🧩

性质	数学表达式	含义
对称性	$R_X(-\tau )=R_X(\tau )$	自相关函数是偶函数
非负定性	相关矩阵为非负定矩阵	保证功率谱密度非负
最大值特性	$|R_X(\tau )|\le R_X(0)$	相关性随时间差增大而衰减

1.3 直观理解 🧠

自相关函数度量了信号与其时移版本的相似度：

• 高自相关：数据点间存在强依赖关系
• 低自相关：数据点近似独立
• 周期性自相关：信号中存在周期性成分

1.4 示例 🌊

对随机相位正弦波 $X(t)=a\cos (\omega t+\Theta )$，其中 $\Theta \sim U(0,2\pi )$：

$R_X(\tau )=\frac{a^2}{2}\cos (\omega \tau )$

其周期性反映信号具有单一频率分量 $\omega$。

📊 二、功率谱密度（Power Spectral Density, PSD）

2.1 数学定义 📐

对平稳过程，功率谱密度是自相关函数的傅里叶变换：

$S_X(\omega )={\int }_{-\infty }^{\infty }{R}_X(\tau )e^{-j\omega \tau }d\tau$

物理意义：表征信号功率在各频率分量上的分布密度。积分全频带得总功率：

${\int }_{-\infty }^{\infty }{S}_X(\omega )d\omega =R_X(0)$

2.2 重要性质 📋

• 非负性：$S_X(\omega )\ge 0$
• 偶对称性：$S_X(-\omega )=S_X(\omega )$
• 维纳-辛钦定理：自相关函数与功率谱密度构成傅里叶变换对

2.3 常见信号的谱密度 📈

信号类型	自相关函数	功率谱密度	特点
白噪声	$R_X(\tau )={\sigma }^2\delta (\tau )$	$S_X(\omega )={\sigma }^2$	功率均匀分布在所有频率上
正弦波	$R_X(\tau )=\frac{a^2}{2}\cos ({\omega }_0\tau )$	$S_X(\omega )=\pi a^2[\delta (\omega -{\omega }_0)+\delta (\omega +{\omega }_0)]$	在特定频率处有脉冲

🔄 三、自相关函数与功率谱密度的关系

3.1 维纳-辛钦定理（Wiener-Khinchin Theorem）🔗

维纳-辛钦定理建立了自相关函数与功率谱密度的互换关系：

S X ( ω ) = F { R X ( τ ) } = ∫ − ∞ ∞ R X ( τ ) e − j ω τ d τ S_X(\omega) = \mathcal{F}\{R_X(\tau)\} = \int_{-\infty}^{\infty} R_X(\tau) e^{-j\omega \tau} d\tau SX​(ω)=F{RX​(τ)}=∫−∞∞​RX​(τ)e−jωτdτ

R X ( τ ) = F − 1 { S X ( ω ) } = 1 2 π ∫ − ∞ ∞ S X ( ω ) e j ω τ d ω R_X(\tau) = \mathcal{F}^{-1}\{S_X(\omega)\} = \frac{1}{2\pi}\int_{-\infty}^{\infty} S_X(\omega) e^{j\omega \tau} d\omega RX​(τ)=F−1{SX​(ω)}=2π1​∫−∞∞​SX​(ω)ejωτdω

3.2 分析意义 💡

• 时域难以观察的周期性和噪声特性，在频域可通过谱峰/平坦度直观识别
• 功率谱密度将复杂信号分解成不同频率成分，便于特征提取
• 在机器学习中，这种转换可以揭示数据中的隐藏模式和周期性

🔄 四、各态历经性（Ergodicity）

4.1 定义与意义 📝

各态历经性：若平稳过程满足各态历经性，则可用单一样本的时间平均替代统计平均。

${\lim }_{T\to \infty }\frac{1}{2T}{\int }_{-T}^{T}x(t)x(t+\tau )dt=R_X(\tau )$

4.2 实际应用价值 🛠️

• 允许从有限长度的单条时间序列估计自相关函数
• 为现实数据分析提供理论支撑
• 机器学习中，意味着可以从单一数据序列中推断整个随机过程的统计特性

💻 五、Python实现与应用

5.1 使用Python计算自相关函数 🐍

import numpy as np
import matplotlib.pyplot as plt
from statsmodels.tsa.stattools import acf# 生成示例数据：AR(1)过程
np.random.seed(42)
n = 1000
phi = 0.8
ar1 = np.zeros(n)
for i in range(1, n):ar1[i] = phi * ar1[i-1] + np.random.normal(0, 1)# 计算自相关函数
lags = 40
acf_values = acf(ar1, nlags=lags)# 可视化
plt.figure(figsize=(10, 6))
plt.stem(range(lags + 1), acf_values)
plt.axhline(y=0, linestyle='--', color='gray')
plt.axhline(y=1.96/np.sqrt(n), linestyle='--', color='red')
plt.axhline(y=-1.96/np.sqrt(n), linestyle='--', color='red')
plt.title('AR(1)过程的自相关函数', fontsize=14)
plt.xlabel('滞后(Lag)', fontsize=12)
plt.ylabel('自相关系数', fontsize=12)
plt.tight_layout()
plt.show()

5.2 功率谱密度计算 📊

from scipy import signal# 生成带噪声的正弦信号
fs = 100  # 采样频率
t = np.arange(0, 10, 1/fs)  # 10秒，采样率100Hz
f1, f2 = 5, 15  # 信号频率
x = np.sin(2*np.pi*f1*t) + 0.5*np.sin(2*np.pi*f2*t) + 0.5*np.random.randn(len(t))# 计算PSD
f, Pxx = signal.welch(x, fs, nperseg=256)# 可视化
plt.figure(figsize=(10, 6))
plt.semilogy(f, Pxx)
plt.title('功率谱密度', fontsize=14)
plt.xlabel('频率 [Hz]', fontsize=12)
plt.ylabel('PSD [V^2/Hz]', fontsize=12)
plt.grid(True)
plt.tight_layout()
plt.show()

🚀 六、机器学习中的应用场景

6.1 信号检测与特征提取 🔍

在接收信号 $Y(t)=aS(t-\tau )+N(t)$ 中，若有用信号 $S(t)$ 与噪声 $N(t)$ 独立，其互相关函数:

$R_{SY}(\tau )=a{R}_S(\tau )$

应用：

• 雷达和声纳系统中的目标检测
• 生物医学信号处理中的异常检测
• IoT传感器数据中的模式识别

6.2 时间序列分析与预测 📈

ARMA模型：自回归滑动平均模型的谱密度具有有理函数形式，便于参数化建模。

例如AR(2)模型：$X_t=0.4X_{t-1}+0.4X_{t-2}+{ϵ}_t$

# AR(2)模型参数估计与频谱分析
from statsmodels.tsa.arima.model import ARIMA
import pandas as pd# 生成AR(2)过程
n = 500
ar2 = np.zeros(n)
ar2[0:2] = np.random.randn(2)
for i in range(2, n):ar2[i] = 0.4*ar2[i-1] + 0.4*ar2[i-2] + np.random.normal(0, 1)# 拟合AR模型
model = ARIMA(ar2, order=(2,0,0))
results = model.fit()
print("AR(2)模型参数估计结果:")
print(results.summary())# 计算理论谱
ar_params = np.r_[1, -results.params[1:]]
w, h = signal.freqz(1, ar_params, worN=1000)
psd = np.abs(h)**2# 可视化理论谱
plt.figure(figsize=(10, 6))
plt.plot(w/np.pi, psd)
plt.title('AR(2)模型的理论功率谱', fontsize=14)
plt.xlabel('归一化频率 (×π rad/sample)', fontsize=12)
plt.ylabel('功率谱密度', fontsize=12)
plt.grid(True)
plt.tight_layout()
plt.show()

6.3 深度学习中的应用 🧠

• 特征工程：将时域数据变换到频域作为神经网络输入特征
• 时频分析：通过短时傅里叶变换(STFT)获取时变频谱图作为2D卷积网络输入
• 异常检测：监测功率谱密度变化识别异常行为

🔬 七、高级谱分析技术

7.1 多元谱分析 🔄

交叉功率谱：研究两个信号之间的频率关系
S X Y ( ω ) = F { R X Y ( τ ) } S_{XY}(\omega) = \mathcal{F}\{R_{XY}(\tau)\} SXY​(ω)=F{RXY​(τ)}

相干函数：衡量两信号在特定频率上的线性相关性
${\gamma }_{XY}^2(\omega )=\frac{|S_{XY}(\omega ){|}^2}{S_X(\omega )S_Y(\omega )}$

7.2 现代谱估计方法 📊

方法	优点	缺点	适用场景
周期图法	计算简单	方差大	长数据序列
Welch方法	降低方差	频率分辨率降低	一般应用
AR模型法	频率分辨率高	需模型阶数选择	短数据序列
多锥方法	方差小、分辨率高	计算复杂	高精度要求

# 比较不同谱估计方法
methods = ['periodogram', 'welch', 'lombscargle']
plt.figure(figsize=(12, 8))# 生成含噪声的多频率信号
np.random.seed(0)
t = np.linspace(0, 10, 1000)
x = np.sin(2*np.pi*5*t) + 0.5*np.sin(2*np.pi*10*t) + 0.3*np.random.randn(len(t))# 周期图法
f1, pxx1 = signal.periodogram(x, fs=100)
plt.subplot(3, 1, 1)
plt.semilogy(f1, pxx1)
plt.title('周期图法', fontsize=12)
plt.ylabel('PSD')# Welch方法
f2, pxx2 = signal.welch(x, fs=100, nperseg=256)
plt.subplot(3, 1, 2)
plt.semilogy(f2, pxx2)
plt.title('Welch方法', fontsize=12)
plt.ylabel('PSD')# Lomb-Scargle方法 (适用于非均匀采样)
f3 = np.linspace(0.1, 50, 500)
pxx3 = signal.lombscargle(t, x, f3)
plt.subplot(3, 1, 3)
plt.semilogy(f3, pxx3)
plt.title('Lomb-Scargle方法', fontsize=12)
plt.xlabel('频率 [Hz]')
plt.ylabel('功率')plt.tight_layout()
plt.show()

🎯 八、总结与展望

8.1 谱分析方法的核心价值 💎

• 时频域转换：通过傅里叶变换在时域和频域间建立桥梁
• 特征提取：揭示数据中的周期性和频率特征
• 降维与去噪：通过频域筛选保留关键信息，去除噪声

8.2 在机器学习中的重要性 🤖

• 为特征工程提供理论基础
• 作为信号预处理的关键步骤
• 在时间序列预测中提供频域视角

8.3 未来发展方向 🚀

• 时频联合分析：小波变换、希尔伯特-黄变换等高级技术
• 深度时频学习：将谱分析与深度学习紧密结合
• 非平稳信号分析：拓展谱分析方法应对非平稳特性

---

🔗 思考与练习

1. 如何从实际数据中估计自相关函数的可靠性？
2. 尝试使用本文提供的Python代码分析一个实际时间序列数据。
3. 研究非平稳信号的谱分析方法，如何处理金融数据等非平稳序列？