接手水文频率分析这个活儿之前,我对皮尔逊三型曲线的全部印象停留在课本公式上。直到有一天,手里攥着一份28年的年最大洪峰流量资料,需要在两天内给出几个设计频率对应的洪峰流量,我才发现事情没那么简单。打开Excel,均值用AVERAGE一算就出来了,Cv也不算难,可到了Cs(偏态系数)这一步,查表、插值、目估适线,再配合那谁都没画对过的“频率格纸”,整个人都麻了。后来我把整套流程搬进Excel,用滚动条调参数、用GAMMA.INV替代查表、用VBA自动寻优,做成了一套能反复用的适线工具。隔壁实验室老王试用完第一句话:“早该这么搞了。”今天把这套东西完整拆开讲一遍,从资料整理、参数初估、频率坐标转换,到理论曲线计算、滚动条调参和网格搜索寻优,手把手带你从一份光秃秃的水文资料,走到一张能直接进报告的设计频率曲线图。
1. 头疼的真正根源:P-III曲线不是画出来的,是“调”出来的
1.1 洪水频率分析为什么绕不开P-III
水文频率分析的核心目标很简单:根据已有的实测流量或降雨资料,推求某个“重现期”对应的设计值。比如百年一遇设计洪水,意思就是每年出现超过这个流量的概率是1%,这个值直接决定堤防高度、桥梁孔径、排水管涵尺寸,属于规划层面的硬指标。
国内水利行业最常用的线型就是皮尔逊三型分布,简称P-III。它本质上是一个带偏态的三参数Gamma分布,能比较好地描述洪水、暴雨这类“下限有限、上限无穷、偏态明显”的变量。几十年来从设计洪水到城市暴雨强度公式,基本都是它在扛大梁。所以哪怕现在有很多新方法,搞水文的人还是绕不开P-III适线。
P-III有三个参数:均值、变差系数Cv、偏态系数Cs。均值好理解,就是平均水平;Cv描述离散程度,相当于相对标准差;Cs描述分布的不对称性,洪水资料几乎都是正偏,Cs大于0。麻烦就麻烦在这个Cs上——它对曲线形状影响极大,却又最不好估,小样本下算出来的值波动非常大,所以才有“矩法初估、目估适线”这一套经验做法。
1.2 让新手当场破防的三个细节
第一个卡点是参数估计。均值、Cv用Excel函数几分钟搞定,Cs需要三阶矩,手算公式不复杂,但牵涉到离差立方求和,数据一多就爱出错。更麻烦的是,Excel自带的SKEW函数算出来的是样本偏度系数,带无偏修正,跟水文上常用的矩法偏度系数有细微差别,直接套用会差一点,虽然不大,但在适线这种“差之毫厘谬以千里”的环节里,最好搞清楚自己用的到底是哪个。
第二个卡点是频率格纸。水文图的横轴不是普通线性刻度,而是“概率格纸”——左边疏右边密,或者反之,目的就是把理论频率曲线从S形拉得接近直线,方便目估适线和外延。过去做图得去买印好的概率格纸,再手工点绘;现在大家都用Excel,但Excel散点图的横轴默认是线性坐标,直接把频率百分比放上去,画出来就是一条奇怪的S形,肉眼根本没法判断是否拟合良好。
第三个卡点是适线。所谓适线,就是调整Cv、Cs(有时也调均值),让理论频率曲线尽量“贴合”经验频率点。传统做法是在透明纸上蒙着曲线反复试画,经验越多调得越快。新手呢?只能一遍遍改参数、重新画图、再改,效率低不说,还容易失去耐心。我当年的崩溃点就是这里:明明Excel就在手边,却硬是把它用成了“高级计算器”,完全没发挥出它动态联动的本事。
1.3 这套Excel工具的整体布局
后来我把整个流程理顺了,做成了一个多Sheet的工作簿,核心就五个区块:
- 资料区:原始水文序列(流量或雨量),以及排序后的经验频率点。
- 参数区:均值、Cv、Cs三个可调参数,放在固定单元格里。
- 计算区:用公式实时计算每个频率下的理论值。
- 图表区:一个散点图,同时展示经验频率点(离散点)和理论频率曲线(连续线)。
- 调参面板:滚动条控件,拖动就能改Cv、Cs,图表立刻刷新。
整套工具的核心逻辑,就是用三个可调参数驱动一整列理论值计算,再由这一列理论值驱动图表更新。参数一变,线条就动,适线从“反复重画”变成了“拖动到手感对了为止”。下面按步骤拆开讲,每一步对应一个Sheet或者一个功能区,你照着做就能复现一版。
2. 资料准备与矩法初估:把原始水文数据变成三个参数
2.1 从原始资料到年最大序列
水文站给的原始资料往往是逐日流量或者逐日雨量,一个月30个值,一年365个值,堆在一起根本没法直接用。做频率分析的第一步,是从长序列里提取“年最大值”——每年挑一个最大日流量(或最大场次降雨量),组成一个n年的极值样本。
如果你手里的资料是整齐的Excel表格,比如一列日期、一列数值,最快的方法是插入数据透视表。把年份拖到行区域,把流量/雨量拖到值区域,值字段设置改成“最大值”,一个透视表就把n年极值全提出来了,比自己拉公式快得多。如果你的资料已经是逐年最大值,那就跳过这步,直接复制一列出来。
拿到极值序列之后,先做两件事:一是按从大到小排序,二是检查有没有明显异常值。排序直接用“数据→排序”,降序排列就行。异常值要特别留意——比如28年资料里突然蹦出一个比其他所有年份大好几倍的数值,这往往是历史特大洪水或者记录错误。历史洪水是宝贵信息,处理方式跟普通实测值不太一样,后面会专门说。
2.2 三个参数的Excel计算公式
假设排序后的流量放在B2:B29(n=28),那么三个矩法参数在Excel里可以这样算:
| 参数 | 公式 | 说明 |
|---|---|---|
| 均值 x̄ | =AVERAGE(B2:B29) | 算术平均值 |
| 标准差 σ | =STDEV.S(B2:B29) | 样本标准差,分母n-1 |
| Cv | =STDEV.S(B2:B29)/AVERAGE(B2:B29) | 变差系数,相对离散程度 |
| Cs | 见下方自定义公式 | 偏态系数,关键参数 |
Cs的矩法公式是:
[ C_s = \frac{n}{(n-1)(n-2)} \cdot \frac{\sum (x_i - \bar{x})^3}{\sigma^3} ]
Excel里没有现成的三阶矩函数,SKEW函数算的是经过无偏修正的版本,跟水文上常用的矩法CS略有差异。为了完全可控,我建议直接在单元格里写自定义公式。比如用辅助列算离差三次方:假设B列是流量,C2填=(B2-AVERAGE($B$2:$B$29))^3,下拉到C29,那么Cs就可以写成:
=COUNT(B2:B29)/((COUNT(B2:B29)-1)*(COUNT(B2:B29)-2)) * SUM(C2:C29) / (STDEV.S(B2:B29)^3)
其中COUNT用来动态统计样本数,以后数据增删不用改公式。这样算出来的就是标准矩法Cs,跟教材对得上。
2.3 初估结果先别急着用,做个合理性检查
矩法初估只是起点,连“半成品”都算不上。拿到均值、Cv、Cs之后,我习惯先做三个快速检查:
第一,Cv是否在合理区间。国内大多数流域的年最大洪峰流量,Cv一般在0.3到1.2之间,干旱半干旱地区可能更大,湿润地区相对较小。如果你算出Cv只有0.05,那大概率资料有问题;如果Cv超过2,得看看是不是混入了特大值。
第二,Cs是否在经验范围内。水文上有个粗略经验:Cs通常取Cv的2到4倍。如果矩法算出来Cs是Cv的10倍,那多半是样本太短或者有异常值,直接用会画出一条形状诡异的曲线。这时候宁可先用Cs=2Cv或3Cv做初值,后面适线再修。
第三,均值是否合理。把均值跟相邻站、多年平均值对比一下,量级不对就查原始资料。我曾经遇到过雨量资料单位从毫米变成米的情况,算出来均值差了1000倍,幸好做了这一步检查。
初估参数合理了,再进下一步。别小看这几分钟的检查,它能帮你省掉后面一晚上的折腾。
3. 频率格纸与坐标转换:Excel散点图里藏着的“非对称坐标”
3.1 概率格纸到底在做什么
先想一个问题:为什么水文频率曲线非得用特殊坐标纸?直接拿普通坐标纸画流量和频率不行吗?其实也可以画,但画出来的理论曲线是一条非常弯曲的S形,两头翘得厉害,中间又特别陡,人眼很难判断离散点是不是“贴”在线上。如果能把横坐标按照正态分布的规律“掰开”——让小概率区域拉宽、中间概率区域压缩——曲线就会被扳得接近一条直线,适线变成“看离散点是否围绕直线分布”,直观得多。
这就是概率格纸的原理。水文上常用的频率格纸,横轴刻度不是均匀的百分比,而是跟正态分布分位数对应的不均匀间距。比如0.01%到0.1%之间的间距,在格纸上是正常坐标的将近2倍宽,而50%附近则被压缩。这样处理之后,理论上如果数据服从正态分布,累积频率曲线就是一条直线。
3.2 用NORM.S.INV把频率转成x坐标
Excel里没有现成的“水文频率格纸”功能,但我们可以用NORM.S.INV(旧版叫NORMSINV)把频率P转换成标准正态分位数,这就是频率格纸的横坐标数值。
转换关系很简单:
[ x = NORM.S.INV(P) ]
比如P=0.5(50%)时,x=0;P=0.01(1%)时,x≈-2.326;P=0.99(99%)时,x≈2.326。下面是一组常用刻度对应的转换值,留作参照:
| 频率P (%) | 0.01 | 0.1 | 1 | 5 | 10 | 20 | 50 | 80 | 90 | 95 | 99 | 99.9 | 99.99 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| NORM.S.INV(P) | -3.719 | -3.090 | -2.326 | -1.645 | -1.282 | -0.842 | 0.000 | 0.842 | 1.282 | 1.645 | 2.326 | 3.090 | 3.719 |
注意,这里P用的是累积概率。在水文频率图上,横轴习惯从左到右频率从小到大(0.01%到99.99%),对应NORM.S.INV值从左到右从-3.72递增到3.72。
实际操作中,我在数据表里加一列辅助列,把每个经验频率P用NORM.S.INV转成x。比如原始频率在A2,辅助单元格写:
=NORM.S.INV(A2)
然后用这一列作为散点图的x轴数据,流量作为y轴数据。这样画出来的图,横轴位置天然就是频率格纸的效果,分布均匀、疏密有致。
3.3 把横轴标签改成“看着像频率格纸”的样子
散点图画出来之后会遇到一个新问题:横轴显示的是-4、-2、0、2、4这些分位数,不是0.01%、1%、50%、99%这种频率值。如果把这样的图贴进报告,那是不合格的,得把轴标签改成水文频率刻度。
有两个思路。一个是做“辅助标签系列”,比较稳妥,步骤如下:
- 先在一列里列出你想要的横轴频率刻度,比如0.01%、0.1%、1%、5%、10%、20%、30%、40%、50%、60%、70%、80%、90%、95%、99%、99.9%、99.99%。
- 旁边一列用NORM.S.INV转成对应的x位置。
- 再在旁边一列给这些刻度配一个固定的y值,这个值要小于图中所有流量数据的最小值,比如0或者比最小值再低一截,目的是让标签系列的点落在图的下边缘。
- 把这个辅助系列加到散点图里,然后右键添加数据标签,标签内容选择“单元格中的值”,范围选频率刻度那一列。
- 最后把原来的x轴隐藏(刻度线、标签都设为无),把辅助系列的标记点设为无填充无线条。一张带标准频率刻度的图就出来了。
另一个思路是用自定义数字格式,但受限于Excel的刻度标签机制,很难做到一一对应,我不推荐。老老实实用辅助系列,一次做好之后存成模板,以后每张图复制过去换数据就行。
3.4 经验频率点怎么算
经验频率其实不复杂,最常用的是魏伯公式(Weibull):
[ P = \frac{m}{n+1} ]
其中m是降序排列后的序号(最大值为1,最小值为n),n是样本数。例如n=28时,最大值的经验频率是1/29≈3.45%,最小值是28/29≈96.55%。
在Excel里,如果你已经按降序排好序,直接用行号减1就是m。比如B列是流量,从B2开始,m就在A2写=ROW()-1,P在C2写=A2/(COUNT($B$2:$B$29)+1),下拉即可。
经验频率点的散点图是离散的,用空心圆标记;理论曲线是连续的,用无标记的平滑线。两类数据放在同一个图表里,目估适线就靠看它们贴不贴。
4. 理论频率曲线的核心:用GAMMA.INV替代查离均系数表
4.1 一张表引发的血案
以前做P-III适线,最烦的一步就是查离均系数Φ值表。那张表按频率P和偏态系数Cs交叉列出Φ值,A4纸打印出来密密麻麻好几页,查的时候先找行、再找列,遇到表里没有的Cs还得插值。更折磨人的是,翻了好几版教材,表的数据还有细微差异,抄错一个数字,整条曲线就歪了。
后来我想明白一件事:P-III分布本质上就是三参数Gamma分布,而Excel的GAMMA.INV函数可以直接求Gamma分布的分位数,为什么还要查表?
4.2 从分布参数到一行Excel公式
P-III分布的概率密度函数是:
[ f(x) = \frac{\beta^\alpha}{\Gamma(\alpha)}(x-a_0)^{\alpha-1}e^{-\beta(x-a_0)} ]
其中α是形状参数,β是尺度参数(注意这里是倒数关系),a0是位置参数。三个参数可以由均值x̄、Cv、Cs换算得到:
[ \alpha = \frac{4}{C_s^2} ]
[ \beta = \frac{2}{\bar{x} \cdot C_v \cdot C_s} ]
[ a_0 = \bar{x}\left(1 - \frac{2C_v}{C_s}\right) ]
而Excel的GAMMA.INV函数,语法是GAMMA.INV(probability, alpha, beta),注意它这里的beta是“尺度参数”的另一种定义,跟上面密度函数里的β互为倒数。也就是说,传给Excel的第三个参数应该是:
[ \text{Excel_beta} = \frac{1}{\beta} = \frac{\bar{x} \cdot C_v \cdot C_s}{2} ]
有了这三个值,任意频率P对应的理论流量x_P就是:
[ x_P = a_0 + GAMMA.INV(1-P, \alpha, \text{Excel_beta}) ]
如果均值x̄在G1单元格,Cv在G2,Cs在G3,频率P在A列,那Excel公式可以写成:
=$G$1*(1-2*$G$2/$G$3) + GAMMA.INV(1-A2, 4/($G$3^2), $G$3*$G$1*$G$2/2)
这里有个值得注意的地方:GAMMA.INV的probability参数是累积概率,也就是“小于等于该值”的概率。但水文频率P表示的是“大于等于该值”的概率(超过概率),所以要用1-P。P越小(越稀遇),1-P越接近1,GAMMA.INV返回的数值越大,对应大流量,逻辑完全正确。
4.3 用实际算例验证公式对不对
拿一个典型例子验证。假设某站年最大洪峰流量均值 x̄=100 m³/s,Cv=0.5,Cs=1.0,推求百年一遇(P=1%)设计流量。
先算参数:
- α = 4/(1²) = 4
- Excel_beta = 100×0.5×1/2 = 25
- a0 = 100×(1-2×0.5/1) = 100×(1-1) = 0
代入公式:
- x_1% = 0 + GAMMA.INV(0.99, 4, 25)
用Excel算一下,GAMMA.INV(0.99, 4, 25)大约等于225.8,也就是说百年一遇设计流量约为226 m³/s,大概是均值的2.26倍。这个倍数关系跟国内很多流域的设计成果是接近的,说明公式没有大的方向性问题。如果把Cv提高到1.0,那设计值会大幅上升,这也符合同等频率下变差系数越大设计值越高的常识。
验证通过之后,就可以放心用了。以前查表插值要花十分钟,现在一列公式拉到尾,不到一秒钟。
4.4 生成理论频率曲线的数据表
理论曲线不能只算几个点,得生成一串平滑的点,才能连成线。我的做法是在工作表里单独放一列频率序列,从P=0.01%开始,按对数或按分位数均匀取点,一直到P=99.99%,总共取50到80个点。具体来说,可以先列出你想要的频率P值,然后用上面的公式生成对应的理论流量。
一个常用的频率点序列大致是:0.01%、0.02%、0.05%、0.1%、0.2%、0.5%、1%、2%、3%、5%、10%、20%、30%、40%、50%、60%、70%、80%、90%、95%、97%、98%、99%、99.5%、99.8%、99.9%、99.95%、99.99%。别看只有28个点,在Excel平滑线散点图里已经能画出非常顺滑的曲线。
画的时候,x坐标同样用NORM.S.INV(P)转换,y坐标用理论流量。这样理论曲线和经验点就在同一套“频率格纸”坐标系里了,贴合程度一目了然。
这里有一个容易踩的坑:P取值不能是0或1。因为GAMMA.INV(0, α, β)=0,GAMMA.INV(1, α, β)=正无穷,Excel直接报错。所以频率序列两端要留出余量,用0.01%和99.99%这种极限值,但永远别取0和1。
5. 目估适线到自动寻优:滚动条、规划求解与VBA怎么选
5.1 为什么矩法初估之后还得“调”
理论上,矩法估计出来的参数已经是对资料的最佳拟合了,为什么还要人工调整?因为矩法是从整体上保证统计特征一致的,并不保证在高频段或者稀遇频率段“贴合经验点”。在小样本情况下,个别特大值或特小值会把矩法参数拉偏,导致曲线在P=1%附近的估计值明显偏离经验点的趋势。这时候就需要适线——在合理范围内微调Cv和Cs,让曲线更贴合经验点的整体走势,尤其是曲线尾部(小频率、大流量)的走势。
传统适线靠手,把透明纸蒙在图上反复试画,半天调一版。现在有了动态图表,完全可以做到“所见即所得”。
5.2 滚动条调参面板
Excel开发工具里的表单控件滚动条,是实现动态调参最简单的手段。前提是先把“开发工具”选项卡调出来:文件→选项→自定义功能区→勾选“开发工具”。
然后按下面步骤做:
- 在参数区预留两个单元格,一个放Cv,一个放Cs。我习惯G2放Cv,G3放Cs。
- 点击“开发工具→插入→表单控件→滚动条(窗体控件)”,在图表旁边拖出一个横向滚动条。
- 右键滚动条→设置控件格式,关键设置:
- 最小值:20(代表0.20)
- 最大值:150(代表1.50)
- 步长:1(每次动0.01)
- 页步长:5(点空白处跳0.05)
- 单元格链接:$G$1(偏好放一个辅助单元格)
- 在G2单元格写
=G1/100,把控件值转换成Cv小数。 - 再放一个滚动条控制Cs。因为Cs一般取Cv的倍数,我习惯用一个倍率系数k来控制,比如k从1.0到6.0,步长0.1,然后G3写
=$G$2*k。
这样拖动滚动条,Cv和Cs实时变化,而理论曲线那一整列因为公式引用了G2、G3,也会跟着刷新,图表立即重绘。调参变成了像调音响旋钮一样的操作,手感非常直观。
这个阶段建议配合一个“残差平方和”指标,不然纯靠眼睛容易走偏。在数据表里加一列,对每个经验频率点P_i,用GAMMA.INV公式算出对应的理论流量,然后跟实测流量求差平方。SS=SUM(差值平方)。理论上SS越小拟合越好,但也不能光追SS小,还要看曲线在关键频率段的走势是否合理。
要注意的是,经验点的频率是离散的,而理论曲线是连续的,计算残差时直接用“经验点频率”去算理论值即可,一一对应,不需要插值。28个点的SS,在Excel里一个SUMSQ就搞定了。
5.3 规划求解:能不能全自动?
滚动条调参毕竟还是手动,如果想自动找最优参数,可以用Excel的规划求解(Solver)。先用加载项把规划求解启用来:文件→选项→加载项→转到→勾选“规划求解加载项”。
设置思路:
- 目标单元格:SS所在的单元格,设为“最小值”。
- 可变单元格:Cv和Cs两个参数。
- 添加约束:
- Cv ≥ 0.1
- Cs ≥ 0.5(防止偏度系数变负)
- Cs/Cv ≥ 1
- Cs/Cv ≤ 6
- 求解方法选择“非线性GRG”。
规划求解在简单情况下效果不错,几分钟之内能算出一个让SS最小的Cv、Cs组合。但它有两个问题:一是非线性优化不一定收敛到全局最优,容易陷在局部极小值;二是如果迭代过程中Cs跑到0附近,GAMMA.INV会报错,整个计算表都亮红灯。所以我的经验是:先用VBA或手动滚动条粗调到一个合理区间,再交给规划求解精修。
5.4 VBA网格搜索:最简单可靠的“穷举法”
规划求解虽然自动化,但有时候不如“穷举”来得放心。所谓网格搜索,就是让Cv和Cs在合理范围内取一系列离散值,逐一计算SS,找出最小值。计算量不大,Cv从0.20到1.50步长0.02,Cs倍率从1.0到6.0步长0.1,总共大约66×51≈3366种组合,每个组合算28个点的GAMMA.INV,合计不到10万次函数调用,VBA跑几秒钟就完事。
下面这段VBA代码可以直接用,假设数据表里A2:A29放的是经验频率P,B2:B29是实测流量,G1放均值,G2放Cv,G3放Cs:
Sub P3GridSearch() Dim ws As Worksheet Dim meanVal As Double, cv As Double, cs As Double Dim ss As Double, minSS As Double Dim bestCv As Double, bestCs As Double Dim p As Double, qObs As Double, qTheory As Double Dim i As Long, j As Long, k As Long Dim n As Long Set ws = ThisWorkbook.Sheets("计算表") n = ws.Cells(Rows.Count, 1).End(xlUp).Row - 1 meanVal = ws.Range("G1").Value minSS = 1E+30 For i = 20 To 150 Step 2 cv = i / 100 For j = 10 To 60 Step 1 cs = (j / 10) * cv ss = 0 For k = 2 To n + 1 p = ws.Cells(k, 1).Value qObs = ws.Cells(k, 2).Value qTheory = meanVal * (1 - 2 * cv / cs) + _ Application.GammaInv(1 - p, 4 / (cs ^ 2), cs * meanVal * cv / 2) ss = ss + (qObs - qTheory) ^ 2 Next k If ss < minSS Then minSS = ss bestCv = cv bestCs = cs End If Next j Next i ws.Range("G2").Value = bestCv ws.Range("G3").Value = bestCs ws.Range("G4").Value = minSS MsgBox "搜索完成,最优Cv=" & bestCv & ",Cs=" & bestCs & ",SS=" & minSS End Sub这段代码思路很清楚:两层循环遍历Cv和Cs,内层对每个经验点算理论流量,累加残差平方。搜索完成后把最优值写回参数区,图表自动更新。要注意Application.GammaInv在Excel 2010及以上版本都可用,如果用的是新版Excel,也可以换成WorksheetFunction.Gamma_Inv。
跑完网格搜索之后,我通常还会用滚动条在最优值附近再微调一下,因为SS最小不一定等于目视最优——有时候稍微牺牲一点点SS,能让曲线尾部更符合物理规律(比如不出现负流量)。这种“机器算一遍、人眼过一遍”的流程,我觉得是最稳妥的。
6. 配套的上层应用:从计算表到规范成果图,以及那些容易翻车的细节
6.1 一张能直接进报告的成果图
适线完成之后,图表还需要“打扮”一下才能见人。我的习惯是这样:
- 经验点用空心圆,大小大概6磅,颜色用黑色,不填充——打印出来也清楚。
- 理论曲线用平滑线,黑色实线,线宽1.5磅。
- 横轴用辅助标签系列伪装成频率刻度,从0.01%到99.99%都要有。
- 纵轴用普通数值轴,标题写“流量 (m³/s)”或者“降雨量 (mm)”。
- 图表标题写“某站年最大洪峰流量频率曲线”。
- 图例可以省略,或者放在右下角。
导出到Word或报告里时,别直接复制图表对象,容易带上一堆格式问题。建议选中图表→复制→在Word里“粘贴为图片”,这样不管对方用什么版本打开都不会乱版。如果报告要求矢量图,那就另存为PDF,再从PDF里裁剪。
6.2 Cv/Cs取值边界与频率分析的地雷
Cv、Cs不是随便取的,水文上有一些硬性约束和常识性边界。
Cs理论上可以小于0,但洪水、暴雨资料几乎都是正偏,Cs小于0的曲线形状是左偏的,物理上不明显,一般直接排除。所以约束条件里我会写Cs ≥ 0.5,防止优化器跑偏。
Cs取太大也要警惕。当Cs很大时,α=4/Cs²会变得很小,Gamma分布的形状会极其偏态,理论最小值a0可能会变成一个明显偏高的正数,导致曲线在低频段出现一个“地板”,所有理论流量都大于等于a0。如果a0超过了0,而实际资料里有接近0的年份,这条曲线就有问题。遇到这种情况,要回头检查资料里是不是混入了异常值,或者样本太短导致Cs失真。
另一个经典问题是历史特大洪水的处理。如果28年实测资料里有一场200年一遇的历史洪水,直接把它当作普通样本参与矩法计算,Cv会被严重拉大,曲线整体偏高。正确做法是采用“不连序”样本计算方法,给历史洪水一个加权经验频率。这个内容展开能写一整篇,这里先提个醒:做设计洪水时,资料里一旦有历史洪水,千万别简单按连序样本处理。
还有一个小坑:经验频率点的绘制要对应降序排位。很多新手把原始序列按时间顺序直接算P=m/(n+1),结果经验点乱成一团,没办法跟理论曲线比较。一定要先降序排列,再编号。
6.3 工作簿的健壮性:防止公式被误改、数据自动扩展
这套工具做出来之后,除了自己用,大概率还会在小组里流传。老王当时就用我的版本,一不小心把GAMMA.INV公式给覆盖了,整张图瞬间崩掉。从那以后我养成了三个习惯:
一是把“计算区”的工作表保护起来。选中允许修改的单元格(参数区、资料区),右键→设置单元格格式→保护→取消锁定;其他公式单元格保持锁定。然后“审阅→保护工作表”,设置一个密码。这样别人只能动参数和数据,碰不了公式。
二是用Excel表格(Table)而不是普通区域。把资料区的数据选中后按Ctrl+T转成表格,然后图表数据范围引用表格列名,这样以后在表格末尾追加数据,图表的范围会自动扩展,不用手动改每一个系列的引用。多少年的资料都能无缝更新。
三是重要成果另存一份“快照”。适线完成、参数定稿之后,把参数区、SS值、计算数据复制到另一个Sheet,粘贴为数值。这样即使以后某个公式意外失效,成果数据还在,不会影响出报告。
6.4 这套方法可以怎么扩展
P-III这套思路打通之后,其他线型其实都是同一个套路。比如耿贝尔极值分布,理论频率值可以用=-LN(-LN(1-P))这种简单公式算,连GAMMA.INV都省了;对数P-III分布,无非是把流量取对数之后再做P-III适线,最后再反变换回来。只要理解了“频率→坐标转换→理论值→拟合”这个循环,换线型只是换一行核心公式的事。
我后来还把这套模板用在了城市暴雨强度公式的推求上。暴雨资料按年度最大值采样,短历时降雨用P-III适线推设计雨强,再按回归拟合强度公式里的A、C、n三个参数,整个过程都是在这个Excel框架里完成的。从洪水频率分析到暴雨强度公式,底层逻辑完全一致。
最后说一个小技巧:在调参面板旁边,我习惯放一个“设计成果表”,直接列出几个规范要求的频率对应的设计值,比如P=0.1%、0.2%、0.5%、1%、2%、5%、10%,引用理论曲线计算公式实时刷新。每次调完参数,不光是图在动,这张表里的设计值也在动,报告里需要的数据随手就抄。你可以在自己复刻这套工具时也加上这一块,等到出报告那天,就知道有多省事了。