洪水、强风与地震的概率尺度

一座桥需要跨过河流,也会遇到强风与地震。水文资料说“百年一遇”,风速表按重现期列值,地震资料则可能给出“某段时间内的超越概率”。这些说法都在连接作用强度与发生机会,但所统计的对象和时间口径并不相同。先分清“什么量超过什么阈值”,才能把统计结果接到桥梁计算中。

完成本专题后,学生能够:

  1. 概率与时间:说明年最大值、年超越概率、重现期和考察期超越概率之间的联系,解释“百年一遇”为何不表示每隔一百年发生一次。
  2. 极值与分位:用给定的Gumbel模型求重现期分位值,辨认GEV分布的参数、支持域和尾部差异,并说明资料外推的不确定性。
  3. 作用与性能:区分洪水、风与地震的危险性以及结构失效概率,把统计输入与水力、气动、动力和抗力模型连接起来。

先确定统计对象

同一场洪水可以记录最大流量、最高水位和洪量;同一场风暴可以记录短时阵风或较长时距的平均风;同一次地震可以用震级描述震源,也可以用峰值地面加速度描述桥位地面运动。名称中都出现“最大”,并不意味着这些量可以共用一条概率曲线。

工程问题 可选统计变量 必须随资料说明的条件
桥孔过水与冲刷 年最大洪峰流量,单位为m³/s 测站、流域、水文年、调蓄条件及流量换算方法
梁底高程与淹没 指定位置的年最高水位,单位为m 高程基准、河道形态、潮汐与上下游回水条件
桥梁抗风 年最大风速,单位为m/s 测高、地表、平均时距、风向和仪器变更
桥梁抗震 桥位地震动强度,如PGA或指定周期谱加速度 场地、分量、谱周期、阻尼及危险性模型

年最大值序列每年只保留一个最大值。若有30个资料完整的年份,便得到30个年最大值,而不是把每年的数万次采样都当作独立的“年极值”。水文年可以采用适合流域的年度划分,但前后必须一致。缺测恰好发生在洪峰或台风期间时,该年的记录最大值可能低于真实年最大值,不能因表格中仍有一个数就视为完整。

另一种取样方式是超阈值序列:保留超过阈值的独立事件峰值,一年可能没有,也可能有数个。同一场洪水的连续高水位、同一场台风的相邻阵风并非许多独立事件,需要按过程分组。阈值太低会混入普通波动,太高又只剩少数样本。事件峰值分布还需配合事件发生率,才能转化为年度危险性;不能把事件样本的分位数直接叫作“百年一遇”。NIST的飓风模型将事件到达过程和事件最大风速分开处理,体现了这种区别。NIST:飓风风速概率模型

基础回顾:年度有无超越

本节可快速回顾;首次推导及次数分区见第2章。本专题主线接着区分年概率与事件年率,再讨论极值分位、有限资料外推及非平稳条件下的结构性能。

令随机变量 XX 表示某桥位一年内的最大作用强度,阈值为 xx,年最大值的分布函数为

FX(x)=P(X≤x).F_X(x)=P(X\le x).

对连续分布,超过阈值的年概率为 p(x)=1−FX(x)p(x)=1-F_X(x)。按年最大值口径定义重现期,若 TT 的数值以年计,则

p(xT)=1T,FX(xT)=1−1T.p(x_T)=\frac1T,\qquad F_X(x_T)=1-\frac1T.

xTx_T 是与 TT 对应的分位值,不是“未来第 TT 年的预测值”。十年一遇对应每年10%的超越概率,百年一遇对应每年1%的超越概率。重现期变大,选取的是分布更靠右的量级;年超越概率则变小。USGS也用年超越概率解释洪水重现期,以免将其误读为固定周期。USGS:百年一遇洪水

在各年独立、分布保持不变的模型下,今年已发生一次超越,并不会使明年“少一次机会”;明年的概率仍为 pp。相邻两年均超越的概率为 p2p^2,虽然小,却不为零。如果长时间观察,以“哪些年份至少超越一次”为事件,平均等待年数为 1/p1/p;某一次实际等待可能只有一年,也可能远长于这个平均值。

例题1 五十年中遇到一次的机会

假定各年独立、统计条件不变。百年一遇阈值对应年超越概率 p=0.01p=0.01;连续 nn 年均不超越的概率是 (1−p)n(1-p)^n,因此

Pn=1−(1−p)n,P50=1−0.9950=39.50%.P_n=1-(1-p)^n,\qquad P_{50}=1-0.99^{50}=39.50\%.

50年只是本例考察窗口;39.50%表示阈值至少被超越一次,不是桥梁倒塌概率。考察时间增加,这一概率应增加但不超过100%;np=0.5np=0.5 则是超越年份数的期望,二者不能混用。完整的10年、50年窗口与零次、恰好一次分区见第2章的次数对照。

年超越率与年超越概率

地震危险性分析常先给出年超越率 λ(x)\lambda(x):每年地震动超过强度 xx 的事件次数期望,单位为年⁻¹。它是计数速率;年超越概率 p(x)p(x) 则是一年内至少出现一次超越的概率,无量纲。一个计数过程可以一年发生多次,因而二者不能直接等同。

若超越事件服从恒定速率的泊松过程,在 nn 年内的事件数 NnN_n 满足 E[Nn]=λn\mathrm E[N_n]=\lambda n,零次事件的概率为 e−λne^{-\lambda n},所以

P(Nn≥1)=1−e−λn,p=1−e−λ,λ=−ln⁡(1−p).P(N_n\ge1)=1-e^{-\lambda n},\qquad p=1-e^{-\lambda},\qquad \lambda=-\ln(1-p).

式中 nn 以年计,λ\lambda 取每年的数值,指数中的 λn\lambda n 无量纲。若将平均事件间隔记为 Tλ=1/λT_\lambda=1/\lambda,它与年最大值口径的 Tp=1/pT_p=1/p 严格不同;当概率很小时,才有 p≈λp\approx\lambda 和 Tp≈TλT_p\approx T_\lambda。例如 λ=0.01年−1\lambda=0.01\ \text{年}^{-1} 对应 p=0.9950%p=0.9950\%,并非恰好1%。

在同一个泊松模型中,以 p=1−e−λp=1-e^{-\lambda} 代入 1−(1−p)n1-(1-p)^n,会恰好得到 1−e−λn1-e^{-\lambda n}。这不是两套相互矛盾的公式,而是同一过程按“年度有无事件”和“事件计数”两种方式表达。USGS的地震危险性说明从震源发生率、给定地震下的场地超越概率出发,汇总成危险性曲线,再进行这种时间概率转换。USGS:地震危险性基础

例题2 由考察期概率求地震动年超越率

某教学危险性曲线中的一个地震动强度值,在恒定泊松模型下具有“50年内10%的超越概率”。求对应年超越率、年超越概率和两种重现期,并说明 1/p1/p 与 1/λ1/\lambda 的统计对象为何不同。

解。 从 0.10=1−e−50λ0.10=1-e^{-50\lambda} 得

λ=−ln⁡0.9050=0.00210721年−1.\lambda=-\frac{\ln0.90}{50}=0.00210721\ \text{年}^{-1}.

因此

p=1−e−0.00210721=0.00210499≈0.2105%,p=1-e^{-0.00210721}=0.00210499\approx0.2105\%,

Tλ=1λ≈474.56年,Tp=1p≈475.06年.T_\lambda=\frac1\lambda\approx474.56\ \text{年},\qquad T_p=\frac1p\approx475.06\ \text{年}.

二者接近但不相等。计算中把 10%/50=0.2%10\%/50=0.2\% 当作年概率只是近似。把求得的 λ\lambda 回代,1−e−50λ=0.101-e^{-50\lambda}=0.10,可核对转换方向。本题没有给出危险性曲线横坐标,因而不能仅凭“50年10%”求出PGA,更不能把它换算成震级。这一概率水平仅为教学输入。

泊松模型要求在采用的尺度上合理描述事件到达。余震成簇、断层更新过程或随时间变化的发生率可能需要其他模型;不能在未说明事件筛选和时间模型时,给所有地震目录直接套用恒定泊松公式。

从概率求强度:Gumbel与GEV

上面的计算从概率求“遇到几次的机会”,并未给出洪峰或风速有多大。要得到 xTx_T,还需要描述年最大值的分布。普通时刻的风速分布与年最大风速分布并不相同;一年取最大值的过程会改变分布形状。

最大值型Gumbel分布

以位置参数 μ\mu、尺度参数 σ>0\sigma>0 表示最大值型Gumbel分布:

FX(x)=exp⁡[−exp(−x−μσ)].F_X(x)=\exp\left[-\exp\left(-\frac{x-\mu}{\sigma}\right)\right].

μ\mu 和 σ\sigma 与 xx 单位相同,σ\sigma 不等于样本标准差。令 FX(xT)=1−1/TF_X(x_T)=1-1/T,对两层指数依次取对数,得到

xT=μ−σln⁡[−ln(1−1T)],T>1.x_T=\mu-\sigma\ln\left[-\ln\left(1-\frac1T\right)\right],\qquad T>1.

注意这是最大值型表达;最小值型Gumbel的符号和尾部方向不同。NIST:Ⅰ型极值分布

例题3 同一风速模型下的三个分位值

沿用抗风部分的教学参数,取 μ=25.0m/s\mu=25.0\ \mathrm{m/s}、σ=2.5m/s\sigma=2.5\ \mathrm{m/s},它们分别对应抗风公式中的 bb、aa。求10年、50年、100年分位风速,并检查重现期翻倍是否使风速翻倍。

解。 代入分位式,得

TT(年) FX(UT)F_X(U_T) UTU_T(m/s)
10 0.90 30.63
50 0.98 34.75
100 0.99 36.50

取一位小数,后两项为34.8 m/s、36.5 m/s,与抗风算例一致。50年至100年,重现期翻倍,风速仅增加约5.02%;在空气密度和气动力系数不变时,速度压增加约10.30%。重现期的比例不能直接当作风速或风力的比例。

数值核对可反向进行:将36.5004 m/s代回分布函数,应得约0.99,因此年超越概率约为0.01。本题参数是假定的,所得风速没有自动包含桥面高度、地形、平均时距与风向转换。

GEV分布怎样改变尾部

广义极值分布(GEV)在位置、尺度之外加入形状参数 ξ\xi。本专题采用如下符号约定:

FX(x)=exp⁡{−[1+ξx−μσ]−1/ξ},σ>0,ξ≠0.F_X(x)=\exp\left\{-\left[1+\xi\frac{x-\mu}{\sigma}\right]^{-1/\xi}\right\},\qquad \sigma>0,\quad\xi\ne0.

括号必须满足支持域条件

1+ξx−μσ>0.1+\xi\frac{x-\mu}{\sigma}>0.

ξ\xi 无量纲。ξ>0\xi>0 时上尾较重,支持域有下端点;ξ<0\xi<0 时存在上端点 μ−σ/ξ\mu-\sigma/\xi;ξ→0\xi\to0 时连续过渡到Gumbel分布。NOAA的极端水位技术报告用三种形状参数比较了这三类尾部。NOAA:极端水位统计方法,技术报告67

在上述支持域和符号约定下,反解可得

xT=μ+σξ{[−ln(1−1T)]−ξ−1}.x_T=\mu+\frac{\sigma}{\xi}\left\{\left[-\ln\left(1-\frac1T\right)\right]^{-\xi}-1\right\}.

利用 a−ξ=exp⁡(−ξln⁡a)≈1−ξln⁡aa^{-\xi}=\exp(-\xi\ln a)\approx1-\xi\ln a,可验证 ξ→0\xi\to0 时回到Gumbel分位式。有的软件以相反符号定义形状参数,抄入参数前应对照分布公式。拟合出的上端点也只是该统计模型的边界,不能据此断言自然界存在已经证实的绝对上限。

图中固定位置参数24 m/s、尺度参数3 m/s,比较形状参数−0.2、0、+0.2,并显示当前可调形状参数对应的分位。三曲线共用坐标,支持域端点与50年、100年分位同步标出。ξ\xi 改变尾部;资料年限影响外推不确定性。该页只演示已给模型,不把参数差异冒充置信区间。

GEV提供一组可比较的尾部模型,并不规定所有水文、风与地震资料都必须采用它。USGS洪峰频率方法中的对数皮尔逊Ⅲ型分布就是另一条技术路线。模型选择要结合资料性质和项目采用的方法,不能因为几个点拟合得更贴近就断定远端尾部更可信。USGS:洪峰频率分析指南17C

短资料怎样支持较长重现期

30年记录可以用于估计100年分位值,但这包含模型外推,不等于已观测到100年的变化。先把各年极值按时间排列,检查仪器、测站、地表和河道变更,再排序观察尾部。课堂上可用 Fi≈i/(N+1)F_i\approx i/(N+1) 给升序样本一个绘图位置,其中 NN 是样本数;这些位置便于作图,并非每个观测值的已知真实概率。

随后才选择分布并估计参数。最大似然、矩法和L矩法利用资料的方式不同,有限样本下可能给出不同尾部。可同时比较经验点、概率图、分位图和目标重现期附近的结果;全体样本上的一个拟合优度指标,不能单独保证远端分位值准确。NIST的极端风速工具列出了多种分布与估计方法,适合用于这种对照。NIST:极端风速统计概览

分位估计应同时给出资料年限、点估计和置信区间。增加记录通常有助于缩小抽样误差,但新增年份可能出现更大极值,使分位点估计上升。历史洪痕、已知某年超过某水位、测站量程之外的洪峰,也包含信息;它们应表示为区间或删失记录,不能随意补成一个精确值。USGS的PeakFQ输出将观测、历史资料、阈值与置信界限放在同一频率图上。USGS:PeakFQ

“100年分位值的95%置信区间”描述估计方法的不确定性,并不表示“未来100年有95%概率不超越该区间上界”,也不等于结构可靠度95%。概率模型描述自然变异,置信区间描述利用有限资料估计模型时的误差;两层不确定性应分别解释。

当统计条件随时间变化

平稳假定意味着所研究变量的分布在考察期间保持不变。上游水库调度、城市化、河床变化、测风环境改变以及气候趋势,都可能使长期序列不再来自同一分布。NOAA关于降水频率更新的说明已将气候趋势作为需要处理的问题。NOAA:降水频率与气候趋势

若每年超越同一固定阈值的概率变为 pjp_j,且各年仍可视为独立,则

Pn=1−∏j=1n(1−pj).P_n=1-\prod_{j=1}^{n}(1-p_j).

若采用发生率随时间变化的非齐次泊松模型,则为

P(考察期内至少一次)=1−exp⁡[−∫0tdλ(t)dt].P(\text{考察期内至少一次})=1-\exp\left[-\int_0^{t_d}\lambda(t)\,\mathrm dt\right].

这时要说明所说的“百年一遇”属于哪个基准时期,或哪个年份的分布。检出趋势后,先辨认仪器与环境变化,再判断怎样建模;既不能机械延用旧常数,也不能仅凭短序列的上升就无限外推。若年际存在明显依赖,上述乘积公式同样需要调整。

从作用危险性走到桥梁性能

水文计算通常先得到相应概率的洪峰流量,再通过河道断面、糙率、上下游水位与桥孔条件计算水位和流速。若这些关系固定且单调,流量分位值可通过模型映射到水位;遇到潮汐顶托、调蓄或河道变化时,单一对应可能不再成立。同名重现期的降雨、洪峰与水位也不保证在同一次过程中同时出现。

风速分位值还要转化为桥面处的输入,由气动力系数、风场频谱与结构模态求响应。地震动危险性曲线给出桥位强度与超越率的关系,结构反应则依赖周期、阻尼、地震动频谱与持续时间、支承及非线性。概率水平相同,只表示选取的超越机会口径相近,不表示三个灾种具有相同的结构后果。

作用超过某个设计阈值,不等于结构失效;作用低于该阈值,也不保证任何受损状态下都不失效。 结构失效需对一个明确极限状态比较抗力 RR 与作用效应 SS,例如 g=R−S<0g=R-S<0。这里的 SS 已经是流量、风速或地面运动经过相应模型产生的需求,不能把m³/s、m/s或地震动加速度直接与构件抗力相减。

若把一年内的灾害情景概括为强度变量 XX,同一时段、同一失效定义下可按全概率关系理解

Pf=∫P(失效∣X=x)dFX(x).P_f=\int P(\text{失效}\mid X=x)\,\mathrm dF_X(x).

其中条件失效概率还受到材料、几何、边界、劣化及模型误差影响。只有在强度变量足以表示所研究情景时,这个一维表达才适用;疲劳、冲刷累积或多次地震损伤可能需要保留过程与历史。水文、抗风和抗震共享的是概率语言,连接结构性能时仍须沿各自的物理过程计算。

阅读与计算任务

选取一份洪水、风或地震资料,写出变量、单位、位置、时间口径与超越阈值,再判断它给的是年概率、年率还是考察期概率。复算三个例题中的一个,并说明哪一项模型假定若改变,答案就需要重算。最后把所得统计量沿“作用强度—结构响应—极限状态”写成三步,指出还缺哪些物理参数。这样,重现期才能成为可解释的设计输入,而不只是一个看起来很大的年数。

迁移核对。 把例题2的情景改为“50年内20%”,仍采用恒定泊松过程。分别求年超越率、年超越概率、平均事件间隔与按年度有无超越定义的重现期;说明两种“平均时间”为什么稍有差别。

核对结果。 年超越率约0.004463 年⁻¹,年超越概率约0.4453%;平均事件间隔约224.07年,年度有无超越的平均等待年数约224.57年。前者允许一年内多次事件,后者只计该年是否发生过超越;这里的差值不是四舍五入造成的。