混合概率分布范文
混合概率分布范文（精选7篇）
混合概率分布 第1篇
近几十年来, 风力发电成本不断下降, 已经具有了与常规电源竞争的实力。风电并网和规划的相关研究[1,2]正日益得到重视。
风速概率分布是体现风电场风能资源的重要指标之一, 也是风电场规划运行的重要参考[3]。常用的风电场风速分布的拟合模型有瑞利分布[4]、对数正态 (Lognormal) 分布[5]和威布尔 (Welbull) 分布[6]等, 其中两参数威布尔分布是在绝大多数情况下拟合效果最好、应用最为广泛的一种, 它在拟合长周期 (如年或部分月份) 风速分布时, 效果良好[6]。
然而, 随着风电场装机容量比重的不断增加, 为保证风电场接入后电网依然保持安全稳定运行, 需要利用短时分布和特殊时段分布信息加强系统运行分析。但当需要研究更短周期或某些特殊时段的风速概率分布特性时, 气候变化等随机因素对分布影响明显增强[7]。经某风电场实际数据统计发现, 同一年内不同月份风速分布随季节气候变化不同, 分布往往存在差别。不同年份的同一月份分布也不相同, 常出现两峰甚至三峰分布情况。其规律已很难用两参数威布尔函数准确逼近。
国内外学者对多峰分布均做过相关研究。文献[8]考虑面积约束后以3次Hermite插值函数对直方图进行逼近, 效果良好, 但缺乏统计意义。其他文献多为求解混合分布。文献[9]介绍了一种用遗传算法求解混合分布的方法, 计算复杂费时;文献[10]介绍了Kececioglu方法;文献[11]采用L-M算法求解参数, 但Kececioglu方法和L-M方法都需要人工作图计算初值, 工作量较大, 有待改进。
针对以上问题, 本文提出采用七参数混合威布尔函数逼近特殊时段风速分布的方法。七参数威布尔函数是两个三参数威布尔函数的叠加, 每个威布尔函数具有各自的尺度、形状和位置参数, 可以适应双峰型风速概率密度函数逼近的需要。本文采用极大似然估计法[12]并根据实测风速数据建立混合威布尔分布的极大似然方程。考虑到特殊时段风速分布的两个峰值间可能具有较大的差异性, 且有时需要更多地关注某一峰值或某些风速段的概率密度, 本文采取对这些风速段进行加权的方式构造似然函数, 以提高这些风速段概率密度函数的拟合优度。在加权似然函数的基础上, 本文采用数论布点法[13]求解待求的7个参数, 该方法只要求目标函数连续, 对初值选择不敏感, 迭代过程只需计算目标函数。
1 风电场短时风速分布的统计特性
取某规划风场2006年和2007年11月份 (冬季大风期) 10 min平均风速采样数据绘制直方图, 如图1所示。分别对两个分布进行两参数威布尔分布的极大似然估计拟合, 可看出二者威布尔函数形状差别较大。这是因为2007年11月大于6 m/s的风速样本个数要大于2006年, 因此该年威布尔函数相对2006年有“右移”趋势。直观看2006年风速分布存在3 m/s～4 m/s和6 m/s～7 m/s两个风频峰值。而2007年风速分布存在6 m/s～7 m/s和10 m/s～11 m/s两个风频峰值, 具有双峰分布的特点。用两参数威布尔分布拟合后, 2006年风频最大误差出现在4 m/s～5 m/s风速段, 2007年风频最大误差出现在8 m/s～9 m/s段, 其相对误差分别达到53.2%和29.4%;而且对两年拟合显著程度进行皮尔森卡方 (Pearson chi-square) 检验 (显著性水平α=0.05) 后, 2006年结果未通过检验, 可见两参数威布尔分布逼近有局限性。在实际应用时, 究竟选用两参数还是混合威布尔分布, 主要根据分布的峰数决定。对分布数据滤波 (去除影响判断的毛刺) 后对风频进行扫描, 若发现几个相邻连续风频先上升后下降, 则可判断出现一个峰。若滤波处理后的分布只有一峰, 则选用两参数威布尔分布进行拟合;若有多峰, 则选用混合威布尔分布进行拟合。图1显然需要选用混合威布尔分布进行拟合。
2 混合威布尔函数拟合双峰分布
三参数威布尔函数是混合威布尔函数的基础, 比两参数威布尔函数多了一个描述分布偏移程度的位置参数。其概率密度函数f (v) 和概率分布函数F (v) 分别为:
F (v) =1-exp-
式中:
是集合 (u, ∞) 的示性函数;k, c, u分别为非负的形状、尺度和位置参数。
在三参数分布基础上, 七参数混合威布尔概率密度函数fmix (v) 和概率分布函数Fmix (v) 可分别表示为:
fmix (v) =rf (v|c1, k1, u1) + (1-r) f (v|c2, k2, u2) (3)
Fmix (v) =rF (v|c1, k1, u1) + (1-r) F (v|c2, k2, u2) (4)
混合威布尔分布函数的参数计算是一个难点。矩估计是统计学常用的点估计方法[14], 它可通过矩估计对样本分布进行大致逼近。以样本低谷风频段的起始风速vdivide为分界风速, 将风速总体样本集V分为V1和V2两个子集, 其中V1为风速不大于vdivide的样本子集, 含N1个样本;V2为风速大于vdivide的样本子集, 含N2个样本。由V1和V2可获得混合威布尔分布的矩估计参数。
以各样本子集风速最小值作为位置参数umi (i=1, 2) 的矩估计值:
形状参数的矩估计值kmi (i=1, 2) 由如下方程求出:
式中:vavi (i=1, 2) 分别为V1和V2的平均风速;σi (i=1, 2) 分别为V1和V2的风速样本方差;Γ为伽马 (Gamma) 函数。
于是可得尺度函数cmi (i=1, 2) 的矩估计值为:
百分比参数rm可通过样本比例得到:
至此得到所有参数的矩估计值。一般来说, 矩估计的精度较差, 因此本文仅采用矩估计值作为优化算法的参考初值, 然后采用数论布点方法求得混合威布尔分布参数的极大似然估计。
另外, 由以上分析可见, 七参数混合威布尔分布就是两个三参数威布尔分布的“叠加”, 同样的思路可拓展到十一参数及更多参数的分布, 以适应三峰或更多峰的情形。本文仅以短时分布最易出现的双峰分布为例说明问题。
3 混合威布尔分布的参数优化
3.1 混合威布尔分布的加权极大似然方程
为求得混合分布参数优化解, 可由混合威布尔概率密度函数构造似然函数:
由此转化为求L的最大值问题。普通极大似然估计的无偏性和一致性使得每个风频段具有相同的“权”, 事实上, 根据不同的研究需要 (如风机额定运行风速要比非额定运行风速重要) , 所有风频段均分配相同权未必合适, 可考虑对不同风频段进行权分配。权分配应根据实际风场情况有所偏重, 为此构造加权极大似然函数如下:
式中:Pi为样本个体vi对应的权, 可根据研究需要设定, 也可通过计算优化调节。
3.2 采用数论布点法计算优化参数
数论布点方法[13]是数论与近似分析交叉的产物, 其实质是在S维单位立方体CS上找到一个均匀散布点集, 由数论方法得到的均匀散布点集通常称为数论网格。
用数论布点方法求函数f (z) , z∈D (D是一个S维矩形a, b=) 的极大值点z*, 使M=f (z*) 的基本思想是:在D上取一个数论网格P={zk, k=1, 2, , n}, 如果f (z) 是连续的, D是有界闭集, 当n充分大时, 有zk*∈P使得f (zk*) 接近M, zk*接近z*。生成数论网格的方法为:
再令:
得到a, b上数论网格{zk= (zk1, zk2, , zkn) , k=1, 2, , n}。为提高求解精度, 使用序贯优化算法, 即在迭代过程中使D的边长收缩, 收缩比为本步与上步求解区域边长比。具体计算步骤如下:
步骤1:初始:m=0, D (0) =D, a (0) =a, b (0) =b;
步骤2:布点:产生数论网格P (m) ={z
步骤3:选优:求P (m) ∪{z
步骤4:迭代:令c (m) = (b (m) -a (m) ) /2, 若max c (m) i <δ, δ为预先给定正数, 取z*=z (m) , M=M (m) , 停止迭代, 否则进入步骤5;
步骤5:收缩区域:令a
D范围的确定对优化结果影响较大。矩估计结果已可大致反映逼近函数, 因此可为D提供参考。一般地, 若矩估计的形状或尺度参数大于1, 相应D的范围下限可设为ki_min=kmi/2或ci_min=cmi/2, 上限可设为ki_max=kmi+3.0或ci_max=cmi+3.0;若形状或尺度参数小于1, 为避免计算溢出, 形状参数下限设为0.05, 上限设为0.95, 位置参数的数论布点优化结果一般与矩估计相差不大, 可设下限ui_min=umi/2, 为避免对数运算真数为0, 上限ui_max=umi-0.01。优化结果r范围变化较大, 若rm<0.8, 则设定上限rmax=rm+0.2, 否则rmax=1;下限统取为rmin=rm/2。
采用文献[13]提供的生成向量产生数论网格。m=1时布点数记为n1, m为其他值时均取为n2。为在第1步就找到较好的M (m) 而之后又不会过分降低计算速度, 取n1>n2。为保证收敛速度合理, 可设定第1步收缩比为r=0.8, 其余各步收缩比为r=0.5。
3.3 加权值优化
特殊风频段的逼近精度与该段权有重要关系, 权过低或过高都会导致逼近精度降低。风频段权上限大于2后一般逼近误差已经较大, 因此权多集中在区间[1,2]内, 变化范围较小, 为简化计算可采用较简单的步进加权方法。每次权增加一个步长, 随着权的增大, 特定风频段的逼近误差一般先减小后增大, 鞍点对应权值即可作为权优化结果。
特殊风频段的优化加权结果有时可能会极大地降低全局逼近精度, 甚至无法通过皮尔森卡方检验。这种情况下就要根据规划或工程需要有所取舍, 可综合考虑逼近精度和拟合效果选取一个略小于优化结果的较优值。
4 算例分析
对某风场2008年1月每天17时至20时期间的10 min平均风速样本进行统计, 并对样本数据进行极大似然估计拟合, 结果如图2所示。
由图2可见, 该特殊时段风频分布为双峰分布, 两峰风频段分别为4 m/s～5 m/s和8 m/s～9 m/s。其中8 m/s～9 m/s风频较大达0.207, 直方图中该风频段变化陡峭。用两参数威布尔函数似然逼近后, 风频为0.127, 相对误差达38.6%, 而该风频段恰是风机额定运行风速区间, 显然两参数威布尔分布不能准确描述该特殊时段的分布, 应考虑采用混合威布尔函数逼近。
4.1 一般极大似然估计拟合
采用本文方法, 分界风速vdivide为7 m/s, 分别以矩估计和一般极大似然估计 (不考虑加权) 对该样本进行混合威布尔函数逼近, 以均方百分比误差 (MSPE) [14]来检验3种方式的拟合误差;以皮尔森卡方检验方法检验不同方式的显著程度, 取显著性水平0.05, 将样本数较少的风频段合并后共划分为11个区间, 3种拟合方式的全局预测误差均方 (MSPE) 、重要风频段 (8 m/s～12 m/s) MSPE比较结果和皮尔森卡方检验结果见表1前3行。
直观看来, 七参数威布尔分布的矩估计和极大似然估计逼近反映出特殊时段的风速分布情况。分析表1结果可知, 两参数分布的似然估计通过了皮尔森卡方检验, 且全局MSPE小于矩估计, 但重要段MSPE较大, 逼近效果欠佳。作为粗略估计, 矩估计全局MSPE较大, 虽未能通过皮尔森卡方检验, 但相比两参数分布, 关键段的MSPE得到有效降低。作为矩估计的优化结果, 七参数分布的一般极大似然估计、全局和重要时段的MSPE均小于两参数分布情况, 且通过皮尔森卡方检验, 逼近效果良好。
4.2 加权极大似然估计拟合
为进一步提高重要风频段的逼近精度, 以重要段 (8 m/s～12 m/s) 的MSPE作为逼近指标, 采用步进加权方法进行优化, 加权步长取为0.1。不同权值对应结果见表1, 逼近效果见图3。
由图3可见, 随权重的增加, 重要风频段的MSPE先减小后增大, 说明特殊段权重较小或较大时逼近精度都会变差;全局MSPE随权重增大而增大, 说明优化过程在提高了局部逼近精度的同时牺牲了全局的逼近精度。在权重为1.7时, 重要风频段的MSPE虽然达到最小值0.114 8, 但对应的全局MSPE达到0.154 6, 误差较大;在权重大于1.3时, 所有结果均拒绝检验, 说明权重较大后全局逼近精度已下降到皮尔森卡方检验不可接受的程度。在接受检验的4组优化结果中, 权重1.3的全局MSPE较权重为1.7时减小了0.52, 而重要段的MSPE仅增加了0.001 9, 可见特殊段权重为1.3的混合威布尔分布兼顾了全局和重要风频段的逼近精度。对应不同的特殊双峰分布时段, 可根据规划的不同要求在不同权重下选取不同的拟合结果, 甚至有时为衡量规划风场额定风速区间对接入电网的冲击也可忽略皮尔森卡方检验。针对本文算例, 采用权重为1.3的优化结果较为合适。
5 结语
针对短时风速概率分布易出现双峰分布的特性提出了混合威布尔逼近方法。其中, 七参数混合威布尔分布能够准确反映短时风速的双峰分布情况, 而更多个三参数威布尔分布的混合可适应更多峰的情形。
极大似然估计模型复杂, 计算时间随峰数和样本数量的增加而延长, 是否可以找到一种逼近精度与计算速度均有所提高的新方法仍有待研究。
概率分布列问题解法解析 第2篇
例1 5封不同的信, 投入三个不同的信箱, 且每封信投入每个信箱的机会均等, 是三个箱子中放有信件数目的最大值.求的分布列.
分析:三个箱子中放有信件数目的最大值ξ取最大值5, 我们容易想到, 因为5封信全部放在一个信箱中.但三个箱子中放有信件数目的最大值ξ取最小值是几呢?颇费思量, 当然经过深入思考后不难得知是1, 2, 2这种情形中的2.
解:ξ的分布列为:
二、对不同情形的发生进行分类
例2甲、乙两队进行一场排球比赛, 根据以往经验, 单局比赛甲队胜乙队的概率为0.6, 本场比赛采用五局三胜制, 即先胜三局的队获胜, 比赛结束.设各局比赛相互间没有影响, 令为本场比赛的局数, 求的概率分布. (精确到0.000 1)
解:当ξ=3时, 表示甲胜3局或乙胜3局.则P (ξ=3) =0.63+0.43=0.28.
当ξ=4时, 表示前3局甲胜2局, 第4局甲胜或前3局乙胜2局, 第4局乙胜.则P (ξ=4) =C320.620.40.6+C320.420.60.4=0.374 4.
当ξ=5时, 表示前4局甲胜2局, 第5局甲胜或前4局乙胜2局, 第5局乙胜.则P (ξ=5) =C420.620.420.6+C420.420.620.4=0.345 6.
则ξ的分布如下:
三、对不同元素的组合进行分类
例3从8个男生5个女生中抽取6个学生参加义务劳动, 其中女生的人数是随机变量, 求ξ的分布列.
解:ξ的分布列如下:
四、对n次独立重复试验发生了k次进行分类
例4某工厂规定, 如果工人在一个季度里有1个月完成生产任务, 可得奖金90元;如果有2个月完成生产任务, 可得奖金210元;如果3个月都完成生产任务, 可得奖金330元;如果3个月都未完成任务, 则没有奖金.假设某工人每月获得奖金与否是等可能的, 求此工人在一个季度里所得奖金的分布列.
可得分布列如下:
五、对n次独立重复试验第k次才发生进行分类
例5某射手有5发子弹, 射击一次命中概率为0.9, 如果命中就停止射击, 否则一直到子弹用尽, 求耗用子弹数ξ的分布列.
解:ξ的取值有1, 2, 3, 4, 5.
当ξ=1时, 即第一枪就中了,
故P (ξ=1) =0.9;
当ξ=2时, 即第一枪未中, 第二枪中了,
基于概率分布模型的暴雨研究 第3篇
黄河中下游地区暴雨主要集中在每年的5~9月这153天的时间里[9]。在此期间每天发生暴雨的概率很低, 而总天数 (153天) 较大, 可见暴雨为稀有事件, 所以其概率分布可用Poisson分布来拟合。本文采用Poisson分布模式来研究黄河中下游若干站点近50年暴雨概率分布特征, 暴雨频数分布规律。同时还用建立的模式计算得到各地每年发生n次以上暴雨的概率, 为各地暴雨预报提供参考。
1 资料与方法
1.1 资料
气象统计中规定, 日降雨量≥50 mm为暴雨日[10]。数据来源与中国气象局中国地面国际交换站。选取黄河中下游若干站点历年暴雨季 (5~9月) 逐日降水量资料, 统计近50年 (1961年至2010年) 暴雨 (日降雨量≥50 mm) 日数及其发生频率。所选各测站站名如表1所示。
1.2 方法
1.2.1 泊松分布
泊松分布是1837年由法国数学家泊松 (PoissonS.D.1781-1840) 首次提出的。泊松分布在各领域的研究相当广泛, 是一种经典的描述稀有事件频率分布的概率模式。
泊松分布可以看作二项分布的极限分布。当n很大时, p很小时, 可以用泊松分布来计算二项分布。参数λ是单位时间 (或单位面积, 体积) 内随机事件的平均发生率.即
服从泊松分布的随机变量X的概率分布为:
该分布仅有一个参数λ (恒为正) 。
递推公式为:
可由 (3) 式得一年中 (5~9月) 发生k次暴雨的概率P。那么, 50年中应有50P年发生k次暴雨, 进而可得50年间每年发生k次暴雨的理论年数。结合各站值和公式 (1) 则可建立使用黄河中下游的暴雨Poisson分布模型。
根据给出的模型, 可以算出各地每年发生n次以上暴雨的概率。
每年发生n次以上暴雨的概率为:
1.2.2 暴雨分布概率模式研究
选用4个代表雨量站50年 (1961年至2010年) 逐日降水资料, 计算一年中暴雨季 (5~9月) , 每天发生暴雨的概率 (表2) 。
由表可以看出, 郑州历年 (5~9月) 共发生了94次暴雨, 延安历年 (5~9月) 共发生了38次暴雨, 西安历年 (5~9月) 共发生了34次暴雨, 太原历年 (5~9月) 共发生了116次暴雨, 且他们每年中 (暴雨季) 每天发生暴雨的概率分别为0.0122, 0.0050, 0.0032, 0.0152。
以郑州为例, 研究模式拟合问题。
由表2知, 郑州 (1951年至2000) 50年中 (5~9月) 共发生了94次暴雨。每年5~9月有15 3天, 在此期间每天发生暴雨的概率为0.0122.此概率值甚小, 而总天数153天较大, 可见暴雨是稀有事件。所以, 其概率分布可用Poisson分布来拟合。
为求郑州每年 (5~9月) 发生k次暴雨的概率, 先计算分布参数, 然后再按照递推公式 (3) 可以求得一年中郑州发生k次暴雨的概率。易得一年中郑州发生k次暴雨的理论年数。
根据中国地面国际交换站提供数据可以查询各站点每年 (暴雨季) 发生k次暴雨的实际观测数据 (表5) 。
各站点Poisson分布的概率模式的理论年数与实测年数的直方图比较如图1, 2, 3, 4所示。
2 检验与分析
2.1 方法
2.1.1 相关系数
定义:设 (X, Y) 是一个二维随机变量, 则称:
为X与Y的相关系数。
Corr (X, Y) 与协方差Cov (X, Y) 是同符号的, 即同为正, 或同为负, 或同为零。这说明, 从相关系数的取值可以反映出X与Y的正相关, 负相关和不相关。
2.1.2 分析
由表6可以看出郑州, 西安, 延安, 太原四个站点实测频数与理论频数的相关系数分别为0.926, 0.993, 0.994, 0.928。由此可见郑州, 西安, 延安, 太原暴雨季 (5~9月) 暴雨日数符合Poisson分布模式。即Poisson分布模式能较好地描述上述站点的暴雨概率分布。
检验结果已表明郑州, 西安, 延安, 太原暴雨概率分布能服从Poisson分布, 那么根据表3给出的模式和参数, 和公式 (5) 可以计算出各地每年暴雨季暴雨发生n次以上暴雨的概率 (表7) 。
3 结论
用Poisson分布理论描述黄河中下游若干站点 (西安, 延安, 郑州, 太原) 暴雨季暴雨分布符合情况较好。同时用建立的Poisson理论模式, 求得各地每年暴雨季中发生n次以上暴雨的概率, 可为各地暴雨预报提供参考。尽管本文选取的站点有限, 但根据降雨的区域连续性可以推测出相近或相邻区域的暴雨情况。本文只能在概率层面上预测出暴雨发生的次数, 未能具体到发生暴雨的具体时间或时间段, 这也为以后得继续研究指明了方向。
摘要：暴雨作为黄河中下游地区主要的气象灾害, 每年都对经济社会产生重要的影响。本文在全面的文献调查基础下选取黄河中下游四个代表雨量站:郑州, 西安, 延安和太原的1960年至2010年每年5月到9月的逐日降水资料, 运用泊松概率分布模式揭示该区域的暴雨统计特征。通过研究我们得出了一个结论, 泊松分布模型能较好地反映该地区暴雨的概率特征。在此基础上我们预测出雨量站一年内发生N次暴雨的概率。研究结果在气象统计和暴雨预报方面有重要的意义。
关键词：暴雨,概率分布,泊松分布,气象统计,黄河中下游
参考文献
[1]丁裕国, 申红艳, 江志红, 等.气候概率分布理论及其应用新进展[J].气象科技, 2009 (3) .
[2]陈艳秋, 袁子鹏, 盛永, 等.基于概率分析的暴雨事件快速评估模型[J].气象与环境学报, 2006 (5) .
[3]于新文, 丁裕国.中国东部地区暴雨的概率特征——基于泊松分布的统计模拟[J].自然灾害学报, 2006 (4) .
[4]王丙参, 魏艳华, 孙春晓.泊松分布参数估计的比较研究[J].四川理工学院学报 (自然科学版) , 2011 (5) .
[5]颜亦琪, 易建军, 孙华安.泊松分布在水文频率计算中的应用[J].人民长江, 2010 (12) .
[6]Park J S, Jung H S.Modeling Ko-rean extreme rainfall using a Kappa distribution and maximum likelihood estimate[J].Theoretical and Applied Climatology, 2002, 7 (2) :55-64.
[7]Khan B, Iqbal MM, Yousufza i.MAK (2011) Flood risk assessment of River Indus of Pakistan4:115–122.Doi:10.1007/s12517-009-0110-9.
[8]雷鸣, 高治定, 宋伟华.黄河中游河龙区间区域性暴雨特性分析[J].人民黄河, 2012 (8) .
[9]刘健, 于兰兰, 翟建青.近50a黄河三角洲极端降水事件特征研究[J].人民黄河, 2013 (2) .
混合概率分布 第4篇
随着国民经济发展和居民生活水平提高, 基于电力电子技术的非线性居民电器日益广泛使用, 且插入式电动汽车和民用新能源获得快速发展。这些分布式谐波源具有类似的大小并遍布在整个电力网络中。由于它们独自产生的谐波并不显著, 过去一般忽略由家用负荷注入的谐波, 但大量的类似负荷产生的集体谐波却是非常可观的, 其对电力系统的影响与危害已不可忽视[1,2,3]。居民非线性负荷产生的谐波同传统工业负荷具有明显的差别, 主要是居民负荷在电力系统中的分散性和工作状态的随机性, 评估和消除它们引起的谐波畸变在理论和工程上非常困难。
国外学者较早关注居民负荷产生的谐波, 并已开展了较多的研究工作。学者Wilsun Xu研究团队发现居民谐波电压畸变主要由系统背景谐波确定, 谐波电流畸变主要由居民负荷和系统背景谐波共同确定[1];研究了紧凑型荧光灯等单个家用电器的谐波负荷模型[2,3,4,5];并建立了基于概率论的自下向上的居民负荷谐波建模体系[6,7]。文献[8-9]定量分析了不同荧光灯负荷渗透率下配电网的谐波危害。文献[10]分析了大量单相电力电子型负荷谐波的分散效应和消减效应。
居民负荷的集体谐波评估本质上是一个概率性谐波求和问题。在概率性谐波求和方面, 国内外学者也已开展了较多研究[11,12,13,14,15,16], 主要有正态分布法、蒙特卡洛方法、中心极限定理法、基于Laguerre多项式法等。正态分布法[11]适用于评估相互之间完全独立、且大量存在的类似负荷所产生的谐波;但是考虑居民负荷类型、用户行为习惯或电力监管政策的影响, 各个居民负荷产生的谐波并不完全独立, 采用正态分布法将产生较大的误差。蒙特卡洛法[12]需要大量复杂的现场数据, 计算耗时多, 难以普遍应用。Laguerre多项式法[13,14,15]计算中会产生截断误差, 且尺度因子、修正系数的选取及样本的标准化还需进一步研究。中心极限定理法[16]联合概率密度函数的确定非常复杂, 而且需要大量的测试数据和复杂的运算, 难以实际应用。文献[17-18]对多谐波源节点的谐波责任进行评估划分, 以确定不同支路的谐波发射责任。现有方法主要从电气设备安全即电磁兼容角度进行谐波水平评估, 一般以集体谐波概率最大值, 即超过某一概率值 (如95%) 谐波电压或电流的大小, 作为谐波水平评估标准。而对于集体谐波的概率性分布关注较少, 而概率性分布体现了谐波来源, 对于确定谐波治理措施和监管政策具有重要指导意义。
云理论是国内李德毅教授提出的一种综合了模糊理论和概率理论的数学方法, 在设备敏感度评估、大坝裂缝监控、图像分割、电能质量数据分析等领域获得成功应用[19,20,21,22]。云理论是一种比概率论正态分布精确, 同时比联合分布简单易算的随机相量求和方法。本文提出一种基于云理论的居民负荷集体谐波概率分布的研究方法, 可以给出精确的居民负荷谐波概率分布, 对于研究居民谐波来源、确定谐波治理措施和监管政策具有重要理论意义。
1 云理论基本概念
1.1 云模型
云模型[23]的定义:给定定量论域U, C是U上的定性概念, 若定量数值x∈U是定性概念C的一次随机实现, x对C的确定度μ (x) ∈[0, 1]是具有稳定倾向的随机数 (不是一个固定的数值) , 即
则x在论域U上的分布称为云, 记为C (X) 。每一个x称为一个云滴。
云的数字特征用期望Ex、熵En和超熵He这3个量表示, 如图1所示, 能够把物理概念的随机性和模糊性很好地集合在一起, 形成了定性与定量之间的映射。
由云滴xi (i=1, 2, …, n) 组成的数据样本期望Ex、熵En和超熵He的计算表达式为:
1.2 云变换
云变换[24,25,26,27]把任意一个不规则的空间数据分布进行数学变换, 生成若干概念的云模型集, 使之成为若干个大小不同的云的叠加。云变换在一定的误差范围内可将任意频率分布函数分解为若干基云的叠加, 一般情况下, 叠加的云的个数越多, 变换误差越小。其数学表达式为:
式中:g (x) 为实际数据的概率分布函数;fi (x) 为云模型的期望函数;ai为幅度系数;n为产生云的个数;ε为误差阈值。
云变换的实质是采用若干个云模型的期望曲线的叠加来拟合实际数据的频率分布。
2 居民负荷集体谐波评估流程
基于云理论的居民负荷集体谐波评估流程主要包括5个部分, 如图2所示。
流程1:测量点电压、电流数据采集
中国居民配电网通常采用“闭环设计、开环运行”的模式, 低压配电网采用220V/380V电压单相或三相供电, 形成“电网—配电站—分支箱—表箱—家庭”的多级辐射电能传输网。谐波测量点可选择家庭、馈线及台区进行电压、电流数据采集, 如附录A图A1所示。
流程2:电压、电流测量数据傅里叶分析
对测量点所测电压、电流数据进行离散傅里叶 (DFT) 分析, 为满足计算精度和速度要求, 一般1s连续采样10个周期, 每个周期256个点。
流程3:电压、电流谐波值预处理
现有研究表明居民供电母线电压畸变主要由系统背景谐波确定, 谐波电流畸变主要由居民负荷和系统背景谐波共同确定[1]。本文主要研究居民负荷所产生的谐波电流。采用系统谐波阻抗法来区分负荷谐波和背景谐波, 即当系统谐波阻抗角在第一象限和第四象限时, 认为谐波电流由系统侧流向负荷侧;当系统谐波阻抗角在第二象限和第三象限时, 谐波电流的流向是由负荷侧流向系统侧。本流程剔除位于第一象限和第四象限的样本点以确保研究数据为负荷所产生的谐波。
流程4:谐波电流幅值云模型分析
谐波电流幅值概率密度函数采用云模型方法表示如式 (8) 所示。
式中:Ih和EIh分别为谐波电流幅值和期望。
其中谐波电流幅值Ih由前述流程得出, 每个谐波电流幅值看作1个云滴, 再根据式 (2) 、式 (3) 、式 (4) 计算得到谐波电流幅值的期望、熵和超熵, 形成谐波电流幅值的云模型表示, 进而计算其概率密度函数, 得出谐波电流幅值的概率分布。
式 (8) 是一个没有明确解析表达式的概率密度函数, 本文通过式 (9) 积分变换后采用高斯—埃尔米特 (Gauss-Hermite) 积分法[20]对任意谐波电流幅值I0积分以求得相应的谐波电流幅值的概率密度f (I0) , 从而进一步获得95%概率最大值。
流程5:云变换法确定谐波电流幅值概率密度函数
基于云变换的曲线拟合算法流程如图3所示。
步骤1:设置初始误差阈值ε, 迭代次数初值k=1。
步骤2:根据谐波电流幅值数集生成其概率分布曲线gk (Ih) 。
步骤3:寻找gk (Ih) 的最大值, 并作为第k个峰的幅度系数ak。
步骤4:在gk (Ih) 的最大值对应的Ih两边[Ih-Δ, Ih+Δ]取点, 根据式 (2) 、式 (3) 计算该峰所对应的云模型期望EIhk和熵Enk, 根据式 (6) 得到云模型期望函数fk (Ih) 。Δ根据需要设置, 例如Δ=0.1。
步骤5:从gk (Ih) 中减去已知云模型的数据分布fk (Ih) , 再迭代寻找最大值, 即重复步骤3, 4, 5, 得到多个基于云的数据分布函数。
步骤6:经过数次迭代后当gk (Ih) 的绝对值的最大值小于ε时迭代终止输出式 (5) 表示的曲线的拟合函数。
3 典型居民负荷谐波电流发射特征
文献[28-30]针对典型家用电器的谐波电流发射特性进行了测试及研究, 考虑国内居民家庭实际拥有情况和品牌型号的差异, 选取普及率较高典型家用电器的测试结果列入表1。
4 应用研究
4.1 数据来源
采用电力录波设备对上海某居民小区的4户家庭进行了现场测量, 如图4所示, 具有5个电流测量点和1个电压测量点。电流测量点分别为4户家庭供电进户线和总线, 电压测量点为公共连接点。对上述测点进行了为期1个月的测试, 采用本文方法对测试数据开展相关研究。
4.2 4户家庭负荷曲线和谐波幅值曲线
4户家庭及总馈线某周日连续24h的15min平均负荷曲线见附录A图A2, 4户家庭及总馈线该日连续24h的3次谐波5s平均幅值曲线见附录A图A3, 下文中图5数据亦为该日数据。
由于居民负荷谐波中主要是较低次频率的谐波, 故分析主要针对3, 5, 7次谐波电流幅值概率分布进行研究。以3次谐波为例, 由附录A图A2和图A3可以看出, 谐波伴随非线性用电负荷的使用同步产生, 电器产生谐波远大于系统背景谐波。
4.3 谐波电流四象限分布
该日中5个测量点的3次谐波电流四象限分布见附录A图A4, 横轴和纵轴分别为谐波电流实部和虚部。由图可知: (1) 3次谐波电流主要集中在第二、三象限, 表明大部分时间居民负荷谐波电流大于系统背景谐波电流, 经矢量叠加后表现为从负荷侧流向系统侧。由于本文研究由居民负荷侧产生的谐波, 因此将仅考虑位于第二、三象限的采样数据。 (2) 3次谐波电流分布没有一定的规律, 不具有典型特定的形状 (如环形、圆形、扇形等) , 不服从特定的概率分布, 这一特征对5次、7次谐波电流同样适用。
4.4 谐波电流概率分布
由测点5的谐波电流测量数据, 可绘制3, 5, 7次谐波电流幅值概率分布曲线图, 如图5所示, 可以看出: (1) 分布曲线和正态分布、瑞利分布等差异明显, 不服从某一特定的分布; (2) 分布曲线具有多峰性。分布曲线的多峰性主要由两方面因素引起, 一方面是由于不同居民电器的谐波电流发射特性的差异性, 例如, 微波炉和荧光灯产生的谐波电流幅值差异明显。另一方面是由于不同居民电器使用概率的差异性, 例如, 微波炉在一天中使用时间一般比荧光灯短, 表现为其概率密度很低。在概率分布曲线上, 微波炉所产生的谐波电流幅值大而概率密度低, 荧光灯所产生的谐波电流幅值小而概率密度高。
4.5 基于云分布和云变换的谐波电流估计
以测点5的3, 5, 7次谐波电流幅值测量数据为例呈现基于云分布和云变换的谐波电流估计方法。由于拉盖尔多项式法和联合概率分布方法不仅要求足够大的谐波样本数, 而且需要准确的网络参数, 在居民配电网中难以满足, 因此本文和较实用的正态分布方法的分析结果进行横向对比。
基于本文第1节和第2节所提理论方法, 可分别绘制测点5的3, 5, 7次谐波电流幅值云分布曲线、正态分布曲线和云变换拟合曲线, 为了同时反应居民用电周末和工作日的用电差异情况, 分析时选取了前述算例日 (周日) 和某工作日 (周二) 的数据进行分析, 见图5、图6。基于云变换的拟合曲线和实际分布曲线的参数见表2。
可见, 云分布和正态分布曲线形态较为接近, 两者均拟合成单峰, 不能反映实际分布曲线的多峰性;基于云变换的拟合曲线能较好地逼近实际分布曲线, 反映出实际分布的多峰性。除3次谐波的峰1外 (相对误差大但绝对误差很小) , 各次谐波电流幅值、概率密度实测值与云算法估计值相对误差均较小, 谐波电流幅值相对误差基本均在10%以内, 在谐波电流幅值估计误差允许的范围内, 云变换法得到的逼近曲线能够完全拟合谐波实测值的频率分布曲线, 准确反映了实际谐波的分布变化趋势, 能拟合实际多峰情况, 清晰反映谐波电流幅值分布, 整体拟合精度与生成云的个数有关, 云的个数越多, 精度越高。
基于分布曲线, 可分别计算出实测分布、正态分布、云分布曲线和云变换拟合曲线所对应的95%概率最大值, 如表3所示。由表3中可以看出, 云变换法在谐波电流幅值95%概率值计算上误差显著小于正态分布法和云分布法。
由图5、图6和表3可以看出, 周末的各次谐波电流幅值的95%概率最大值均显著高于工作日, 即周末居民用电谐波发射水平显著高于工作日。
由表2可以看出, 周末和工作日居民负荷谐波电流幅值分布中各峰值的谐波电流幅值基本相当, 如:3次谐波最大电流峰值中周末曲线峰5谐波电流幅值 (3.02A) 和工作日曲线峰4谐波电流幅值 (2.90 A) 相当;周末曲线峰2谐波电流幅值 (0.42A) 和工作日曲线峰3谐波电流幅值 (0.45A) 相当;周末曲线峰1谐波电流幅值 (0.10A) 和工作日曲线峰1谐波电流幅值 (0.08A) 相当, 是相同负荷电器谐波源的直观反映。
根据各次谐波电流幅值多峰性概率分布曲线, 可用于判别产生该谐波电流的电器类型或电器类型组合。如图5 (a) 所示3次谐波电流分布曲线具有5个峰, 各个峰所对应的谐波电流幅值和概率密度见表2。谐波是伴随负荷电器的使用产生的, 每个电器谐波发射电流幅值累积为测试点谐波电流幅值, 电器的使用时长决定了该类型电器谐波发射电流幅值的概率密度, 故谐波电流幅值概率分布曲线的横纵坐标分别对应电器谐波电流幅值和其使用时长。和表1进行对比发现, 峰1所对应的谐波其幅值约为0.10A、概率密度约1.18, 幅值较小但概率密度最大, 由使用时间较长但谐波电流发射水平低的荧光灯、电视机等电器所产生;峰5所对应的谐波其幅值约3.02A、概率密度约0.02, 幅值大但概率密度很低, 由使用时间较短但谐波电流发射水平高的微波炉所产生。谐波电流幅值概率分布曲线峰值和居民电器间具有很强的匹配性, 这对识别谐波来源和制定谐波抑制措施具有重要意义。
由于居民负荷电器种类的多样性和用电习惯的随机性, 导致居民负荷谐波不服从特定分布且具有多峰性, 现有居民谐波水平评估一般以集体谐波概率最大值, 即超过某一概率值 (如95%) 谐波电压或电流的大小, 作为谐波水平评估标准, 对谐波的概率性分布关注较少。本文所提方法在谐波电流幅值概率密度曲线形态、拟合精度和95%概率最大值计算精度三方面均显著优于现有方法, 即为居民负荷功率特性、使用时长特征及负荷类型的反映。由图5可见, 基于云变换的拟合曲线能较好地逼近实际分布曲线, 和实际曲线形态一致, 反映出实际分布的多峰性, 显著优于只能模拟单峰的正态分布法;由表2可见, 和实际值相比, 云算法谐波电流幅值及概率密度值相对误差水平较小, 在工程允许误差范围内;由表3可以看出, 95%概率最大值计算精度显著优于正态分布法, 相对误差最大仅4.27%。证明了所提方法应用于居民负荷谐波发射概率性分布分析中具有可行性, 但由于实际居民用电的随机性太强, 要对分析结果进行严格校验在工程上很难实现。
注:“—”表示无此数据。
5 结语
本文提出一种基于云理论的居民负荷谐波集体影响的评估方法。研究表明, 居民负荷谐波电流分布没有一定的规律, 不服从特定的概率分布, 具有多峰性, 采用正态分布方法进行评估存在较大的误差。本文所提云变换方法可以精确拟合居民负荷谐波的概率分布曲线。根据居民负荷概率分布曲线特征, 有可能判别产生该谐波电流的电器或者电器类型或不同组合类型, 对识别谐波来源和制定谐波抑制措施具有重要意义。
附录见本刊网络版 (http://www.aeps-info.com/aeps/ch/index.aspx) 。
摘要：电力电子型居民负荷和新能源装置在配电网中应用越来越广泛, 其产生的集体谐波非常可观, 对电力系统的影响与危害已不可忽视。文中提出一种基于云理论的居民负荷集体谐波评估方法, 介绍了云模型和云变换的基本理论, 阐述了基于云理论的居民负荷集体谐波评估流程。通过实际案例证明, 采用云变换方法可精确拟合谐波电流概率分布曲线, 能较好地呈现居民负荷谐波电流幅值分布曲线的多峰性, 其95%概率最大值估计精度显著优于正态分布法。研究结果对分析居民谐波来源、确定谐波治理措施和监管政策具有重要意义。
混合概率分布 第5篇
高压输电网络中,注入转移分布因子(ISF)体现了支路有功潮流与各节点有功注入之间近似的线性映射关系[1,2],表示了系统内某一节点有功注入的变化对支路有功潮流变化的影响。电力网络中几种重要的线性化分布因子,如功率传输分布因子(PTDF),支路开断分布因子(LODF)以及开断传输分布因子(OTDF),都可以由注入转移分布因子推导获得[2]。因此,注入转移分布因子在电力系统最优潮流[3]、阻塞管理[4,5,6]和安全评估[7,8]等应用中具有重要的作用。显然,研究如何获取更为准确、合理的注入转移分布因子具有重要的现实意义。
传统的注入转移分布因子估计方法是以直流潮流模型为基础的一种确定性估计方法[2,9,10,11]。在该方法中,注入转移分布因子可直接由支路电抗矩阵推导得到。然而,基于直流潮流模型的估计方法可能存在元件参数不准确、平衡节点设置与系统实际平衡策略不符、信息更新不及时、难以追踪系统运行点及拓扑变化等缺陷[12,13,14],从而存在较大的估计误差。
随着电力系统监测技术的发展[15,16],系统运行人员能够通过数据采集与监控(SCADA)或广域测量系统(WAMS)获取大量以固定时间间隔采集的量测数据,支路有功潮流与节点有功注入数据均可从中筛选出来,从而使利用量测数据进行注入转移分布因子的估计成为可能。根据支路有功潮流与节点有功注入之间的近似线性关系,文献[13-14,17] 提出了多种利用量测数据基于最小二乘原理的注入转移分布因子的估计方法。此类确定性最小二乘估计方法具有如下特点:1仅需使用量测数据进行注入转移分布因子的估计,避免了由于支路参数不准确、拓扑结构更新不及时所引起的估计误差[14]; 2因子求取过程中,不再预置潮流平衡节点[18,19],避免了由于平衡节点设置与实际系统平衡策略不符所造成的估计误差;3量测数据蕴含电力系统运行状态及拓扑结构变化的信息[20,21],利用量测数据进行注入转移分布因子估计的方法,可以自动、实时地适应系统运行状态与拓扑结构的变化调整节点注入与线路潮流间的映射关系[14]。以上特点使得基于量测的最小二乘估计方法的估计精度优于基于物理模型的直流潮流方法[13,14,17]。
然而,需要注意的是,即使基于量测的注入转移分布因子估计方法更加准确,但其估计结果仍可能受到测量噪声与线性化假设的影响而存在误差。注入转移分布因子的估计误差会造成对支路潮流估计的偏差,进而影响到电力系统安全分析与调度决策结果的有效性。文献[22]研究了注入转移分布因子不确定情况下的鲁棒性调度决策方法,阐明了注入转移分布因子估计误差对于调度决策结果的显著影响,但能够量化给出估计误差分布情况的注入转移分布因子的不确定性估计方法目前尚未见诸文献。
针对上述研究现状,本文提出一种依据量测信息的注入转移分布因子的概率估计方法。方法首先在输电网络线性化假设条件下,以注入转移分布因子为回归系数,建立了支路有功潮流关于各节点有功注入的线性回归模型;进而,根据贝叶斯线性回归理论[23,24,25]建立注入转移分布因子的后验概率分布模型;同时,为适应大规模电力系统中高维注入转移分布因子的概率估计需求,采用吉布斯采样数值解法[26,27]对线性贝叶斯回归模型进行求解,最终得出注入转移分布因子的后验概率分布。
1注入转移分布因子的概率估计原理
1.1含有注入转移分布因子的线性回归模型
对于输电网络,以下线性化条件可近似满足[2]: 1支路电阻值远小于电抗值;2支路两端电压相角差很小,可以忽略;3系统中各节点电压标幺值约等于1。因而在忽略线路对地导纳的前提下,支路有功潮流与各节点有功注入之间存在着近似的线性关系,即系统内某一支路有功潮流可以近似表达为各节点有功注入的线性组合[13]:
式中:PBranch,k为流过支路k的有功潮流;PNode为对应系统内N个节点有功注入的N维列向量;Mk为对应于支路k的N维注入转移分布因子向量,表示系统内N个节点的有功注入对支路k上有功潮流的影响。
随着电力系统量测技术的发展,SCADA系统或WAMS可以提供大量的时序量测数据,从中容易筛选出支路有功潮流与节点有功注入值,充足的数据使采用回归方式来估计注入转移分布因子成为可能。考虑到回归估计中残差的存在,在式(1)的基础上,可以建立包含注入转移分布因子向量的线性回归模型:
式中:εk为服从正态分布N(0,σk2)的残差变量。
1.2注入转移分布因子的贝叶斯概率估计
贝叶斯线性回归[23]是贝叶斯推理中针对线性回归模型的一种方法,它可以用于对线性回归模型中未知系数的概率估计,以获取系数的概率分布。
式(2)中待估计的变量为关于支路k的注入转移分布因子向量Mk和回归残差的方差标量σk2,为不失一般性,视各待估计变量为相互独立的随机变量。根据贝叶斯公式,所有随机变量的联合后验概率密度函数可以表示为:
式中:P (Mk,σk2|PBranch,k,PNode)为给定观测值PBranch,k与PNode下Mk与σk2的联合后验概率密度函数;P(Mk,σk2)为关于Mk与σk2的联合先验概率密度函数;P (PBranch,k|Mk,σk2,PNode)为似然函数; P(PBranch,k)为通过统计获得的支路有功潮流PBranch,k的边缘概率密度函数。
易知在已知式(3)所示联合概率分布的情况下, 支路k注入转移分布因子向量Mk的联合后验概率分布可以通过如下积分获得:
式中:P(Mk|PBranch,k,PNode)为注入转移分布因子向量Mk的联合后验概率分布。
式(3)、式(4)提供了一种注入转移分布因子向量联合后验概率分布的估计方式。然而,由于当前电力系统规模庞大,注入转移分布因子向量维数通常较高,式(4)难以通过解析的方式求解出来。因此,本文采用吉布斯采样数值算法来近似求解各注入转移分布因子的后验概率分布。
2注入转移分布因子后验概率分布的数值解法
吉布斯采样是马尔可夫链蒙特卡洛算法中的一种,它通过从各变量已知的条件分布中采样,统计采样值,近似得出全部(部分)变量的联合分布或单一变量的边缘分布[26,27]。该方法可以避免复杂的解析推导,且输入、输出明确,易于编程快速实现对注入转移分布因子后验概率分布的估计。
吉布斯采样算法的思路是根据式(3)所示的联合概率分布,在Mk与σk2其一为已知的条件下,对另一部分进行采样,如此迭代,最终统计形成各自的概率分布。
如图1所示,吉布斯采样过程主要包括如下4个步骤。
步骤1:初始化。为σk2赋初值并设置迭代计数变量,开始采样过程。此处,σk2的初值由普通最小二乘估计获得。
步骤2:对Mk的联合后验概率分布采样。当 σk2已知时,式(3)中关于Mk与σk2的联合后验分布将成为仅关于向量Mk的联合后验分布。根据贝叶斯公式,为了得到Mk的联合后验概率分布表达式, 必须首先设定其先验分布与似然函数。
在贝叶斯线性回归中,Mk的先验分布通常被设置为多维正态分布,可表达为:
式中:Zk为N维均值列向量;Σk为N ×N阶对角协方差阵。
本文中,Zk被设置为零向量,Σk的对角线元素赋值为1 000σk2。由于式(2)中残差变量εk服从正态分布N(0,σk2),因而式(3)中的似然函数为:
式中:σk2的当前值由其最新采样值来更新(初次迭代为所赋初值)。
由先验分布与似然函数表达式,根据式(3)所示关系,Mk的后验概率分布将是一个联合正态分布, 该分布的均值向量Zk*与协方差矩阵Σk*分别为:
从而,可通过标准多维随机变量采样程序[24], 获得向量Mk的采样值。存储采样值,并用它们更新Mk的当前值。
步骤3:对σk2的后验概率分布采样。同样,当Mk由其采样值替代时,式(3)将成为关于σk2的边缘分布。在贝叶斯线性回归中,σk2通常被设定为服从逆伽马先验分布:
式中:Tk为先验分布的形状参数;θk为先验分布的规模参数;Γ(Tk)为伽马函数在Tk处的取值。本文中,Tk与θk的值都设置为0.5。
将式(9)与式(6)代入式(3)中化简后发现,σ2k的后验概率分布也是逆伽马分布,该分布的形状参数为Tk+0.5,规模参数为θk+0.5(PBranch,kPNodeMk)T(PBranch,k-PNodeMk)。
根据σk2的边缘后验概率分布函数,通过单变量采样程序便可以获取σk2的采样值,存储采样值。
步骤4:收敛检验。检验采样次数是否达到设定值,如果已经达到设定次数,结束采样程序,输出存储的Mk与σk2的采样值;如果未达到,以最新采样值更新σk2的当前值并返回步骤2。
步骤1至步骤4描述了对于注入转移分布因子向量Mk的吉布斯采样过程。最后,通过对向量Mk中各元素采样值的统计,便可得到各注入转移分布因子的边缘概率分布。
3算例分析
本文利用截取于河南省220kV输电网络的21节点41支路系统中的实际SCADA量测数据,验证本文方法的有效性。系统的拓扑结构如图2所示。
3.1注入转移分布因子概率估计
选取301组时序量测数据(采集间隔为5 min) 对注入转移分布因子进行概率估计,另301组数据用于交叉检验。为简洁起见,将这两组数据分别记为数据集A与数据集B。
从数据集A中筛选出PBranch,k与PNode的观察值并代入注入转移分布因子的贝叶斯后验概率分布模型,利用吉布斯采样方法获得注入转移分布因子向量Mk的多组采样值。为保证估计效果,采样次数设置为18 000次。
通过统计注入转移分布因子采样值,可获得其概率估计结果。表1给出了关于支路21的6个注入转移分布因子部分分位点的估计结果。以表1中M2115与M2116为例(其中上标15,16为节点号),其采样轨迹如图3所示。
轨迹图反映了吉布斯采样器的采样过程,图中所示各样本值经统计即可获得图4中的概率密度曲线,两图从不同角度反映出注入转移分布因子一致的不确定分布情况。
理论上来说,单个注入转移分布因子的边缘后验概率分布服从t分布[23],该分布具有样本容量越大越趋近于正态分布的性质。算例中注入转移分布因子的采样次数为18 000次,样本容量较大,因此, 图4中注入转移分布因子的概率密度分布形状趋近于正态分布,这直观地表明吉布斯采样值获得的后验概率分布良好地收敛到了其实际分布,显示了本文方法的有效性。
综上,注入转移分布因子的概率估计结果可以总结为表2。
由表2可见,概率估计可以获得注入转移分布因子的中位数与期望值,在2.5%与97.5%分位点之间的注入转移分布因子值界定了其95%置信区间。 例如,M2115的95% 置信区间为[-0.042 3, -0.018 2],这意味着因子M2115取该区间以外的值为小概率事件,实际运行中几乎不会发生。关于支路21的各个注入转移分布因子的95%置信区间长度较短,因子波动较小,这表明高压输电网络确实具有较高的线性化程度。另外值得注意的是,以上求得的注入转移分布因子的概率分布或不确定性区间体现了本文所提出的注入转移分布因子概率估计方法与已有研究中确定性估计方法的本质区别。
3.2与注入转移分布因子确定性估计方法的对比
基于直流潮流模型的估计方法[2,9,10,11]及利用量测数据的最小二乘相关的估计方法[13,14,17]是两类主要的注入转移分布因子的确定性估计方法。通过直流估计方法,一组固定的注入转移分布因子值可以从支路导纳矩阵中直接推导得出。而基于数据的最小二乘相关估计方法随着采集得到新数据,在线更新注入转移分布因子估计值。其中,文献[13]使用潮流模拟实验证明了基于测量的最小二乘估计方法相比于基于直流潮流模型的估计方法优势明显。此处,本文利用河南电网某输电网络的实际参数与运行数据,将两种确定性方法与本文方法进行对比,以说明本文方法的有效性。
仍以支路21相关的注入转移分布因子为例进行说明。3种方法获得的注入转移分布因子估计结果列于表3。其中,直流潮流估计方法中选取节点16为参考节点。
直流估计方法估计结果依赖于直流潮流模型, 基于物理模型估计方法的固有缺陷将降低其估计的精度。例如:在直流潮流估计方法获得的结果中,与参考节点相关的注入转移分布因子均为零值[10],当网络实际平衡策略与计算中参考节点设置不一致时,直流估计获得的注入转移分布因子将是不准确的。相反,最小二乘估计与本文提出的概率估计方法均使用量测数据进行估计,使这一问题得以避免。 观察表3中最小二乘估计结果及概率估计期望值可见,节点11的注入功率对支路21的有功潮流影响最弱,相应注入转移分布因子值接近于零,然而直流估计结果显示,节点11的功率注入对支路21的有功潮流影响最强,这显然与基于数据的估计方法结果相冲突,显示了直流潮流估计方法的不足。
此外,由表3可见,注入转移分布因子的最小二乘估计结果接近概率估计期望值,这表明概率估计的期望值完全可以作为注入转移分布因子的点估计结果加以使用。
3.3支路有功潮流预估结果对比
注入转移分布因子的估计效果可从其对支路潮流估计的准确程度反映出来。在系统任意注入情况下,根据式(2),利用Mk与σk2的采样值,即可获得支路k有功潮流样本并统计得到其概率分布。由注入转移分布因子概率估计结果获得的支路有功潮流概率分布包含更加丰富的信息,可以为调度人员预防线路过载情况发生提供更为全面的决策依据。
从数据集B中取出一组节点注入数据进行说明,其中节点11~16的有功注入分别为-204.68, -94.48,141.07,-104.25,76.84,876.70 MW。仍以支路21为例,由已知注入转移分布因子概率估计结果,可计算得出支路21有功潮流的期望值、2.5% 分位点、97.5% 分位点的概率预测结果分别为152.051 8,149.745 4,154.386 2MW,而其量测值为151.537 0 MW。可见,支路21所传输有功功率的真实量测值在概率估计值的95%置信区间内,且预测期望值接近量测值,是量测值的良好点估计结果。
按照这种方式,应用数据集B中301组量测数据进行交叉检验,对比由概率估计、直流估计与最小二乘估计方法获得的注入转移分布因子应用于支路有功潮流预估的有效性。检验结果如图5所示。
由图5可见,依据3种注入转移分布因子估计结果得到的支路有功潮流均可追踪支路潮流量测值的变化趋势。但同时可见,本文提出的概率估计方法与最小二乘估计方法和直流估计方法相比,追踪效果更好,3种方法对应的平均绝对误差(MAE)指标分别为1.023 1,1.133 7,2.352 0,也验证了这一直观感觉。另外,图5中的黄色分布带表示概率估计注入转移分布因子计算所得支路有功潮流预测值的95%置信区间,在测试时间段内,实际有96.01%的功率点落于该分布带中,体现了由本文注入转移分布因子概率估计结果所得支路有功潮流概率分布的合理性。
另外,本文方法估计所得注入转移分布因子计算支路有功潮流与交流潮流计算结果的对比见附录A,系统拓扑结构变化后本文方法的适用性测试见附录B。
4结语
本文提出一种基于量测数据的注入转移分布因子的概率估计方法。方法依据贝叶斯线性回归理论建立了注入转移分布因子的概率分布估计模型,并采用吉布斯采样数值解法对模型进行了有效求解。 文中以河南省局部输电系统的实际运行数据,对方法的有效性进行了测试,测试结果表明:相比较于基于直流潮流、最小二乘的注入转移分布因子估计方法,本文方法能够提供更为准确的注入转移分布因子的点估计结果,同时,估计得到的注入转移分布因子的概率分布结果可以真实地反映注入转移分布因子估计结果的偏差范围,为调度人员提供更为全面、 可靠的电网功率分布信息,以预防输电线路过载情况的发生。
混合概率分布 第6篇
关键词：结构可靠性,敏感性分析,混合不确定性,概率,区间
0 引言
传统的结构可靠性分析[1,2,3,4]基于概率方法, 需要大量实验样本构造参数的精确概率分布。为弥补概率方法对大样本量的依赖性, 近年来国内外学者提出了两类基于概率与非概率区间[5,6]的混合可靠性模型, 用于信息量缺乏的复杂结构的可靠性分析。第一种混合模型将所有参数分成两类, 一类具有足够的信息, 另一类信息不充分, 并分别使用随机分布和区间来描述两类参数。如文献[7]基于两阶极限状态方程, 提出了一种处理线性问题的混合可靠性分析方法;文献[8]通过综合概率与非概率方法, 得出一种针对不确定结构的鲁棒性设计技术;文献[9]提出了一种求解概率-区间混合问题的高效可靠性分析方法;文献[10]构建了一种概率-凸模型混合可靠性模型, 并提出了高效求解算法。第二种混合模型中的不确定变量使用带有区间参数的概率分布来表示。该模型首先由Elishakoff等[11]提出, Elishakoff等探讨了随机振动结构最差响应的反优化问题;文献[12]将混合不确定模型应用于结构系统, 并计算了由区间参数导致的失效概率区间。
目前概率-区间混合可靠性分析的研究总体上仍处于起步阶段, 而对其敏感性分析的研究更为缺乏。对于实际结构而言, 可靠度敏感性分析具有非常重要的工程意义。通过敏感性分析, 可定量获得相关参数对于结构可靠性的影响程度, 从而有效指导实际结构的设计。针对第一种混合模型, 文献[13]初步构造了几类敏感性分析指标, 用以分析结构可靠度对于区间参数的影响程度。本文借鉴了该方法的研究思路, 提出了一种针对第二类混合模型的可靠度敏感性分析方法, 定量获得了结构失效概率区间对于重要区间分布参数的一阶梯度。
1 传统的一次二阶矩可靠度分析方法
对于只含有独立随机变量X的结构, 可靠性定义如下:
式中, Pr{}为求解概率的函数;X为由独立随机变量构成的向量;g为功能函数。
使用一次二阶矩方法 (first order reliability method, FORM) 计算可靠性时, 需要将原始空间中的随机变量X= (x1, x2, , xn) 映射到标准正态空间U= (u1, u2, , un) 中[14], 即
式中, 为随机变量xi的累积概率函数;ф-1为标准正态分布累积概率函数的逆函数。
由此可以得到功能函数的映射关系:
其中, T为随机变量从原始空间转换到标准空间的转换函数。可进一步转换为如下优化问题进行分析:
通过上式获得的最优点U*即为最大可能点 (most probable point, MPP) 。β=‖U*‖为可靠度指标, 通过公式, 可近似获得可靠度。
2 概率-区间混合可靠性分析
对于许多实际工程问题, 难以获得足够的信息来构造某些不确定参数的精确概率分布, 而假设的概率分布又可能导致很大的分析误差。本文使用混合模型对结构中的不确定性进行度量, 即使用概率分布描述不确定参数, 而其中的关键分布参数因为信息的缺失只能给定变化区间而非精确值。
2.1 带区间分布参数的随机变量可靠性模型
对于一带区间分布参数向量Y的随机向量X, 从原始空间转换到标准正态空间的过程可表示如下:
式中, Y为区间分布参数构成的m维向量。
通过上述转换, 可以得到标准正态空间中的极限状态方程:
由于区间参数Y的影响, G (U, Y) =0所表示的极限状态在参数空间中不再是一个曲面, 而是由两个边界面构成的带, 如图1所示, 图中距离原点最近与最远的极限状态曲面具有不同的Y值。通过分析两边界面的可靠性, 可以得出混合可靠性指标β:
式中, βU和βL分别为上下边界面的可靠度指标, ML和MU分别为上下边界面上的最大可能点, 对应的区间向量值分别为YU和YL。
混合可靠度指标β不再是确定值, 而是一波动区间, 这种可靠度指标的不确定性是由于引入区间分布参数而导致的。由此, 可获得结构失效概率区间:
为计算边界βL与βU, 可构造如下两个优化问题:
与
显然, βL与βU的求解为典型的嵌套优化问题。
2.2 混合可靠性的求解
式 (7) 所表示的等概率转换过程中, 累积概率关于区间分布参数的单调性对可靠性区间的分析和计算至关重要, 文献[15]已详细分析了区间分布参数与极限状态边界面之间的关系。对于极限状态方程中仅含单个区间分布参数的问题, 有如下性质[15]:
性质1如果某随机变量的累积概率函数对其区间分布参数单调, 那么极限状态面映射到标准正态空间所得区域的两条边界面必然分别对应区间分布参数的两个边界, 故M点 (最大可能点, 包括ML与MU) 必在分布参数区间的边界上达到。
对于多区间变量问题有如下性质:
性质2如果极限状态方程中的多个随机变量的累积概率函数分别对各自的区间分布参数单调, 则两个边界面上的M点的区间分布参数取值一定是其边界值的组合。
在上述性质中, 不确定变量的累计概率关于其分布参数的单调性具有非常重要的作用。对工程中常用的随机分布单调性进行分析可发现[15], 所有常用随机变量的累计概率都为其分布参数的单调函数, 故性质1和性质2对一般问题是成立的。
对于式 (10) 和式 (11) 中的优化问题, 本文基于文献的解耦策略[9], 构造一高效的序列迭代格式进行求解, 即内外层相对独立地进行迭代, 内层为以区间分布参数为变量的极值求解, 外层为以随机参数为变量的可靠度计算, 最终使内外层同时达到稳定解。基于性质1和性质2, 可将内层优化转换为如下线性规划问题而不失精度:
式中, Uk、Yk分别为U和Y的第k步迭代的初始值;Yi是向量Y中的第i个值, 外层优化采用HL-RF迭代:
算法流程如下:
(1) 输入初始点U0与Y0;
(2) 建立线性规划问题如式 (12) 所示, 并利用单纯形法进行求解, 获得Yk;
(3) 利用式 (13) 计算Uk;
(4) 检验是否收敛。如果|βk-βk-1|/|βk-1|ε1与|G (Uk) |ε2同时成立, 则进入下一步;否则, 令k←k+1, 回到步骤 (2) 。其中ε1与ε2为给定的较小值;
(5) 由可靠性指标βL得出失效概率上边界PfU。
上述由式 (12) 和式 (13) 组成的迭代过程仅求解了PfU。对于失效概率下边界PfL可完全按照类似的过程进行求解, 仅需要将式 (12) 改为最大化问题即可。
3 混合可靠度敏感性分析
传统的可靠度敏感性分析通常指可靠度的值对某个量的一阶梯度。本文的问题由于受区间参数的影响, 分析得出的可靠度指标不再是特定值, 而是一区间, 为此分析中将研究结构可靠度区间对于相关参数的敏感程度。具体分析之前, 先作如下规定:
式中, δP、分别为失效概率区间的宽度与中心值;δi、分别为分布参数区间的宽度与中心值。
下面将推导出如表1所示的6类敏感性指标的计算格式。
由上文中的单调性分析已知, 对于工程中常用的随机分布, 失效概率的最大值与最小值只会出现在分布参数区间的边界上, 所以6类敏感性指标的推导过程中都分两种情况进行分析:①PfU出现在YiU处, PfL出现在YiL处;②PfU出现在YiL处, PfL出现在YiU处。并作如下定义:
(1) 类型一:为失效概率区间宽度对分布参数区间宽度的敏感性。
①PfU出现在YiU处, PfL出现在YiL处:
下一步计算。此时, M点在区间分布参数的边界上, 令Pfb表示PfL或PfU, Yib表示YiL或YiU, 可得
其中, u为参数Yi所属随机变量在标准正态空间中对应的坐标值。可见, 问题最终转化为u/Yib的求解。基于等概率变换, 可知
针对工程中常用的随机分布, 式 (19) 可给出解析表达式。
②PfU出现在YiL处, PfL出现在YiU处:
问题最终转化为求解。
(2) 类型二:为失效概率中心值对分布参数区间宽度的敏感性。
①PfU出现在YiU处, PfL出现在YiL处:
②PfU出现在YiL处, PfL出现在YiU处:
(3) 类型三:为失效概率的区间宽度对分布参数区间中心值的敏感性。
①PfU出现在YiU处, PfL出现在YiL处:
②PfU出现在YiL处, PfL出现在YiU处:
(4) 类型四:为失效概率的区间中心值对分布参数区间中心值的敏感性。
①PfU出现在YiU处, PfL出现在YiL处:
②PfU出现在YiL处, PfL出现在YiU处:
(5) 类型五:为失效概率的区间宽度对确定性分布参数的敏感性:
(6) 类型六:为失效概率的区间中心值对确定性分布参数的敏感性:
4 算例
4.1 管状悬臂梁
图2所示的管状悬臂梁[9]受到外力F1、F2、P, 以及扭矩T的作用。功能函数定义为强度Sy与固定端圆周下表面处最大应力σmax之差:
模型中4个固定参数取值为:L1=120mm, L2=60mm, θ1=5°, θ2=10°。
σmax表达式如下:
式中, τzx为剪切力。
其中σx按照下式计算:
式中, A为横截面的面积;I为横截面对圆心的极惯性矩;M为弯矩;c为横截面到圆心的距离。
式 (32) 中第一项是由轴力引起的应力, 第二项是由弯矩M引起的应力。
不确定参数及其分布如表2所示。该问题包含两种分布类型:正态分布与极值Ⅰ型分布, 共7个随机变量;每个随机变量包含一个区间参数, 共7个区间变量。
敏感性分析结果如表3所示, 表中每一项结果都同时给出了中心差分法的计算值。可以发现, 本文方法获得结果与中心差分法获得的结果非常接近, 表明本文方法具有较高的计算精度。另外, 利用中心差分法计算时, 需要重新进行结构可靠性的求解, 而采用本文方法时, 只需要在原可靠性分析的基础上即可快速获得敏感性结果, 效率较高, 为复杂结构的分析及可靠性优化设计提供了高效计算工具。
另外, 为定量比较各参数对结构可靠性的敏感程度, 对表3中的结果进行了量纲一化处理:
式中, s表示变量a对变量b的敏感性指标进行量纲一化之后的结果;b0为b的名义值。
量纲一化后敏感性分析结果如表4所示, 通过量纲一分析, 任何一种敏感性类型中参数的敏感程度可进行定量比较。如类型Ⅰ中, 随机变量d的均值所对应敏感性指标最大, 表明随着d区间宽度变化, 失效概率区间的宽度相对其他几个区间分布参数以最快速度变化。所以, 为保证结构的失效概率波动幅度较小, 可考虑将管径d均值的波动控制在一定范围之内。
4.2 车门结构
汽车车身是整车最重要的组成部分之一, 而车门是车身的重要部件。车门是由薄板冲压成型并通过焊接连成一个整体的受力结构, 承担载荷的部件主要有外门板、内门板、上加强板、下加强板等。刚性好、不易下沉是衡量车门性能的一项重要指标。图3所示为某款轿车的前门结构, 考虑车门在门锁位置承受900N垂直载荷的工况, 因车门变形后垂直方向的最大位移可表征其刚度大小, 所以选择该方向最大挠度作为功能函数, 建立极限状态方程:
式中, T为板料厚度;F (T) 为车门挠度值;f为车门挠度的上限, 这里取4.169mm。
建立有限元模型计算挠度值, 采用壳单元, 单元数目为18 500。
选择车门上五个重要部件进行分析, 分别是:门外板, 门内板, 上防撞梁, 下防撞梁, 窗框下方加强板。由于受板料制造精度以及冲压过程的影响, 板材厚度存在波动性。本算例中以正态分布度量料厚的波动, 由于缺乏实验数据, 给定均值和标准差变动区间而非精确值, 各个随机变量的分布参数取值情况如表5所示。每一个随机变量均有一个分布参数具有不确定性, 前三个随机变量的均值是区间变量, 后两个变量的标准差为区间变量。由于可靠性分析过程需重复调用有限元计算, 将导致计算量过大的问题。所以, 前期通过拉丁超立方采样得出30个样本点, 在此基础上利用径向基函数[16]构造响应面, 并通过多个验证点校核响应面的精度。
mm
基于该响应面, 使用本文方法计算了四类敏感性指标, 得到量纲一化之后的结果, 如表6所示。针对T1的均值得出的四个敏感性指标, 类型Ⅳ的绝对值最大, 而且为负值, 表明:随着正态分布随机变量T1取值区间中心值逐渐增大, 平均失效概率以较快的速度减小, 即增大车门外板的料厚, 会使整体的平均失效概率显著减小, 结构可靠性增强。针对类型Ⅰ的所有结果中, 随机变量T2的均值所对应的敏感性指标最大, 表明随着T2的均值取值区间宽度变大, 失效概率区间的宽度增大速度较快, 而且增速要高于其他区间参数由于宽度增大带来的影响。所以, 如需保证车门质量的稳定性, 很重要的一点就是将车门内板料厚度T2的均值波动幅度控制在一定范围之内。整个结果中, 绝对值最大的指标出现在类型Ⅳ中T2的均值项中, 其值为负。另外, 较大的敏感性指标大部分出现在内外板料厚度T1、T2处, 表明内外板厚度对车门垂向刚度贡献要高于其余几个部件, 所以在考虑车门的垂向刚度时, 需要重点关注。
5 结语
本文构建了一种基于概率-区间混合不确定模型的可靠度敏感性分析方法。给出了6种可靠度敏感性指标的计算方法, 定量分析了区间分布参数对于结构失效概率区间的影响程度。基于概率转化过程的单调性, 分两种情况给出了各个敏感性指标的计算格式。数值算例表明, 本文方法具有较高的计算效率和精度。该可靠度敏感性分析方法可用于不确定信息缺乏的实际工程结构的分析, 并有望应用于工程结构的可靠性优化设计之中。
混合概率分布 第7篇
中国经济的快速发展引发众多交通问题, 而优先发展公共交通是解决我国城市交通问题的关键之一。为落实城市公共交通的主体地位, 就必须保证公交线网的优化调整和公交营运调度质量, 为此公交运营数据的实时采集必不可少。现状公交数据的采集多依靠人工调查法, 但数据处理工作量大, 耗费成本高。而基于公交IC卡的数据采集技术受外界因素影响较小, 所得数据库信息量大、具体全面且真实可靠, 是1种快速、经济、有效的客流调查方法。与此同时, 如何利用公交IC卡数据挖掘出更多、范围更广的信息也成为了1个亟待解决的课题。
目前, 不少学者在单条公交线路客流OD矩阵反推的课题上作了大量的相关研究, 模型一般集中于概率论模型[1]、最大熵法[2]、双层网络规划模型[3]、结构化模型[4]、重力模型和最小二乘法等[5]。这些方法主要分析了公交乘客出行行为的规律, 对公交站点影响乘客下车概率的研究较少[6,7,8,9]。其中, 赵锦焕等[1]在研究公交乘客下车概率时, 提出1个重要的影响因素-公交站点换乘功能, 但其量化表示只考虑该线路站点的公交换乘条数, 未考虑周边站点对其站点的辐射影响。胡郁葱等[10]在IC卡数据挖掘和统计基础上, 对线路站点OD矩阵、区域OD矩阵进行了推算。陈素平等[11]在基于小票法调查数据的基础上提出1种公交线路OD分析方法, 可准确得出公交线路的客流OD, 同时还可以得到线路沿途各站点的上下客人数、断面流量等信息。
公交线路客流OD矩阵的分布是基于公交站点建立的, 不同的站点对公交客流量的吸引强度不同, 因此, 站点吸引是反推公交线路客流OD矩阵的首要考虑问题。而站点吸引是公交站点对公交客流量吸引强度的简称, 其影响因子有很多。
笔者在以上研究的基础上, 全面分析了公交乘客的出行行为特性以及公交站点对乘客出行的影响, 提出了公交出行距离、站点用地性质以及站点换乘能力等3个影响乘客下车概率的主要因素, 构建了基于站点吸引的公交客流OD分布概率模型, 并采用C++语言编程得以实现, 以此分析公交IC卡乘客的空间出行特性。
1 模型建立
公交IC卡乘车记录的起讫点判别, 即判断刷卡记录对应的上车站点以及该上车站点对应的下车站点。起讫点的判别是公交IC数据处理的关键环节, 也是模型建立的核心。
大多数的城市公交采用全线一票制, 就是上车1次性刷卡完毕, 下车不需要再次刷卡, 因此刷卡记录中只有乘客上车站点, 不记录下车站点。由于公交定线定站的运营特征以及城市居民的出行特点, 决定了公交出行在路线的选择和客流的分布上具有一定的规律性和稳定性。利用这些特性, 笔者提出基于站点吸引的公交OD分布概率模型。
1) 下车站点的判定模型。以1条公交线路单向运行1趟的数据分析为例。
假设此线路有m个站点 (包含首末站) ;Di为i站点下车人数 (i=1, 2, …, m) ;Si为i站点上车人数 (i=1, 2, …, m) ;P′ij为某乘客在公交站点i上车, 在站点j下车的概率, (i=1, 2, …, m;j=1, 2, …, m) 。
其中, 上车站点已进行过判断, 因此, 可利用Access软件对IC卡数据进行统计分析得到Si。
根据公交单向运行的特性, 首站没有下车乘客, 因此D1=0。
在第2个站点下车的乘客来自于首站上车的乘客, 因此D2=S1×P′12。
在第3个站点下车的乘客来自于首站和第2个站点上车乘客, 因此D3=S1×P′13+S2×P′23。
依次类推得到各站点下车人数计算公式
2) 下车概率的确定方法。在分析已有的城市公交客流调查数据的基础上, 得出了影响乘客下车概率的因素主要有3个:公交出行距离;下车站点附近的土地利用性质;下车站点的换乘能力。
(1) 居民公交出行受到出行距离的影响。出行者选择公交出行时, 其出行站数主要集中在某个范围内, 属于中长距离的出行。也就是说当乘坐站数到达此范围内, 出行者下车的概率较大;反之, 当出行距离过长或过短时, 出行者很少采用公交出行的方式。可以看出, 下车概率随途经站点数量服从泊松分布。因此, 只考虑途经站点数目而得到的下车概率为
式中:Fij为乘客在i站点上车j站点下车的概率;λ为平均公交出行途经站点数量, 当i站点以后的站点数量小于平均公交出行途经站点数量时λ=m-i, m为线路单向站点个数。
(2) 公交出行受到各站点周围用地性质的影响。附近有大型商业区、火车站、汽车站等的一些站点发生与吸引的客流量较多。对于公交站点来说, 若某站点上车的人数越多, 就说明该站点的公交客流发生量越大;若某站点下车的人数越多, 说明该站点的公交客流吸引量越大。在i站发生量一定时, j站的吸引量越大, 那么从i站到j站的OD分布量越大;反之, 在j站吸引量一定时, i站的发生量越大, 那么从i站到j站的OD分布量也越大。因此, 可以利用各站点上、下车人数来量化站点周围的用地性质, 那么定义Wi为公交站点i用地吸引强度系数, 即
(3) 公交出行受到站点换乘能力的影响。由于下车站点附近的可换乘公交线路越多, 站点的换乘能力也相应越强, 同时该站点下车换乘的人也越多。因此, 可以利用各站点附近的可换乘公交线路条数来量化站点的服务水平, 即站点的换乘能力是用以站点为中心, 方圆300~500 m区域内, 可利用换乘的公交线路条数来表示。定义Ki为公交站点i的换乘能力系数, 即
式中:li为i站点的可换乘公交线路数。
以上可知, 下车概率Pij与居民公交出行距离、站点用地吸引强度和站点换乘能力正相关, 即
计算公式为
再将所得概率Pij归一化, 其公式为
将式 (6) 代入下车站点的判定模型公式 (1) , 即可得出各站点的下车人数。
3) 该模型的特点。基于站点吸引的公交客流OD分布概率模型, 是根据IC卡数据站点的上车人数, 利用下车概率理论来计算各站点的下车人数。该方法约束条件较少, 运算简单, 可得到较为准确的结果, 但是此法不能判断单个乘客的下车站点, 也不能得出单个乘客完整的出行过程。
2 案例分析
以济南市的公交IC卡运营数据为基础, 结合Microsoft Access 2010数据库技术和Microsoft Visual Studio 2010中的C++语言编程, 采用基于站点吸引的公交OD分布概率模型, 对济南市1条公交线路的公交客流情况进行分析。
2.1 公交的出行距离分布
文献[6]以某市1条公交线路为跟车调查对象, 对该线路乘客的出行站数进行统计分析, 绘制了乘客出行站数频率分布曲线图, 并通过检验发现在显著性水平0.05下该调查数据服从泊松分布。同时对该市的其他9条公交线路调查数据进行统计校验, 发现在显著性水平0.05下, 绝大多数公交线路都服从泊松分布。
笔者以济南市83路公交车为跟车调查对象, 获取该线路出行者的乘车站数统计数据, 并利用文献[6]中的结论, 采用泊松分布拟合得到相应的泊松分布参数, 这里特指83路公交线路的平均途经站点个数, 即λ=10。
2.2 公交站点可换乘线路分布
根据济南公交线路一览表, 可统计出83路公交线路中, 以各站点为中心, 方圆300~500m内所有可以换乘的公交线路条数, 统计结果见表1。
2.3 站点用地吸引强度
由于公交出行具有很强的往返性, 1个站点的发生量和吸引量基本保持平衡, 也就是该站的发生量同时可以反映站点的吸引量, 因此可以通过下行方向的上车发生率近似得到上行方向的下车吸引强度。其中, 上车发生率即是1个站点的上车人数与该线路所有站点的上车人数的比值。
以济南市83路公交出行的数据为例, 借助Access软件对其数据进行处理, 得出公交线路的各站上车人数, 统计结果分为上行方向和下行方向, 可见表3。
本文以济南市83路公交出行的数据为例, 使用Access将其数据库打开, 应用筛选、交叉查询等功能得出公交线路的各站上车人数, 统计结果分上行方向与下行方向, 见表2。
2.4 生成单条公交线路OD的数据处理
基于以上影响公交出行的3个因素的调查统计数据, 利用基于站点吸引的公交客流OD分布概率模型, 采用C++语言编程进行数据的处理, 得出83路公交线路OD分布, 见表3。
2.5 模型精度校验
该案例采用核查线法进行精度校验。具体实施方案是选用83路公交线路中的3个路段断面为基准断面, 调查基准断面的小时客流量 (包括上行方向和下行方向) 。同时, 由所得OD分布表中统计相应断面的理论日客流量, 并折算成断面小时客流量。对各断面调查的断面小时客流量与理论的断面小时客流量进行比较, 计算相对误差, 见表4。
由表4可见, 选定的3个基准断面的调查断面小时客流量与理论断面小时客流量的相对误差基本在10%以内。考虑到二者都只包括使用IC卡的乘客, 因此误差在10%以内是符合公交OD反推要求的。将同样的数据运用赵锦焕等[1]提出的概率论模型反推, 得出结果的相对误差大于10%。由此可见, 改进后的下车站点服务水平这一影响因子更加贴合实际, 大大提高了模型的精度。总体看来, 笔者所提出的基于站点吸引的公交客流OD分布概率模型能够在前人的基础上更准确地推算出下车站点的信息, 是目前比较适用的1种公交客流OD推算方法。
3 结束语
笔者从公交出行距离、下车站点的用地性质以及下车站点的服务水平等3个影响因素出发, 构建了推算单条公交线路客流OD分布的概率模型。实例分析的结果表明基于站点吸引的公交客流OD分布概率模型方便简单、易于实现, 且3个影响因素的综合考虑大大提高了公交客流OD反推的精度。
公交线路客流OD是使用IC卡乘客出行和使用纸币的乘客出行的总和。然而, 在本文中, 相关参数的获取仅仅依靠公交IC卡数据, 只包括使用IC卡的乘客出行, 而忽略了使用纸币的乘客出行。因此如何科学全面计算各参数, 从而进一步提高公交客流OD矩阵的推算精度和适用性, 是下一步研究重点。
参考文献
[1]赵锦焕, 李文权, 张健, 等.公交客流OD矩阵推算[J].交通运输工程与信息学报, 2009, 7 (4) :124-128.
[2]刘翠, 陈洪仁.公交线路客流OD矩阵推算方法研究[J].城市交通, 2007, 5 (4) :81-83.
[3]周晶, 张伦珂.利用IC卡数据估计公交OD矩阵的模型及算法[J].系统工程理论与实践, 2006 (4) :131-135.
[4]朱从坤, 丁建霆, 陈瑜.公交线路OD反推的结构化模型研究[J].哈尔滨工业大学学报, 2005, 37 (6) :851-853.
[5]戴霄.基于公交IC信息的公交数据分析方法研究[D].南京:东南大学, 2006.
[6]Doblas J, Benitez F G.An approach to estimating and updating origin-destination matrices based upon traffic counts preserving the prior structure of a survey matrix[J].Transportation Research Part B:Method Logical, 2005, 39 (1) :565-591.
[7]Park Jin Young, Kim Dong-Jun.Use of smart card data to define public transit Use in seoul, South Korea[J].Transportation Research Record:Journal of the Transportation Research Board, 2008 (263) :3-9.
[8]窦慧丽, 刘好德, 杨晓光.基于站点上下客人数的公交客流OD反推方法研究[J].交通与计算机, 2007, 25 (2) :79-81.
[9]周涛, 翟长旭, 高志刚.基于公交IC卡数据的OD推算技术研究[J].城市交通, 2007, 5 (3) :48-52.
[10]胡郁葱, 梁杰荣, 梁枫明.基于IC卡数据挖掘获取公交OD矩阵的方法[J].交通信息与安全, 2012, 30 (4) :66-70.
混合概率分布范文
声明：除非特别标注，否则均为本站原创文章，转载时请以链接形式注明文章出处。如若本站内容侵犯了原著者的合法权益，可联系本站删除。
