当前位置：首页 > news >正文

EMD-SVM 太阳能功率预测

news 2025/10/29 7:21:24

含数据读取→EMD 分解→SVM 回归→多步预测→结果保存/绘图，完成 EMD-SVM 太阳能功率预测
基于 15 min 采样光伏功率 + 4 维气象数据

一、文件结构

EMD_SVM_Solar/  
├─ main.m                % 主脚本  
├─ emd_decomp.m          % 调用 EMD 分解  
├─ build_train_svm.m     % 训练 SVM 子函数  
├─ forecast_svm.m        % 滚动预测子函数  
├─ calc_metrics.m        % 计算 MAPE/RMSE/R²  
├─ data/  
│  ├─ solar_power.xlsx   % 时间序列功率  
│  └─ weather.xlsx       % 同期气象数据  
└─ results/              % 输出图、指标表

二、数据格式

solar_power.xlsx
列 A：时间（yyyy-mm-dd HH:MM:SS）
列 B：功率（MW）
weather.xlsx
列 A-D：辐照度(W/m²)、温度(°C)、风速(m/s)、湿度(%)

（数据须对齐、无缺失；如缺失请先插值。）

三、核心代码

main.m —— 一键运行

%% 1. 读数 & 归一化
clear; clc; close all;
power = readmatrix('data/solar_power.xlsx','Range','B:B');
meteo = readmatrix('data/weather.xlsx');          % 4×N
N = length(power);% 归一到 [0,1]
power = (power - min(power)) / (max(power) - min(power));%% 2. EMD 分解
[IMFs, res] = emd_decomp(power);          % 返回 cell 数组
K = numel(IMFs);                          % IMF 个数%% 3. 构造监督数据集
lag = 24;   % 用过去 6h 预测下一步
[X, Y] = deal(cell(K,1));                 % 每个 IMF 独立
for k = 1:K[X{k}, Y{k}] = buildData(IMFs{k}, meteo, lag);
end%% 4. 训练/验证/测试划分
ratio = [0.7 0.15 0.15];
idx = splitData(N-lag, ratio);mdl = cell(K,1);
for k = 1:Kmdl{k} = build_train_svm(X{k}(idx.train,:), Y{k}(idx.train));
end%% 5. 测试集滚动预测
horizon = 96;   % 预测未来 24h
pred = zeros(horizon, K);
lastWin = cellfun(@(x) x(end-lag+1:end), IMFs, 'uni', 0);
meteoFut = meteo(end-horizon+1:end, :);   % 假设已知未来气象for t = 1:horizonfor k = 1:Kfeat = [lastWin{k}'; meteoFut(t,:)'];pred(t,k) = predict(mdl{k}, feat');lastWin{k}(1:end-1) = lastWin{k}(2:end);lastWin{k}(end) = pred(t,k);end
end
yhat = sum(pred,2);                       % 反归一化前
yhat = yhat*(max(power)-min(power))+min(power);%% 6. 评估 & 绘图
ytrue = readmatrix('data/solar_power.xlsx','Range', sprintf('B%d:B%d', N-horizon+2, N+1));
figure; plot(ytrue,'b'); hold on; plot(yhat,'r--');
legend('真实','预测'); xlabel('采样点'); ylabel('功率 MW');
title(sprintf('EMD-SVM 未来 24h 预测\\nMAPE=%.2f%%',calc_metrics(ytrue,yhat)));writematrix([ytrue, yhat],'results/pred24h.csv');

emd_decomp.m —— 调用 MATLAB 内置 emd

function [IMFs, res] = emd_decomp(sig)[IMFs, res] = emd(sig,'MaxNumIMF',8,'Display',0);IMFs = mat2cell(IMFs, size(IMFs,1), ones(1,size(IMFs,2)));
end

buildData.m —— 生成监督样本

function [X, y] = buildData(sig, meteo, lag)N = numel(sig);X = []; y = [];for t = lag:N-1X = [X; [sig(t-lag+1:t)'  meteo(t+1,:)]]; %#ok<AGROW>y = [y; sig(t+1)];end
end

splitData.m —— 时序划分

function idx = splitData(N, ratio)trainN = floor(N*ratio(1));valN   = floor(N*ratio(2));idx.train = 1:trainN;idx.val   = trainN+1:trainN+valN;idx.test  = trainN+valN+1:N;
end

build_train_svm.m —— 训练 ε-SVR（RBF）

function mdl = build_train_svm(X, y)% 网格搜索 [C,gamma]（可换成贝叶斯优化）C_set   = 2.^(-2:2:8);g_set   = 2.^(-5:2:3);best = struct('mse',Inf);for C = C_set, for g = g_setmdl = fitrsvm(X, y, 'KernelFunction','rbf', ...'KernelScale',1/g, 'BoxConstraint',C, ...'Standardize',true);yhat = predict(mdl, X);mse = mean((y-yhat).^2);if mse < best.mse, best = struct('mdl',mdl,'mse',mse); endend, endmdl = best.mdl;
end

calc_metrics.m —— 指标计算

function mape = calc_metrics(ytrue, yhat)mape = mean(abs((ytrue-yhat)./ytrue))*100;fprintf('RMSE=%.3f MW, R²=%.3f\n', ...sqrt(mean((ytrue-yhat).^2)), 1-var(ytrue-yhat)/var(ytrue));
end