做风电和光伏出力的随机性建模时绕不开两个经典分布表征风速随机性的Weibull分布以及描述光伏出力占比特征的Beta分布。很多论文里都有“风速服从Weibull分布、光伏功率服从Beta分布”这句话但真正想在Matlab里把两套模型拟合出来、组合起来、再评估效果并没那么顺手。这篇博文我从实际项目角度出发把从数据生成、参数估计、组合建模到结果可视化的完整流程拆开讲清楚代码直接可跑。适合正在写新能源方向论文、做微电网容量配置或者搞功率预测的读者。不需要你有多深的数学功底但最好用过Matlab基础语法至少知道histogram和plot的区别。看完之后你能复现风电和光电概率分布的组合分析也能根据自己场站的数据替换出结果。1. 组合建模思路为什么是Weibull与Beta1.1 风电随机性风速Weibull模型的来龙去脉风电的随机性本质上来自风速的随机性。风速不是一个稳定值它受气压、地形、温度等多种因素影响实测数据通常表现为右偏、厚尾的形态。两参数Weibull分布的概率密度函数是f(v) (k/λ) * (v/λ)^(k-1) * exp(-(v/λ)^k)其中v是风速k是形状参数λ是尺度参数。k控制密度曲线的形状k小于1时曲线在零点附近陡峭上升风速大多集中在极小区间k等于2附近时曲线接近瑞利分布这是很多风资源报告里默认的情况k大于3时峰形变得尖锐风速集中在某个值附近。λ则大致决定风速的数值级别直接和平均风速挂钩。我见过很多初学者纠结要不要用三参数Weibull就是多一个位置参数μ把分布起点从0平移到某个风速。工程上两参数基本够用三参数虽然拟合偏差小一些但参数辨识的稳定性变差数据量不够时经常出现不收敛或者拟合参数明显不合理。除非你手上风速数据存在很明显的“零风速时段较多”且地点特殊否则优先用两参数。不过这里有个关键点需要说清楚标题里说的是“风电的Weibull分布”。实际建模有两种做法一种是直接对风速样本做Weibull拟合然后用功率曲线把风速分布映射为风电功率分布另一种是干脆把风电功率样本当随机变量直接拟合成Weibull。我个人更推荐第一种因为风速的Weibull分布物理意义清晰参数稳定而功率样本里大量0值和满发值会让Weibull拟合变得很难看。1.2 光伏随机性Beta分布为什么刚好合适光伏出力的随机性来自光照辐照度、温度、云层遮挡等因素。与风速不同光伏出力有一个天然边界功率只能从0到额定容量。用标幺值表示就是0到1之间的数。Beta分布恰好定义在[0, 1]区间内概率密度函数是f(x) x^(α-1) * (1-x)^(β-1) / B(α, β)其中B(α, β)是Beta函数α和β都是大于0的形状参数。α和β的取值直接影响分布形态αβ1时是均匀分布α和β都大于1时分布呈单峰形态且峰值位置由两者比值决定α小于β时密度曲线偏向左侧表示光伏出力偏小的时间更多这对很多多云地区很符合。从数学上看Beta分布几乎是描述光伏出力的“天选之子”——有界、单峰、左右偏态灵活可调。实际处理时光伏出力数据经常包含夜间零出力时段如果全年逐小时数据直接送进拟合会把分布拉向0.1以下的极端区间。所以工程上常规做法是分离零出力和正出力对正出力部分做归一化然后拟合Beta分布或者使用零膨胀Beta模型。后面对这个坑会有详细说明。1.3 两条组合路线怎么选“组合研究”这四个字在不同论文里含义差别挺大。我归纳为两条主流路线第一条是加权混合分布。把风电功率分布和光伏功率分布按权重叠在一起形成总出力的概率密度函数。设风电装机容量占比为w则组合密度函数是f_total(x) w * f_wind(x) (1-w) * f_solar(x)。这个形式简单适合做解析推导但它的物理含义不是“风电和光伏同时出力的总和”而是“随机抽取一个时刻该时刻出力来自风电或光伏的概率密度”。换句话说它描述的是整体电力出力的分布特征适合宏观分析。第二条是联合抽样或卷积。利用风电、光伏各自的分布抽样然后把两者出力相加可按装机容量比例缩放得到总出力样本再估计总出力分布。这个方法更贴近工程实际微电网容量配置、可靠性评估、储能容量优化都是基于这种“两个随机源叠加”的逻辑。要做选择关键看你的目标。如果论文方向偏概率建模、解析表达式走混合模型路线如果偏规划运行、需要模拟总出力序列走联合抽样路线。本文两种都会给出实现你按需取用。2. 数据准备与分布参数估计2.1 数据准备没有场站数据怎么验证算法很多人卡在第一步不是因为没有方法而是没有数据。手上有真实风电站和光伏电站的历史数据自然没问题直接读进来用。如果没有可以用Matlab自带随机数函数合成一套带已知真值的数据用来验证拟合程序是否写对。% 固定随机种子保证结果可复现 rng(42); % 一年小时级数据共8760个点 n 8760; % 用已知参数的Weibull分布生成风速样本 v_true_k 2.1; % 真实形状参数 v_true_lambda 6.5; % 真实尺度参数 (m/s) v_sample wblrnd(v_true_lambda, v_true_k, n, 1); % 用已知参数Beta分布生成光伏出力样本标幺值 pv_true_a 2.2; % 真实alpha pv_true_b 3.8; % 真实beta pv_sample betarnd(pv_true_a, pv_true_b, n, 1);这里用wblrnd生成Weibull样本用betarnd生成Beta样本两个函数都是Matlab Statistics and Machine Learning Toolbox里的。注意wblrnd的前两个参数顺序是(scale, shape)也就是(λ, k)和很多论文里的书写习惯不一致我当年第一次用就因为这个顺序搞反拟合出的参数完全对不上。后面所有用到的地方都要留意这个顺序。合成数据的最大好处是你在拟合后可以把估计值和真值对比直接判断算法实现是否正确。比如上面生成的v_sample用wblfit拟合后得到的λ应该接近6.5k接近2.1对pv_sample用betafit拟合后得到的α接近2.2β接近3.8。如果对不上就要回头检查代码了。2.2 Weibull参数估计的Matlab实现细节Matlab里做Weibull参数估计最简单的办法是直接用wblfit它基于极大似然估计。一行代码就能得到参数结果% 对风速样本拟合Weibull分布 phat wblfit(v_sample); lambda_fit phat(1); % 尺度参数 k_fit phat(2); % 形状参数跑完以后对照一下lambda_fit约6.52k_fit约2.08左右与真值接近。之所以不完全相等是因为有限样本估计本身就存在抽样误差这很正常。如果你想让结果更稳定可以加大样本量到87600甚至更多误差会进一步缩小。除了极大似然工程上还有矩估计和经验公式。矩估计的思路是利用Weibull分布的均值和方差与参数的解析关系反解参数。Weibull分布的均值μ与λ、k的关系是μ λ * Γ(1 1/k)其中Γ是伽马函数。如果已知样本均值和标准差σ可以用一个近似公式估算kk ≈ (σ/μ)^(-1.086)然后通过λ μ / Γ(11/k)反算λ。这个公式在风速资源评估中很流行手算或者在没有统计工具箱的环境下都能用。在Matlab里也可以直接用gamrnd配合迭代去解但既然有wblfit我建议你直接用内置函数把重点放在后续组合分析上不要在基础拟合上重复造轮子。有一点必须提醒wblfit对输入数据的范围很敏感。如果风速向量里有NaN、0或者负值拟合结果可能异常。风速为0时Weibull分布的概率密度在某些k值下趋近无穷极大似然迭代可能会出问题。后面专门讲零值处理。2.3 Beta参数估计的Matlab实现细节Beta分布的参数估计同样可以直接用Matlab的内置函数% 对光伏出力样本拟合Beta分布 abhat betafit(pv_sample); alpha_fit abhat(1); beta_fit abhat(2);betafit同样是极大似然估计内部用迭代算法求解。跑完以后alpha_fit差不多在2.2附近beta_fit在3.8附近。这个函数用起来很简单但有几个坑是文档里不会写的第一输入数据必须严格处于[0, 1]区间。如果你用光伏出力(kW) / 装机容量(kW)得到标幺值理论上是0到1但因为测量误差或者被四舍五入可能出现1.000000001这样的值betafit直接报错或者警告。处理办法是把数据做一次clipx min(max(x, eps), 1-eps)。第二数据中的0和1会让Beta分布的边界出现奇异。光伏出力为0的夜间时段如果全扔进去Beta拟合会试图用一个极端的参数组合去匹配这个零质量结果就是α被压到很小曲线变得特别难看。建议做法是先用逻辑索引把pv_sample拆成pv_zero等于0的部分和pv_pos大于0的部分只对pv_pos做归一化后拟合Beta。后面评估总体分布时再带上零出力的概率质量。第三betafit在数据量太少时可能不收敛。如果你只有几十个点拟合结果可能出现NaN。一般样本量最好在100以上做新能源出力分析时最少也得一个月逐小时数据也就是720个点通常没问题。3. 完整Matlab实现拟合、组合与可视化3.1 从风速Weibull到风电功率分布的转换现在进入核心实操。对风速做Weibull拟合拿到参数以后下一步是把风速分布转换为风电功率分布。这需要一条功率曲线。常规风机的功率曲线近似如下切入风速v_ci额定风速v_r切出风速v_co。低于切入风速或者高于切出风速时出力为0切入到额定之间近似线性爬升额定到切出之间保持满发。% 典型风机功率曲线参数 v_ci 3; % 切入风速 (m/s) v_r 12; % 额定风速 (m/s) v_co 25; % 切出风速 (m/s) p_rated 1; % 额定功率标幺化 % 将风速样本映射为风电功率样本 function p windPowerCurve(v) p zeros(size(v)); idx_ramp (v v_ci) (v v_r); idx_full (v v_r) (v v_co); p(idx_ramp) (v(idx_ramp) - v_ci) / (v_r - v_ci); p(idx_full) p_rated; end p_wind_sample windPowerCurve(v_sample);这段代码把风速样本映射成了风电功率样本量纲已经统一到标幺值。接下来你可以把p_wind_sample的分布画出来大概率会看到两个尖峰一个在0附近对应风速低于切入风速一个在1附近对应满发时段。中间爬升段的分布比较平坦。这种现象在实际风电功率数据里非常常见所以直接对功率样本强行做Weibull拟合时效果普遍不好这也印证了我前面说的“先拟合风速再映射功率”的思路。如果你想从风速Weibull分布解析推导风电功率的密度函数也不是不行但涉及到分段变换和雅可比行列式公式比较繁琐。工程上直接用Monte Carlo抽样加核密度估计ksdensity就能得到功率分布曲线简单又够用。后面组合阶段也是基于抽样样本所以这里不需要强行写出解析表达式。3.2 加权混合模型与蒙特卡洛联合抽样的实现到了最核心的组合环节。我按两条路线分别给代码。先看加权混合模型。假设风电装机占比为w光伏装机占比为1-w组合分布密度为两者概率密度的加权和。风电功率密度用ksdensity从p_wind_sample估计光伏功率密度直接用Beta拟合后的理论密度函数。% 路由1加权混合模型 x linspace(0, 1, 500); w 0.6; % 风电装机占比 % 风电功率的经验密度估计 [f_wind, xi] ksdensity(p_wind_sample, x, Support, [0, 1]); % 光伏功率的理论Beta密度 f_solar betapdf(x, alpha_fit, beta_fit); % 组合密度 f_mix w * f_wind (1 - w) * f_solar; % 画图对比 figure; plot(x, f_wind, r-, LineWidth, 1.5); hold on; plot(x, f_solar, b-, LineWidth, 1.5); plot(x, f_mix, k--, LineWidth, 2); legend(风电功率密度, 光伏功率密度, 组合密度); xlabel(标幺化功率); ylabel(概率密度);注意ksdensity的Support参数要设为[0, 1]否则风电功率分布会在0附近估算出负区间的密度这在物理上没有意义。另一种更贴近工程的做法是蒙特卡洛联合抽样。既然已经得到了风速Weibull分布和光伏Beta分布的参数就可以直接从这两个分布抽取一批样本映射到功率并叠加成总出力。% 路由2蒙特卡洛联合抽样 N 10000; v_sim wblrnd(lambda_fit, k_fit, N, 1); pv_sim betarnd(alpha_fit, beta_fit, N, 1); % 风速模拟样本映射为风电功率 p_wind_sim windPowerCurve(v_sim); % 总出力 风电出力(按容量占比) 光伏出力(按容量占比) p_total w * p_wind_sim (1 - w) * pv_sim; % 估计总出力的概率密度 [f_total, x_total] ksdensity(p_total, Support, [0, 1]); figure; histogram(p_total, 50, Normalization, pdf, FaceAlpha, 0.3); hold on; plot(x_total, f_total, k-, LineWidth, 2); xlabel(总出力标幺值); ylabel(概率密度);两条路线得到的结果在含义上略有不同。混合模型得到的是“整体出力的解析密度”联合抽样得到的是“两个电源真正叠加后的总出力分布”。如果你做储能容量配置或可靠性评估用后者如果你写综述性文章、只描述分布形态用前者。顺带提一个参数问题w的取值直接决定结果形态。按装机容量比取值是最常规的比如风电场50MW、光伏电站50MWw0.5。但有的论文会按“保证率”或“置信度”来优化w常见做法是让组合分布尽可能贴近历史总出力分布用KL散度最小化来搜索w。具体方法我放到第4部分讨论。3.3 结果可视化的输出要点做概率分布分析图比数字重要。审稿人和导师第一眼看的是分布曲线是否平滑、拟合是否贴合、组合形态是否合理。我的经验是多输出四类图基本就能覆盖需求。第一张图风速直方图与Weibull拟合曲线叠加。用histogram的Normalization设为pdf密度直方图才能和理论密度函数在同一个尺度上对比。figure; histogram(v_sample, 50, Normalization, pdf, FaceAlpha, 0.4); hold on; v_line linspace(min(v_sample), max(v_sample), 300); plot(v_line, wblpdf(v_line, lambda_fit, k_fit), r-, LineWidth, 2); xlabel(风速 (m/s)); ylabel(概率密度); legend(实测直方图, Weibull拟合);第二张图光伏出力直方图与Beta拟合曲线叠加。如果对正出力部分单独拟合记得把零出力的概率单独标注在图上否则看到直方图在0附近有一根很高的柱子会让人误以为拟合失败。第三张图组合模型对比。把风电功率密度、光伏功率密度、组合密度画在同一张图标出装机占比。这张图信息量很大直接说明组合权重的影响。第四张图总出力样本的直方图和核密度曲线。这是联合抽样路线的结果输出配合累计分布函数ecdf绘制CDF曲线方便后续对比不同装机配比下的出力特性。所有图建议统一设置坐标系范围标幺值功率的x轴固定为0到1风速图的坐标按实际范围调整。字体大小在出版投稿场景建议统一设12磅以上这里不再赘述。4. 常见问题与实战避坑指南4.1 零值、边界值与参数不收敛的处理这是我在实际数据处理里踩过最多的坑。风速数据经常有0值光伏数据到了晚上几乎全是0这两种情况如果直接塞进wblfit和betafit结果可能完全偏离你的预期。先处理风速零值。两参数Weibull的支撑集理论上是[0, ∞)但概率密度在v0处的行为取决于kk小于1时f(0)趋向无穷k大于1时f(0)为0。如果样本中0风速比例较高比如静风频繁的内陆地区单靠Weibull无法同时拟合“零值频率”和“正风速形态”。一个相对实用的办法是使用混合模型把零风速的离散概率单独建模对正风速部分用Weibull拟合。类似的做法在风资源领域叫“离散-连续混合分布”。如果只是为了论文分析且零值比例低于5%直接忽略零值、只对正风速拟合通常不会有太大问题。再看光伏的0和1边界。betafit要求输入严格在0到1之间0和1本身会让似然函数出现退化。处理逻辑分两步第一步把数据拆成零出力和正出力两组记下零出力比例p_zero第二步对正出力数据用pv_pos min(max(pv_pos, eps), 1-eps)做边界收缩再进入betafit。画总体分布时要用“p_zero * δ(0) (1-p_zero) * Beta密度”的叠加形式其中δ(0)表示0点的冲激质量。4.2 拟合优度评估不要只信p值参数拟合完当然要评估做得好不好。Matlab有kstest可以做Kolmogorov-Smirnov检验但带拟合参数的KS检验有个坑因为参数本身是从样本里估计出来的直接使用标准临界值会让检验变得保守大概率拒绝原假设。尤其是样本量在几千以上时KS检验对微小偏差非常敏感几乎必然拒绝“样本来自该分布”的原假设。我实际项目里的做法是多重验证。首先眼见图把直方图和拟合曲线叠加肉眼判断形态是否吻合。这听起来不够“科学”但在处理海量新能源数据时视觉判断其实非常有效。其次算拟合误差RMSE在直方图密度值和理论密度值之间计算均方根误差数值越小越好。然后用KS距离作为参考指标但不要只看p值而是看统计量本身的大小。最后可以做一个Q-Q图横轴为理论分位数、纵轴为样本分位数如果散点贴近yx直线说明拟合良好。如果RMSE偏离明显优先怀疑参数估计是否正确再看数据预处理是否有问题比如0值是否处理、边界是否收缩、样本量是否足够。我见过一个案例某同学拟合光伏Beta分布时没有去掉夜间零出力结果α被拟合到0.6密度曲线在0附近翘得很高RMSE大得离谱去掉零值以后一切恢复正常。4.3 组合权重的选择策略与优化思路组合权重w的选择直接决定组合分布长什么样。最朴素的办法是按装机容量比简单透明可复现。但如果风电和光伏的实际利用率差异很大比如当地弃风严重、风电实际出力远低于额定容量直接用装机比会让组合分布偏高估风电贡献。这种情况下建议按“可用容量”或“平均出力比”来设定权重。更精细的做法是让组合分布逼近历史总出力分布。你可以取一个候选w的网格比如0.1到0.9步长0.01对每个w计算组合分布与历史总出力经验分布的KL散度选择KL散度最小的w作为最优权重。KL散度计算的核心代码如下% 假设 f_hist 是历史总出力的核密度估计f_mix 是当前权重下的组合密度 kl_div sum(f_hist .* log(f_hist ./ (f_mix 1e-12))) * mean(diff(x));加上一个极小值1e-12是为了防止除零和出现Inf。这种基于数据驱动的权重优化思路比拍脑袋定权重要有说服力得多写在论文里也是一个加分项。我最后还想多一句以上所有方法本质上都是对随机性的“概率近似”。真实的风电光伏出力还受时间相关性、季节变化、天气过程等因素影响单靠静态分布并不能完全刻画。分布建模更像是给你提供一把标尺让你在做规划和调度时心里有底。如果你后续想做更精细的时序模拟可以考虑把Weibull参数按季节分段或者给Beta分布加上随时间变化的参数。这个方向我最近也在试后续有结果再分享。