使用周期图检测时间序列的季节性

Question

我正在使用周期图捕获时间序列的季节性，我想使用前十个频率分量来创建季节性时间序列，到目前为止，我绘制了周期图：

data=elec_price[:48*365]
from scipy.signal import periodogram
f, Pxx_den = periodogram(data)

数据是价格的一年子集，间隔为 30 分钟 , 并在按幅值排序后从fourier_coefficients列表中捕获前十个频率分量

fourier_coefficients=Pxx_den.tolist()
fourier_coefficients=pd.DataFrame(fourier_coefficients,columns=['amplitude'])

最高振幅的最高频率是： 365,2,730,22,52,5,729,8 ，接下来我需要做的是使用这些最高频率分量来获得时间序列的季节性，我生成了每个频率分量的正弦波，并将它们加在一起绘制时间序列，我想我不确定这是不是正确的方法，因为我记得频率分量应该有一个虚部但我无法从周期图数据中找到它，其次我只是假设这些分量中的每一个都是正弦波。

sin_waves=fourier_coefficients.drop(columns=['frequency coneficient','amplitude']).copy()
sin_waves[' sin_wave_1']=0
sin_waves[' sin_wave_2']=0
sin_waves[' sin_wave_3']=0
sin_waves[' sin_wave_4']=0
sin_waves[' sin_wave_5']=0
sin_waves[' sin_wave_6']=0
sin_waves[' sin_wave_7']=0
sin_waves[' sin_wave_8']=0
for i in range(8761):
    sin_waves[' sin_wave_1'][i]= fourier_coefficients['amplitude'][365]*math.sin(math.pi*i*365/8761)
    sin_waves[' sin_wave_2'][i]= fourier_coefficients['amplitude'][2]*math.sin(math.pi*i*2/8761)
    sin_waves[' sin_wave_3'][i]= fourier_coefficients['amplitude'][730]*math.sin(math.pi*i*730/8761)
    sin_waves[' sin_wave_4'][i]= fourier_coefficients['amplitude'][22]*math.sin(math.pi*i*22/8761)
    sin_waves[' sin_wave_5'][i]= fourier_coefficients['amplitude'][52]*math.sin(math.pi*i*52/8761)
    sin_waves[' sin_wave_6'][i]= fourier_coefficients['amplitude'][5]*math.sin(math.pi*i*5/8761)
    sin_waves[' sin_wave_7'][i]= fourier_coefficients['amplitude'][729]*math.sin(math.pi*i*729/8761)
    sin_waves[' sin_wave_8'][i]= fourier_coefficients['amplitude'][8]*math.sin(math.pi*i*8/8761)
sin_waves['accumulated_sin_wave']=(sin_waves[' sin_wave_1']+sin_waves[' sin_wave_2']
                + sin_waves[' sin_wave_3']+sin_waves[' sin_wave_4']+sin_waves[' sin_wave_5']+
                sin_waves[' sin_wave_6']+sin_waves[' sin_wave_7']+sin_waves[' sin_wave_8'])

8761 是数据中的样本数。
然后我绘制了sinaves['accumulated_sin_wave']

这是第一年的季节性成分图：

所以，我想知道我应该怎么做，而不是仅仅绘制这些正弦波来捕获最高频率的真实分量。

Answer 1

嗯，对我来说更好更简单的方法就是使用 FFT 过滤 IFFT 模型，它看起来很简单，首先我使用以下方法为我的数据绘制 FFt：

one_year_values=48*365
data=elec_price[:one_year_values]
from numpy.fft import fft, fftfreq, ifft

n=len(data)
freqs=fftfreq(n)
mask = freqs>0


fft_values=fft(data)
fft_list=np.copy(fft_values).tolist()
fft_filtered=np.copy(fft_list).tolist() # copying from the original fft list in order to filter it

power= 2*(np.abs(fft_values/one_year_values))**2 # calculating the power make sit easy to compare between frequency components, because they are complex numbers, and the power of them gets the absolute values of them.
power=power.tolist()

在此之后，我按幅度对功率列表进行了排序，并使用了前 10 个频率，我的数据是： 0,365,17155,1,17519,11,17509,366,17154,2,17518 请注意，0 是数据中的 DC 分量，因此您可以使用或不使用它，具体取决于您是否希望在季节性中进行偏移。 另外，请注意，这些频率中的每一个都是对的，每个对是彼此的镜像，因为 FFT 是一个图，其中每个频率分量及其负频率具有相同的功率。 在选择了我的前 5 个频率后，我过滤了 FFT_filtered 列表：

tuple(fft_filtered)
for i in range(len(fft_filtered)):
    fft_filtered[i]=0
for i in (0
,365
,17155
,1
,17519
,11
,17509
,366
,17154
,2
,17518):
    fft_filtered[i]=fft_list[i]

然后我使用 IFFT 反转 FFt 以获取使用前 5 个频率的数据图：

reverse=np.real(ifft(fft_filtered))
#%%
plt.plot(reverse)
plt.legend(loc='best')
plt.ylabel('seasonal factor')
plt.xlabel('time stamps')
plt.title('1st  year seasonality-top 5 freq')
plt.show(block=False)

这是反向列表的情节： 1 year seasonality using top 5 frequencies