Page 112 - 2022年第53卷第5期

Gibbs 抽样和 Metropolis 算法 ［30］ 结合起来采用混合 Gibbs 抽样的方法来构造马尔科夫链。其中 Gibbs 抽
样方法 ［31］ 研究已较为成熟，本文不再赘述，Metropolis算法 ［32-33］ 的实现方法如下：
（1）初始化：t=0，选择一个对称提议分布 q（x），满足 q（x|y）=q（y|x），给定一个起始样本点 x ，迭
t
代终止值为T；
（2）令t=t+1，从q（x|x ）中生成候选样本x′；
t-1
（3）计算接收概率α，α=min{1，x′/ x }；
t-1
（4）从均匀分布中抽取随机 α ，若 α ≤α，则接受候选样本，x = x′； 否则，拒绝候选样本，并令
t
t
t
x = x ；
t
t-1
（5）重复（2）—（4）步直到迭代终止。
采用混合Gibbs抽样同时产生多条Markov链，若这几条链稳定下来则认为混合Gibbs抽样收敛。
2.4 可靠度评价 在实际工程中，往往需要给出在一定可靠度条件下参数的估计值，设参数估计的
可靠度为R，分布函数与可靠度的关系为：
ê ( ) ξ ù
é
x
R( ) x = 1 - F( ) x =expê- η ú ú （9）
ë û
对上式取对数可得：
( ) 1 ξ
1
x R = η ln R （10）
式中 x 为两参数威布尔分布在可靠度 R 条件下的估计值。将基于贝叶斯理论更新得到的两参数后验
R
威布尔分布估计结果代入式（10），即可计算出在可靠度 R条件下的砂砾料级配特征参数的估计值 x 。
R
此外根据贝叶斯方法更新后的威布尔分布作为总体分布，产生的随机数可作为压实质量评估模型中
坝料级配特征参数的输入指标。

3 算例验证

假设砂砾石料 P5 含量服从形状参数 ξ=3、尺度参数 η=2 的威布尔分布，从中产生 25 个随机数作
为先验样本数据集，5 个随机数作为现场试验数据集，采用本文方法编制的程序对上述小样本数据
进行威布尔分布的贝叶斯估计求解。在对后验分布进行模拟仿真求解时，本文构建了迭代次数为
10 000、20 000、30 000 的三条 Markov 链，仿真过程中监控形状参数 ξ 和尺度参数 η 输出的模拟估计
值，经过统计分析用于判断迭代过程是否收敛。三条马尔科夫链输出结果的统计分析结果如表 1 所
示。从表 1 和图 2 可以看出，随着迭代次数的增加，MC 误差在逐步的减小，且三条链在 2.5%、50%
和 97.5% 这三个分位数对应的估值基本相同，这表明在迭代过程中 Markov 链逐渐收敛并趋于稳定。
当 Markov 链达到平稳分布时，将迭代计算结果的均值作为后验形状参数和尺度参数的估计值，后验
形状参数 ξ=3.1950 和尺度参数 η=2.0380 与真实值的相对误差分别为 6.5% 和 1.9%，在误差允许范围
内，由此验证了本文小样本条件下贝叶斯估计方法和编制程序的准确性。
上述通过一次蒙特卡罗模拟验证了本文方法的有效性，为了消除随机性的影响，先后进行了100
次蒙特卡罗模拟，并分别采用了本文方法与传统 Bootstrap 方法进行贝叶斯估计对比分析。后验分布

表1 后验分布参数的模拟估计值
待估参数 均值 标准差 MC误差 2.5% 50% 97.5% 样本数
3.1927 0.5200 0.0322 2.1963 3.2014 4.1699 10000
形状参数ξ 3.1720 0.5068 0.0222 2.0956 3.1815 4.1183 20000
3.1950 0.5016 0.0180 2.1518 3.2113 4.1291 30000
2.0385 0.1510 0.0094 1.7489 2.0312 2.3411 10000
尺度参数η 2.0364 0.1468 0.0064 1.7581 2.0317 2.3415 20000
2.0380 0.1500 0.0054 1.7496 2.0347 2.3408 30000

—612 —