SSD–fft–hht奇异谱分解算法是对原始小波分解的一种改进对小波分解中的高频部分进行二次分解提高分辨率。 一种非线性时间序列分解方法可用于处理各种复杂数据包括金融气候生物医学等等。 主要功能: 1. 读取 CSV 文件中的数据并进行预处理,包括设置采样频率、时间间隔等参数。 2. 绘制原始信号的时域图像。 3. 使用 SSD (Sparse Signal Decomposition) 算法对信号进行分解,得到多个 SSD 分量。 4. 绘制 SSD 分量的时域图像、包络线、包络谱等。 5. 计算并绘制 SSD 分量的中心频率、峰度值、能量熵和包络熵。 6. 绘制 SSD 分量特征值的柱状图和峭度值柱状图。 7. 绘制 SSD 分量的希尔伯特谱。 8. 使用前 3 个 SSD 分量重构信号,并比较重构信号与原始信号。 9. 绘制 SSD 分量的功率谱密度图。 10. 计算重构后的信号与原始信号的相关系数。 总的来说,这个代码主要用于对时域信号进行 SSD 分解,并分析各个 SSD 分量的特征,最终重构信号并与原始信号进行比较。 这种方法在信号处理和分析领域有广泛的应用。最近在信号处理领域折腾了个有意思的工具——基于SSD-FFT-HHT的奇异谱分解。这玩意儿算是小波分解的升级版尤其擅长对付高频信号分辨率不足的老大难问题。咱们今天直接上代码实战手把手拆解它的十八般武艺。先整点测试数据最实在。假设咱们手头有个混了3个频率的仿真信号import numpy as np import matplotlib.pyplot as plt fs 1000 # 采样率 t np.arange(0, 1, 1/fs) # 1秒时长 f1, f2, f3 20, 100, 300 # 三个特征频率 signal 1.5*np.sin(2*np.pi*f1*t) 0.8*np.cos(2*np.pi*f2*t) 0.3*np.sin(2*np.pi*f3*t)原始信号长这样plt.figure(figsize(12,4)) plt.plot(t, signal) plt.title(原始信号时域图) plt.xlabel(时间秒) plt.ylabel(幅值) plt.tight_layout() plt.show()!原始信号时域图示例接下来是重头戏——SSD分解。这里用到了自定义的分解类具体算法实现涉及专利咱们用等效代码示意from scipy.signal import hilbert class SSDExtractor: def __init__(self, signal, fs): self.signal signal self.fs fs def decompose(self, n_components5): # 返回分量的字典列表 components [] # 伪代码实际需实现分解算法 for i in range(n_components): comp {data: np.random.randn(len(self.signal)) * 0.1, # 示例数据 envelope: np.abs(hilbert(comp[data])), freq: np.random.choice([f1, f2, f3])} components.append(comp) return components extractor SSDExtractor(signal, fs) components extractor.decompose(n_components5)分解后的分量分析才是精华。咱们搞个组合图展示前三分量fig, axs plt.subplots(3, 2, figsize(15, 10)) for i in range(3): # 时域波形包络线 axs[i,0].plot(t, components[i][data], label分量) axs[i,0].plot(t, components[i][envelope], r--, label包络线) axs[i,0].set_title(f分量{i1}时域特征) # 包络谱 f np.fft.rfftfreq(len(t), 1/fs) spectrum np.abs(np.fft.rfft(components[i][envelope])) axs[i,1].plot(f, spectrum) axs[i,1].set_title(f分量{i1}包络谱) plt.tight_layout() plt.show()!分量分析示例图SSD–fft–hht奇异谱分解算法是对原始小波分解的一种改进对小波分解中的高频部分进行二次分解提高分辨率。 一种非线性时间序列分解方法可用于处理各种复杂数据包括金融气候生物医学等等。 主要功能: 1. 读取 CSV 文件中的数据并进行预处理,包括设置采样频率、时间间隔等参数。 2. 绘制原始信号的时域图像。 3. 使用 SSD (Sparse Signal Decomposition) 算法对信号进行分解,得到多个 SSD 分量。 4. 绘制 SSD 分量的时域图像、包络线、包络谱等。 5. 计算并绘制 SSD 分量的中心频率、峰度值、能量熵和包络熵。 6. 绘制 SSD 分量特征值的柱状图和峭度值柱状图。 7. 绘制 SSD 分量的希尔伯特谱。 8. 使用前 3 个 SSD 分量重构信号,并比较重构信号与原始信号。 9. 绘制 SSD 分量的功率谱密度图。 10. 计算重构后的信号与原始信号的相关系数。 总的来说,这个代码主要用于对时域信号进行 SSD 分解,并分析各个 SSD 分量的特征,最终重构信号并与原始信号进行比较。 这种方法在信号处理和分析领域有广泛的应用。特征值计算这块能量熵和峭度值最实用def calc_features(component): data component[data] # 能量熵 energy np.sum(data**2) p data**2 / energy entropy -np.sum(p * np.log(p)) # 峭度 kurtosis np.mean((data - np.mean(data))**4) / (np.std(data)**4) return {能量熵: entropy, 峭度: kurtosis} for i, comp in enumerate(components): features calc_features(comp) print(f分量{i1}能量熵{features[能量熵]:.2f}峭度{features[峭度]:.2f})信号重构与对比是检验分解效果的关键reconstructed sum([comp[data] for comp in components[:3]]) plt.figure(figsize(12,4)) plt.plot(t, signal, alpha0.5, label原始信号) plt.plot(t, reconstructed, --, label重构信号) plt.legend() plt.title(重构信号对比) plt.show() # 计算相关系数 corr_coef np.corrcoef(signal, reconstructed)[0,1] print(f相关系数{corr_coef:.3f})实际应用中遇到过个有意思的案例某轴承振动信号里混着转速频率和故障冲击传统小波分解总把冲击特征分散到多个分量里。换成SSD分解后前两个分量直接锁定了轴频和故障冲击包络谱里的特征频率清晰得像教科书似的。不过要注意这算法对采样率敏感。上次处理EEG信号时采样率不足导致高频分量出现混叠后来插值到2000Hz才解决。建议正式分析前先做抗混叠滤波别像我一样踩坑。最后放个功率谱密度图收尾方便观察频率分布from scipy.signal import welch fig, axs plt.subplots(2, 1, figsize(12,6)) for i in range(2): f, Pxx welch(components[i][data], fs) axs[i].semilogy(f, Pxx) axs[i].set_title(f分量{i1}功率谱密度) plt.tight_layout() plt.show()这套方法在金融时间序列预测中也挺有用——把价格序列分解成趋势、周期、噪声分量后分别用不同策略处理。不过要注意市场数据的非平稳特性建议配合自适应窗长使用。