滑带或滑面是坡体内部岩土体受到揉皱及剪切逐渐形成的软弱带或剪切面，影响滑坡的演化过程与变形稳定性. Jiang等^[1-2]研究不同部位滑带土的颗粒级配、矿物组成和化学成分的空间差异特征. Lu等^[3]探究滑带土颗粒级配的分形维数对黄土坡滑坡微观结构和抗剪强度参数空间变异性的影响. 以上研究表明，在滑带的形成过程中，不同位置岩土体的物质成分、颗粒级配与微观结构是空间变化的，抗剪强度参数存在空间变异性. 目前，滑带力学参数的获取可以通过工程地质类比法、室内外力学试验与强度参数反演等手段. 由于测量误差、岩土体固有变异与转换模型误差等因素的存在^[4-5]以及取样数量与取样位置的限制，通过试验手段获取的抗剪强度参数存在误差，不能代表整个研究区域. 对于强度参数反演方法，通常结合力学试验数据，利用滑坡的位移监测数据^[6-7]和滑坡实际稳定性状况^[8-9]等，对输入强度参数进行反演. 为了建立精确的滑带或滑面抗剪强度参数模型，既需要考虑滑带土力学强度参数的空间变异性与不确定性，也需要利用滑坡的后验信息.

对于土性参数的空间变异性与不确定性，利用随机场理论、地统计学理论与可靠度理论可以研究计算模型的不确定性^[10-16]. 这些理论方法依赖于参数的基本统计特征，这些参数信息的获取较困难. 为了解决参数先验信息获取困难或不准确的问题，贝叶斯方法被应用到岩土工程的不确定性建模中，Wang等^[12]介绍了岩土体参数表征过程中多种不确定性影响的贝叶斯反分析框架. Zhang等^[13]采用马尔可夫链蒙特卡洛（MCMC）的贝叶斯方法，利用边坡的失稳观测信息对边坡的多个参数进行随机反演. Wang等^[14]利用MCMC和最大似然的贝叶斯方法，对公路边坡的摩擦角和锚固力进行随机反演. Gong等^[15]结合随机场理论与MCMC贝叶斯方法，利用测点信息对场地的岩土参数特征进行表征. 蒋水华等^[16]以简单边坡模型为例，基于可靠度的贝叶斯方法利用试验数据对随机场参数进行反演. 以上研究表明，贝叶斯方法可以有效地利用外源信息对参数先验信息进行反演修正，在随机变量反演方面有较多的应用.

综上所述，本文构建基于随机场-贝叶斯方法的滑面抗剪强度参数反演框架，研究更新前、后各滑面抗剪强度参数统计特征的变化规律. 基于PAWN方法，研究不同位置条块滑面抗剪强度参数变量对稳定性的敏感度. 在获取更新前、后滑面强度参数样本的基础上，探究库水位下降工况下滑面抗剪强度参数空间变异性对滑坡变形的影响规律.

贝叶斯理论能够利用多源信息数据，将岩土体参数变量（如渗透系数、黏聚力或摩擦角等）的先验信息更新为未知参数的后验信息分布^[12-19]. 多源信息数据可以分为直接信息和间接信息数据. 直接信息通常为通过地质勘察、室内试验和原位勘测等手段获得岩土体参数向量 ${{X}}$在特定空间位置的输入信息 ${g_i}$，间接信息数据包括滑坡的地表位移、深部位移和地下水位等监测信息以及其他工程输出信息 ${m_i}$，其中 $i$为岩土体采样点、探测点或监测点的空间位置. 后验概率密度函数为

(1) $f_{x\left| {{m_1},{m_2}, \cdots ,{m_{{l_1}}}} \right.}^{''}\left( x \right) = {{{C}}_1}f_x'\left( x \right)L\left( x \right).$

式中： $f_{x\left| {{m_1},{m_2}, \cdots ,{m_{{l_1}}}} \right.}^{''}\left( x \right)$为多组间接信息数据 ${m_i}$条件下的后验概率密度函数； ${l_1}$为间接信息的总数目； $f_x'\left( x \right)$为 ${{X}}$分量x的先验概率密度函数； ${{{C}}_1}$为比例系数； $L\left( x \right)$为似然函数，用来描述信息的不确定性，能够将直接信息 ${g_i}$或间接信息数据 ${m_i}$与预测信息 ${h_i}$考虑到似然函数中，以间接信息数据对变量进行更新，似然函数如下：

(2) $L\left( {{x_i}} \right) \propto {\rm{Pr}}\left( {{h_i}\left| {x = {m_i}} \right.} \right)\left( {i = 1,2,3, \cdots ,{l_1}} \right),$

(3) $L\left( x \right) = \prod {L\left( {{x_i}} \right) = {{{C}}_2}\prod {{\varphi _\zeta }\left[ {{h_i} - {m_i}} \right]} } .$

式中： ${h_i}$为计算模型对应间接信息数据 ${m_i}$的预测值，可以通过有限元与离散元数值模型或物理数学模型计算获得； ${\varphi _\zeta }$为偏差 $\zeta $的概率密度函数，通常假设每次的偏差 ${\zeta _i}$相互独立，表征为均值为0、标准差为某一值的正态分布. ${m_i}$与 ${h_i}$间的偏差 ${\zeta _i}$来源于测量设备、手段、模型计算等误差.

$f_{x\left| {{m_1},{m_2}, \cdots ,{m_{{l_1}}}} \right.}^{''}\left( x \right)$的求解至关重要，解析解的获取较困难，尤其是对于高维度变量，如随机场模型. 基于结构可靠度的贝叶斯更新方法可以解决该问题^[16,19]，该方法定义极限功能函数 ${{H}}$，处于功能函数失效区域 $\varOmega $上的样本即为更新后的变量样本：

(4) $ \begin{array}{*{20}{c}} {{{H}}\left( {d,{{X}}} \right) = d - {{C}}L\left( {{X}} \right)};&{\varOmega = \left[ {{{H}}\left( {d,{{X}}} \right) \leqslant 0} \right]} \end{array}.$

式中： $d$服从均匀分布， ${{C}}$为常数. 当采用结构可靠度方法时，通常需要将变量 $d$和 ${{X}}$转化为独立标准正态变量 ${{U}} = \left[ {{U_0};{U_1}; \cdots ;{U_{{l_2}}}} \right] \in {{\bf{R}}^{{l_2} + 1}}$，其中 ${l_2}$为随机场变量 ${{X}}$的数目，变换公式如下：

(5) $d = \varPhi \left( {{U_0}} \right).$

式中： $\varPhi $为标准正态累计概率密度函数.

(6) ${{X}} = {{T}}\left( {{U_1}; \cdots ;{U_{{l_2}}}} \right).$

式中： ${{T}}$为变量 ${{X}}$转换为标准正态变量的函数，可以通过Rosenblatt或Nataf变换来实现. 式（4）可以转化为下式，破坏区域为 ${\varOmega _U}$.

(7) $ \left.\begin{array}{c} {{H}}\left( {{U}} \right) = {U_0} - {\varPhi ^{ - 1}}\left\{ {{{C}}L\left[ {{{T}}\left( {{U_1}; \cdots ;{U_{{l_2}}}} \right)} \right]} \right\};\\ {\varOmega _U} = \left[ {{{H}}\left( {{U}} \right) \leqslant 0} \right]. \end{array}\right\}$

采用子集模拟的方法产生属于破坏区域中的样本，定义 ${{Z}}_{\rm{e}}=\left[{{U}}\in {\bf{R}}^{{l}_{2}+1}:{H}\left({{U}}\right) \leqslant 0\right]$为极限功能函数 ${{H}}$的失效事件. 假设事件 ${Z_{\rm{e}}}$由一系列中间事件 ${Z_1},{Z_2}, \cdots ,{Z_{{n_1}}}$组成，这些中间事件满足 $0 < {Z_1} < $ $ {Z_2} < \cdots < {Z_{{n_1}}}$关系，则 ${{{Z}}_{\rm{e}}}$事情的破坏概率 ${\rm{Pr}}\left( {{{{Z}}_{\rm{e}}}} \right)$可以表示为

(8) ${\rm{Pr}}\left( {{{{Z}}_{\rm{e}}}} \right) = {\rm{Pr}}\left( {\bigcap\limits_{w = 1}^{{n_1}} {{Z_w}} } \right) = \prod\limits_{w = 1}^{{n_1}} {{\rm{Pr}}\left( {{Z_w}\left| {{Z_{w - 1}}} \right.} \right)} .$

式中：中间事件 $ {Z}_{w}=\left[{{U}}\in {\bf{R}}^{{l}_{2}+1}:{H}\left({{U}}\right) \leqslant {b}_{w}\right]$, ${b_1} > $ $ {b_2} > \cdots > {b_w} = 0$.

子集模拟的简要步骤如下.

1）基于随机变量 $d$和 ${{X}}$的概率密度函数，产生 $J$组样本 ${{{U}}^j}(j = 1,2, \cdots ,J)$代入功能函数 ${{H}}$中计算得到J个系统的响应值. 按照升序排列，第 $\left( {1 - {p_0}} \right) J$个响应值是 ${b_1}$，符合 $ {{H}}\left({{U}}\right) \geqslant {b}_{1}$区域的事件概率 ${\rm{Pr}}\left( {{Z_1}} \right) = {p_0}$，其中 ${p_0} \in \left[ {0.1,0.3} \right]$，通常取 ${p_0} = 0.1$.

2）在第1）步产生 ${p_0} J$样本的基础上，根据自适应条件抽样方法产生 $\left( {1 - {p_0}} \right) J$个满足上述区域的样本，重新构成 $J$个样本，再次进行升序排列，则第 $(1 - $ $ {p_0} ) J$个响应值是 ${b_2}$. 这种情况下 ${\rm{Pr}}\left( {\left. {Z > {Z_2}} \right|Z > {Z_1}} \right)$的概率为 ${p_0}$，会有 ${p_0} J$个样本的系统响应值在 $ {{H}}\left({{U}}\right) \geqslant {b}_{2}$的区域内，用这 ${p_0} J$个样本作为抽样方法产生下一层 $\left( {1 - {p_0}} \right) J$个的“种子”.

3）不断重复步骤2），若模拟进行到第 ${n_1}$层后， ${p_0} J$个更新样本全部处于失效区域内，则停止后续模拟，整个模拟过程共产生 $J + \left( {{n_1} - 1} \right)\left( {1 - {p_0}} \right) J$个样本.

自适应条件抽样方法通过构建候选样本与抽取样本的联合分布函数，抽取符合破坏区域的样本^[20]. 该方法的优势如下：基于每次模拟过程中接受率参数 $p$，不断地调整联合分布函数，达到高效率抽样. 其中 $p$为马尔科夫链产生条件样本的效率，即为基于“种子”产生的条件样本符合破坏区域的概率. 自适应条件抽样方法的简要步骤如下.

在第k+1层子集模拟中，将k层模拟中的样本 $\left\{ {{{U}}_k^{\left( j \right)}:j = 1, \cdots ,{N_{\rm{s}}}} \right\},{N_{\rm{s}}} = {p_0} J$作为“种子”. 自适应过程需要分成多步 $i_{\rm{t}} = 1,2, \cdots ,{N_{\rm{a}}}$进行抽样，其中 ${N_{\rm{a}}}$为自适应步数，要满足 ${p_{\rm{a}}} = {{{N_{\rm{a}}}}/ {{N_{\rm{s}}}}}$. 在每一迭代i_t步中，从 $\left\{ {{{U}}_k^{\left( j \right)}:j = 1, \cdots ,{N_{\rm{s}}}} \right\}$无放回地选择 ${N_{\rm{a}}}$个样本作为条件抽样的种子，因此在该层子集模拟中的条件抽样可以分为 ${v_k} = \left( {i_{\rm{t}} - 1} \right){N_{\rm{a}}} + 1, \cdots , i_{\rm{t}} {N_{\rm{a}}}$次进行.

第 ${v_k}$次条件抽样基于已有的种子样本均值 $\;{{\mu }}$和互相关参数 $\;{{\rho }} = \sqrt {1 - {{{\sigma }}^2}} $组成的协方差矩阵，从多维正态分布中抽取 $\left( {1 - {p_0}} \right) J$个样本. $\;{{\mu }}$和 $\;{{\rho }}$可以通过计算 $\left\{ {{{U}}_k^{\left( j \right)}:j = 1, \cdots ,{N_{\rm{s}}}} \right\}$样本的均值和标准差来获得，公式如下.

(9) ${{\hat \mu }} = \frac{1}{{{N_{\rm{s}}}}}\sum\limits_{j = 1}^{{N_{\rm{s}}}} {{{U}}_k^j} .$

(10) ${{{\hat \sigma }}^2} = \frac{1}{{{N_{\rm{s}}} - 1}}\sum\limits_{j = 1}^{{N_{\rm{s}}}} {{{\left( {{{U}}_k^{\left( j \right)} - {{\hat \mu }}} \right)}^2}} .$

(11) ${{\sigma }} = \min \;\left\{ {{\lambda _{i_{\rm{t}}}}{{{\sigma }}_0},1.0} \right\}.$

(12) $\lg\; {\lambda _{i_{\rm{t}} + 1}} = \lg {\lambda _{i_{\rm{t}}}} + {\zeta _{i_{\rm{t}}}}\left[ {{{\hat p}_{i_{\rm{t}}}} - {p^*}} \right].$

(13) ${\hat p_{i_{\rm{t}}}} = \frac{1}{{{N_{\rm{a}}}}}\sum\limits_{{l_k} = 1}^{{N_{\rm{a}}}} {{{\hat E}}\left[ {{p_{{v_k}}}} \right]} .$

式中： ${\lambda _{i_{\rm{t}}}}$为缩放系数； ${\lambda _1} \in \left( {0,1.0} \right)$为初始缩放系数，可以随机生成； ${{{\sigma }}_0}$为基于式(10)计算每一层的初始样本得到的初始标准差； ${\hat p_{i_{\rm{t}}}}$为该步的接受率； ${p^*}$为最大接受率0.44； ${{\hat E}}\left[ {{p_{{v_k}}}} \right]$为每一 ${v_k}$步产生种子的平均接受率. 通过每一步接受率与最大接受率的比较，不断调整 ${\lambda _{i_{\rm{t}}}}$，实现自适应更新.

根据 ${{X}}$的统计特征（均值、变异系数、相关系数与波动范围），建立相应区域的随机场模型.

由于斜坡的稳定性通常由滑面上的平均强度决定，而不是某一点的强度，采用二维局部平均法^[21]进行随机场的离散. 局部平均法通过局部平均过程将点强度离散到一定的区域，具有更好的适应性. 为了建立滑面黏聚力 ${{c}}$与摩擦角 ${{\varphi}} $的二维平稳相关随机场模型 ${{X}}\left( {{z_x}{\text{、}}{z_y}} \right)$，其中 ${z_x}{\text{、}}{z_y}$为滑面位置点坐标. 将 ${{c}}$模拟为对数正态分布随机场 ${{{X}}_{{c}}}\left( {{z_x},{z_y}} \right)$，将 ${{\varphi}} $模拟为正态分布随机场 ${{{X}}_{{\varphi }}}\left( {{z_x},{z_y}} \right)$，计算公式如下.

(14) ${{\varPsi }} = \left[ {\begin{array}{*{20}{c}} 1&{{\rho _{{{c\varphi}} }}} \\ {{\rho _{{{c\varphi}} }}}&1 \end{array}} \right]{\rm{,}}\;{{{L}}_{{1}}}{{{L}}_{{1}}}^{\rm{T}} = {{\varPsi }},\;{{\varSigma }} = {{{L}}_{{2}}}{{{L}}_{{2}}}^{\rm{T}}.$

(15) $\begin{array}{l} {{\xi }} = \left[ {{{{\xi }}_{{c}}},{{{\xi }}_{{\varphi }}}} \right] = \left[ {{{\left[ {{\xi _{{c_1}}},{\xi _{c_2}}, \cdots ,{\xi _{c_q}}} \right]}^{\rm{T}}},{{\left[ {{\xi _{{\varphi _1}}},{\xi _{{\varphi _2}}}, \cdots ,{\xi _{{\varphi _q}}}} \right]}^{\rm{T}}}} \right]. \end{array} $

(16) $\begin{array}{l} {\tilde{ X}}\left( {{z_x},{z_y}} \right) = \left[ {{{{\tilde{ X}}}_{{c}}}\left( {{z_x},{z_y}} \right),{{{\tilde{ X}}}_{{\varphi }}}\left( {{z_x},{z_y}} \right)} \right] = {{{L}}_{{2}}}{{\xi }}{{{L}}_1}^{\rm{T}}. \end{array} $

(17) ${{{X}}_{{c}}}\left( {{z_x},{z_y}} \right) = \exp \;\left[ {{{{\tilde{ X}}}_{{c}}}\left( {{z_x},{z_y}} \right)} \right].$

(18) ${{{X}}_{{\varphi }}}\left( {{z_x},{z_y}} \right) = {{F}}_x^{ - 1}\left\{ {\varPhi \left[ {{{{\tilde{ X}}}_{{\varphi }}}\left( {{z_x},{z_y}} \right)} \right]} \right\}.$

式中： ${{\varPsi }}$为 ${{c}}$与摩擦角间互相关系数 $\;{\rho _{{{c\varphi}} }}$组合的互相关矩阵； ${{{L}}_{{1}}}$为 ${{\varPsi }}$经过乔列斯基分解得到的下三角矩阵； ${{{L}}_{{2}}}$为自相关矩阵 ${{\varSigma }}$经过乔列斯基分解得到的下三角矩阵； ${{\xi }}$为超立方拉丁抽样得到的独立标准正态随机变量；q为变量数目，等于条块总数； ${\tilde{ X}}\left( {{z_x},{z_y}} \right)$为相关的标准高斯随机场； ${{F}}_x^{ - 1}\left\{ {\varPhi \left[ \;{ } \right]} \right\}$为正态分布逆函数.

自相关矩阵 ${{\varSigma }}$由不同条块间 $\tilde i{\text{、}}\tilde j$的自相关系数 $\;\rho \left( {l{_{\tilde i}},l{_{\tilde j}}} \right)$组成，其中沿滑面分布任意2个条块的信息如图1所示，计算公式如下：

(19) $ \begin{split} \rho \left( {l{_{\tilde i}},l{_{\tilde j}}} \right) =& \frac{1}{{l{_{\tilde i}}l{_{\tilde j}}\cos \alpha '\cos \beta '}} \times \\ &\int_{{x_p}\left( {\tilde i} \right)}^{{x_k}\left( {\tilde i} \right)} {\left[ {\int_{{x_p}\left( {\tilde j} \right)}^{{x_k}\left( {\tilde j} \right)} \rho \left( {z_x,z_y} \right){\rm{d}}x\left( {\tilde j} \right)} \right]} {\rm{d}}x\left( {\tilde i} \right). \end{split} $

(20) $ \left. \begin{array}{l} {z_x} = \left| {{z_x}\left( {\tilde j} \right) - {z_x}\left( {\tilde i} \right)} \right|,\\ {z_y} = \left| {\tan \beta '\left[ {{z_x}\left( {\tilde j} \right) - {x_p}\left( {\tilde j} \right)} \right] - \tan \alpha '\left[ {{z_x}\left( {\tilde i} \right) - {x_p}\left( {\tilde i} \right)} \right] + } \right.\\ \quad\;\;\; \left. {{z_y}\left( {\tilde i} \right) - {y_p}\left( {\tilde j} \right)} \right|. \end{array} \right\} $

(21) $\rho \left( {z_x,z_y} \right) = \exp \;\left[ { - \left( {\frac{{2z_x}}{{{\delta _{\rm{h}}}}} + \frac{{2z_y}}{{{\delta _{\rm{v}}}}}} \right)} \right].$

式中： $\alpha '$和 $\;\beta '$为任意2个条块 $\tilde i{\text{、}}\tilde j$滑面的倾角， $l{_{\tilde i}}{\text{、}} $ $ l{_{\tilde j}}$为2个条块滑面长度， ${x_p}$、 ${y_p}$为条块滑面起点坐标， ${x_k}$、 ${y_k}$为条块滑面终点坐标， $z_x{\text{、}}z_y$分别为条块 $\tilde i{\text{、}}\tilde j$滑面内任意两点的水平和竖直距离， ${\delta _{\rm{h}}}$、 ${\delta _{\rm{v}}}$为水平波动范围和竖直波动范围.

为了探究滑面不同位置的抗剪强度参数对滑坡稳定性的贡献程度，采用新的全局敏感性分析方法PAWN^[22-23]来量化不同条块滑面的黏聚力与摩擦角变量对滑坡稳定系数变量的敏感性. 对于研究随机场变量对输出值的敏感性影响，该方法相比于常规的敏感性分析方法，不受输入变量和输出变量计算过程的限制，计算步骤相对简单，结果有较高的可靠性.

不同于基于方差的全局敏感性分析方法，PAWN通过计算输出变量 $y'$的非条件累计概率密度函数 ${F_{y'}}\left( {y'} \right)$和某一变量 $x'$条件下输出变量 $y'$的条件累计概率密度函数 ${F_{y'|x'}}\left( {y'|x'} \right)$最大距离的均值、最大值或中位数来表征某一变量对输出变量的贡献度 $S$，最大距离可以通过Kolmogorov–Smirnov test来获取. Pianosi等^[22]对PAWN进行一些改进，改进后的PAWN对产生输入变量的抽样方法不再限制，提高了计算效率，须设置的参数数目变少. 改进后的PAWN方法如下.

将每一个输入变量 $x'$等间距分成 $\tilde n$个子区间 ${\mathcal{I}_{\tilde w}}$，产生条件样本 ${\rm{Y}}{{\rm{C}}_{\tilde w}}$和非条件样本 ${{{Y}}_{\tilde U}}$. 改进后的敏感性指标由下式表示：

(22) $\hat S = {\rm{sta}}{{\rm{t}}_{\tilde w = 1,\cdots,\tilde n}}\mathop {\max }\limits_{y'} \left| {{{\hat F}_{y'}}\left( {y'} \right) - {{\hat F}_{y'|x'}}\left( {y'|x' \in {\mathcal{I}_{\tilde w}}} \right)} \right|.$

式中：stat 表示取均值、最大值或中位数，由用户自定义^[22]，结合实例对各种状态的结果进行比较； $\hat S$为0~1.0， $\hat S$越小表明变量 $x'$对 $y'$的影响越小； ${\hat F_{y'}}\left( {y'} \right)$为样本 ${{{Y}}_{\tilde U}}$的经验累积分布函数； ${\hat F_{y'|x'}} \left( {y'|x' \in {\mathcal{I}_{\tilde w}}} \right)$为样本 ${\rm{Y}}{{\rm{C}}_{\tilde w}}$的经验累积分布函数.

基于式（22）计算得到敏感性指标具有一定的近似误差，尤其是当样本容量较小时，评估输入变量敏感性指标的可靠性是很重要的. 伪参数 ${\hat S_{{\rm{dummy}}}}$的引入可以用于评估一定样本容量下敏感性指标的精度，通过产生的伪参数 ${\hat S_{{\rm{dummy}}}}$与 $\hat S$进行比较，识别输入变量对结果的影响程度. 当 $\hat S$> ${\hat S_{{\rm{dummy}}}}$时，说明输入因素对输出结果重要. 对于 ${\hat S_{{\rm{dummy}}}}$的计算，利用Bootstrap抽样方法对非条件样本 ${{{Y}}_{\tilde U}}$进行一定次数的重抽样. 利用式（22）计算重抽样子集的敏感性指标分布，得到置信区间，判断敏感性指标的差异能否足够区分相关的输入变量. 将条块的黏聚力和摩擦角视为输入因素，共有2q个变量，通过随机场方法产生，每一输入变量样本数为10 000个， $\tilde w$通常取默认设置值10.

白水河滑坡位于长江主干道的南岸，距三峡大坝坝址56 km，滑坡前缘临空面较陡（高度为30~50 m，坡度为20°~30°），滑面平均倾角为25°. 白水河滑坡划分为预警区和非预警区两部分（见图2）. 滑坡地质勘察和滑坡监测信息表明，变形破坏区域主要分布在预警区内，后缘和边界处出现剧烈变形，预警区外没有明显的宏观变形迹象，钻孔取出的岩芯没有明显的滑带.

自三峡水库蓄水到135 m后，滑坡内出现强烈的宏观变形. 2003年6月到2006年9月，库水位在135 m到138 m间波动，坡体变形稳定. 2006年9月到2007年2月，库水位逐步升到155 m；到2007年6月，库水位首次下降至145 m，伴随着季节性强降雨（每年6月~9月），滑坡产生了较大的一次阶跃变形，具有加速滑动的迹象，位于3-3′剖面处发现多处裂缝^[24-25]. 如图2所示，滑坡后缘沿江公路南侧（山上方向）约10⁴ m³的岩土体发生塌滑，滑坡东侧和后缘边界基本贯通，西侧边界裂缝呈羽状断续展布. 卢书强等^[24]根据专业监测成果与勘查报告分析认为，白水河滑坡在库水位快速下降和强降雨共同作用下是不稳定的. 根据滑坡不同发育阶段的稳定系数评估方法^[9-10]可知，滑坡3-3′剖面在当时自重、强降雨与库水位首次下降的工况下处于不稳定状态^[25]，稳定系数为0.95~1.00，为该工况下的反演分析提供了条件.

滑坡体物质由碎石、块石、角砾和粉质黏土、黏土组成，厚度为7.75~38.50 m，总体趋势是滑坡后部滑体厚度小，前缘滑体厚度大. 滑带以含碎石或者含角砾粉质黏土为主，部分滑带岩土的物质成分为角砾土和黏土，滑面可见明显的磨光面和擦痕. 滑带厚度为0.20~1.30 m，平均厚度为0.7 m.

滑坡稳定系数计算和滑坡渗流变形模拟所需的参数有容重 $\gamma $、弹性模量 $E$、泊松比 $\;\mu$和饱和渗透系数 ${K_{\rm{s}}}$等，如表1所示. 为了确定滑带土黏聚力和摩擦角的初始统计特征，根据李远耀等^[26]对三峡库区400处堆积层滑坡滑带土黏聚力和摩擦角的统计结果，得出黏聚力主要服从对数正态分布，摩擦角服从正态分布. 分别给出I-V类易滑地层滑坡的滑带土黏聚力和摩擦角的统计特征，其中白水河滑坡所处地层为侏罗系香溪组，滑带土的基本特征与I类易滑地层滑坡描述的滑带土特征相同. 经统计可得，I类易滑地层滑坡的滑带土天然峰值黏聚力和摩擦角的均值分别为32.9 kPa、19.6°，残余黏聚力和摩擦角的均值分别为30.4 kPa、13.7°. 基于白水河滑坡勘察报告^[25]中滑带土样的天然固结快剪试验结果，综合确定黏聚力和摩擦角的均值分别为32 kPa、15°，标准差分别为8 kPa、2°. 经统计可知，岩土体黏聚力和摩擦角的相关系数为−0.2~−0.7^[27]，本文取−0.5；不同成因土性参数的 ${\delta _{\rm{h}}}$一般为20~80 m^[28], ${\delta _{\rm{v}}}$一般为2~6 m，且 ${\delta _{\rm{v}}}$的波动范围通常为 ${\delta _{\rm{h}}}$波动范围的1/10. 取抗剪强度参数的 ${\delta _{\rm{h}}}$为40 m， ${\delta _{\rm{v}}}$为4 m. 若获取的先验信息越接近实际，则反演更新的效率越高；反之，则反演更新的效率越低. 基于贝叶斯方法的更新机制可知，反演后的样本会不断地接近真实情况.

基于自适应条件抽样算法的贝叶斯方法的具体步骤如下.

1）建立白水河滑坡3-3′剖面的二维地质模型，条块划分见图3，共划分了55个条块，导出所有条块的滑面位置坐标信息. 图中，H_i为高程，D_i为距离.

2）基于黏聚力与摩擦角的先验统计信息特征值，如均值、变异系数、波动范围与变量间的相关系数等，建立滑面的初始随机场模型，生成随机变量样本 ${{X}}$.

3）采用剩余推力法计算白水河滑坡3-3′剖面的稳定系数，其中地下水位由有限元软件通过施加相应的时间段的工况而获取.

4）基于白水河滑坡3-3′剖面的实际稳定性状况（欠稳定），评估该剖面的实际稳定系数为0.95~1.00. 结合正态分布95%的置信区间，确定 ${m_i}$为0.975，偏差 ${\zeta _i}$为0.01.

5）在结构可靠度的贝叶斯更新框架下，利用白水河滑坡3-3′剖面的实际稳定性状况与预测状况的差异，对滑面随机场模型的变量样本进行筛选，获得更新后的变量.

6）基于更新后的变量，采用非侵入式随机有限元法^[11]开展白水河滑坡3-3′剖面的渗流变形分析，探究滑面抗剪强度参数空间变异下的滑坡变形以及反演更新前、后样本下的变形差异.

依据3-3′剖面的地质模型，建立如图4所示的有限元模型，滑带厚度为0.5 m. 模型采用孔压/位移耦合的CPE4RP平面应变四边形单元，该减缩积分单元在满足位移精度要求的条件下计算代价相对较低. 坡体共划分747个四边形单元、884个节点，滑带共划分170个四边形单元、342个节点. 采用理想弹塑性本构和摩尔库仑破坏准则，模型底部约束竖直位移和水平位移，滑坡前缘库水位以下部分为水头边界，水位以上部分为降雨边界.

参数反演工况为2007年2月1日到2007年7月1日，库水位从155 m降到145 m，期间并施加强降雨^[25]. 为了探究库水位下降条件下滑面抗剪强度参数空间变异对滑坡变形的影响，分析反演更新前、后样本下的滑坡变形差异，采用的研究工况为库水位从175 m下降到145 m.

为了研究式（22）中stat 取不同状态（均值、最大值或中位数）对敏感性结果的影响，结合PAWN方法，分别计算不同状态下55个条块的黏聚力与摩擦角共110个变量对稳定系数的敏感度指标，如图5所示. 图中，1~55号为对应条块号的黏聚力变量，56~110号为对应条块号的摩擦角变量. 从图5可以看出，stat 取均值和中位数的敏感度指标结果基本相同，stat 取最大值的敏感度指标略大于其他2种结果，且变量间敏感度指标的波动幅度相对较大，但从整体来看，3种状态的敏感度指标变化趋势较相似，均能够反映出不同条块敏感度指标的变化趋势. 总之，由于最大值为数据出现的极端情况，不能有效地反映数据的整体性，建议stat 取均值和中位数来进行分析.

基于stat 取均值的敏感度计算结果进行分析. 结果表明， ${S_{{\rm{dummy}}}}$为0.02，所有变量的 $S$均大于 ${S_{{\rm{dummy}}}}$，说明各条块滑面的抗剪强度参数变量对滑坡稳定系数的影响不可忽略. 如图6所示为条块滑面抗剪强度参数对稳定系数的敏感度. 可以看出，黏聚力变量和摩擦角变量的敏感度指标从滑坡后缘到前缘的变化趋势基本相同. 敏感度指标较小的条块号为28~40，敏感度指标相对较大的区域主要分布在滑坡前缘和滑坡后缘道路附近. 滑面抗剪强度参数的敏感性分析表明，滑面不同位置的参数变化对稳定系数有着不同的影响作用，敏感度较大的区域参数设置较为重要，对滑坡的稳定性起着较大的作用，在滑坡加固设计时应优选这些区域^[28].

基于结构可靠度的贝叶斯自适应反演流程，共更新10 000组抗剪强度参数随机场样本. 在该过程中，共经历了17层子集模拟. 如图7所示为10 000组稳定系数F_S的先验与后验概率密度函数（PDF）分布曲线. 可见，更新前的稳定系数为1.119~1.233，稳定系数的均值和标准差分别为1.176和0.034；更新后的稳定系数为0.95~1.00，均值和标准差分别为0.986和0.011. 可见，更新后的稳定系数区间与预定区间吻合，说明该方法的有效性.

如图8、9所示分别为3个不同位置处条块黏聚力与内摩擦角稳定系数的先验PDF分布曲线和后验PDF分布曲线，条块分别处于滑坡前、中与后部位置. 可以看出，与先验PDF分布曲线相比，更新后PDF分布曲线变瘦、变高，即更新后条块滑面的抗剪强度参数变量的变异性变小.

如图10、11所示分别为反演前、后各条块黏聚力与内摩擦角的均值 $\bar c {\text{、}}\bar\varphi$的变化. 反演前各条块的黏聚力均值约为32.75 kPa，这是因为建立的初始随机场是平稳的，即变量均值不随位置而变化. 从图10可知，更新后条块黏聚力均值远大于先验均值的数目占总条块数的56%，小于或接近于先验均值的数目占总条块数的44%，更新后所有条块黏聚力均值的平均值为33.77 kPa，总体上更新后条块黏聚力均值略呈增大趋势. 对于摩擦角均值，反演前各条块的摩擦角均值约为14.5°，更新后所有条块的摩擦角均值都小于该值，总体上有从后缘到前缘的下降趋势，更新后所有条块摩擦角均值的平均值为11.64°. 由此可知，平稳随机场不一定适用于所有的岩土体，参数具有平稳特性的岩土体在一定范围内的形成过程和沉积环境相对稳定，没有经历过地质体的局部运动，岩土体的物质成分与结构随空间位置的变化不大. 说明一些岩土工程对非平稳随机场模拟的必要性，对滑坡来讲，滑坡的演化伴随着滑带内部岩土体的物理化学成分、颗粒组成、受力状态、地质环境与力学性质的不断改变，滑坡的定向滑动可能使得岩土体的结构与物理力学性质存在一定的变化趋势. 以黄土坡临江1号滑坡为例，滑带岩土体的颗粒级配、物质组成与力学参数沿着主滑方向呈现一定的变化趋势^[3].

如图12、13所示为反演前、后各条块黏聚力与内摩擦角的标准差变化曲线. 可知，更新后大多数条块的黏聚力与内摩擦角的标准差均小于反演前的标准差，表明经过后验信息的更新，输入参数的不确定性减小. 反演前、后从滑坡后缘到前缘条块滑面参数的标准差变化幅度发生一定的改变. 由随机场理论可知，各条块滑面参数的标准差会受波动范围的影响，波动范围越大，一定范围内条块间参数的变异性越小，反之则越大. 可知，反演后参数的波动范围发生了改变.

为了探究更新前、后随机场参数变化对滑坡变形的影响，对库水位从175 m下降到145 m时的滑坡位移进行分析. 如图14~16所示为参数更新前、后滑坡在A、B剖面处的位移曲线，其中图16为位移变异系数曲线. 图中，h_s为剖面深度，d_s为剖面位移. 可以看出，反演前、后位移沿深度的变化趋势一致，但反演后A、B剖面处各深度位移的最值、均值与变化范围均比反演前的大，这符合反演后滑坡稳定性较差的状况. 从图16可知，滑面附近的位移变异系数较其他深度部位的大，说明抗剪强度参数的变化对滑面的位移影响相对较大. 对于A剖面，更新后全部深度范围的变异系数比更新前小；对于B剖面，0~27 m深度处的变异系数比更新前大，大于27 m深度处的变异系数比更新前小. 这由于B剖面位于滑坡前缘，该部位受库水波动的作用明显，位移变形容易出现极端的情况. 比较A、B剖面的位移可知，A剖面的位移最值与均值大于B剖面，说明滑坡在A剖面处变形更明显，位移变化特征符合白水河滑坡的变形规律.

（1）参数反演结果表明，反演前建立的平稳随机场经过更新后变为非平稳随机场，尤其是反演后的摩擦角均值变化趋势更明显. 结合滑带的形成环境分析可知，非平稳随机场的模拟相对适合于滑带. 对于条块滑面参数的标准差，反演后多数条块的抗剪强度参数标准差减小，变化幅度受到影响.

（2）对于白水河滑坡3-3′剖面，全局敏感性分析表明，各条块滑面的抗剪强度参数变量对滑坡稳定系数的影响不可忽略，变量对滑坡稳定系数影响较大的区域分布在滑坡前缘和滑坡后缘道路附近，滑坡防治工程应优选这些影响较大的区域.

（3）经过贝叶斯反演后，抗剪强度参数的不确定性减小，滑坡的位移变形增大，滑坡前缘坡表的位移波动变化较大，能够为滑坡分析与防治提供较危险的工况. 与传统的确定性分析相比，本文既考虑了滑带抗剪强度参数的空间变异性，也考虑了滑坡的后验稳定性，可以获得滑坡变形的极端情况及位移主要分布区间.