如何去除时间序列中的特定谐波？

Posted 2023-02-25

【中文标题】如何去除时间序列中的特定谐波？

【英文标题】：How to remove particular harmonics in time series?

【发布时间】：2020-05-27 19:27:03

【问题描述】：

我有一个时间序列，其中包含 2014 年到 2019 年的每小时（移动平均）数据。

我使用fft() 函数查找谐波，并使用quantmod::findPeaks() 查找最高值的位置，结果得到以下位置：8 56 2193 4384 6575 8766 306737 308928 311119 313310 315447 315495（时间序列中有 315499 个元素)。

现在我需要将其转换回时域，但用零代替这些位置的数据值。我该怎么做？

另外，我注意到当我绘制fft() 函数（FFT 的绝对值）时，x 轴仍然显示 2014 年到 2019 年。它不应该变成 Hz 轴（频率）吗？

z 是我的时间序列

FFT <- fft(z)
FFT <- abs(FFT/FFT[1])
plot(FFT)

PK <- findPeaks(FFT, thresh = T)
PK

PK 的结果是我提到的 12 个峰值，代表最高的谐波失真。 我需要删除那些。我相信我可以简单地对 12 个山峰执行FFT[8]=0，但我不确定。 （我试过了，但它不会让我回到时域）

我还需要证明生成的 FFT 图确实是傅里叶频谱，但 x 轴不是频率。关于如何进行更正的任何想法？

【问题讨论】：

你能给我们看看minimal reproducible example好吗？请注意，zz <- Re(fft(fft(z),inverse=TRUE))/length(z) 给出的结果等于原始 z ...

谢谢你，Ben，代码是正确的，谢谢，我只需要先去除谐波。我编辑了问题以包含示例。

findPeaks() 函数来自哪个包？ quantmod? （我不认为它是基础 R 的一部分）

是的，来自quantmod

【参考方案1】：

模拟方波+（低）噪声

set.seed(101)
ncycles <- 20
len <- 50
z <- rep(rep(0:1,each=len),ncycles) + rnorm(len*ncycles,mean=0, sd=0.02)
plot(z,type="l")

FFT 并找到峰值

f <- fft(z)
F <- abs(f/f[1])
PK <- quantmod::findPeaks(F,thresh=0.01)

绘制 sqrt（功率谱）：

plot(F,type="l")
abline(v=PK,col=2)

将峰值归零：

fz <- f
## position of 'peaks' is actually *one before* PK (not sure why?)
fz[PK-1] <- 0

逆变换，与原来的比较：

zz <- Re(fft(fz,inverse=TRUE))/length(f)
plot(z,type="l")
lines(zz,type="l",col=2)

【讨论】：

在功率谱中，x轴仍然为“时间”，有没有办法将其转换为频率？我的意思是，它是频域，对吧？所以它不应该是 x 轴上的时间。

我怀疑您正在将数据制作成时间序列，这不是必需的。如果你想取消时间序列，你可以使用c()