混沌时间序列预测：基于时空RBF神经网络的MATLAB实现

简介：这是一份面向本科、硕士阶段教研学习的Matlab实现资料，聚焦时空RBF神经网络（时空RBF-NN）在混沌时间序列预测中的应用。压缩包内共10个文件，包含3个.m脚本、3个.mat数据集和4个png结果图，整体仅1.32MB，轻量易用；其中脚本涵盖标准RBF与时空RBF两种模型的Mackey-Glass时间序列预测实现，数据集用于直接加载实验数据，结果图分别对应训练误差、测试误差与预测对比，便于逐项分析模型性能。目前已有281人学习下载，适合正在开展神经网络预测、混沌时间序列方向课题研究或课程实验的本科与硕士读者。借助该资源可快速复现时空RBF-NN的完整预测流程，节省代码调试与数据准备时间；同时项目结构简洁、注释清晰，为后续修改网络结构、替换其他混沌数据集或对比不同预测算法提供了扎实的Matlab代码基础。

1. 混沌时间序列预测：为什么用时空 RBF-NN

Mackey-Glass 方程生成的序列，是评估非线性预测模型的经典基准。它的特点在于确定性的微分方程却能呈现非常接近随机的长期状态，一步预测容易做，多步预测时误差会随混沌发散迅速放大。用普通 RBF 网络可以捕捉短期关联，但面对混沌序列的时间延迟依赖，隐含层的静态映射很容易丢掉历史状态。时空 RBF-NN（ST-RBF）把时间延迟作为输入通道重新组织，让网络在学习下一时刻的值时显式看到前面几个时刻的状态向量，这是它在多步预测上更稳的关键。下面从 Mackey-Glass 序列构造开始，逐步拆解普通 RBF 与时空 RBF-NN 的 MATLAB 实现，适合准备做混沌预测、神经网络预测的本科和研究生在 MATLAB 2014/2019a 上复现。

2. Mackey-Glass 序列的构造与数据预处理

2.1 从微分方程到离散时间序列

Mackey-Glass 方程的形式是：

dx/dt = β·x(t−τ) / (1 + x(t−τ)^n) − γ·x(t)

当n > 17且τ足够大时，系统进入混沌区。工程里最常用的参数组合是β=0.2、γ=0.1、n=10、τ=17，生成的序列既混沌又不至于发散太快。MATLAB 里可以直接用dde23求解时滞微分方程，也可以像下面这样用欧拉法递推，后者更容易嵌入到后续预测脚本：

% Mackey_Glass_Time_Series_Prediction_RBF.m 中常见的数据生成段 dt = 0.01; tau = 17; x = zeros(1200, 1); x(1) = 1.2; history = ones(tau, 1) * 1.2; % 历史状态缓存 for i = 2:length(x) % dx/dt = 0.2*history(end)/(1+history(end)^10) - 0.1*x(i-1) x(i) = x(i-1) + dt * (0.2*history(end)/(1+history(end)^10) - 0.1*x(i-1)); history = [history(2:end); x(i-1)]; end x = x(200:end); % 去掉前 200 个瞬态点，只剩稳定混沌段 x = x(:); % 统一为列向量

这段代码用history保存最近tau个状态，history(end)就是当前步的时滞项；循环内先计算欧拉增量再更新历史数组。需要注意x(1)是标量，而history是tau维向量，二者类型要一致。去掉前 200 点的原因是延迟微分方程初值附近有一段瞬态，不剔除会影响后续构造样本的质量。资源包里的Data.mat保存的就是这样一个已经清洗好的序列，你不需要每次重新生成，但理解它来自哪里，后面调窗口长度时才知道怎么对应。

2.2 用时间延迟重构样本，让 RBF 看到历史

RBF 网络做的是输入到输出的映射，如果你把某个时刻的值x(t)直接作为输入去预测x(t+1)，它学到的几乎是恒等映射，对混沌序列没有任何泛化能力。常用做法是把过去m个时刻的值组装成输入向量：

[x(t-m+1), x(t-m+2), ..., x(t)] -> x(t+1)

这里的m就是嵌入维数。对 Mackey-Glass 序列，m=18到m=30是比较合理的区间，取太小会丢信息，取太大则 RBF 的隐层节点数会快速膨胀。构造样本时，先把数据切成训练集和测试集，再滑动生成矩阵：

% 构造样本：X 每行是一个长度为 m 的历史窗口，Y 是下一时刻的目标值 m = 24; N = length(x); X = zeros(N-m, m); Y = zeros(N-m, 1); for k = 1:N-m X(k, :) = x(k:k+m-1); Y(k) = x(k+m); end % 划分：前 70% 训练，后 30% 测试 trainNum = round(0.7 * (N-m)); Xt = X(1:trainNum, :); Yt = Y(1:trainNum); Xv = X(trainNum+1:end, :); Yv = Y(trainNum+1:end);

这里X的每一行是一个长度为 24 的滑动窗口，滑动步长为 1，所以相邻窗口之间高度相关。很多人在这一步直接用X(k, :) = x(k:k+m-1)就完了，却忽略了一句：如果序列没有先转成列向量，x(k:k+m-1)拿回来的是行向量还是列向量取决于x的形状。代码里先x = x(:)就是为了避免这个隐晦 bug。数据划分也很关键，窗口滑过的时间轴有先后顺序，不能随机打乱，否则会把未来的信息泄漏进训练集，测试结果会虚高。资源里的Time_SeriesTesting.png和Time_SeriesTraining.png用的就是这种按时间先后划分的样本组织方式。

2.3 归一化与范围选择

RBF 神经网络里，隐含层的高斯核函数exp(-||x-c||^2 / (2*spread^2))对输入尺度非常敏感。由于 Mackey-Glass 序列跑起来后幅值基本落在 0.4 到 1.4 之间，如果不归一化，协方差矩阵的距离计算会被数值较大的维度主导，网络中心也会偏向幅值大的区域。我一般会先做 min-max 归一化，把数据压到[0.1, 0.9]：

% min-max 归一化到 [0.1, 0.9] mn = min(x); mx = max(x); x_norm = (x - mn) / (mx - mn) * 0.8 + 0.1;

归一化之后训练出的网络权重，在预测时要先把预测结果反归一化回原始量纲：

x_pred = (y_pred - 0.1) / 0.8 * (mx - mn) + mn

这里把上下界设为 0.1 和 0.9 而不是 0 和 1，是为了避开高斯函数输出在边界处过于接近 0 的问题。如果输入里恰好有训练集之外的极值，比如测试段出现新的更大峰值，min-max 归一化会因为固定了mx而产生轻微偏移。更稳妥的做法是用训练集的mn和mx来归一化训练、测试和后续预测，不要让测试数据的信息参与归一化参数计算。

表 2-1 列出 Mackey-Glass 序列建模时最常见的参数组合：

这套参数可以保证你在对比普通 RBF 和时空 RBF-NN 时，至少不是数据集构造上的差异，而是模型结构本身的差异。

3. 普通 RBF 网络预测：从原理到 MATLAB 实现

3.1 RBF 网络结构与中心选择

RBF（径向基函数）网络是一种三层前馈结构：输入层把向量原样送入，隐含层每个节点对应一个中心向量c_i，用径向基函数计算输入到中心的距离，输出层把各个基函数的响应做加权线性组合。最常用的是高斯核：

φ(x, c_i) = exp(-||x - c_i||^2 / (2*sigma_i^2))

这里的sigma_i在 MATLAB 的newrb中对应 spread 参数。如果 spread 太小，每个中心只管到非常近的点，网络会在输入空间留下大量空隙；如果 spread 太大，所有基函数都重叠，预测会趋于平均，混沌序列中的波峰波谷会被抹平。常见做法是先固定 spread=1.0，用newrb的误差目标参数自动增加中心数。

由于 Mackey-Glass 序列的轨道在相空间里位于一个带状流形上，不是全部输入空间均匀分布，所以 RBF 中心的选择更看重样本覆盖度而不是数量。newrb的实现方式是反复从训练集里选点作为中心，每轮添加一个中心，直到满足均方误差目标或达到最大节点数。代码里RBF.mat保存的就是这样训练好的一个net。

3.2 用 newrb 训练普通 RBF 预测器

newrb是最直接的调用方式，不需要手写梯度下降：

% 普通 RBF 网络训练 rng(0); % 固定随机种子，便于复现 spread = 1.0; % 基函数扩展常数 goal_mse = 0.001; % 期望均方误差 max_neurons = 100; % 隐含层最大神经元数 net_rbf = newrb(Xt', Yt', goal_mse, spread, max_neurons); % 测试集预测 Yp = sim(net_rbf, Xv'); % 反归一化并计算误差 Yp_denorm = (Yp' * 0.8 + 0.1) .* (mx - mn) + mn; Yv_denorm = (Yv * 0.8 + 0.1) .* (mx - mn) + mn; test_mse = mean((Yp_denorm - Yv_denorm).^2); fprintf('RBF 测试 MSE: %.4f\n', test_mse);

这里Xt和Yt是第 2 章构造出的训练样本，newrb期望输入是矩阵，但每一列是一个样本，所以要用转置Xt'。sim函数返回的Yp是一个行向量，和Yv的形状不一致，写代码时最容易在这一步出错。反归一化时，Yp是网络的原始输出，必须按第 2.3 节的反变换公式还原；如果你直接拿归一化的误差去和原始数据对比，得到的 MSE 会看起来很小，但没有任何现实意义。

运行完这段代码后，资源包里的Time_SeriesTraining.png就是普通 RBF 在训练集上的预测曲线叠加，Time_SeriesTrainingMSE.png是训练过程的迭代 MSE。你在自己的机器上复现时，可能看到newrb的训练速度很快——这是因为 RBF 网络只需要解一个线性方程组得到输出层权重，不存在 BP 的反复迭代。

newrb常用参数如表 3-1 所示：

3.3 隐含层数量失控时的处理

newrb的max_neurons设得太大时，训练误差可能降到机器精度，但测试误差反而上升。我一般会把goal_mse设成 0.005 而不是 0.001，并观察max_neurons是否被触顶。如果训练结束时节点数已经到了max_neurons，说明 spread 选得太大或goal_mse过于苛刻，需要先调 spread。

另一个问题是spread对输出层权重的条件数影响很大。高斯核响应接近 0 的节点，其输出权重会在最小二乘解里被放得很大，导致测试时轻微扰动就引起剧烈波动。常见做法是优先检查sim(net_rbf, Xt')在训练集上的误差，如果训练集上的预测和目标值之间有规律性偏移，先不要怀疑网络结构，回到数据归一化那段检查训练集和测试集是否用了各自的mn/mx。

4. 时空 RBF-NN 的改进结构与训练流程

4.1 普通 RBF 在混沌序列上的瓶颈

普通 RBF 把整个历史窗口一次性映射到输出，隐含层并没有显式地表达时间关系。它学到的是“窗口内值的组合方式”，而不是“这些值随时间演化的规律”。Mackey-Glass 序列的混沌特性意味着相邻状态在相空间里的距离很接近，但演化方向可能完全不同。普通 RBF 对这类邻域的判别力不足，多步递归预测时误差会沿着轨道指数放大，所以单步预测不错，多步预测就发散了。

时空 RBF-NN 的基本思路是把时序信息拆成两个维度：空间状态（当前和邻近时刻的采样点）和时间演化（这些点之间的先后关系）。实现上，最直接的一种是做输入端的延迟增强：在原来的m维窗口之外，把窗口内的每个采样点再附上它的时间索引和一步差分。相当于把每个样本从一维序列展开成一张时间-空间状态表。

4.2 构造时空输入特征

这里以Mackey_Glass_Time_Series_Prediction_ST_RBF.m为例，常见做法是先给每个历史点附加相对时间戳，再构造差分项：

% 构造时空 RBF 输入：原始窗口 + 一阶差分 + 时间位置 m = 24; N = length(x_norm); % 一阶差分 dx = diff(x_norm); % 长度 N-1 % 对每个 k，取 x(k:k+m-1) 和 dx(k:k+m-2)，再补一段时间索引 len = N - m - 1; Xst = zeros(len, m + (m-1) + 1); Yst = zeros(len, 1); for k = 1:len win = x_norm(k:k+m-1); dwin = dx(k:k+m-2); time_idx = (0:m-1)' / m; % 时间位置归一化 Xst(k, :) = [win(1:end-1), dwin, win(end), time_idx(1)]; Yst(k) = x_norm(k+m); end

这段代码把输入从原来的m维扩展到(m-1) + (m-2) + 1 + 1 = 2m-1维，特点是保留了窗口内相邻点的差分信息，同时额外加了一个时间定位值。为什么这么做？因为混沌序列在时间轴上是一段反复折返的轨道，窗口形态相同不代表下一步相同，而差分项能区分上升段和下降段。time_idx相当于告诉网络当前的步长位置，避免不同相位但相同形态的窗口混在一起。

需要注意的是，不是所有“时空 RBF”都采用这种差分展开。文献中另一种常见做法是把 RBF 的每个中心c_i从空间向量扩展为[c_i, t_i]，在高斯距离中同时计算空间距离和时间距离。那样需要修改核函数，newrb不支持，得自己写核矩阵。资源里的ST_RBF.mat对应的是前一种扩展输入的做法，因为它仍然可以调用newrb，只是在输入矩阵上做了工程改造。

4.3 训练顺序与结果保存

时空 RBF-NN 的输入维度变高后，训练样本数量如果不够，高维空间会变得稀疏。所以训练集占比可以提高到 80%，同时把max_neurons从 100 降到 80，防止中心过多造成过拟合。训练代码和普通 RBF 几乎一样，区别只是在输入矩阵上：

% ST-RBF 训练 spread = 1.2; % 维度变高，适当增大 spread net_st = newrb(Xst(1:trainNum, :)', Yst(1:trainNum)', 0.005, spread, 80); % 保存结果，方便 Results_graphs.m 加载绘图 save('ST_RBF.mat', 'net_st', 'Xst', 'Yst', 'mn', 'mx', 'm');

这里把spread提高到 1.2，是因为扩展输入之后的欧几里得距离变大了，如果还用 1.0，高斯响应大部分接近 0，网络学不到东西。save保存的mn和mx是第 2.3 节归一化用的参数，Results_graphs.m里要反归一化时必须用到。如果不存归一化参数，生成的Time_SeriesTesting.png会因为量纲不一致而出现明显偏差。

表 4-1 给出普通 RBF 与时空 RBF-NN 在输入构成上的差异：

这个表可以帮你判断为什么Time_SeriesTestingMSE.png里 ST-RBF 的测试误差通常比普通 RBF 低：差分项在混沌序列的局部趋势上提供了更强的判别信息。

5. 训练验证与结果分析：从 MSE 曲线看模型行为

5.1 如何从 MAT 文件恢复并预测

资源里给出了RBF.mat和ST_RBF.mat，这两份文件分别保存了训练好的网络对象。你可以直接加载它们，不需要重新训练。但注意它们是在x_norm这种归一化数据上训练出来的，加载后预测也要先按同样的方式归一化，否则结果会完全偏离。典型的验证脚本是这样：

clear; close all; load('RBF.mat'); % net_rbf load('ST_RBF.mat'); % net_st load('Data.mat'); % 原始数据 x, 以及 mn, mx 等 % 构造测试窗口 m = 24; N = length(x); Xv = zeros(N-m, m); Yv = zeros(N-m, 1); for k = 1:N-m Xv(k, :) = x(k:k+m-1); Yv(k) = x(k+m); end % 普通 RBF 预测 yr = sim(net_rbf, Xv'); % 时空 RBF 预测需要重新构造差分输入，这里略去相同部分 yerr_r = abs(yr' - Yv); fprintf('普通 RBF 最大绝对误差: %.4f\n', max(yerr_r));

这里把Yv设为真实下一时刻值，yr'和Yv做逐点差，取绝对值后找最大值，能直观看到混沌序列在峰值附近的预测偏差。Data.mat里的x可能已经是原始数据，也可能已经归一化，你在加载后先执行plot(x)看一眼幅值范围，再决定要不要反归一化。

5.2 读 MSE 图时重点看什么

Time_SeriesTrainingMSE.png和Time_SeriesTestingMSE.png是两份不同数据上的误差曲线。不要只看最终 MSE 数值，要看曲线形态：

• 如果训练 MSE 单调下降并出现平台期，说明网络容量足够；
• 如果训练 MSE 已经很低，但测试 MSE 在某个步数后突然跳高，典型过拟合信号；
• 如果两条曲线都有周期性抖动，说明输入窗口没有对齐混沌序列的轨道周期，需要调嵌入维数。

资源里的Results_graphs.m脚本会把这些图重新绘制出来。我一般会把它改造成逐点误差折线图，这样比 MSE 曲线更能定位误差集中在序列的哪些区段。注意Time_SeriesTestingMSE.png呈现的只是 30% 测试段上的统计值，真正的混沌预测能力要看多步递归预测，而不是单步统计。

5.3 对比两种模型的验证方法

要公平对比 RBF 和 ST-RBF，必须保证训练/测试集完全相同。最简单的方式是先加载Data.mat，用第 2 章的窗口构造逻辑生成一份固定的Xv/Yv，然后分别用net_rbf和net_st预测。不要分两次重新滑窗，因为如果序列长度是奇数，第二次滑窗的偏移会让两边的测试样本错开一个相位，对比就失去意义了。

可以用下面这种方式快速得到对比表：

% 计算两种模型在相同测试集上的 MSE mse_rbf = mean((yr' - Yv).^2); mse_st = mean((y_st' - Yv).^2); ratio = mse_rbf / mse_st; fprintf('MSE 比值 RBF/ST-RBF: %.2f\n', ratio);

ratio大于 1 说明 ST-RBF 在整体误差上更小；但注意这个指标只适用于相同步长预测。如果你做的是递归多步预测，需要把预测值送回输入窗口，一步步迭代往前推，这时误差会指数放大，两种模型的差距也会拉大。资源里的单步预测图能看出短期好坏，真正的混沌预测能力要看多步迭代后的对比。

表 5-1 建议的记录指标：

这种表在论文里可以当实验设置，不用额外造图。

6. 调参与复现技巧：让预测结果更稳定

6.1 spread 与神经元数量的联动调整

复现这类项目最常见的翻车点，是把newrb的goal_mse设成 0.0001、max_neurons设成 500，结果训练集预测接近完美，测试集一团糟。RBF 网络的两个超参数是联动的：spread决定每个基函数的覆盖半径，max_neurons决定能容纳多少覆盖区域。如果输入维度是 24，spread=0.5 时每个中心只覆盖大约 2 个特征跨度的区域，整个输入空间可能需要几百个中心；spread=2.0 时中心数量明显减少，但边界模糊。我习惯先固定max_neurons=100，从 0.5 开始按 0.1 步长扫到 2.0，选测试 MSE 最低的点。

% 手动扫描 spread，而不是用默认值 for s = 0.5:0.1:2.0 net = newrb(Xt', Yt', 0.005, s, 100); Yp = sim(net, Xv'); mseList(round(s*10)) = mean((Yp' - Yv).^2); end [bestMSE, idx] = min(mseList); bestSpread = 0.5 + (idx-1)*0.1;

注意newrb每次训练都会重新选中心，训练结果带随机性。即使固定随机种子，不同 MATLAB 版本也可能给出不同中心。做一个十次重复实验，取平均 MSE 比只看单次结果更有说服力。

6.2 多步预测的验证技巧

混沌时间序列预测不能只做一步验证。把测试集上的单步预测结果当作输入继续迭代，代码上很简单：先预测t+1，把这个值接在窗口末尾、丢掉窗口最前面的值，再预测t+2。

% 递归多步预测示意 win = Xv(1, :); % 初始窗口 steps = 20; preds = zeros(steps, 1); for j = 1:steps yhat = sim(net_st, win'); preds(j) = yhat; win = [win(2:end), yhat]; % 窗口前移一位 end

如果 20 步之后预测序列和真实序列仍然保持同样的周期形态，说明模型学到了轨道结构；如果很快就衰减成一个常数，说明它只是记住了最近邻，没有学到动态。资源里的Time_SeriesTesting.png就是这种单步预测结果的叠加图，你可以自己把递归多步的线画上去，对比会非常直观。

6.3 用 Results_graphs.m 输出论文级图表

拿到项目后先别急着跑主程序，打开Results_graphs.m看它加载了哪些变量。我一般会做三处修改：一是把坐标轴字号调大到 11，二是用ylim固定纵轴范围，三是把普通 RBF 和 ST-RBF 的预测线画在同一个小图里。这样出的图，投稿或放进毕业论文时不用再重新画。

figure; plot(1:length(Yv), Yv, 'k-', 'LineWidth', 1); hold on; plot(1:length(Yv), yr, 'r--', 'LineWidth', 1); legend('真实值', 'RBF预测', 'Location', 'northwest'); xlabel('时间步'); ylabel('x(t)');

修改后注意保存路径，Results_graphs.m默认生成的图片会覆盖资源里已有的 PNG，建议先saveas(gcf, 'my_result.png')另存为别的文件名，再打开生成的 PNG 与初始图片做像素级对比，就可以判断复现代码是否落在同一组模型参数上。

本文还有配套的精品资源，点击获取