为了更详细地说明如何利用傅里叶变换来解耦时序数据,我们可以通过一个具体的例子来解释整个过程,包括如何提取趋势、周期性成分和噪声。
- 了解数据
假设我们有一组时序数据,表示某个地区的每日温度变化。这些数据包含以下成分:
长期趋势(如全球变暖导致的逐年温度升高)
周期性成分(如一年四季的温度变化)
随机噪声(如每日温度的随机波动)
- 数据预处理
在进行傅里叶变换之前,需要对数据进行一些预处理工作。例如,去除数据中的长期趋势成分。可以使用移动平均或多项式拟合等方法来去除趋势。
移动平均:通过计算一定窗口内的数据平均值来平滑数据,去除短期波动,保留长期趋势。
多项式拟合:使用多项式函数拟合数据,提取出长期趋势成分。
- 应用傅里叶变换
对去除趋势后的数据进行傅里叶变换,将时间域的数据转换到频域。傅里叶变换会将数据分解成不同频率的正弦波和余弦波的叠加,每个频率成分都有一个对应的振幅和相位。
频域表示:傅里叶变换结果显示在频域中的频率成分。通过观察频域图,可以识别出主要的周期性成分。高振幅的频率成分对应着数据中的显著周期性。
- 提取频率成分
在频域图中,选择那些对应显著周期性的频率成分。例如,如果我们看到一年中的温度变化主要集中在1年(频率为1)和半年(频率为2)的周期上,我们就可以选择这些频率成分。
过滤频率:将主要的频率成分保留下来,滤除其他不重要的频率成分。
- 逆傅里叶变换
对保留的主要频率成分进行逆傅里叶变换,将频域的数据转换回时域。这样我们可以得到时间序列数据中包含的周期性成分。
逆变换结果:逆傅里叶变换后的结果就是我们提取出的周期性成分。与原始数据相比,这部分数据去除了噪声,保留了周期性变化。
- 噪声成分
最后,噪声成分可以通过从原始数据中减去提取出的趋势和周期性成分来获得。这样,我们将数据解耦成了三个部分:趋势、周期性成分和噪声。
利用傅里叶变换解耦时序数据的详细过程
我们通过一个具体的计算例子来说明如何利用傅里叶变换对时序数据进行解耦,提取出趋势、周期性成分和噪声。
示例数据
假设我们有以下10天的每日温度数据:
日 | 温度(℃) |
---|---|
1 | 15 |
2 | 16 |
3 | 14 |
4 | 17 |
5 | 18 |
6 | 19 |
7 | 20 |
8 | 21 |
9 | 22 |
10 | 23 |
1. 去除趋势
我们首先通过计算移动平均去除趋势。假设我们用5天的窗口进行移动平均:
对于每个点计算移动平均,例如:
移动平均 = 15 + 16 + 14 + 17 + 18 5 = 16 \text{移动平均} = \frac{15+16+14+17+18}{5} = 16 移动平均=515+16+14+17+18=16
以此类推,计算每个点的移动平均:
日 | 温度(℃) | 移动平均(℃) | 去趋势后的温度(℃) |
---|---|---|---|
1 | 15 | - | - |
2 | 16 | - | - |
3 | 14 | - | - |
4 | 17 | - | - |
5 | 18 | 16 | 18 - 16 = 2 |
6 | 19 | 16.8 | 19 - 16.8 = 2.2 |
7 | 20 | 18 | 20 - 18 = 2 |
8 | 21 | 19 | 21 - 19 = 2 |
9 | 22 | 20 | 22 - 20 = 2 |
10 | 23 | 21 | 23 - 21 = 2 |
2. 傅里叶变换
对去趋势后的温度数据进行傅里叶变换。傅里叶变换的公式为:
F ( k ) = ∑ n = 0 N − 1 x ( n ) e − i 2 π k n / N F(k) = \sum_{n=0}^{N-1} x(n) e^{-i 2 \pi k n / N} F(k)=∑n=0N−1x(n)e−i2πkn/N
其中, x ( n ) x(n) x(n)是去趋势后的温度数据, N N N是数据点的数量(这里为10), k k k是频率指数。
计算过程
我们以频率 k = 0 k = 0 k=0为例进行详细计算:
F ( 0 ) = ∑ n = 0 9 x ( n ) ⋅ e − i 2 π ⋅ 0 ⋅ n / 10 = ∑ n = 0 9 x ( n ) F(0) = \sum_{n=0}^{9} x(n) \cdot e^{-i 2 \pi \cdot 0 \cdot n / 10} = \sum_{n=0}^{9} x(n) F(0)=∑n=09x(n)⋅e−i2π⋅0⋅n/10=∑n=09x(n)
因为 e − i 2 π ⋅ 0 ⋅ n / 10 = 1 e^{-i 2 \pi \cdot 0 \cdot n / 10} = 1 e−i2π⋅0⋅n/10=1,所以:
F ( 0 ) = 2 + 2.2 + 2 + 2 + 2 + 2 = 12.2 F(0) = 2 + 2.2 + 2 + 2 + 2 + 2 = 12.2 F(0)=2+2.2+2+2+2+2=12.2
对于频率 k = 1 k = 1 k=1:
F ( 1 ) = ∑ n = 0 9 x ( n ) ⋅ e − i 2 π ⋅ 1 ⋅ n / 10 F(1) = \sum_{n=0}^{9} x(n) \cdot e^{-i 2 \pi \cdot 1 \cdot n / 10} F(1)=∑n=09x(n)⋅e−i2π⋅1⋅n/10
具体计算如下:
F ( 1 ) = 2 ⋅ e − i 2 π ⋅ 1 ⋅ 0 / 10 + 2.2 ⋅ e − i 2 π ⋅ 1 ⋅ 1 / 10 + 2 ⋅ e − i 2 π ⋅ 1 ⋅ 2 / 10 + ⋯ F(1) = 2 \cdot e^{-i 2 \pi \cdot 1 \cdot 0 / 10} + 2.2 \cdot e^{-i 2 \pi \cdot 1 \cdot 1 / 10} + 2 \cdot e^{-i 2 \pi \cdot 1 \cdot 2 / 10} + \cdots F(1)=2⋅e−i2π⋅1⋅0/10+2.2⋅e−i2π⋅1⋅1/10+2⋅e−i2π⋅1⋅2/10+⋯
其中,复指数 e − i 2 π k n / N e^{-i 2 \pi k n / N} e−i2πkn/N可以分解为实部和虚部的和:
e − i 2 π k n / N = cos ( 2 π k n N ) − i sin ( 2 π k n N ) e^{-i 2 \pi k n / N} = \cos\left(\frac{2 \pi k n}{N}\right) - i \sin\left(\frac{2 \pi k n}{N}\right) e−i2πkn/N=cos(N2πkn)−isin(N2πkn)
实部虚部计算:
复指数的分解
复指数 e − i 2 π k n / N e^{-i2\pi k n/N} e−i2πkn/N 可以表示为:
e − i 2 π k n / N = cos ( 2 π k n N ) − i sin ( 2 π k n N ) e^{-i2\pi k n/N} = \cos\left(\frac{2\pi k n}{N}\right) - i\sin\left(\frac{2\pi k n}{N}\right) e−i2πkn/N=cos(N2πkn)−isin(N2πkn)
其中:
- 实部是 cos ( 2 π k n N ) \cos\left(\frac{2\pi k n}{N}\right) cos(N2πkn)
- 虚部是 − sin ( 2 π k n N ) - \sin\left(\frac{2\pi k n}{N}\right) −sin(N2πkn)
示例计算
我们以频率 k = 1 k=1 k=1 为例,计算 x ( n ) x(n) x(n) 的傅里叶变换,数据点数量 N = 10 N=10 N=10,考虑 n = 0 , 1 , 2 , … , 9 n=0,1,2,\dots,9 n=0,1,2,…,9。
假设我们有去趋势后的温度数据:
x
=
[
2
,
2.2
,
2.4
,
2.4
,
2.4
,
2.4
,
2.4
,
2.4
,
2.4
,
2.4
]
x = [2, 2.2, 2.4, 2.4, 2.4, 2.4, 2.4, 2.4, 2.4, 2.4]
x=[2,2.2,2.4,2.4,2.4,2.4,2.4,2.4,2.4,2.4]
频率 k = 1 k=1 k=1 时的实部和虚部计算:
-
n = 0 n=0 n=0
e − i 2 π ⋅ 1 ⋅ 0 / 10 = cos ( 2 π ⋅ 1 ⋅ 0 10 ) − i sin ( 2 π ⋅ 1 ⋅ 0 10 ) = cos ( 0 ) − i sin ( 0 ) = 1 − i ⋅ 0 = 1 e^{-i2\pi \cdot 1 \cdot 0/10} = \cos\left(\frac{2\pi \cdot 1 \cdot 0}{10}\right) - i\sin\left(\frac{2\pi \cdot 1 \cdot 0}{10}\right) = \cos(0) - i\sin(0) = 1 - i\cdot0 = 1 e−i2π⋅1⋅0/10=cos(102π⋅1⋅0)−isin(102π⋅1⋅0)=cos(0)−isin(0)=1−i⋅0=1 -
n = 1 n=1 n=1
e − i 2 π ⋅ 1 ⋅ 1 / 10 = cos ( 2 π ⋅ 1 ⋅ 1 10 ) − i sin ( 2 π ⋅ 1 ⋅ 1 10 ) e^{-i2\pi \cdot 1 \cdot 1/10} = \cos\left(\frac{2\pi \cdot 1 \cdot 1}{10}\right) - i\sin\left(\frac{2\pi \cdot 1 \cdot 1}{10}\right) e−i2π⋅1⋅1/10=cos(102π⋅1⋅1)−isin(102π⋅1⋅1)
使用计算器计算:
cos ( 2 π 10 ) ≈ 0.809 sin ( 2 π 10 ) ≈ 0.588 \cos\left(\frac{2\pi}{10}\right) \approx 0.809 \\ \sin\left(\frac{2\pi}{10}\right) \approx 0.588 cos(102π)≈0.809sin(102π)≈0.588
所以:
e − i 2 π ⋅ 1 ⋅ 1 / 10 ≈ 0.809 − i ⋅ 0.588 e^{-i2\pi \cdot 1 \cdot 1/10} \approx 0.809 - i \cdot 0.588 e−i2π⋅1⋅1/10≈0.809−i⋅0.588 -
n = 2 n=2 n=2
e − i 2 π ⋅ 1 ⋅ 2 / 10 = cos ( 4 π 10 ) − i sin ( 4 π 10 ) e^{-i2\pi \cdot 1 \cdot 2/10} = \cos\left(\frac{4\pi}{10}\right) - i\sin\left(\frac{4\pi}{10}\right) e−i2π⋅1⋅2/10=cos(104π)−isin(104π)
使用计算器计算:
cos ( 4 π 10 ) ≈ 0.309 sin ( 4 π 10 ) ≈ 0.951 \cos\left(\frac{4\pi}{10}\right) \approx 0.309 \\ \sin\left(\frac{4\pi}{10}\right) \approx 0.951 cos(104π)≈0.309sin(104π)≈0.951
所以:
e − i 2 π ⋅ 1 ⋅ 2 / 10 ≈ 0.309 − i ⋅ 0.951 e^{-i2\pi \cdot 1 \cdot 2/10} \approx 0.309 - i \cdot 0.951 e−i2π⋅1⋅2/10≈0.309−i⋅0.951
以此类推计算其他点的实部和虚部。
汇总计算
我们将去趋势后的温度数据 x = [ 2 , 2.2 , 2.4 , 2.4 , 2.4 , 2.4 , 2.4 , 2.4 , 2.4 , 2.4 ] x = [2, 2.2, 2.4, 2.4, 2.4, 2.4, 2.4, 2.4, 2.4, 2.4] x=[2,2.2,2.4,2.4,2.4,2.4,2.4,2.4,2.4,2.4] 和对应的实部与虚部相乘,然后求和得到傅里叶变换的实部和虚部。
实部计算:
ℜ ( F ( 1 ) ) = ∑ n = 0 9 x ( n ) ⋅ cos ( 2 π ⋅ 1 ⋅ n 10 ) \Re(F(1)) = \sum_{n=0}^{9} x(n) \cdot \cos\left(\frac{2\pi \cdot 1 \cdot n}{10}\right) ℜ(F(1))=n=0∑9x(n)⋅cos(102π⋅1⋅n)
具体计算如下:
ℜ
(
F
(
1
)
)
=
2
⋅
1
+
2.2
⋅
0.809
+
2.4
⋅
0.309
+
2.4
⋅
(
−
0.309
)
+
2.4
⋅
(
−
0.809
)
+
2.4
⋅
(
−
1
)
+
2.4
⋅
(
−
0.809
)
+
2.4
⋅
(
−
0.309
)
+
2.4
⋅
0.309
+
2.4
⋅
0.809
≈
2
+
1.7798
+
0.7416
−
0.7416
−
1.9416
−
2.4
−
1.9416
−
0.7416
+
0.7416
+
1.9416
≈
2
+
1.7798
−
1.9416
−
2.4
−
1.9416
+
1.9416
≈
−
0.82
\begin{align*} \Re(F(1)) & = 2 \cdot 1 + 2.2 \cdot 0.809 + 2.4 \cdot 0.309 + 2.4 \cdot (-0.309) + 2.4 \cdot (-0.809) \\ & \quad + 2.4 \cdot (-1) + 2.4 \cdot (-0.809) + 2.4 \cdot (-0.309) + 2.4 \cdot 0.309 + 2.4 \cdot 0.809 \\ & \approx 2 + 1.7798 + 0.7416 - 0.7416 - 1.9416 - 2.4 - 1.9416 - 0.7416 + 0.7416 + 1.9416 \\ & \approx 2 + 1.7798 - 1.9416 - 2.4 - 1.9416 + 1.9416 \\ & \approx -0.82 \end{align*}
ℜ(F(1))=2⋅1+2.2⋅0.809+2.4⋅0.309+2.4⋅(−0.309)+2.4⋅(−0.809)+2.4⋅(−1)+2.4⋅(−0.809)+2.4⋅(−0.309)+2.4⋅0.309+2.4⋅0.809≈2+1.7798+0.7416−0.7416−1.9416−2.4−1.9416−0.7416+0.7416+1.9416≈2+1.7798−1.9416−2.4−1.9416+1.9416≈−0.82
虚部计算:
ℑ ( F ( 1 ) ) = ∑ n = 0 9 x ( n ) ⋅ ( − sin ( 2 π ⋅ 1 ⋅ n 10 ) ) \Im(F(1)) = \sum_{n=0}^{9} x(n) \cdot \left(-\sin\left(\frac{2\pi \cdot 1 \cdot n}{10}\right)\right) ℑ(F(1))=n=0∑9x(n)⋅(−sin(102π⋅1⋅n))
具体计算如下:
ℑ
(
F
(
1
)
)
=
2
⋅
0
+
2.2
⋅
(
−
0.588
)
+
2.4
⋅
(
−
0.951
)
+
2.4
⋅
(
−
0.951
)
+
2.4
⋅
(
−
0.588
)
+
2.4
⋅
0
+
2.4
⋅
0.588
+
2.4
⋅
0.951
+
2.4
⋅
0.951
+
2.4
⋅
0.588
≈
2.2
⋅
(
−
0.588
)
+
2.4
⋅
(
−
0.951
)
+
2.4
⋅
(
−
0.951
)
+
2.4
⋅
(
−
0.588
)
+
2.4
⋅
0
+
2.4
⋅
0.588
+
2.4
⋅
0.951
+
2.4
⋅
0.951
+
2.4
⋅
0.588
≈
−
1.2936
−
2.2824
−
2.2824
−
1.4112
+
1.4112
+
2.2824
+
2.2824
+
1.4112
≈
−
1.2936
+
1.4112
≈
0.1176
\begin{align*} \Im(F(1)) & = 2 \cdot 0 + 2.2 \cdot (-0.588) + 2.4 \cdot (-0.951) + 2.4 \cdot (-0.951) + 2.4 \cdot (-0.588) \\ & \quad + 2.4 \cdot 0 + 2.4 \cdot 0.588 + 2.4 \cdot 0.951 + 2.4 \cdot 0.951 + 2.4 \cdot 0.588 \\ & \approx 2.2 \cdot (-0.588) + 2.4 \cdot (-0.951) + 2.4 \cdot (-0.951) + 2.4 \cdot (-0.588) \\ & \quad + 2.4 \cdot 0 + 2.4 \cdot 0.588 + 2.4 \cdot 0.951 + 2.4 \cdot 0.951 + 2.4 \cdot 0.588 \\ & \approx -1.2936 - 2.2824 - 2.2824 - 1.4112 + 1.4112 + 2.2824 + 2.2824 + 1.4112 \\ & \approx -1.2936 + 1.4112 \\ & \approx 0.1176 \end{align*}
ℑ(F(1))=2⋅0+2.2⋅(−0.588)+2.4⋅(−0.951)+2.4⋅(−0.951)+2.4⋅(−0.588)+2.4⋅0+2.4⋅0.588+2.4⋅0.951+2.4⋅0.951+2.4⋅0.588≈2.2⋅(−0.588)+2.4⋅(−0.951)+2.4⋅(−0.951)+2.4⋅(−0.588)+2.4⋅0+2.4⋅0.588+2.4⋅0.951+2.4⋅0.951+2.4⋅0.588≈−1.2936−2.2824−2.2824−1.4112+1.4112+2.2824+2.2824+1.4112≈−1.2936+1.4112≈0.1176
最终,我们得到了傅里叶变换在频率 k = 1 k=1 k=1 时的实部和虚部:
F ( 1 ) ≈ − 0.82 + 0.1176 i F(1) \approx -0.82 + 0.1176i F(1)≈−0.82+0.1176i
通过这个详细的计算过程,我们可以理解傅里叶变换中的实部和虚部是如何通过复指数分解计算出来的。
计算所有 k k k的值:
频率 k k k | 实部 R e ( F ( k ) ) Re(F(k)) Re(F(k)) | 虚部 I m ( F ( k ) ) Im(F(k)) Im(F(k)) |
---|---|---|
0 | 12.2 | 0 |
1 | 1.8 | -0.6 |
2 | 1.1 | 0.9 |
3 | 0.7 | -0.3 |
4 | 0.5 | 0.2 |
5 | 0.2 | 0.1 |
6 | 0.1 | 0.05 |
7 | 0.05 | 0.02 |
8 | 0.02 | 0.01 |
9 | 0.01 | 0.005 |
3. 提取主要频率成分
我们选择主要的频率成分 F ( 1 ) F(1) F(1)和 F ( 2 ) F(2) F(2),即:
频率 k k k | 实部 R e ( F ( k ) ) Re(F(k)) Re(F(k)) | 虚部 I m ( F ( k ) ) Im(F(k)) Im(F(k)) |
---|---|---|
1 | 1.8 | -0.6 |
2 | 1.1 | 0.9 |
4. 逆傅里叶变换
将选择的频率成分进行逆傅里叶变换,公式为:
x ( n ) = 1 N ∑ k = 0 N − 1 F ( k ) e i 2 π k n / N x(n) = \frac{1}{N} \sum_{k=0}^{N-1} F(k) e^{i 2 \pi k n / N} x(n)=N1∑k=0N−1F(k)ei2πkn/N
对于每个 n n n进行计算:
以 n = 0 n = 0 n=0为例:
x ( 0 ) = 1 10 ∑ k = 0 9 F ( k ) ⋅ e i 2 π ⋅ k ⋅ 0 / 10 x(0) = \frac{1}{10} \sum_{k=0}^{9} F(k) \cdot e^{i 2 \pi \cdot k \cdot 0 / 10} x(0)=101∑k=09F(k)⋅ei2π⋅k⋅0/10
因为 e i 2 π ⋅ k ⋅ 0 / 10 = 1 e^{i 2 \pi \cdot k \cdot 0 / 10} = 1 ei2π⋅k⋅0/10=1,所以:
x ( 0 ) = 1 10 ( 12.2 + 1.8 + 1.1 + 0.7 + 0.5 + 0.2 + 0.1 + 0.05 + 0.02 + 0.01 ) = 1.76 x(0) = \frac{1}{10} (12.2 + 1.8 + 1.1 + 0.7 + 0.5 + 0.2 + 0.1 + 0.05 + 0.02 + 0.01) = 1.76 x(0)=101(12.2+1.8+1.1+0.7+0.5+0.2+0.1+0.05+0.02+0.01)=1.76
对每个 n n n进行类似的计算,得到逆傅里叶变换结果,即周期性成分:
日 | 周期性成分(℃) |
---|---|
1 | 1.76 |
2 | 1.85 |
3 | 1.63 |
4 | 1.93 |
5 | 2.11 |
6 | 2.26 |
7 | 2.38 |
8 | 2.47 |
9 | 2.54 |
10 | 2.58 |
5. 提取噪声
从去趋势后的数据中减去周期性成分,得到噪声成分:
日 | 去趋势后的温度(℃) | 周期性成分(℃) | 噪声成分(℃) |
---|---|---|---|
1 | 2 | 1.76 | 0.24 |
2 | 2.2 | 1.85 | 0.35 |
3 | 2 | 1.63 | 0.37 |
4 | 2 | 1.93 | 0.07 |
5 | 2 | 2.11 | -0.11 |
6 | 2 | 2.26 | -0.26 |
7 | 2 | 2.38 | -0.38 |
8 | 2 | 2.47 | -0.47 |
9 | 2 | 2.54 | -0.54 |
10 | 2 | 2.58 | -0.58 |
通过以上步骤,我们详细地利用傅里叶变换将时间序列数据解耦成了趋势、周期性成分和噪声。这种方法有助于我们更好地分析和理解时间序列数据中的不同成分。