从时序数据中提取特征分量
原始数据
import matplotlib.pyplot as plt from matplotlib import font_manager fname="/usr/local/python3.6/lib/python3.6/site-packages/matplotlib/mpl-data/fonts/ttf/simhei.ttf" zhfont = font_manager.FontProperties(fname=fname) plt.figure(figsize=(14, 4)) plt.plot(dataset) plt.title("原始数据", fontproperties=zhfont, fontsize=12) plt.grid() plt.show()
方法一:傅里叶变换
import datetime import numpy as np import pandas as pd from scipy.fft import fft, ifft import matplotlib.pyplot as plt from matplotlib import font_manager fname="/usr/local/python3.6/lib/python3.6/site-packages/matplotlib/mpl-data/fonts/ttf/simhei.ttf" zhfont = font_manager.FontProperties(fname=fname) # 样本:288个5分钟点,即一天 N = len(dataset) # 采样间隔(秒) T = 300 # 每5分钟一个数据点 # 采样频率(Hz) Fs = 1 / T # 即每300秒一个样本 # 总观察时间(秒) D = N * T # 一整天的时间长度 # 生成时间向量 t = np.arange(0, D, T) # 生成示例时序数据(两个正弦波的叠加) # 假设每天有一个主要周期为24小时(低频)和一个次要周期为12小时(高频) # y = 1.0 * np.sin(2.0 * np.pi * (1 / (24*3600)) * t) + 0.5 * np.sin(2.0 * np.pi * (1 / (12*3600)) * t) y = dataset # 计算傅里叶变换 yf = fft(y) # 频率向量 xf = np.fft.fftfreq(N, T) # 绘制频谱图 plt.figure(figsize=(14, 4)) plt.plot(xf, np.abs(yf)) plt.title("将原始数据从时域转成频域", fontproperties=zhfont, fontsize=12) plt.xlabel("Frequency (Hz)", fontsize=12) # plt.ylabel("Amplitude") plt.grid() plt.show() # 提取各个正弦分量 def extract_sine_component(yf, freq_index): yf_component = np.zeros_like(yf) yf_component[freq_index] = yf[freq_index] if freq_index != 0: yf_component[-freq_index] = yf[-freq_index] # 对称频谱 y_component = ifft(yf_component) return np.real(y_component) # 找出幅值最大的10个频率分量及其索引(只处理前N//2部分) amplitudes = np.abs(yf[:N//2]) indices = np.argsort(amplitudes)[-10:] plt.figure(figsize=(14, 4)) y_reconstructed = np.zeros_like(y) for idx in indices: y_component = extract_sine_component(yf, idx) y_reconstructed += y_component plt.plot(t, y_component, label=f'Frequency: {xf[idx]:.5f} Hz') plt.title("基于傅里叶变换将原始数据拆分成多个正弦波,取Top10", fontproperties=zhfont, fontsize=12) plt.legend() plt.grid() plt.show() # 绘制重建后的时序数据与原始数据的对比 plt.figure(figsize=(14, 4)) plt.plot(t, y, label='原始数据') plt.plot(t, y_reconstructed, label='Top10正弦波叠加', linestyle='--') plt.title("原始数据 vs Top10正弦波叠加", fontproperties=zhfont, fontsize=12) plt.legend() plt.grid() plt.show()
作者:Standby — 一生热爱名山大川、草原沙漠,还有我们小郭宝贝!
出处:http://www.cnblogs.com/standby/
本文版权归作者和博客园共有,欢迎转载,但未经作者同意必须保留此段声明,且在文章页面明显位置给出原文连接,否则保留追究法律责任的权利。
出处:http://www.cnblogs.com/standby/
本文版权归作者和博客园共有,欢迎转载,但未经作者同意必须保留此段声明,且在文章页面明显位置给出原文连接,否则保留追究法律责任的权利。