本补充材料是对数字信号处理(Digital Signal Processing)中傅里叶变换相关内容的初步介绍与探讨,以更好理解TimesNet中快速傅里叶变换应用的有效性与局限性(以及日后进一步研究),以及模型如何自动发现时间序列中的周期性模式。仅供参考。
1 傅里叶分析
傅里叶分析的核心思想是:任何周期函数都可以表示为不同频率的正弦和余弦函数的线性组合。对于非周期函数,可以通过傅里叶变换将其分解为连续频率的谐波分量。
1.1 傅里叶级数
对于周期为T
的周期函数f(t)
,其傅里叶级数(Fourier Series)展开为
f(t) = \frac{a_0}{2} + \sum\limits_{n=1}^{\infty} \left[a_n \cos\left(\frac{2\pi n t}{T}\right) + b_n \sin\left(\frac{2\pi n t}{T}\right)\right]
其中系数计算公式为
a_0 = \frac{2}{T} \int_{0}^{T} f(t) \text{d}t
a_n = \frac{2}{T} \int_{0}^{T} f(t) \cos\left(\frac{2\pi n t}{T}\right) \text{d}t
b_n = \frac{2}{T} \int_{0}^{T} f(t) \sin\left(\frac{2\pi n t}{T}\right) \text{d}t
1.2 复数形式的傅里叶级数
使用欧拉公式\text{e}^{\text{i}x} = \cos x + \text{i}\sin x
,可以得到更简洁的复数形式:
f(t) = \sum\limits_{n=-\infty}^{\infty} c_n \text{e}^{\text{i} \frac{2\pi n t}{T}}
其中系数
c_n = \frac{1}{T} \int_{0}^{T} f(t) \text{e}^{-\text{i} \frac{2\pi n t}{T}} \text{d}t
2 傅里叶变换
对于非周期函数,引入傅里叶变换(Fourier Transform)将时域(Time Domain)信号变换到频域(Frequency Domain)。傅里叶变换的理论极为复杂,此处仅介绍以下几种。
2.1 连续傅里叶变换
连续傅里叶变换(Continuous Fourier Transform,CFT)的定义为
F(\omega) = \int_{-\infty}^{\infty} f(t) \text{e}^{-\text{i}\omega t} \text{d}t
逆傅里叶变换则为
f(t) = \frac{1}{2\pi} \int_{-\infty}^{\infty} F(\omega) \text{e}^{\text{i}\omega t} \text{d}\omega
其中\omega = 2\pi f
是角频率,相关定义如下:
- 普通频率(Ordinary Frequency)
f
:单位时间内完成的周期数,单位是赫兹(\text{Hz}
)。例如f = 50\text{Hz}
表示每秒振动50次。 - 角频率(Angular Frequency)
\omega
:单位时间内转过的角度(以弧度计),单位是弧度/秒(\text{rad/s}
)
2.2 离散傅里叶变换
实际应用中通常处理的是离散时间信号。设x_0, x_1, \cdots, x_{N-1}
是N
个等间隔采样的数据点,离散傅里叶变换(Discrete Fourier Transform,DFT)定义为
X_k = \sum\limits_{n=0}^{N-1} x_n \text{e}^{-\text{i} \frac{2\pi}{N} k n},\ k = 0, 1, \cdots, N-1
逆离散傅里叶变换为
x_n = \frac{1}{N} \sum\limits_{k=0}^{N-1} X_k \text{e}^{\text{i} \frac{2\pi}{N} k n},\ n = 0, 1, \cdots, N-1
此处X_k
表示频率为\frac{k}{N}
倍采样频率的复振幅。
3 快速傅里叶变换
快速傅里叶变换(Fast Fourier Transform,FTT)是计算离散傅里叶变换的高效算法,由Cooley和Tukey于1965年提出。
3.1 基本原理
FTT的核心思想是分治策略,将DFT分解为更小的DFT。
当N = 2^m
(即N
是2的幂)时,可以将DFT分解为
X_k = \sum\limits_{n=0}^{N-1} x_n \text{e}^{-\text{i} \frac{2\pi}{N} k n} = \sum\limits_{r=0}^{\frac{N}2-1} x_{2r} \text{e}^{-\text{i} \frac{2\pi}{N} k (2r)} + \sum\limits_{r=0}^{\frac{N}2-1} x_{2r+1} \text{e}^{-\text{i} \frac{2\pi}{N} k (2r+1)}
令
E_k = \sum\limits_{r=0}^{\frac{N}2-1} x_{2r} \text{e}^{-\text{i} \frac{2\pi}{\frac{N}2} k r} \quad \text{(偶数点DFT)}
O_k = \sum\limits_{r=0}^{\frac{N}2-1} x_{2r+1} \text{e}^{-\text{i} \frac{2\pi}{\frac{N}2} k r} \quad \text{(奇数点DFT)}
则原DFT可表示为
X_k = E_k + \text{e}^{-\text{i} \frac{2\pi}{N} k} O_k
由于E_k
和O_k
的周期都是\frac{N}2
,我们还有:
X_{k+\frac{N}2} = E_k - \text{e}^{-\text{i} \frac{2\pi}{N} k} O_k
由此,一个N
点的DFT被分解为两个N/2
点的DFT。
3.2 算法复杂度
直接计算DFT:O(N^2)
次复数乘法和加法
FFT算法:O(N \log_2 N)
次复数乘法和加法
当N = 1024
时,FFT比直接DFT快约100倍;当N = 32768
时,快约2000倍。
3.3 在TimesNet中的应用
在TimesNet中,FFT用于提取时间序列的周期性,具体步骤对应论文中的公式
\mathbf{A}=\text{Avg}(\text{Amp}(\text{FFT}(\mathbf{X}_{\text{1D}})))
其中相关运算的详细解释如下:
\text{FFT}(\mathbf{X}_{\text{1D}})
:对长度为T
的时间序列进行FFT,得到频域表示\text{Amp}(\cdot)
:计算复数的振幅(模长),即|X_k| = \sqrt{\text{Re}(X_k)^2 + \text{Im}(X_k)^2}
,其中\text{Re}(X_k)
是实部(Real Part),\text{Im}(X_k)
是虚部(Imaginary Part)\text{Avg}(\cdot)
:在通道维度C
上取平均值
在此之后,选择前k
个最显著的频率
f_1,\cdots,f_k= \underset{f_*\in \{1,\cdots,[\frac{T}{2}]\}}{\arg\text{Topk}} (\mathbf{A})
对应的周期长度为
p_1,\cdots,p_k=\left \lceil \frac {T}{f_1} \right \rceil ,\cdots,\left \lceil \frac {T}{f_k} \right \rceil
3.4 频域的对称性
由于实值信号的傅里叶变换具有共轭对称性:
X_{N-k} = X_k^*
其中“^*
”表示复共轭。因此振幅谱是对称的:
|X_{N-k}| = |X_k|
这就是为何在TimesNet中只需要考虑频率范围\{1,\cdots,[\frac{T}{2}]\}
,避免了冗余计算。
4 实际计算考虑
4.1 振幅与功率谱
在实际应用中,除了振幅谱外,还常用功率谱:
P_k = |X_k|^2
功率谱更能突出主导频率成分。
4.2 窗函数效应
有限长度的采样相当于对无限长信号加矩形窗,这会引入频谱泄漏。为减少此效应,可在FFT前加窗函数(如汉宁窗、汉明窗)。
TimesNet中直接使用矩形窗,这对于周期性检测通常是足够的。
《从傅里叶变换到快速傅里叶变换》有2条评论