ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

CEEMDAN-ISOS-VMD-GRU-ARIMA:非平稳时间序列预测的分工流水线

CEEMDAN-ISOS-VMD-GRU-ARIMA:非平稳时间序列预测的分工流水线 简介这是一套面向时间序列预测学习与毕业设计的完整Python实现融合CEEMDAN、ISOS、VMD、GRU与ARIMA等算法适合计算机、电子信息、数学等专业学生用于课程设计、期末大作业或毕业设计。压缩包共3个文件包含1个Python主程序与2个CSV数据集整体仅52KB代码采用参数化编程并配有保姆级注释几乎一行一注释基于Anaconda、PyCharm与TensorFlow环境即可运行方便初学者从数据读取到模型训练逐步理解完整流程。资源已收获387人次学习浏览可作为CEEMDAN分解、ISOS优化、VMD模态分解与GRU-ARIMA组合预测的参考实现实现从数据预处理、分解优化到组合预测与结果对比的完整链路。作者为资深算法工程师程序结构清晰、参数易改附带可直接运行的实际城市CSV数据既能支撑算法对比实验也可作为论文图表或项目预研的快速起点。1. CEEMDAN-ISOS-VMD-GRU-ARIMA把时间序列预测做成一条分工明确的流水线做负荷预测或者量化策略的人应该都遇到过这种尴尬用ARIMA预测非线性序列残差大到没法看换成GRU趋势又拟合不稳末端直接发散。这个标题里的组合模型核心不是堆模型而是把一条非平稳序列拆成不同复杂度成分再按成分分配模型。CEEMDAN负责粗分解ISOS优化VMD参数做二次分解GRU抓非线性高频分量ARIMA收趋势与残差最后重构得到预测。下面按可复现的顺序把这套链路拆开讲清楚适合拿到源码但不知道怎么调参的Python从业者。2. 四个模块的分工CEEMDAN 拆粗、VMD 拆细、GRU 抓非线性、ARIMA 收残差2.1 为什么 CEEMDAN 之后还要 VMD一次分解不干净这里先说清楚一个容易被忽略的点CEEMDAN和VMD不是竞争关系是接力关系。EMD类方法的目标是把非平稳序列展开成一组本征模态函数IMF每个IMF在同一时刻尽量只包含一个频率成分。但真实序列里突变、噪声和有效信号经常挤在同一频段EMD会模态混叠。EEMD用加噪声的方式缓解混叠代价是分解不完备重构有残差噪声。CEEMDAN在每一层分解时加入自适应白噪声边分解边抵消噪声得到的IMF集合再相加可以近似无损重构原始序列这一步是整条流水线的地基。但CEEMDAN对高频段的分离并不彻底。实际跑数据时会发现前几个IMF的频谱仍然严重重叠特别是当采样点上有强噪声时IMF1里往往同时包含真实事件和噪声。如果直接拿这种IMF训练GRU模型会去拟合噪声换到测试集上就是翻车。VMD则不同它把信号分解问题放进变分框架约束每个模态在频域内围绕自己的中心频率紧支撑对非平稳高频串更有效。所以常见做法是CEEMDAN全局拆再把过零率高、波动剧烈的那部分IMF相加成一条高频复合序列交给VMD二次分解。顺带提一个参数上的认知VMD一次分解的效果完全被K和alpha两个参数决定。K是模态数K给小了多个频率成分挤在同一个模态里K给大了同样的频率被拆成好几个虚假模态后面逐个预测会浪费计算量。alpha是惩罚因子控制模态带宽alpha太小模态会碎成脉冲alpha太大模态之间重叠。这两个参数没有解析解只能靠搜索这就是ISOS进场的原因。注意CEEMDAN重构是线性叠加VMD重构也是线性叠加所以整条流水线最后预测重构时可以直接把所有分量预测累加不需要额外解码网络。这个性质是混合分解类模型能成立的前提。2.2 GRU 和 ARIMA 的分工边界高频交给 GRU低频与趋势交给 ARIMA很多混合模型做得不理想问题往往出在分工没做对。ARIMA本质是线性模型它对低频趋势、慢变周期这类低复杂度成分非常擅长参数少、可解释、外推稳定。但把ARIMA扔到高频段上它会把一切波动都当成随机噪声预测曲线被过度平滑。GRU作为门控循环网络能学习高频段里反复出现的局部模式和非线性状态转换但它的趋势外推能力很差让GRU预测单调趋势经常出现预测值持续走平甚至反向末端漂移没法看。所以分工边界要按序列成分划VMD二次分解出来的子模态归GRU逐一预测CEEMDAN分解得到的最后几个低频IMF和趋势残差归ARIMA。这里有一个实践经验中频IMF怎么分我一般再算一次过零率过零率高于0.05的分量进GRU低于0.05的进ARIMA。这个阈值不是唯一标准但比凭肉眼挑IMF靠谱。如果某个中频IMF方差贡献很小、过零率又高说明它基本是噪声成分可以直接不预测重构时也不会带来明显误差。GRU在整条链路里的角色更像“细节补全”。VMD已经帮它把频谱分开了GRU只需要学每个子模态在滑动窗口内的短时演化不需要从头理解整条序列的趋势。训练难度大幅下降这也是为什么分解后再用GRU通常比直接用GRU预测原始序列效果稳定得多。ARIMA的角色则相反它吸收的是去掉高频后的“骨架”保证最终预测不会像纯神经网络那样失去均值回归能力。2.3 ISOS 在这个链路里到底优化谁VMD 参数与 GRU 超参数ISOS这个缩写在不同实现里给出过不同展开一种叫Improved Sparrow Optimization Algorithm另一种叫Improved Swarm Optimization Strategy。拿到任何源码第一件事不是纠结缩写全称而是看它的改进策略落在哪里。常见的改进点有三个发现者位置更新里加自适应t分布扰动、追随者更新用莱维飞行、初始化用混沌映射代替均匀随机。这些改进都是为了同一个目标在有限的迭代次数内让种群既能快速收缩到好解区域又能偶尔跳出局部最优。ISOS在这条链路里负责两个搜索任务而且必须分阶段。第一个任务是在VMD分解之前搜索K和alpha。目标函数最常用包络熵也就是所有模态包络熵之和。包络熵越低说明每个模态的包络谱越集中频率混叠越少。第二个任务是在GRU训练之前搜索隐藏层单元数、学习率和滑动窗口长度。目标函数用验证集的MAPE因为GRU是随机初始化训练同一个超参组合跑两次结果也有波动常见做法是每个组合固定随机种子后训练一次用验证集MAPE做排序。为什么不推荐网格搜索VMD目标函数虽然单次便宜但ISOS每评估一个个体都要完整做一次VMD分解GRU更贵一次训练可能要几十秒。网格搜索在三维超参空间里要跑几百次ISOS用十来个种群迭代几十次虽然也有重复计算但因为它会保留历史最优并朝最优区域收缩同样预算下更容易找到可用参数。这个区别在超参空间大时尤其明显。注意ISOS优化VMD参数时搜索到的K要保持整数alpha保持连续。种群初始化时K用随机整数alpha用均匀分布。目标函数里不要忘记给K过大的情况加一个小惩罚否则优化器会发现模态数越多包络熵总和越小然后给你一个K15的荒唐结果。3. 先跑通最小闭环CEEMDAN 分解与 ISOS-VMD 参数寻优3.1 环境准备与数据切分从 Python 安装到虚拟环境先说环境。整个链路用到的库比较多强烈建议在虚拟环境里安装直接往系统Python里塞容易跟别的项目打架。我用VSCode做Python开发一般先建一个venv再装依赖PyCharm里也可以用conda环境。Python版本选3.8到3.10之间比较稳3.10以下安装TensorFlow都流畅3.11以上个别旧版本vmdpy会有编译问题。主要依赖如下pip install EMD-signal vmdpy pandas numpy scikit-learn statsmodels pmdarima tensorflow说明PyEMD的pip包名是EMD-signal导入时用import PyEMD。vmdpy的导入名和包名一致。pmdarima提供auto_arima内部依赖statsmodels。如果只是为了先跑通分解可以只装前四个库后面建模再补TensorFlow。数据切分是时间序列预测里最容易被忽略的步骤普通机器学习里的随机打乱在这里绝不能做。时序数据必须按时间顺序切import pandas as pd import numpy as np df pd.read_csv(series.csv, parse_dates[date]) raw df[value].values.astype(float) # 先处理缺失值时间序列没有删除这一说只能填充 raw pd.Series(raw).ffill().bfill().values # 按时间顺序切分前70%训练后30%测试 split int(len(raw) * 0.7) train_data, test_data raw[:split], raw[split:]说明缺失值用ffill前向填充和bfill后向填充是时序数据的默认做法不要用均值填充否则会引入未来信息。切分比例一般取70/30或80/20如果数据有强周期性最好保证测试集里包含至少一个完整周期否则验证结果没有代表性。3.2 CEEMDAN 分解PyEMD 的 trials 和 epsilon 怎么设分解代码非常短难在参数怎么定from PyEMD import CEEMDAN import numpy as np def ceemdan_decompose(series, trials100, epsilon0.005): ceemdan CEEMDAN(trialstrials, epsilonepsilon) imfs ceemdan(series) return np.asarray(imfs) # shape: (n_imfs, n_points)说明trials是添加白噪声的次数每次加不同的噪声实现再平均噪声会相互抵消。trials越大分解越稳定但计算时间线性增长我一般先用100如果IMF数量波动明显就调到200。epsilon是噪声幅值系数相对原序列标准差的比例0.005是稳妥起点。太小退化成普通EMD模态混叠重现太大IMF1会变成纯噪声占比很高的伪分量。拿到IMF之后不要急着建模先做一次“体检”观察每个分量的过零率和方差def zero_cross_rate(series): sign np.sign(series) return np.sum(np.diff(sign) ! 0) / len(series) for i, imf in enumerate(imfs): print(fIMF{i}: zero_cross_rate{zero_cross_rate(imf):.4f}, fvar{np.var(imf):.4f})说明过零率反映振荡快慢方差反映能量大小。过零率高且方差大的分量是高频段的主要贡献者要相加起来做VMD二次分解。过零率接近0.5的分量基本是白噪声可以不进模型。最后一个IMF通常是趋势项方差大、过零率低直接走ARIMA。这个体检结果也决定了第4章的分量分配不要跳过。3.3 ISOS-VMD 参数寻优包络熵目标函数与搜索范围把高频分量相加成一条序列后交给VMD。VMD本身是确定性的但K和alpha得调from vmdpy import VMD from scipy.signal import hilbert def envelope_entropy(imf): envelope np.abs(hilbert(imf)) p envelope / (np.sum(envelope) 1e-12) return -np.sum(p * np.log(p 1e-12)) def vmd_cost(K, alpha, high_freq): u, _, _ VMD(high_freq, alpha, 0, int(K), 1, 1, 1e-7) return np.sum([envelope_entropy(u[i, :]) for i in range(int(K))])说明vmdpy的VMD函数签名是VMD(f, alpha, tau, K, DC, init, tol)。tau是噪声容忍参数取0表示不允许噪声项DC1表示保留直流分量init1表示用均匀分布初始化中心频率tol是收敛阈值。envelope_entropy计算时加1e-12是为了防log(0)。包络熵越小模态越纯净。然后跑ISOS搜索。这里给一个教学用的简化骨架重点在流程def isos_vmd_optimize(high_freq, pop10, max_iter30, seed42): rng np.random.default_rng(seed) pop_K rng.integers(3, 11, pop) # K 搜索范围 [3, 10] pop_alpha rng.uniform(200, 3000, pop) # alpha 搜索范围 best_cost float(inf) best_params None for it in range(max_iter): for i in range(pop): cost vmd_cost(pop_K[i], pop_alpha[i], high_freq) if cost best_cost: best_cost cost best_params (int(pop_K[i]), pop_alpha[i]) # 发现者位置更新让 alpha 向当前最优个体靠拢 best_alpha best_params[1] pop_alpha pop_alpha 0.5 * (best_alpha - pop_alpha) rng.normal(0, 20, pop) # 防止越界 pop_alpha np.clip(pop_alpha, 200, 3000) # 定期用 t 分布扰动防止种群早熟 if it % 5 0: idx rng.integers(0, pop, max(1, pop // 3)) pop_alpha[idx] rng.standard_t(df3, sizelen(idx)) * 50 pop_alpha np.clip(pop_alpha, 200, 3000) return best_params, best_cost说明这是麻雀搜索里发现者位置更新的简化版省略了追随者和警戒者。完整ISOS的改进点一般还包括自适应步长、混沌初始化、莱维飞行。搜索范围方面K取3到10是因为时间序列分解中K超过10后容易出现虚假模态alpha取200到3000覆盖了从窄带到较宽带的范围。收缩系数0.5让个体快速向当前最优靠拢所以每5代做一次t分布扰动来防止早熟。t分布的自由度df3时尾部更厚跳出局部最优的能力比正态扰动强。注意VMD结果有个参数敏感性玄学同一组K和alpha输入序列长度变化100个点最优参数可能就变了。所以不要在训练集上优化完参数后拿去硬套测试集正确做法是在训练集上优化预测时保持最优参数不变只在分解新数据时重新调用VMD。4. GRU 与 ARIMA 建模分量分配、训练与重构4.1 分量分配策略什么给 GRU、什么给 ARIMACEEMDAN分解出来的IMF配合VMD二次分解后的子模态数量可能达到十几个。每一个都训练一个GRU不现实所以要先分类高频复合序列经VMD得到的子模态全部给GRU逐个建模。CEEMDAN原始IMF中过零率中等、方差贡献大于1%的给GRU但可以共享一个模型结构仅用不同数据训练。过零率极低、方差占比大的趋势IMF和残差给ARIMA。方差贡献极小且过零率接近0.5的认为是噪声分量直接忽略。这套分配规则的核心是让每个模型处理它擅长的序列。GRU处理的是局部波动模式ARIMA处理的是全局趋势骨架。如果趋势分量里混进一点高频波动ARIMA也能容忍但GRU无法容忍输入序列里同时存在两种时间尺度。这也是为什么VMD必须先做否则GRU的输入维度会失控。预测长度对齐是所有分量最终相加的前提。测试集长度固定每个GRU要输出同样长度的预测序列ARIMA也要预测同样长度。如果某个分量因为数据窗不足导致预测长度偏短重构时就要对齐尾部这个操作最容易引入误差两个长度不一的numpy数组直接相加会直接报错。4.2 GRU 训练代码滑动窗口、学习率与递归预测对每一个分配到GRU的分量先构造滑动窗口数据集def make_dataset(series, window12): X, y [], [] for i in range(len(series) - window): X.append(series[i:iwindow]) y.append(series[iwindow]) return np.array(X).reshape(-1, window, 1), np.array(y) X, y make_dataset(imf, window12) split int(len(X) * 0.8) X_train, X_val, y_train, y_val X[:split], X[split:], y[:split], y[split:]说明窗口长度window是GRU最敏感的超参数之一常见取值范围8到24需要ISOS搜索或经验设定。窗口太长训练样本量减少模型学不到新变化窗口太短模型看不到一个完整周期滞后明显。如果数据有明显周期窗口至少要覆盖一个周期。reshape成三维是因为Keras的GRU要求输入形状为(batch, timesteps, features)。模型定义与训练import tensorflow as tf def build_gru(window, units32, lr0.001): model tf.keras.Sequential([ tf.keras.layers.GRU(unitsunits, input_shape(window, 1), return_sequencesFalse), tf.keras.layers.Dropout(0.2), tf.keras.layers.Dense(1) ]) model.compile(optimizertf.keras.optimizers.Adam(learning_ratelr), lossmse) return model model build_gru(12, units32, lr0.001) callback tf.keras.callbacks.EarlyStopping( monitorval_loss, patience10, restore_best_weightsTrue ) model.fit(X_train, y_train, validation_data(X_val, y_val), epochs100, batch_size32, callbacks[callback], verbose0)说明units取32还是64取决于分量复杂度方差大、包络熵高的分量给64平滑分量给32。学习率0.001是Adam的默认值如果训练loss震荡优先降到0.0005如果loss下降太慢再考虑升到0.002。Dropout0.2用来抑制过拟合分量多的时候逐个训练容易记忆噪声。EarlyStopping的patience10表示验证集loss连续10轮不降就停止并恢复最优权重。多步预测用递归滑窗def recursive_forecast(model, last_window, n_steps): forecasts [] cur last_window.astype(np.float32).copy() for _ in range(n_steps): pred model.predict(cur.reshape(1, -1, 1), verbose0)[0, 0] forecasts.append(pred) cur np.roll(cur, -1) cur[-1] pred return np.array(forecasts)说明last_window是测试集起点前的一段历史值。递归预测会把上一步预测值当作下一步输入所以误差会累积。如果训练时窗口输入全部来自真实值预测时输入全部来自预测值分布不一致这就是GRU滞后问题的根源之一。缓解手段是在训练数据里混入少量预测噪声或直接用scheduled sampling策略。这个坑在第5章还会有详细展开。4.3 ARIMA 建模与最终重构auto_arima 的 d 怎么定低频分量和趋势项交给ARIMA最简单的方式是用auto_arima自动定阶from pmdarima import auto_arima import numpy as np def arima_forecast(series, n_steps): model auto_arima( series, start_p0, max_p5, start_q0, max_q5, d1, seasonalFalse, stepwiseTrue, traceFalse, error_actionignore, suppress_warningsTrue ) return np.asarray(model.predict(n_periodsn_steps))说明d1表示对非平稳序列做一阶差分。趋势项通常是非平稳的一阶差分后变平稳如果序列本身已经平稳d0更合适过差分会把长期趋势信息削弱。stepwiseTrue让auto_arima用贪心搜索速度比穷举快很多。error_actionignore保证某些阶数组合失败时不中断程序。低频序列的样本量往往不大p和q超过5后模型参数过多容易过拟合因此max_p和max_q锁在5。最终重构是所有预测分量求和final_forecast np.zeros(len(test_data)) for pred in all_gru_predictions: final_forecast pred final_forecast arima_prediction # 如果之前做过归一化这里必须反归一化 # final_forecast final_forecast * std mean说明因为CEEMDAN和VMD都是线性重构预测阶段的分量求和就是最终预测。注意归一化的反转必须在误差计算之前完成否则RMSE和MAPE都失去意义。这里还有个细节被判定为噪声而忽略的分量在重构时不要加它的预测值因为它的真实值贡献本来就接近0加上去只会把噪声误差带进最终结果。5. 避坑指南分解、优化、预测三个阶段最容易翻车的地方5.1 CEEMDAN 分解末尾出现震荡现象分解得到的IMF在序列两端出现大幅摆动重构后两端误差明显大于中间。原因CEEMDAN处理有限长序列时边界处缺乏数据包络拟合不准确导致IMF在两端发散。trials设置过小或epsilon过大会放大这种边界效应。解决最简单的办法是把两端各截掉5%的数据再分解但预测时需要还原对应位置。更稳妥的做法是在分解前做镜像延拓把原始序列首尾各扩展一段对称数据分解后截掉延拓部分。此外把epsilon从0.005降到0.002trials提升到150边界震荡会明显减轻。这个坑很隐蔽因为分解的中间段看起来完全正常只有把重构误差按时间位置画出来才看得到两端翘起。5.2 ISOS-VMD 每次跑出的参数都不一样现象同一份训练集同一套ISOS代码两次运行得到的最优K和alpha不同有时K差两个数。原因ISOS是随机优化器种群初始化、t分布扰动都是随机的。更麻烦的是包络熵目标函数在部分参数区域比较平坦几个不同的K和alpha组合包络熵差值很小优化器随机性占主导。解决第一固定随机种子numpy的seed和python的random都要固定第二把K的搜索范围缩到合理区间比如已知数据周期性强就设K为4到7不要从1到15第三在目标函数中加过分解惩罚比如K每超过8惩罚项为0.05乘包络熵让优化器对“多拆出虚假模态”有明确排斥。如果条件允许同一组参数跑三次取目标函数最小的结果作为最终配置比一味加大迭代次数更有效。5.3 GRU 预测曲线整体滞后一个时间步现象在测试集上画预测曲线发现预测值比真实值晚了一拍曲线整体右移但MAPE数值居然不高。原因递归预测时输入窗口里最后一个值其实是上一步的预测值模型倾向于输出跟输入相近的值导致预测序列比真实序列延迟。训练时窗口全部是真实值预测时窗口逐步变成预测值输入分布不一致模型没见过这种输入。解决一是训练阶段加入scheduled sampling以一定概率用预测值替代窗口末尾的真实值再训练二是在recursive_forecast里对预测值做一次轻量校正比如用ARIMA对相邻两步预测残差建模并补偿三是直接改成多步输出把GRU的Dense(1)换成Dense(n_steps)一次预测未来多步而不是逐步递归。我一般优先做scheduled sampling因为它不改变模型结构只改训练数据构造方式。5.4 auto_arima 报错或 d 阶数选得离谱现象auto_arima运行时报“Data must be stationary”或者给出d2/3的高阶差分预测结果是一条水平线。原因趋势项序列不是平稳序列auto_arima的stepwise搜索可能找不到稳定解。d过高会过度差分把长期趋势信息全部抹掉预测自然退化成均值。解决在喂给ARIMA之前先做ADF检验确定d值from statsmodels.tsa.stattools import adfuller adf_stat, p_value, *_ adfuller(series) # p_value 0.05 表示序列平稳d 取 0否则 d 取 1如果ADF检验得到的p值大于0.05序列一阶差分后再测一次通常d1就够。不要轻易让auto_arima自动选d直接传d1或d0让它只搜索p和q。这条经验避免了大量报错。5.5 重构后的误差反而比单个模型更大现象把GRU和ARIMA的预测叠加后RMSE比直接用GRU预测原始序列还高。原因分量预测的误差之间存在正相关尤其高频子模态之间特征相似GRU对它们的预测误差会叠加。另一个原因是某个高方差子模态预测崩了把整体重构结果带崩。解决对每个分量预测结果按方差贡献做加权而不是简单求和。虽然理论上线性分解重构应该等权相加但那是分解与重构都精确时成立的预测阶段每个分量都有误差权重应该与分量的可预测性相关。简单的做法是用验证集上每个分量的预测误差倒数做权重更工程的做法是加一个线性回归层把所有分量预测作为特征拟合最终预测权重由回归学习得到。这样重构结果通常不会比单模型更差。6. 验证与封装让每个模块都经得起消融6.1 消融实验证明每一步都不是白加的拿到完整源码后先不要急着用全部模块跑结果按顺序做消融ARIMA单独跑、GRU单独跑、CEEMDAN-GRU、CEEMDAN-VMD-GRU、最后是完整链路。模型组合RMSEMAPE (%)单一 ARIMA----单一 GRU----CEEMDAN GRU----CEEMDAN VMD GRU----CEEMDAN ISOS-VMD GRU ARIMA----表格里的指标位留给自己跑完实验后填重点是每次改动只加一个模块。如果发现加了VMD后误差没下降先看ISOS搜索到的K是不是过分解如果加了ARIMA后误差反而上升大概率是趋势分量和GRU分量之间有重叠分配策略需要调整。消融实验是判断“这个模块到底有没有用”的唯一标准否则整条链路就是一个无法解释的黑匣子。6.2 完整 Pipeline 的封装要点这套流程涉及分解、优化、训练、预测、重构五个阶段不封装的话代码很快失控。我习惯用一个配置字典驱动全流程把重要参数集中管理。关键设计是CEEMDAN的分解结果和ISOS-VMD的最优参数要缓存到本地文件npz和json因为分解和参数寻优是最耗时的环节重新跑数据时不应该重复计算。GRU训练结果保存成模型文件预测时直接加载。接口只暴露fit和predict两个方法内部流程对调用方不可见。封装时另一个容易忽视的问题是随机种子管理。ISOS、GRU初始化都要固定随机种子否则每次运行结果都不同没法复现论文里的指标。我会把所有seed放进配置字典里每跑一次实验自动记录当次的配置和指标这样调参失败时还能回溯是哪一步出了问题。我现在的习惯是拿到任何一条时间序列先做分解再谈模型这个动作已经变成默认流程。希望帮到你。本文还有配套的精品资源点击获取
返回列表