简介一份融合时间卷积网络TCN、双向门控循环单元BiGRU与注意力机制的多变量时间序列预测Matlab完整源码面向自动化、电气、计算机等专业高年级本科生及研究生适用于课程设计、期末大作业和毕业设计等需要完成预测建模与对比分析的场景。压缩包共八个文件其中六个m脚本分别承担数据整理、模型构建、误差计算和雷达图展示等功能一个mat文件保存训练后的网络结构一个xlsx为风电场运行特征数据集整体大小仅4.26MB轻量易用且便于修改调试。已有八十八人学习下载运行环境为Matlab2023b主程序一键执行后会在命令窗口输出R2、MSE、MAE、MAPE四项回归指标同时生成可视化图表便于快速评估模型效果。代码采用参数化编程模型超参数集中可改注释明细、思路清晰数据结构已按多输入单输出格式整理读者可直接替换自有数据并调整历史步长与隐藏单元数量快速完成不同变量组合下的对比实验也可作为毕业设计中的基线模型进一步扩展。1. 把 TCN 塞进 BiGRU 之前先想清楚你要解决什么问题多变量时间序列预测的麻烦从来不在“模型不够深”而在“趋势、周期、突变混在一起时同一个模型既要盯得住长期依赖又不能在局部抖动上反应过度”。纯 BiGRU 擅长捕捉双向上下文但序列一长GRU 的隐状态会逐渐“忘记”较早的信息纯 TCN 靠膨胀卷积把感受野撑大却缺乏对关键时间步的选择能力。TCN-BiGRU-Attention 这个组合之所以在负荷预测、股价、交通流、设备剩余寿命这类场景里被反复使用是因为它把三件事拆给了三个模块TCN 先用不同膨胀率的卷积核把多变量输入转成高维特征BiGRU 再沿着时间维正向反向各扫一遍捕捉依赖最后 Attention 决定“到底哪几个时间步的哪几个特征对预测结果影响最大”。对新手来说最大的误区是拿到源码就改参数跑数据却不理解每一层张量形状的变化对熟手而言真正的价值在于搞清楚这套结构在 Matlab 里如何高效实现、验证集上哪些指标能说明“Attention 确实在工作”。本文按“原理 → 数据准备 → Matlab 实现 → 训练与调参 → 验证技巧”的顺序展开所有代码均可直接复制运行。2. TCN、BiGRU、Attention 在时序预测里的角色划分与 Matlab 选型理由2.1 为什么先用 TCN 做“特征抽取器”而不是直接堆 BiGRUTCN 的结构核心是因果膨胀卷积causal dilated convolution。因果意味着 t 时刻的输出只依赖 t 及之前的输入不会“看到未来”这是时间序列预测的硬约束膨胀率dilation rate让卷积核在不增加参数量的情况下扩大感受野。例如两层 kernel_size3、dilation[1,2] 的 TCN感受野是 1(3-1)*(12)7 个时间步而参数量只有两个 3×1 卷积核。Matlab 用dlconv实现这个操作非常直接不需要像 Python 那样引入额外的包。% TCN 单层实现膨胀因果卷积 dilation 2; % 膨胀率控制感受野扩张速度 filterSize 3; % 卷积核长度 numFilters 64; % 输出通道数 % 在 dlarray 上做因果卷积裁剪掉右侧多余输出以保持时间对齐 convOut dlconv(X, W, B, Stride, 1, DilationFactor, dilation, Padding, causal);Padding, causal是 Matlab 2019b 及以后版本支持的关键参数它自动在序列左侧补零保证因果性。如果你用的是 2018a 或更老版本就得手动在序列前面补(filterSize-1)*dilation个零再裁剪掉末尾等量数据。TCN 输出张量形状为[numFilters, seqLen, numObservations]在喂给 BiGRU 前需要用permute调整维度顺序以符合lstmLayer/gruLayer的输入约定。2.2 BiGRU 的双向扫描如何补上 TCN 的“上下文盲区”TCN 的感受野是单向扩张的虽然大但对“前文已经出现、后文才显形”的模式不敏感。BiGRU 在这里的作用不是再抽一遍特征而是把 TCN 输出的每个时间步嵌入放到正反两个 GRU 里各自跑一遍得到两个隐藏状态并拼接最终形成带上下文信息的特征序列。在 Matlab 中用bilstmLayer或手动封装两个gruLayer都可以但如果你需要显式拿到前向和后向的隐藏状态分别做 Attention 拼接建议手动拆% 手动搭建 BiGRU前向、反向各一层 forwardGRU gruLayer(32, OutputMode, sequence); backwardGRU gruLayer(32, OutputMode, sequence); % 注意反向输入需要先 flip 时间维经过 gruLayer 后再 flip 回来 X_reverse flip(X, 2); [Y_forward, ~] predict(forwardGRU, X); [Y_backward, ~] predict(backwardGRU, X_reverse); Y_backward flip(Y_backward, 2); Y_bigru cat(1, Y_forward, Y_backward); % 拼接后形状 [64, seqLen, N]cat(1, ...)在特征维拼接是 Attention 的输入准备。这里隐藏单元选 32拼接后是 64与 TCN 输出通道数对齐后续 Attention 打分时维度一致更好调参。一个值得注意的细节gruLayer的OutputMode必须设置为sequence因为 Attention 需要每个时间步的隐藏状态如果设成last就只剩最后一个时间步的输出Attention 无从计算。2.3 Attention 的权重到底关注的是什么Attention 层接收 BiGRU 的输出矩阵 H ∈ R^{64×T}T 为时间步数通过一个可学习的打分函数计算每个时间步的权重然后按权重加权求和得到上下文向量 c% 加性 Attention 打分与加权求和 % H: [64, T, N] attentionScore tanh(FC1(H)); % FC1 将 64 维打分为 32 维 attentionWeight softmax(FC2(attentionScore), 1); % FC2 打分为 1 维按时间维归一化 % 通过点乘广播并相加得到上下文向量 context sum(H .* attentionWeight, 2); % 形状 [64, 1, N]这里softmax(..., 1)的维度 1 是时间维因为 dlarray 的默认维度排列是[C, T, N]通道、时间、观测。打分函数选择加性 Attention 而不是点积 Attention是因为 BiGRU 的隐藏状态是高维向量点积容易导致梯度消失加性 Attention 多一层全连接做非线性变换效果更稳定。观察attentionWeight的分布是一件很有意思的事——训练好的模型通常会对“突变点”或“周期性拐点”施加更高权重这也是在论文里画注意力热力图的依据。Matlab 里可以用heatmap(squeeze(extractdata(attentionWeight)))直接可视化。3. Matlab 完整源码的数据准备从 CSV 到归一化与滑动窗口3.1 多变量数据集的导入与缺失值处理原始数据通常是 CSV 格式每一列是一个变量每一行是一个时间点。用readtable导入后将表格转为数值矩阵并检查是否存在NaN。多变量序列最忌讳的是逐列单独插值——这会破坏变量间的相关性我一般用fillmissing的movmedian方法它用一个滑动窗口的中位数填充能保留局部趋势% 导入多变量时间序列 dataTable readtable(multivariate_series.csv); data table2array(dataTable); % [T, numFeatures] % 滑动窗口中位数填充缺失值 filledData fillmissing(data, movmedian, 15);窗口大小 15 表示取当前点前后 7 个点共 15 个值的中位数。之后做异常值替换一个实用技巧是用 3σ 准则某列的值超出均值±3倍标准差时用该列上下两个正常采样点的均值替换。这一操作必须逐列做且要记录替换位置——如果测试集里也出现同类异常要用训练集统计出的均值和标准差处理而不是重新计算否则会造成信息泄漏。3.2 归一化的“Fit 到训练集Transform 到全部”原则归一化方法选择zscore零均值单位方差还是mapminmax映射到 [0,1]取决于后续激活函数。Attention 里用到了softmax它的输入绝对值大小不影响归一化结果但 GRU 的tanh激活对输入范围敏感所以对 BiGRU 的输入用 zscore 更合适。关键在实现顺序先用训练集拟合mu和sigma再用同一参数转换验证集和测试集——这也是热词检索里“bilstm代码matlab soc”这类提问最常见的坑。% 按特征维度拟合训练集统计量 mu mean(trainData, 1); sigma std(trainData, 0, 1); sigma(sigma 0) 1; % 防止某列全常数导致除零 % 统一转换三组数据 trainNorm (trainData - mu) ./ sigma; valNorm (valData - mu) ./ sigma; testNorm (testData - mu) ./ sigma;注意mean(trainData, 1)的维度参数 1表示按每列求均值得到形状[1, numFeatures]的行向量。如果省略1mean默认对整个矩阵求标量归一化就错了。3.3 滑动窗口切样本时的顺序陷阱滑动窗口的典型做法是用过去lookback个时间步预测未来horizon个时间步。假设lookback24、horizon6样本 i 的输入是data[i:i23, :]标签是data[i24:i29, 目标列]。一个容易被忽略的问题相邻样本高度重叠如果直接划分训练/验证集会造成验证集泄漏到训练集。常见做法是按时间顺序切——前 70% 时间段的样本归训练集中间 15% 归验证集最后 15% 归测试集而不是随机打乱。function [XTrain, YTrain] createSlidingWindows(data, targetCol, lookback, horizon) numSamples size(data, 1) - lookback - horizon 1; X zeros(lookback, size(data, 2), numSamples); Y zeros(numSamples, 1); for i 1:numSamples X(:, :, i) data(i:ilookback-1, :); Y(i, 1) data(ilookbackhorizon-1, targetCol); end % 转为 dlarray通道维在第二维 XTrain dlarray(X, TCB);我在dlarray(X, TCB)里用了TCB的维度标签——T 是时间C 是通道B 是批次。这是 Matlab 深度学习自定义训练循环最常用的排列方式与dlconv、gruLayer的接口都匹配。切片层数size(data, 2)是所有特征列数取目标列标签时用的是data(...)索引而不是Y X(:, targetCol, i)——后者拿的是输入窗口内的值而不是未来值这个错误在初学代码里出现频率极高。4. TCN-BiGRU-Attention 模型的 Matlab 实现与训练循环4.1 网络结构定义的完整代码使用dlnetwork定义自定义网络。整体流程TCN 特征抽取 → BiGRU 上下文编码 → Attention 加权求和 → 全连接输出。% 定义 TCN-BiGRU-Attention 网络Matlab R2021a inputSize size(XTrain, 1); % 时间步数 lookback numFeatures size(XTrain, 2); % 输入变量数 numHidden 32; numOutputs 1; % 初始化权重 wConv1 initializeGlorot([3, numFeatures, 64]); % [filterSize, inChannels, outChannels] bConv1 zeros(64, 1, single); wGRU_F initializeGlorot([numHidden, 64numHidden]); % GRU 权重矩阵 bGRU_F zeros(numHidden, 1, single); % 构建 dlnetwork 结构 layers [ featureInputLayer(numFeatures, Normalization, none) convolution1dLayer(3, 64, Padding, causal, DilationFactor, 1, Name, tcn_1) reluLayer convolution1dLayer(3, 64, Padding, causal, DilationFactor, 2, Name, tcn_2) reluLayer % 这里不是标准层序列需要自定义训练循环 ];这里有个实现分歧点Matlab 的convolution1dLayer在 2020a 之后才支持DilationFactor参数且dlnetwork对自定义的 GRU 加 Attention 结构支持并不友好。更稳妥的路线是用自定义训练循环即把所有模块写成函数前向传播走自定义modelPredictions梯度用dlgradient自动推导% 自定义前向传播函数 function [Y, context, attWeight] modelPredictions(parameters, X) % TCN 两层膨胀卷积 Z dlconv(X, parameters.conv1.W, parameters.conv1.B, Padding, causal, DilationFactor, 1); Z relu(Z); Z dlconv(Z, parameters.conv2.W, parameters.conv2.B, Padding, causal, DilationFactor, 2); Z relu(Z); % BiGRU手动定义 GRU 单元循环前向传播 [Y_f, ~] gruForward(parameters.gruF, Z); Z_flip flip(Z, 1); [Y_b, ~] gruForward(parameters.gruB, Z_flip); Y_b flip(Y_b, 1); H cat(2, Y_f, Y_b); % [T, 2*numHidden, N] % Attention attScore tanh(dlconv(H, parameters.attW1, parameters.attB1, Padding, same)); attLogit dlconv(attScore, parameters.attW2, parameters.attB2, Padding, same); attWeight softmax(attLogit, 1); context sum(H .* attWeight, 1); % [1, 2*numHidden, N] Y fullyconnect(context, parameters.fcW, parameters.fcB); end4.2 gruForward 的实现细节Matlab 2019b 之后提供了一个较隐晦的接口gruLayer但要在自定义循环里拿到逐时间步的隐状态最常见的做法还是自己实现门控计算。GRU 的前向公式是标准的重置门 r、更新门 z、候选隐状态 h~function [H, hT] gruForward(params, X) % X: [T, C, N] [T, C, N] size(X); h zeros(params.numHidden, 1, N); H zeros(params.numHidden, T, N); for t 1:T xt X(t, :, :); r sigmoid(params.Wr * xt params.Ur * h params.br); z sigmoid(params.Wz * xt params.Uz * h params.bz); hCandidate tanh(params.Wh * xt params.Uh * (r .* h) params.bh); h (1 - z) .* h z .* hCandidate; H(:, t, :) h; end hT h; endparams.Wr、params.Ur等六组权重矩阵构成了 GRU 的全部可学习参数。矩阵乘法在循环里执行T 个时间步跑 T 次循环这是纯 Matlab 实现的性能瓶颈如果样本数很多、lookback 很长可以考虑用dlarray的C维做小批量并行或者直接改用gruLayer后通过dlnetwork提取隐藏状态。但自定义实现的好处是能显式拿到 h方便 Attention 接入——官方层封装后想取中间隐状态要用activations或forward的 Outputs 选项不如自己来得透明。4.3 训练循环、损失函数与学习率调度训练采用 Adam 优化器损失函数用均方误差MSE。多变量预测里如果只预测一个目标变量MSE 足够如果预测多个变量需要对每个输出通道做加权 MSE权重可以按变量量纲反比设置。学习率调度使用余弦退火epoch 从 0 到 maxEpochs学习率从 lr0 衰减到 lr0×0.01% 自定义训练循环 numEpochs 100; learnRate 0.001; trailingAvg []; trailingAvgSq []; gradDecay 0.9; gradDecaySq 0.999; for epoch 1:numEpochs % 按小批量打乱数据 shuffledIdx randperm(size(XTrain, 3)); for i 1:numBatches idx shuffledIdx((i-1)*batchSize1 : min(i*batchSize, end)); XBatch XTrain(:, :, idx); YBatch YTrain(:, :, idx); [loss, grads] dlfeval(modelLoss, parameters, XBatch, YBatch); % 余弦退火学习率 lr learnRate * 0.5 * (1 cos(pi * (epoch-1) / numEpochs)); [parameters, trailingAvg, trailingAvgSq] adamupdate(parameters, grads, ... trailingAvg, trailingAvgSq, epoch, lr, gradDecay, gradDecaySq); end % 验证损失记录用于早停 valLoss computeLoss(parameters, XVal, YVal); if valLoss bestValLoss bestValLoss valLoss; bestParams parameters; end enddlfeval和dlgradient是 Matlab 自动微分的核心接口modelLoss里要先跑前向传播得到预测值再对预测值和真实值算mse最后调dlgradient(loss, parameters)。注意dlgradient只能算标量对结构体/元胞数组的梯度如果 parameters 是包含多个字段的 struct梯度也是同构 struct。批量大小batchSize选 64 或 128——Attention 的 softmax 对 batch 内所有样本独立归一化批量稍微大点不影响收敛方向但显存内存占用会线性增长。早停的 patience 设 10 个 epoch即连续 10 个 epoch 验证损失不下降就恢复最优参数并停止。5. 预测结果评估不只是 RMSE还要看残差和时延5.1 反归一化与多步预测的累积误差控制模型输出是归一化后的值评估前必须用训练集的mu(targetCol)和sigma(targetCol)反变换。多步预测horizon1有两种做法直接多输出模型一次输出 horizon 个值和递归预测用上一步的预测值作为下一步输入。直接多输出的误差更可控因为每一步都是独立输出分支不累积递归预测会把第一小步的误差逐步放大但只要你的训练数据本身噪声不小递归预测的“真实轨迹贴合度”反而更高。下面的代码实现反归一化% 反归一化预测值 predRaw extractdata(YPred); predDenorm predRaw .* sigma(targetCol) mu(targetCol); trueDenorm extractdata(YTest) .* sigma(targetCol) mu(targetCol); % 计算 RMSE、MAE、MAPE rmse sqrt(mean((predDenorm - trueDenorm).^2, all)); mae mean(abs(predDenorm - trueDenorm), all); mape mean(abs((predDenorm - trueDenorm) ./ trueDenorm), all) * 100;MAPE 在真实值接近零时会爆炸如果你的目标变量含有零值建议改用 SMAPE对称平均绝对百分比误差公式是2*|pred-true| / (|pred||true|)。5.2 与 BiGRU、TCN 单模型对比的消融实验模板判断 Attention 是否有效的最直接方法是跑消融TCN-BiGRU-Attention vs. TCN-BiGRU去掉 Attention直接用最后时间步的隐状态接全连接vs. 纯 BiGRU。在同一个数据集、同一个窗口和归一化方案下如果加了 Attention 的模型在 RMSE 上没有明显下降问题通常出在 Attention 打分没有收敛——权重全部变成均匀分布等价于对所有时间步取平均。验证 Attention 是否真的“学会了”关注关键步观察训练过程中attWeight的熵% 计算注意力权重的信息熵诊断是否退化 attEntropy -sum(attWeight .* log(attWeight 1e-12), 1); % 如果熵接近 log(T)说明接近均匀分布Attention 失效 fprintf(Attention entropy: %.4f (max possible: %.4f)\n, ... mean(extractdata(attEntropy), all), log(size(XTrain, 1)));如果熵值始终接近理论最大值常见原因有三个Attention 层的全连接权重初始值过小梯度更新慢学习率偏大导致权重在 tanh 饱和区震荡或者 BiGRU 的隐藏状态本身就没有区分度——可以把numHidden从 32 加到 64增加特征表达能力。5.3 预测时延检验一个容易被忽略的定性指标定量指标之外用交叉相关函数检查预测序列与真实序列的时延偏移。如果 RMSE 很低但预测曲线整体滞后于真实值半步说明模型在“抄近道”——学着把上一个时间步的值搬过来而不是学动态规律。检验方法是用finddelay计算两组序列的滞后% 计算预测与真实的滞后步数 d finddelay(trueDenorm, predDenorm); fprintf(Lag between prediction and ground truth: %d steps\n, d);滞后步数为 0 是理想情况若滞后为负预测超前或正预测滞后且绝对值超过 2 步就要警惕。此时可以尝试缩小lookback或者把horizon的目标从“第 h 步的值”改为“第 1 到 h 步的累计值”让模型学会预测增量而不是绝对值。滚动预测的多步误差是逐步放大的所以报告指标时请分别给出 horizon1、3、6 的结果并说明每一步用了哪些输入——是递归预测还是直接多输出这个细节直接决定了读者能否复现你的数字。本文还有配套的精品资源点击获取
