多层线性模型在碳钢土壤腐蚀规律中的应用研究

D0L:10.13374/.issn1001-053x.2013.11.008 第35卷第11期 北京科技大学学报 Vol.35 No.11 2013年11月 Journal of University of Science and Technology Beijing Nov.2013 多层线性模型在碳钢土壤腐蚀规律中的应用研究 鲁庆，穆志纯 北京科技大学自动化学院，北京100083 ☒通信作者，E-mail:mu@ies.ustb.edu.cm 摘要针对材料在自然士土壤环境中的腐蚀预测问题，提出应用多层线性模型对材料的腐蚀规律进行研究.多层线性模 型是分析具有层次结构数据的新型统计技术，对于地区间土壤腐蚀规律的差异性特点和单地区实验观测样本不足的问 题，可以为不同的区域分别建立腐蚀率模型，使模型假设与实际更为吻合.以碳钢在土壤中的腐蚀数据为研究对象，建 立了腐蚀率的多层线性模型.实验验证了该模型可以准确地拟合和预测碳钢在土壤的腐蚀率变化，优于指数平滑算法和 差分自回归移动平均算法. 关键词碳钢：钢腐蚀：土壤：腐蚀率：多层线性模型 分类号TG172.4:TP391.9 Application study of hierarchical linear models in carbon steel corrosion in soils LU Qing,MU Zhi-chun凶 School of Automation and Electrical Engineering,University of Science and Technology Beijing,Beijing 100083,China XCorresponding author,E-mail:mu@ies.ustb.edu.cn ABSTRACT Aiming at the problem of materials corrosion prediction in the natural soil environment,hierarchical linear modeling was used to study the corrosion rules of materials.The hierarchical linear model is a new statistical technique for hierarchical data analysis.Different models for every region can be established to fit regional soil corrosion characteristics and solve the problem of small samples in a specific region,which makes modeling assumptions more consistent with actual situations.A hierarchical linear corrosion model was built with carbon steel corrosion data in soils being the object of study.The test demonstrates that the model can precisely describe and predict the corrosion rate of carbon steel in soils,and it is superior to exponential smoothing algorithm and differential autoregressive integrated moving average (ARIMA)algorithm. KEY WORDS carbon steel:steel corrosion:soils:corrosion rate:hierarchical linear model 材料在土壤环境的腐蚀规律是自然环境腐蚀研 是为每个地区建立特有的土壤腐蚀性模型.但是， 究的一个重要内容，对地下工程的施工选材防护和 由于自然环境实验的困难性，实验观测数据往往存 材料寿命评估等具有重要意义.近年来，国内外研 在高维度、小样本等问题，具体到本文研究的所用 究机构和学者利用多元回归叫、时间序列分析②、 数据，每个地区只有23组观测数据，很难满足多 神经网络B-4、灰色系统问等方法建立各种腐蚀 元回归、时间序列等方法的建模数据需求.神经网 规律模型，取得了很多成果和进展. 络方法在解释土壤腐蚀规律影响因素方面也存在局 由于受到各环境因素相互影响，不同地区的土 限性. 壤腐蚀规律具有很大的相异性.要准确描述并预 本文应用多层线性模型的原理，以碳钢为例， 测材料在土壤中的腐蚀变化，一个比较合理的策略 把材料在土壤中的腐蚀变化过程用一个两层分层模 收稿日期：2012-10-21

第 35 卷 第 11 期 北 京 科 技 大 学 学 报 Vol. 35 No. 11 2013 年 11 月 Journal of University of Science and Technology Beijing Nov. 2013 多层线性模型在碳钢土壤腐蚀规律中的应用研究 鲁 庆，穆志纯 北京科技大学自动化学院，北京 100083 通信作者，E-mail: mu@ies.ustb.edu.cn 摘 要 针对材料在自然土壤环境中的腐蚀预测问题，提出应用多层线性模型对材料的腐蚀规律进行研究. 多层线性模 型是分析具有层次结构数据的新型统计技术，对于地区间土壤腐蚀规律的差异性特点和单地区实验观测样本不足的问 题，可以为不同的区域分别建立腐蚀率模型，使模型假设与实际更为吻合. 以碳钢在土壤中的腐蚀数据为研究对象，建 立了腐蚀率的多层线性模型. 实验验证了该模型可以准确地拟合和预测碳钢在土壤的腐蚀率变化，优于指数平滑算法和 差分自回归移动平均算法. 关键词 碳钢；钢腐蚀；土壤；腐蚀率；多层线性模型 分类号 TG172.4; TP391.9 Application study of hierarchical linear models in carbon steel corrosion in soils LU Qing, MU Zhi-chun School of Automation and Electrical Engineering, University of Science and Technology Beijing, Beijing 100083, China Corresponding author, E-mail: mu@ies.ustb.edu.cn ABSTRACT Aiming at the problem of materials corrosion prediction in the natural soil environment, hierarchical linear modeling was used to study the corrosion rules of materials. The hierarchical linear model is a new statistical technique for hierarchical data analysis. Different models for every region can be established to fit regional soil corrosion characteristics and solve the problem of small samples in a specific region, which makes modeling assumptions more consistent with actual situations. A hierarchical linear corrosion model was built with carbon steel corrosion data in soils being the object of study. The test demonstrates that the model can precisely describe and predict the corrosion rate of carbon steel in soils, and it is superior to exponential smoothing algorithm and differential autoregressive integrated moving average (ARIMA) algorithm. KEY WORDS carbon steel; steel corrosion; soils; corrosion rate; hierarchical linear model 材料在土壤环境的腐蚀规律是自然环境腐蚀研 究的一个重要内容，对地下工程的施工选材防护和 材料寿命评估等具有重要意义. 近年来，国内外研 究机构和学者利用多元回归[1]、时间序列分析 [2]、 神经网络 [3−4]、灰色系统 [5] 等方法建立各种腐蚀 规律模型，取得了很多成果和进展. 由于受到各环境因素相互影响，不同地区的土 壤腐蚀规律具有很大的相异性 [6] . 要准确描述并预 测材料在土壤中的腐蚀变化，一个比较合理的策略 是为每个地区建立特有的土壤腐蚀性模型. 但是， 由于自然环境实验的困难性，实验观测数据往往存 在高维度、小样本等问题，具体到本文研究的所用 数据，每个地区只有 2∼3 组观测数据，很难满足多 元回归、时间序列等方法的建模数据需求. 神经网 络方法在解释土壤腐蚀规律影响因素方面也存在局 限性. 本文应用多层线性模型的原理，以碳钢为例， 把材料在土壤中的腐蚀变化过程用一个两层分层模 收稿日期：2012-10-21 DOI:10.13374/j.issn1001-053x.2013.11.008

第11期 鲁庆等：多层线性模型在碳钢土壤腐蚀规律中的应用研究 .1459· 型进行描述，将对每个地区的累次观察视为归属于 看出，对于具有层级结构的数据，如果使用传统 该地区，尝试建立单独的地区腐蚀率模型，较好地 的回归分析参数估计算法，得到的参数估计将是有 拟合土壤腐蚀规律的地区性特征，并在一定程度上 偏的 解决小样本等问题 1.2参数估计 多层线性模型利用经验贝叶斯方法，对不同层 1多层线性模型 次的参数估值进行加权，充分利用数据中包含的信 多层线性模型最早应用在社会科学领域，由 息，求得最大似然估值.以式(1)和式(2)所示模型 Lindley和Smith可在1972年提出.与传统的回归 为例，其参数估计过程如下 模型比起来，它不受方差齐性和样本独立性等传统 首先，对式(1)进行最小二乘(OLS)回归方法 统计假设的严格限制，在分析具有层级结构数据方 估值，得到B的第一组估值： 面有很大的优越性母，在心理、经济、教育等多个 学科领域9-1山得到广泛应用，但在自然科学研究 了ioi=卫-a1id; 领域还没有得到应有的重视 、ai:=∑(ai-aa)Ya-a)/∑(a4-a). (4) 1.1数学模型 以二层模型为例，当利用多层线性模型描述个 然后，利用广义最小二乘法(GLS)估计固定效 体变化时，数学模型可以表示如下 应Y. 第一层模型： 定义精度△： Yui =Boi Biati+rti. (1) 46=Var(o:+(:-Bo:）=0+oi(⑤) △=Var(u1+(间：-：》=T1+ 第二层模型： 式中：To0=Var(o),是随机误差4o:的方差：T1= B0:=Y00+W01+0: (2) Var(1i),是随机误差h1i的方差：koi=Var(oi- o),1:=Var(1i-),分别表示由于估计值误 B1i=Y10+Wh1+41i 差引起的差异.可得到y的估值： 第一层模型描述了个体的增长变化，具体到本 00= ∑Ao:7i/∑Ao: 文的研究对象，为特定地区的土壤腐蚀变化规律， 701= 其中Yi表示i地区在t时点的土壤腐蚀率：ai 4o:W:/∑Ao:W: (6) 为时间变量，表示在时点t的材料埋存时间：0 10= △1iY/∑A1 和3：表示i地区腐蚀率模型的截距和斜率：T为 11= 随机误差，r~N(0,2).第一层模型看上去与线 由式(2)可得到参数B的第二组估值： 性回归模型相似，但截距0：和斜率3：不再是常 结：=o0+01W; 数，而是依赖于第二层模型的定义.在第二层模型 (7) i:-10+1W. 中，W:为第二层预测变量，表示不同地区的环境参 多层线性模型用经验贝叶斯(EB)算法对参数B的 数；Yo0、o1、Y10和11分别表示对应二层模型的 两组估值进行加权，最终估值为 截距和斜率：4gi~N(0,T),为地区间偏差；T为 方差协方差矩阵：通过对Y、B和各层方差协方差 6=A0:o:+(1-0)6; (⑧) 6、T进行参数估计，就可以求得单独的地区腐蚀 31=入131i+(1-入1)31 模型 其中入为估值信度，定义为： 将式(2)合并式(1)，可得到混合方程： Aoi T00/(To0 +koi); Yu =Y00+WiY01+Qu1o+auWiY1+Hoi+atHi+7i. 入1i=T1/（T1+ (3) 在本文的实验中，土壤腐蚀台站的单站数据只 在式(3)中，o+a1i+Yi是残差项，其中偏 有23组观测数据，按照传统方法，无法进行可靠 差o:和1：对于同一地区是分别相等的，在不同 的统计分析.如将各台站数据作为一个整体样本进 的地区间则存在差异，只有当不存在组间偏差，例 行建模，又会因为无法满足样本独立同分布假定而 如0和山1：为0时，方程(3)才能满足最小二乘 造成估值有偏.从式（⑧）可以看出，当估值可靠性 (OL$)回归方法关于方差齐次性的假定.在此可以 入较高时，o:在最终结果中将占据较大的比重：如

第 11 期 鲁 庆等：多层线性模型在碳钢土壤腐蚀规律中的应用研究 1459 ·· 型进行描述，将对每个地区的累次观察视为归属于 该地区，尝试建立单独的地区腐蚀率模型，较好地 拟合土壤腐蚀规律的地区性特征，并在一定程度上 解决小样本等问题. 1 多层线性模型 多层线性模型最早应用在社会科学领域，由 Lindley 和 Smith[7] 在 1972 年提出. 与传统的回归 模型比起来，它不受方差齐性和样本独立性等传统 统计假设的严格限制，在分析具有层级结构数据方 面有很大的优越性 [8]，在心理、经济、教育等多个 学科领域 [9−11] 得到广泛应用，但在自然科学研究 领域还没有得到应有的重视. 1.1 数学模型 以二层模型为例，当利用多层线性模型描述个 体变化时，数学模型可以表示如下. 第一层模型： Yti = β0i + β1iαti + rti. (1) 第二层模型： β0i = γ00 + Wiγ01 + µ0i ; (2) β1i = γ10 + Wiγ11 + µ1i . 第一层模型描述了个体的增长变化，具体到本 文的研究对象，为特定地区的土壤腐蚀变化规律， 其中 Yti 表示 i 地区在 t 时点的土壤腐蚀率；αti 为时间变量，表示在时点 t 的材料埋存时间；β0i 和 β1i 表示 i 地区腐蚀率模型的截距和斜率；rti 为 随机误差，rti ∼ N(0, δ2 ). 第一层模型看上去与线 性回归模型相似，但截距 β0i 和斜率 βli 不再是常 数，而是依赖于第二层模型的定义. 在第二层模型 中，Wi 为第二层预测变量，表示不同地区的环境参 数；γ00 、γ01 、γ10 和 γ11 分别表示对应二层模型的 截距和斜率；µqi ∼ N(0, T ) ，为地区间偏差；T 为 方差协方差矩阵；通过对 γ 、β 和各层方差协方差 δ 、T 进行参数估计，就可以求得单独的地区腐蚀 模型. 将式 (2) 合并式 (1)，可得到混合方程： Yti = γ00+Wiγ01+αtiγ10+αtiWiγ11+µ0i+αtiµ1i+γti. (3) 在式 (3) 中，µ0i + αtiµ1i + γti 是残差项，其中偏 差 µ0i 和 µ1i 对于同一地区是分别相等的，在不同 的地区间则存在差异，只有当不存在组间偏差，例 如 µ0i 和 µ1i 为 0 时，方程 (3) 才能满足最小二乘 (OLS) 回归方法关于方差齐次性的假定. 在此可以 看出，对于具有层 级结构的数据，如果使用传统 的回归分析参数估计算法，得到的参数估计将是有 偏的. 1.2 参数估计 多层线性模型利用经验贝叶斯方法，对不同层 次的参数估值进行加权，充分利用数据中包含的信 息，求得最大似然估值. 以式 (1) 和式 (2) 所示模型 为例，其参数估计过程如下. 首先，对式 (1) 进行最小二乘 (OLS) 回归方法 估值，得到 β 的第一组估值： ( βˆ 0i = Y¯ − βˆ 1iα¯ti; βˆ 1i = X(αti − α¯ti)(Yti − Y¯ ti)/ X(αti − α¯ti) 2 . (4) 然后，利用广义最小二乘法 (GLS) 估计固定效 应 γ. 定义精度 ∆： ( ∆ −1 0i = Var(µ0i + (βˆ 0i − β0i)) = τ00 + k0i ; ∆ −1 1i = Var(µ1i + (βˆ 1i − β1i)) = τ11 + k1i . (5) 式中：τ00 = Var(µ0i)，是随机误差 µ0i 的方差；τ11 = Var(µ1i)，是随机误差 µ1i 的方差；k0i = Var(βˆ 0i − β0i)，k1i = Var(βˆ 1i − β1i)，分别表示由于估计值误 差引起的差异. 可得到 γ 的估值：    γˆ00 = X∆0iY¯ ti/ X∆0i ; γˆ01 = X∆0iWiβˆ 0i/ X∆0iW2 i ; γˆ10 = X∆1iY¯ ti/ X∆1i ; γˆ11 = X∆1iWiβˆ 1i/ X∆1iW2 i . (6) 由式 (2) 可得到参数 β 的第二组估值： ( βˆ∗ 0i = ˆγ00 + ˆγ01Wi ; βˆ∗ 1i = ˆγ10 + ˆγ11Wi . (7) 多层线性模型用经验贝叶斯 (EB) 算法对参数 β 的 两组估值进行加权，最终估值为 ( βˆ∗∗ 0i = λ0iβˆ 0i + (1 − λ0i)βˆ∗ 0i ; βˆ∗∗ 1i = λ1iβˆ 1i + (1 − λ1i)βˆ∗ 1i . (8) 其中 λ 为估值信度，定义为： ( λ0i = τ00/(τ00 + k0i); λ1i = τ11/(τ11 + k1i). 在本文的实验中，土壤腐蚀台站的单站数据只 有 2∼3 组观测数据，按照传统方法，无法进行可靠 的统计分析. 如将各台站数据作为一个整体样本进 行建模，又会因为无法满足样本独立同分布假定而 造成估值有偏. 从式 (8) 可以看出，当估值可靠性 λ 较高时，βˆ 0i 在最终结果中将占据较大的比重；如

.1460 北京科技大学学报 第35卷 果可靠性入较低，：的估值将在更大程度上影响 从自由度J-N-1的T分布，其中Bg是B的估 最终结果.多层线性模型通过这种加权的方式，可 值，(W,)P2是，的标准误差，J为第一层模型 以借助其他二层单位的样本数据中包含的信息，对 的样本数，N为第一层模型的自变量个数.如假设 取样很少的单位也可以进行分析建模，相对解决了 成立，则表示对应的因子与因变量Y无关系 小样本问题 为便于说明，上述的计算过程假设各层方差协 2结果与讨论 方差参数已知，在实际计算中，需要对方差、协方差 2.1研究样本和资料库 参数6和T进行估计.通常使用最大似然估计(ML)、 本文以国家材料腐蚀站网中碳钢在土壤中的 约束最大似然估计(REML)、贝叶斯估计等估计方 腐蚀数据进行分析和建模，对不同土壤环境因素下 法.以约束最大似然估计为例，其原理就是通过选 碳钢平均腐蚀率进行预测.按照多层模型分析的需 择参数6和T的估计值，应用EM算法、费舍尔得 要，将原始数据整理为两层数据表，第一层数据如 分算法等，在迭代过程中使观测到这一实际数据Y 表1所示，包含了实样站点、埋样时间和实验材料 的似然值最大并趋于稳定.限于篇幅，在这里就不 的腐蚀率记录共37条数据：第二层数据如表2所 再展开讨论，具体讨论可参见文献[12-13). 示，描述了各试样站的土壤环境数据，包括土壤pH 1.3模型检验 值、有机质含量、全氮含量等理化参数，共21条数 可以通过单参数检验对固定效应参数Y进行 据，涵盖了21个地区试验台站.考虑到数据的可比 假设检验，以考察Y是否具有显著效应14.对于任 性，对表2数据进行了Z-SCORE标准化处理 意Ygs,单参数检验的假设为Ho:4s=0,检验统 计量t=as/(W)/2,服从自由度J-N-1的T分 表1不同区域碳钢的腐蚀率数据 Table 1 Corrosion rate of carbon steel at different sites 布，其中gs是Ygs的估值，(Vg)/2是as的标准 站点 埋存时间/a 腐蚀率/(mma1) 误差，J为第二层模型的样本数，N为第二层模型 沈阳站 3.719 的自变量个数.如假设成立，则表示Ygs对应的因 3 3.215 子对预测因变量B无意义. 4.155 成都站 随机效应参数B的假设检验方式与固定效应 3 2.343 4.363 检验方式类似.对于任意，单参数检验的假设 西安站 2.912 为Ho:月=0，检验统计量t=月/W,)/2,服 5 2.120 表2站点土壤理化性质数据示例 Table 2 Physical and chemical properties of soils at test sites 站点 pH值有机质质量分数/%全氮含量（质量分数/%- 可溶性盐离子的质量分数/% Na+ 土壤电阻/2 HCO3 s0- Mg+ K+ 大庆站10.0500 1.0250 0.0525 0.0669 0.0142 0.0008 0.0033 0.0516 5.1000 华南站5.4500 0.3250 0.0140 0.0031 0.0049 0.0006 0.0002 0.0007 460.0000 沈阳站6.7500 2.3800 0.1090 0.0126 0.0118 0.0017 0.0004 0.0049 32.9000 2.2建立模型 在着显著变异：方差=5.563，X2=179.7,P<0.001. 2.2.1零模型 这说明如果为所有地区建立一个统一的模型，其参 为考察不同地区土壤腐蚀率是否存在明显差 数估计结果是有偏的，即利用多层线性模型进行建 异，首先建立零模型如下 模是必要的.可以通过计算跨级相关系数度量变异 第一层模型： 分布， Yi Boi rt. (9) 第二层模型： ICC Var(uoi)/[Var(uoi)+Var(rti)]= 5.563/(5.563+1.012)=0.861, 00i=Y00+40i, 参数估计结果如表3所示.考察随机部分参 即由不同地区间特性导致的变异占总变异的比例为 数估计结果可以看出，不同地区站点的腐蚀规律存 86.1%

· 1460 · 北 京 科 技 大 学 学 报 第 35 卷 果可靠性 λ 较低，βˆ∗ 0i 的估值将在更大程度上影响 最终结果. 多层线性模型通过这种加权的方式，可 以借助其他二层单位的样本数据中包含的信息，对 取样很少的单位也可以进行分析建模，相对解决了 小样本问题. 为便于说明，上述的计算过程假设各层方差协 方差参数已知，在实际计算中，需要对方差、协方差 参数 δ 和T进行估计. 通常使用最大似然估计 (ML)、 约束最大似然估计 (REML)、贝叶斯估计等估计方 法. 以约束最大似然估计为例，其原理就是通过选 择参数 δ 和T的估计值，应用 EM 算法、费舍尔得 分算法等，在迭代过程中使观测到这一实际数据 Y 的似然值最大并趋于稳定. 限于篇幅，在这里就不 再展开讨论，具体讨论可参见文献 [12-13]. 1.3 模型检验 可以通过单参数检验对固定效应参数 γ 进行 假设检验，以考察 γ 是否具有显著效应 [14] . 对于任 意 γqs，单参数检验的假设为 H0：γˆqs = 0，检验统 计量 t = ˆγqs/(Vγˆqs ) 1/2，服从自由度 J-N-1 的 T 分 布，其中 γˆqs 是 γqs 的估值，(Vγˆqs ) 1/2 是 γˆqs 的标准 误差，J 为第二层模型的样本数，N 为第二层模型 的自变量个数. 如假设成立，则表示 γqs 对应的因 子对预测因变量 β 无意义. 随机效应参数 β 的假设检验方式与固定效应 检验方式类似. 对于任意 βqj ，单参数检验的假设 为 H0：βˆ qj = 0 ，检验统计量 t = βˆ qj/(Vβˆqj ) 1/2，服 从自由度 J-N-1 的 T 分布，其中 βˆ qj 是 βqj 的估 值，(Vβˆqj ) 1/2 是 βˆ qj 的标准误差，J 为第一层模型 的样本数，N 为第一层模型的自变量个数. 如假设 成立，则表示 βqj 对应的因子与因变量 Y 无关系. 2 结果与讨论 2.1 研究样本和资料库 本文以国家材料腐蚀站网中碳钢在土壤中的 腐蚀数据进行分析和建模，对不同土壤环境因素下 碳钢平均腐蚀率进行预测. 按照多层模型分析的需 要，将原始数据整理为两层数据表，第一层数据如 表 1 所示，包含了实样站点、埋样时间和实验材料 的腐蚀率记录共 37 条数据；第二层数据如表 2 所 示，描述了各试样站的土壤环境数据，包括土壤 pH 值、有机质含量、全氮含量等理化参数，共 21 条数 据，涵盖了 21 个地区试验台站. 考虑到数据的可比 性，对表 2 数据进行了 Z-SCORE 标准化处理. 表 1 不同区域碳钢的腐蚀率数据 Table 1 Corrosion rate of carbon steel at different sites 站点 埋存时间/a 腐蚀率/(mm·a –1 ) 沈阳站 1 3.719 3 3.215 成都站 1 4.155 3 2.343 西安站 1 4.363 3 2.912 5 2.120 表 2 站点土壤理化性质数据示例 Table 2 Physical and chemical properties of soils at test sites 站点 pH 值 有机质质量分数/%全氮含量(质量分数)/% 可溶性盐离子的质量分数/% 土壤电阻/ Ω HCO− 3 SO2− 4 Mg+ K+ Na+ 大庆站 10.0500 1.0250 0.0525 0.0669 0.0142 0.0008 0.0033 0.0516 5.1000 华南站 5.4500 0.3250 0.0140 0.0031 0.0049 0.0006 0.0002 0.0007 460.0000 沈阳站 6.7500 2.3800 0.1090 0.0126 0.0118 0.0017 0.0004 0.0049 32.9000 2.2 建立模型 2.2.1 零模型 为考察不同地区土壤腐蚀率是否存在明显差 异，首先建立零模型如下. 第一层模型： Yti = β0i + rti. (9) 第二层模型： β0i = γ00 + µ0i . 参数估计结果如表 3 所示. 考察随机部分参 数估计结果可以看出，不同地区站点的腐蚀规律存 在着显著变异：方差 =5.563，χ 2 = 179.7, P < 0.001. 这说明如果为所有地区建立一个统一的模型，其参 数估计结果是有偏的，即利用多层线性模型进行建 模是必要的. 可以通过计算跨级相关系数度量变异 分布， ICC = Var(µ0i)/[Var(µ0i) + V ar(rti)] = 5.563/ (5.563 + 1.012) = 0.861, 即由不同地区间特性导致的变异占总变异的比例为 86.1%

第11期 鲁庆等：多层线性模型在碳钢土壤腐蚀规律中的应用研究 .1461· 表3碳钢腐蚀率零模型参数估计结果 Table 3 Parameter estimation results of the zero model of carbon steel corrosion rate 固定效应 系数 标准差 T值 P值 Y00 3.211313 0.531361 6.044 0.000 随机效应 标准误差 方差 x2 P值 Hoi 2.35862 5.56307 179.71340 0.000 ru 1.00615 1.01233 2.2.2基本模型 性因素，保持模型的稳定性，经过对比多种模型， 为了解碳钢材料在土壤中的腐蚀规律是否具有 本文采用了约束随机效应1：为0的策略.为确定 线性变化趋势，并且不同地区间其变化趋势是否存 约束是否恰当，对两种模型进行了比较.3o和31： 在差异，在零模型的基础上，建立基本模型如下. 均为随机变量时，模型的偏差度为139.93，自由度 第一层模型： 为4：相应模型（式(10）)的偏差度为142.18，自由 度为2.两个模型偏差度统计量差值为2.25，这个 Yu Boi Bui(ati)+rti (10) 差值服从自由度为2的X2分布，检验结果不显著 第二层模型： (P=0.325),即模型的简化是可行的，因此最终确定 30i=00+40i; 了目前第二层模型的形式. 不同地区土壤环境碳钢腐蚀率基本模型的参 31i=Y10 数估计如表4所示.固定效应的估计结果表明，腐 基本模型的第一层模型定义了碳钢在不同地 蚀率初始状态平均值为4.766(T=7.716),随着埋存 区腐蚀率随时间的变化趋势，第二层模型定义了不 时间的增加，腐蚀率有显著的线性下降趋势，下降 同地区腐蚀率变化参数截距和斜率的变化.在第 速度为0.134(010=10=0.134，T=-3.942).由于 二层模型，当定义：为随机变量时（即3：=10+ 截距和增长率的t检验值均显著，表明这两个参数 山1)，3o和31：的相关系数为0.925，为消除共线 在描述平均增长轨迹时是必要的， 表4碳钢腐蚀率基本模型参数估计结果 Table 4 Parameter estimation results of basic model of carbon steel corrosion rate 固定效应 系数 标准差 T值 P值 00 4.765820 0.617623 7.716 0.000 Y10 -0.134149 0.034030 -3.942 0.001 随机效应 标准误差 方差 x2 P值 Hoi 2.07276 4.29633 217.40681 0.000 rti 0.81737 0.66809 在目前的实验数据中，由于每个站点观测数据 行检验，删除不显著的因子，得到最终模型估计结 较少，不能支撑复杂模型的建模，因此考虑建立线 果如表5所示. 性模型.从上述固定效应的估计结果也可以看出， 首先考察固定效应参数估计结果.对于第一层 线性模型结构是可以接受的 腐蚀率模型的截距，可以看出，pH值（系数=-2.797， 进一步考察随机部分参数估计结果可以看出， T=-3.955,P=0.002)、全氮含量（系数=-5.881，T= 第一层模型截距3o:的随机误差o:的方差为 -3.875,P=0.002)、S02(系数=-21.61，T=-4.010, 4.296,X2值为217.4，P<0.001,即碳钢在不同 P=0.001)、电阻（系数=-1.612，T=-3.823,P=0.002) 地区的土壤腐蚀率初始状态存在着明显差异 四个变量越小，则第一层腐蚀率模型的截距越大， 2.2.3全模型 对于初始腐蚀率具有显著负效应：而Mg+(系 在基本模型的基础上，在第二层模型中逐步增 数=7.844976，T=2.793,P=0.015)和K+(系数 加与土壤腐蚀特性相关的pH值、有机物含量、全 =13.401239,T=3.638,P=0.003)则对于初始腐蚀率 氮含量、HCO3、SO、Mg+、K+、Na+、土壤电 具有显著正效应.在控制上述因素后，不同地区站 阻等理化参数，采用1.3节提出的方法对各变量进 点的初始腐蚀率的均值为3.705mma-1

第 11 期 鲁 庆等：多层线性模型在碳钢土壤腐蚀规律中的应用研究 1461 ·· 表 3 碳钢腐蚀率零模型参数估计结果 Table 3 Parameter estimation results of the zero model of carbon steel corrosion rate 固定效应 系数 标准差 T 值 P 值 γ00 3.211313 0.531361 6.044 0.000 随机效应 标准误差 方差 χ 2 P 值 µ0i 2.35862 5.56307 179.71340 0.000 rti 1.00615 1.01233 — — 2.2.2 基本模型 为了解碳钢材料在土壤中的腐蚀规律是否具有 线性变化趋势，并且不同地区间其变化趋势是否存 在差异，在零模型的基础上，建立基本模型如下. 第一层模型： Yti = β0i + β1i(αti) + rti. (10) 第二层模型： β0i = γ00 + µ0i ; β1i = γ10. 基本模型的第一层模型定义了碳钢在不同地 区腐蚀率随时间的变化趋势，第二层模型定义了不 同地区腐蚀率变化参数截距和斜率的变化. 在第 二层模型，当定义 β1i 为随机变量时(即β1i = γ10 + µ1i)，β0i 和 β1i 的相关系数为 0.925，为消除共线 性因素，保持模型的稳定性，经过对比多种模型， 本文采用了约束随机效应 γ1i 为 0 的策略. 为确定 约束是否恰当，对两种模型进行了比较. β0i 和 β1i 均为随机变量时，模型的偏差度为 139.93，自由度 为 4；相应模型 (式 (10)) 的偏差度为 142.18，自由 度为 2. 两个模型偏差度统计量差值为 2.25，这个 差值服从自由度为 2 的 χ 2 分布，检验结果不显著 (P=0.325)，即模型的简化是可行的，因此最终确定 了目前第二层模型的形式. 不同地区土壤环境碳钢腐蚀率基本模型的参 数估计如表 4 所示. 固定效应的估计结果表明，腐 蚀率初始状态平均值为 4.766 (T=7.716)，随着埋存 时间的增加，腐蚀率有显著的线性下降趋势，下降 速度为 0.134 (β10 = γ10= –0.134，T= –3.942). 由于 截距和增长率的 t 检验值均显著，表明这两个参数 在描述平均增长轨迹时是必要的. 表 4 碳钢腐蚀率基本模型参数估计结果 Table 4 Parameter estimation results of basic model of carbon steel corrosion rate 固定效应 系数 标准差 T 值 P 值 γ00 4.765820 0.617623 7.716 0.000 γ10 –0.134149 0.034030 –3.942 0.001 随机效应 标准误差 方差 χ 2 P 值 µ0i 2.07276 4.29633 217.40681 0.000 rti 0.81737 0.66809 — — 在目前的实验数据中，由于每个站点观测数据 较少，不能支撑复杂模型的建模，因此考虑建立线 性模型. 从上述固定效应的估计结果也可以看出， 线性模型结构是可以接受的. 进一步考察随机部分参数估计结果可以看出， 第一层模型截距 β0i 的随机误差 µ0i 的方差为 4.296，χ 2 值为 217.4，P < 0.001，即碳钢在不同 地区的土壤腐蚀率初始状态存在着明显差异. 2.2.3 全模型 在基本模型的基础上，在第二层模型中逐步增 加与土壤腐蚀特性相关的 pH 值、有机物含量、全 氮含量、HCO− 3 、SO2− 4 、Mg+、K+ 、Na+、土壤电 阻等理化参数，采用 1.3 节提出的方法对各变量进 行检验，删除不显著的因子，得到最终模型估计结 果如表 5 所示. 首先考察固定效应参数估计结果. 对于第一层 腐蚀率模型的截距，可以看出，pH 值 (系数=–2.797， T= –3.955，P=0.002)、全氮含量 (系数 = –5.881，T= –3.875，P=0.002)、SO2− 4 (系数 = –21.61，T= –4.010， P=0.001)、电阻 (系数= –1.612，T= –3.823，P=0.002) 四个变量越小，则第一层腐蚀率模型的截距越大， 对于初始腐蚀率具有显著负效应； 而 Mg+(系 数 =7.844976， T=2.793， P=0.015) 和 K+(系数 =13.401239，T=3.638，P=0.003) 则对于初始腐蚀率 具有显著正效应. 在控制上述因素后，不同地区站 点的初始腐蚀率的均值为 3.705 mm·a −1

.1462 北京科技大学学报 第35卷 表5碳钢腐蚀率全模型参数估计结果 Table 5 Parameter estimation results of the full model of carbon steel corrosion rate 因子 固定效应 系数 标准差 T值 P值 初始腐蚀率/(mma-1) Yoo 3.705291 0.386479 9.587 0.000 pH值 01 -2.797140 0.707255 -3.955 0.002 全氮含量/% Y02 -5.880776 1.517765 -3.875 0.002 S0含量/% Y03 -21.615691 5.390560 -4.010 0.001 Mg+含量/% 04 7.844976 2.808501 2.793 0.015 K+含量/% y05 13.401239 3.683997 3.638 0.003 电阻/ 06 -1.612091 0.421720 -3.823 0.002 腐蚀率变化速度 Y10 -0.697246 0.150994 -4.618 0.000 pH值 y11 0.173946 0.075356 2.308 0.031 全氮含量/% Y12 0.427550 0.147249 2.904 0.009 HC03含量/% Y13 0.427550 0.147249 2.904 0.009 S0含量/% Y14 2.845542 0.740702 3.842 0.001 Mg+含量/% 15 -1.048672 0.286515 -3.660 0.002 K+含量/% Y16 -0.804762 0.353150 -2.279 0.033 Na+含量/% Y17 -2.179878 0.557934 -3.907 0.001 随机效应 标准误差 方差 P值 Hoi 1.63887 2.68589 86.87952 0.000 Tui 0.73006 0.53298 对于第一层腐蚀率模型的斜率，可以看出，pH 31i=Y10+11(pH)+12(Nitrogeni)+13[(HC03-)]+ 值（系数=0.174，T=2.308,P=0.031)、全氮含 1a[(SO2)l+h5[(Mg+)月+Y16lKt)月+m[(Na+)l 量（系数=0.428，T=2.904,P=0.009)、HC03(系 数=0.428，T=2.904,P=0.009)和S02-(系数 式中，Nitrogen表示全氮含量，Resistance表示电 =2.846,T=3.842,P=0.001)四个变量越大，则 阻.建立全模型后，根据土壤理化参数等自变量数 第一层腐蚀率模型的斜率越大，即对于腐蚀率变化 据，并利用经验贝叶斯公式（式（⑧）)，就可以为任一 地区i建立独立的土壤腐蚀率模型.例如对于大庆 速度具有显著正效应：而Mg+(系数=-1.049，T= -3.660,P=0.002)、K+(系数=-0.805，T= 中心站，可计算得出，3oDQ=2.921,3DQ=-1.066, -2.279,P=0.033)和Na+(系数=-2.180，T= 其腐蚀率模型为YDQ=2.921-1.066a(tDQ) -3.907,P=0.001)则对于腐蚀率变化速度具有显 为了检验模型，本文计算了全部37组实验数 著负效应.在控制上述因素后，腐蚀率变化速度为 据的模型估值，并与实际测量数据对比，结果如图 -0.6972,仍具有显著的线性下降趋势(T=-4.618,P 1所示.其中横坐标为实验样本编号：纵坐标为腐 <0.001). 蚀率，单位为mm-1.方形图标表示利用本文求得 继续考察随机效应估计结果，不同地区的碳钢 的二层线性模型计算出的，不同地区各埋样时间的 腐蚀率仍存在着显著的差异（方差=2.68589，X2= 碳钢估计腐蚀率：圆圈图标表示实际测量的实验数 86.88,P<0.001),表明模型可能还存在进一步解 据.可以看出，模型估值与观测数据实际值非常接 释的可能.但是，对比基本模型（方差=4.29633）， 近，实证了依据多层线性模型原理建立的模型能够 全模型方差为2.68589有了显著的下降，下降比例 较好地描述腐蚀率的规律 w=(4.29633-2.68589)/4.29633×100%=37.48%,即实 2.3模型预测能力 验收集的士壤理化参数解释了初始情况的37.48% 为考察依据多层线性模型原理建立模型的预 方差. 测能力，本文采用将任一地区最后一次腐蚀率观测 根据估值结果，全模型表示为如下 数据删除，然后使用剩余数据建模，并检验模型预 第一层模型： 测值与该次观测实际值差异的方式.由于实验数据 Yui Boi+Bui(at)+rti. (11) 少，为保证模型稳定，在地区选择时，确定规则为 该地区至少有两次以上的观测数据，且每次只删一 第二层模型： 组数据.在本文使用的数据中，全部21个土壤腐蚀 Boi =Y00+701(pH;)+702(Nitrogen;)+703[(SO)]+ 试验台站共有13个台站满足条件.删除数据后模型 Yo4[(Mg+)]+705[(K+):]+Y06(Resistance;)+oi: 预测值与删除值的偏差结果如图2所示.尽管由于

· 1462 · 北 京 科 技 大 学 学 报 第 35 卷 表 5 碳钢腐蚀率全模型参数估计结果 Table 5 Parameter estimation results of the full model of carbon steel corrosion rate 因子 固定效应 系数 标准差 T 值 P 值 初始腐蚀率/ (mm·a−1 ) γ00 3.705291 0.386479 9.587 0.000 pH 值 γ01 –2.797140 0.707255 –3.955 0.002 全氮含量/% γ02 –5.880776 1.517765 –3.875 0.002 SO2− 4 含量/% γ03 –21.615691 5.390560 –4.010 0.001 Mg+ 含量/% γ04 7.844976 2.808501 2.793 0.015 K+ 含量/% γ05 13.401239 3.683997 3.638 0.003 电阻/Ω γ06 –1.612091 0.421720 –3.823 0.002 腐蚀率变化速度 γ10 –0.697246 0.150994 –4.618 0.000 pH 值 γ11 0.173946 0.075356 2.308 0.031 全氮含量/% γ12 0.427550 0.147249 2.904 0.009 HCO− 3 含量/% γ13 0.427550 0.147249 2.904 0.009 SO2− 4 含量/ % γ14 2.845542 0.740702 3.842 0.001 Mg+ 含量/% γ15 –1.048672 0.286515 –3.660 0.002 K+ 含量/% γ16 –0.804762 0.353150 –2.279 0.033 Na+ 含量/% γ17 –2.179878 0.557934 –3.907 0.001 随机效应 标准误差 方差 χ 2 P 值 µ0i 1.63887 2.68589 86.87952 0.000 rti 0.73006 0.53298 — — 对于第一层腐蚀率模型的斜率，可以看出，pH 值 (系数 = 0.174, T = 2.308, P = 0.031)、全氮含 量 (系数 = 0.428, T = 2.904, P = 0.009)、HCO− 3 (系 数 = 0.428, T = 2.904,P = 0.009) 和 SO2− 4 (系数 = 2.846,T = 3.842,P = 0.001) 四个变量越大，则 第一层腐蚀率模型的斜率越大，即对于腐蚀率变化 速度具有显著正效应；而 Mg+(系数 = −1.049,T = −3.660,P = 0.002)、 K+(系数 = −0.805,T = −2.279, P = 0.033) 和 Na+(系数 = −2.180，T = −3.907,P = 0.001) 则对于腐蚀率变化速度具有显 著负效应. 在控制上述因素后，腐蚀率变化速度为 –0.6972，仍具有显著的线性下降趋势 (T = −4.618，P <0.001). 继续考察随机效应估计结果，不同地区的碳钢 腐蚀率仍存在着显著的差异 (方差 =2.68589，χ 2= 86.88，P <0.001)，表明模型可能还存在进一步解 释的可能. 但是，对比基本模型 (方差 =4.29633)， 全模型方差为2.68589 有了显著的下降，下降比例 ω=(4.29633–2.68589)/4.29633×100%=37.48%，即实 验收集的土壤理化参数解释了初始情况的 37.48% 方差. 根据估值结果，全模型表示为如下. 第一层模型： Yti = β0i + β1i(αti) + rti. (11) 第二层模型： β0i = γ00+γ01(pHi )+γ02(Nitrogeni )+γ03[(SO2− 4 )i ]+ γ04[(Mg+)i ] + γ05[(K+)i ] + γ06(Resistancei) + µ0i ; β1i = γ10+γ11(pHi )+γ12(Nitrogeni)+γ13[(HCO3 −)i ]+ γ14[(SO2− 4 )i ]+γ15[(Mg+)i ]+γ16[(K+ )i ]+γ17[(Na+)i ]. 式中，Nitrogen 表示全氮含量，Resistance 表示电 阻. 建立全模型后，根据土壤理化参数等自变量数 据，并利用经验贝叶斯公式 (式 (8))，就可以为任一 地区 i 建立独立的土壤腐蚀率模型. 例如对于大庆 中心站，可计算得出，β0DQ =2.921，β1DQ =–1.066， 其腐蚀率模型为 YtDQ = 2.921 − 1.066a(tDQ). 为了检验模型，本文计算了全部 37 组实验数 据的模型估值，并与实际测量数据对比，结果如图 1 所示. 其中横坐标为实验样本编号；纵坐标为腐 蚀率，单位为 mm· −1 . 方形图标表示利用本文求得 的二层线性模型计算出的，不同地区各埋样时间的 碳钢估计腐蚀率；圆圈图标表示实际测量的实验数 据. 可以看出，模型估值与观测数据实际值非常接 近，实证了依据多层线性模型原理建立的模型能够 较好地描述腐蚀率的规律. 2.3 模型预测能力 为考察依据多层线性模型原理建立模型的预 测能力，本文采用将任一地区最后一次腐蚀率观测 数据删除，然后使用剩余数据建模，并检验模型预 测值与该次观测实际值差异的方式. 由于实验数据 少，为保证模型稳定，在地区选择时，确定规则为 该地区至少有两次以上的观测数据，且每次只删一 组数据. 在本文使用的数据中，全部 21 个土壤腐蚀 试验台站共有 13 个台站满足条件. 删除数据后模型 预测值与删除值的偏差结果如图 2 所示. 尽管由于