MATLAB实现基于季节性趋势分解（STL）进行时间序列预测
MATLAB 实现基于季节性趋势分解（STL, Seasonal and Trend decomposition using Loess）进行时间序列预测 的项目实例。该方法适用于具有明显季节性和趋势成分的时间序列（如电力负荷、气温、销售数据等），通过 STL 分解后对各成分分别建模，再组合实现预测。
🎯 项目目标
对时间序列进行 STL 分解：
( y_t = T_t + S_t + R_t )
其中：(T_t) = 趋势项，(S_t) = 季节项，(R_t) = 余项（残差）
对趋势项和余项分别建立预测模型（如 ARIMA、线性回归、MLP 等）
利用已知的季节模式（周期重复）外推季节项
重构预测结果：(\hat{y}{t+h} = \hat{T}{t+h} + \hat{S}{t+h} + \hat{R}{t+h})
⚠️ 注意：MATLAB 没有内置 STL 函数，但可通过 movmean + 自定义 Loess 或使用 替代方案（如 islocalmax + 滤波），或调用 Python 的 statsmodels（需配置引擎）。
为保持纯 MATLAB 实现，我们采用 近似 STL 的经典三步滤波法（Cleveland et al. 方法简化版）。
✅ 推荐方案：使用 MATLAB 官方 stl 替代函数（自定义实现）
由于 R2023a 起 MATLAB 仍未提供 stl，我们将使用一个高效、开源的 MATLAB STL 实现：
🔧 步骤 1：下载或定义 STL 函数
从 MATLAB File Exchange 获取可靠实现：
STL Decomposition by Danil (File ID: #82794)
或使用以下简化版本（适合教学/快速原型）
✅ 简化版 STL（基于滑动平均 + 季节周期对齐）
matlab
function [trend, seasonal, remainder] = simple_stl(y, period, nIter)
% SIMPLE_STL: 简化版 STL 分解（非 Loess，但保留核心思想）
% 输入：
% y - 时间序列 (列向量)
% period - 季节周期（如 12 表示月度数据年周期）
% nIter - 内部迭代次数（默认 2）
% 输出：
% trend, seasonal, remainder
if nargin < 3, nIter = 2; end
n = length(y);
seasonal = zeros(n, 1);
trend = zeros(n, 1);
% 初始化
detrended = y;
for iter = 1:nIter
% Step 1: 提取季节项（按周期分组取中位数/均值）
seasonal = zeros(n, 1);
for i = 1:period
idx = i:period:n;
if ~isempty(idx)
val = detrended(idx);
seasonal(idx) = mean(val); % 可改用 median 更鲁棒
end
end
% 去除季节项
deseason = y - seasonal;
% Step 2: 提取趋势项（用移动平均平滑）
win = min(2floor(period/2)+1, floor(n/5)); % 窗口大小
if mod(win,2)==0, win = win+1; end
trend = smoothdata(deseason, ‘movmean’, win);
% 更新去趋势序列用于下一轮季节提取
detrended = y - trend;
end
remainder = y - trend - seasonal;
end
💡 对于高精度需求，建议使用完整 Loess 版本（可提供）。
📊 Step 2：加载或生成时间序列数据
matlab
% 示例：合成带趋势+季节+噪声的数据
t = (1:200)';
season = 10 sin(2pit/12); % 周期=12（月度）
trend = 0.1 t; % 线性增长趋势
noise = 2 randn(size(t)); % 随机噪声
y = trend + season + noise + 50; % 基准值50
% 可视化
figure;
plot(t, y, ‘k’); xlabel(‘Time’); ylabel(‘Value’);
title(‘Original Time Series with Trend + Seasonality’);
🔍 Step 3：执行 STL 分解
matlab
period = 12; % 假设年度季节性（月度数据）
[trend, seasonal, remainder] = simple_stl(y, period);
% 可视化分解结果
figure;
subplot(4,1,1); plot(y); title(‘Original’); grid on;
subplot(4,1,2); plot(trend); title(‘Trend’); grid on;
subplot(4,1,3); plot(seasonal); title(‘Seasonal’); grid on;
subplot(4,1,4); plot(remainder); title(‘Remainder’); grid on;
📈 Step 4：分别预测各成分
4.1 预测季节项 (S_{t+h})
季节项是周期性的，直接重复最后一个完整周期：
matlab
h = 24; % 预测未来24步
lastCycle = seasonal(end-period+1:end);
future_seasonal = repmat(lastCycle, ceil(h/period), 1);
future_seasonal = future_seasonal(1:h);
4.2 预测趋势项 (T_{t+h})
使用简单线性回归（或 ARIMA、MLP）：
matlab
% 线性拟合趋势
p = polyfit(t, trend, 1); % 一次多项式
t_future = (length(y)+1 : length(y)+h)';
future_trend = polyval(p, t_future);
更高级：可用 arima 模型拟合趋势：
matlab
mdl = arima(‘Constant’, 0, ‘D’, 1, ‘Seasonality’, 0, ‘MALags’, 1);
estMdl = estimate(mdl, trend);
[future_trend, ~] = forecast(estMdl, h, ‘Y0’, trend);
4.3 预测余项 (R_{t+h})
余项通常视为白噪声，预测为0；若存在自相关，可用 AR 模型：
matlab
% 简单处理：设为0（或用AR(1)）
future_remainder = zeros(h, 1);
% 可选：AR(1) 拟合余项
if std(remainder) > 1e-6
armdl = arima(1,0,0);
armdl = estimate(armdl, remainder, ‘Display’, ‘off’);
[future_remainder, ~] = forecast(armdl, h, ‘Y0’, remainder);
else
future_remainder = zeros(h, 1);
end
🔗 Step 5：重构预测结果
matlab
y_pred = future_trend + future_seasonal + future_remainder;
% 可视化完整预测
t_all = [t; t_future];
y_all = [y; NaN(h,1)];
figure;
plot(t, y, ‘b’, ‘LineWidth’, 1.5); hold on;
plot(t_future, y_pred, ‘r–’, ‘LineWidth’, 2);
xlabel(‘Time’); ylabel(‘Value’);
legend(‘Historical’, ‘STL Forecast’);
title(‘STL-Based Time Series Forecasting’);
grid on;
📉 Step 6：评估预测性能（若有真实未来值）
matlab
% 假设有真实未来值 y_true（实际中可能没有）
% y_true = … ;
% rmse = sqrt(mean((y_pred - y_true).^2));
% fprintf(‘Forecast RMSE: %.2f\n’, rmse);
🧠 改进建议
成分 改进方法
趋势项 使用 Prophet、LSTM、SVR、或局部加权回归（Loess）
季节项 多重季节性（如周+年）、傅里叶级数拟合
余项 GARCH 模型（波动聚集）、ResNet 残差学习
整体框架 结合 EMD/VMD + STL，或使用 TBATS、Facebook Prophet
📦 实际应用建议
数据频率：确保 period 设置正确（日数据：7/365；小时数据：24/168）
缺失值处理：STL 对缺失敏感，需插值（fillmissing）
非固定周期：考虑使用 动态时间规整（DTW） 或 小波分解
✅ 总结
本项目展示了如何在 纯 MATLAB 环境 下：
- 实现 STL 分解（简化版）；
- 分别建模趋势、季节、余项；
- 重构未来预测；
- 可视化与评估。
若你有 真实数据集（如电力负荷、零售销量），只需替换 y 和 period 即可直接应用。
如需以下内容，请告知：
完整 .m 脚本文件；
调用 Python statsmodels.tsa.seasonal.STL 的混合方案；
基于 STL + LSTM 的混合预测模型；
多变量 STL（MSTL）扩展。
我会为你定制代码！
