留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

附不等式约束的地表沉降时间序列自回归EIV模型

陶叶青 王坚 刘超

陶叶青, 王坚, 刘超. 附不等式约束的地表沉降时间序列自回归EIV模型[J]. 武汉大学学报 ● 信息科学版, 2020, 45(9): 1455-1460. doi: 10.13203/j.whugis20180268
引用本文: 陶叶青, 王坚, 刘超. 附不等式约束的地表沉降时间序列自回归EIV模型[J]. 武汉大学学报 ● 信息科学版, 2020, 45(9): 1455-1460. doi: 10.13203/j.whugis20180268
TAO Yeqing, WANG Jian, LIU Chao. A Solution for Ground Subsidence Prediction of Time Series Based on Autoregression EIV Model with Inequality Constraints[J]. Geomatics and Information Science of Wuhan University, 2020, 45(9): 1455-1460. doi: 10.13203/j.whugis20180268
Citation: TAO Yeqing, WANG Jian, LIU Chao. A Solution for Ground Subsidence Prediction of Time Series Based on Autoregression EIV Model with Inequality Constraints[J]. Geomatics and Information Science of Wuhan University, 2020, 45(9): 1455-1460. doi: 10.13203/j.whugis20180268

附不等式约束的地表沉降时间序列自回归EIV模型

doi: 10.13203/j.whugis20180268
基金项目: 

国家自然科学基金 41601501

国家自然科学基金 41874029

国家自然科学基金 41404004

江苏省高校自然科学基金 16KJD420001

江苏省建设厅科技规划项目 2017ZD259

详细信息
    作者简介:

    陶叶青, 博士, 副教授, 主要从事测量平差理论与方法研究。yenneytao@163.com

  • 中图分类号: P207

A Solution for Ground Subsidence Prediction of Time Series Based on Autoregression EIV Model with Inequality Constraints

Funds: 

The National Natural Science Foundation of China 41601501

The National Natural Science Foundation of China 41874029

The National Natural Science Foundation of China 41404004

the Natural Science Foundation for Colleges and Universities of Jiangsu Province 16KJD420001

Science Foundation of Jiangsu Provincial Department of Construction 2017ZD259

More Information
    Author Bio:

    TAO Yeqing, PhD, associate professor, specializes in the theories and methods of surveying adjustment. E-mail: yenneytao@163.com

  • 摘要: 地表沉降监测是预防地质灾害发生的有效方法, 时间序列模型是进行地表沉降预测的主要模型。传统AR(autoregression)模型参数估计算法没有顾及模型系数矩阵中的元素含有观测误差的情况。为同时顾及模型观测向量与系数矩阵中的元素均含有观测误差的情况, 应用变量中含有误差(errors-in-variables, EIV)的参数估计模型改进基于时间序列的地表沉降预测模型。为提高AR模型参数解的收敛速度与精度, 应用先验信息构建具有不等式约束的EIV模型, 将建立的附有不等式约束EIV模型参数估计问题转化为非线性模型的二次规划问题, 结合中位函数建立参数估计迭代算法。为论证所建立算法的有效性与可行性, 通过实验分别对AR模型参数估计的最小二乘算法、EIV模型参数估计算法与该算法进行比较。实验证明, 所建立算法具有较高的精度和效率, 是一种可行的方法。
  • 表  1  算例1应用的模拟数据

    Table  1.   Simulation Data of Example 1

    序号 模拟值
    1 15.9189
    2 14.7974
    3 23.2462
    4 7.3572
    5 -10.8052
    6 -18.0692
    7 -10.5677
    8 2.7923
    9 14.6792
    10 13.1268
    11 1.1775
    12 -10.7078
    13 -13.8582
    14 -6.0010
    15 5.4601
    16 9.6419
    17 9.1258
    18 -1.9343
    19 -6.5519
    20 -9.4017
    21 -2.8574
    22 4.8899
    23 7.6215
    24 5.0390
    25 0.6078
    26 -6.3572
    27 -5.7925
    28 -1.3012
    29 2.1084
    30 5.8308
    下载: 导出CSV

    表  2  4种算法得到的参数估值与其精度

    Table  2.   Model Parameters Estimation and Its Accuracy Based on the Four Methods

    算法 a1 a2 a3 精度
    LS -0.3246 -0.4667 0.7709 1.4160
    SVD -0.2725 -0.5461 0.8363 1.4455
    TLS-V -0.3068 -0.5084 0.8043 0.9042
    TLS-M -0.2903 -0.5077 0.8136 0.9203
    下载: 导出CSV

    表  3  监测点形变观测数据

    Table  3.   Observation Data for Settlement of Monitoring Point

    周期 高差/mm
    1 -1.5
    2 9.6
    3 -2.4
    4 -1
    5 -5.6
    6 -1.5
    7 -8.5
    15 -3.6
    16 -2.9
    17 -10.4
    18 -14.3
    19 -3.1
    20 -9.4
    21 -0.6
    8 -3.3
    9 -4.1
    10 -5.3
    11 -1.2
    12 -4.4
    13 -1.6
    14 2
    22 -4.6
    23 -2.3
    24 -3.1
    25 -2.3
    26 -3.9
    27 -2.6
    28 -2.7
    下载: 导出CSV

    表  4  3种算法求解的模型参数估值与检核精度

    Table  4.   Model Parameters Estimation and Accuracy by the Three Different Algorithms

    算法 a1 a2 a3 RMSE/mm
    LS 0.0277 0.3714 0.4470 0.8235
    ICLS 0.3234 0.4131 0.2401 0.7773
    ICTLS 0.2218 0.3747 0.4147 0.7173
    下载: 导出CSV
  • [1] 李达, 邓喀中, 高晓雄, 等.基于SBAS-InSAR的矿区地表沉降监测与分析[J].武汉大学学报·信息科学版, 2018, 43(10): 1 531-1 537 doi:  10.13203/j.whugis20160566

    Li Da, Deng Kazhong, Gao Xiaoxiong, et al. Monitoring and Analysis of Surface Subsidence in Mining Area Based on SBAS-InSAR[J].Geomatics and Information Science of Wuhan University, 2018, 43(10): 1 531-1 537 doi:  10.13203/j.whugis20160566
    [2] Zahmatkesh A, Choobbasti A J. Evaluation of Wall Deflections and Ground Surface Settlements in Deep Excavations[J]. Arabian Journal of Geosciences, 2015, 8 (5): 3 055-3 063 doi:  10.1007/s12517-014-1419-6
    [3] Moeinossadat S R, Ahangari K, Shahriar K. Control of Ground Settlements Caused by EPBS Tunneling Using an Intelligent Predictive Model[J]. Indian Geotechnical Journal, 2017(5):1-10 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=568fc128a04881fa01eb87670e610051
    [4] 向巍, 郭际明, 傅露.基于垂直距离最小二乘拟合的双曲线沉降模型[J].武汉大学学报·信息科学版, 2013, 38(5):571-574 http://ch.whu.edu.cn/article/id/2645

    Xiang Wei, Guo Jiming, Fu Lu. Hyperbolic Settlement Model Based on Least-Squares Orthogonal Distances Fitting[J]. Geomatics and Information Science of Wuhan University, 2013, 38(5):571-574 http://ch.whu.edu.cn/article/id/2645
    [5] 王建民, 张锦.基于高斯过程回归的变形智能预测模型及应用[J].武汉大学学报·信息科学版, 2018, 43(2):248-254 doi:  10.13203/j.whugis20160075

    Wang Jianmin, Zhang Jin. Deformation Intelligent Prediction Model Based on Gaussian Process Regression and Application[J]. Geomatics and Information Science of Wuhan University, 2018, 43(2):248-254 doi:  10.13203/j.whugis20160075
    [6] 段光耀, 刘欢欢, 宫辉力, 等.京津城际铁路沿线不均匀地面沉降演化特征[J].武汉大学学报·信息科学版, 2017, 42(12): 1 847-1 853 doi:  10.13203/j.whugis20150537

    Duan Guangyao, Liu Huanhuan, Gong Huili, et al. Evolution Characteristics of Uneven Land Subsidence Along Beijing-Tianjin Inter-city Railway[J]. Geomatics and Information Science of Wuhan University, 2017, 42(12): 1 847-1 853 doi:  10.13203/j.whugis20150537
    [7] Amiri-Simkooei A R. Application of Least Squares Variance Component Estimation to Errors-in-Variables Models[J]. Journal of Geodesy, 2013, 87(12):935-944 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=6125d49648f95a639a2e805c2d7c295d
    [8] Zhang S L, Tong X H, Zhang K L, et al. A Solution to EIV Model with Inequality Constraints and Its Geodetic Applications[J]. Journal of Geodesy, 2013, 87(1): 23-28 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=786335433327bf83c48572ded36e0eb0
    [9] Fang Xing, Li Bofeng, Alkhatib H, et al.Bayesian Inference for the Errors-in-Variables Model[J].Studia Geophysica et Geodaetica, 2016, 61(1): 1-18 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=102db1aba60e81971988fd72bf7c253a
    [10] 王乐洋, 李海燕, 温扬茂, 等.地震同震滑动分布反演的总体最小二乘方法[J].测绘学报, 2017, 46(3):307-315 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb201703006

    Wang Leyang, Li Haiyan, Wen Yangmao, et al. Total Least Squares Method Inversion for Coseismic Slip Distribution[J]. Acta Geodaetica et Cartographica Sinica, 2017, 46(3):307-315 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb201703006
    [11] Schaffrin B, Wieser A. On Weighted Total Least-Squares Adjustment for Linear Regression[J]. Journal of Geodesy, 2008, 82(7):415-421 doi:  10.1007/s00190-007-0190-9
    [12] 姚宜斌, 熊朝晖, 张豹, 等.顾及设计矩阵误差的AR模型新解法[J].测绘学报, 2017, 46(11):1 795-1 801 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb201711001

    Yao Yibin, Xiong Zhaohui, Zhang Bao, et al. A New Method to Solving AR Model Parameters Considering Random Errors of Design Matrix[J]. Acta Geodaetica et Cartographica Sinica, 2017, 46(11): 1 795-1 801 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb201711001
    [13] Rousseeuw P, Wagner J. Robust Regression with a Distribution Intercept Using Least Median of Squares[J]. Computational Statistics & Data Analysis, 1994, 17:65-76
    [14] Fang Xing. On Non-combinatorial Weighted Total Least Squares with Inequality Constraints[J]. Journal of Geodesy, 2014, 88(8):805-816 doi:  10.1007/s00190-014-0723-y
    [15] Schaffrin B. A Note on Constrained Total Least-Squares Estimation[J]. Linear Algebra & Its Applications, 2006, 417(1):245-258 http://cn.bing.com/academic/profile?id=ba532e7ef9374a84fe38b7b5fe362c35&encoded=0&v=paper_preview&mkt=zh-cn
    [16] 许文源, 王东谦.带有线性不等式约束的最小二乘[J].系统科学与数学, 1984, 4 (1):55-62 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=QK000001694665

    Xu Wenyuan, Wang Dongqian. Least Square Estimation with Linear Inequality Constraints[J]. Journal of Systems Science and Mathematical Sciences, 1984, 4 (1):55-62 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=QK000001694665
  • 加载中
计量
  • 文章访问数:  62
  • HTML全文浏览量:  27
  • PDF下载量:  18
  • 被引次数: 0
出版历程
  • 收稿日期:  2018-07-13
  • 刊出日期:  2020-09-05

附不等式约束的地表沉降时间序列自回归EIV模型

doi: 10.13203/j.whugis20180268
    基金项目:

    国家自然科学基金 41601501

    国家自然科学基金 41874029

    国家自然科学基金 41404004

    江苏省高校自然科学基金 16KJD420001

    江苏省建设厅科技规划项目 2017ZD259

    作者简介:

    陶叶青, 博士, 副教授, 主要从事测量平差理论与方法研究。yenneytao@163.com

  • 中图分类号: P207

摘要: 地表沉降监测是预防地质灾害发生的有效方法, 时间序列模型是进行地表沉降预测的主要模型。传统AR(autoregression)模型参数估计算法没有顾及模型系数矩阵中的元素含有观测误差的情况。为同时顾及模型观测向量与系数矩阵中的元素均含有观测误差的情况, 应用变量中含有误差(errors-in-variables, EIV)的参数估计模型改进基于时间序列的地表沉降预测模型。为提高AR模型参数解的收敛速度与精度, 应用先验信息构建具有不等式约束的EIV模型, 将建立的附有不等式约束EIV模型参数估计问题转化为非线性模型的二次规划问题, 结合中位函数建立参数估计迭代算法。为论证所建立算法的有效性与可行性, 通过实验分别对AR模型参数估计的最小二乘算法、EIV模型参数估计算法与该算法进行比较。实验证明, 所建立算法具有较高的精度和效率, 是一种可行的方法。

English Abstract

陶叶青, 王坚, 刘超. 附不等式约束的地表沉降时间序列自回归EIV模型[J]. 武汉大学学报 ● 信息科学版, 2020, 45(9): 1455-1460. doi: 10.13203/j.whugis20180268
引用本文: 陶叶青, 王坚, 刘超. 附不等式约束的地表沉降时间序列自回归EIV模型[J]. 武汉大学学报 ● 信息科学版, 2020, 45(9): 1455-1460. doi: 10.13203/j.whugis20180268
TAO Yeqing, WANG Jian, LIU Chao. A Solution for Ground Subsidence Prediction of Time Series Based on Autoregression EIV Model with Inequality Constraints[J]. Geomatics and Information Science of Wuhan University, 2020, 45(9): 1455-1460. doi: 10.13203/j.whugis20180268
Citation: TAO Yeqing, WANG Jian, LIU Chao. A Solution for Ground Subsidence Prediction of Time Series Based on Autoregression EIV Model with Inequality Constraints[J]. Geomatics and Information Science of Wuhan University, 2020, 45(9): 1455-1460. doi: 10.13203/j.whugis20180268
  • 地表沉降主要是由于自然和人为因素的作用, 使得地下松散土层及岩层受到压缩而导致地面高程逐渐降低。地表沉降可能会发展成为一种严重的地质灾害, 对道路、房屋建筑、地下管线等造成破坏, 并可能引发次生地质灾害, 受到研究人员的关注[1-3]。地表沉降是随时间演变的地质演化过程, 应用时间序列的回归模型能够描述其演化特征, AR(autoregression)模型是一种广泛应用的时间序列模型。

    应用观测技术对地表进行监测是预测地表沉降发生, 避免引发重大地质灾害的有效方法。雷达干涉测量技术、卫星导航定位技术、三维激光扫描技术等应用于沉降监测的同时, 结合观测数据, 研究人员可应用数学模型对形变进行推估[4-6]。应用数学模型对地表沉降进行推估, 其精度受到两个因素的影响:一是模型的符合程度; 二是模型参数的求解精度。数学模型参数传统解算方法是应用高斯-马尔可夫(Gauss-Markov model, G-M)模型在最小二乘准则(least squares, LS)下求解参数。建立G-M模型的前提是仅模型中的观测向量中的元素含有误差, 描述观测向量与模型参数函数关系的系数矩阵中的元素不含有误差, 这并不符合AR模型的实际特征。AR模型的系数矩阵与观测向量中的元素均由观测元素组成, 因此, 其系数矩阵与观测向量中均含有观测误差, 这一类模型称为变量中含有误差(errors-in-variables EIV)参数估计模型。EIV模型是一种非线性参数估计模型, 其模型参数的估计算法受到关注的同时, EIV模型在地学专业的应用也被关注[7-10]

    地表沉降满足一定的形变特征, 在对其监测数据进行处理时, 具有有效的先验信息可以使用, 根据先验信息建立具有不等式约束的参数估计模型能够提高模型参数解的精度, 而现有成果缺乏对具有约束条件的地表沉降时间序列模型的研究。根据建立的地表沉降时间序列不等式约束EIV模型, 本文将不等式约束问题转化为非线性最优化问题。同时, 为使得在模型中重复出现的同一元素得到相同的改正数, 应用中位函数计算其改正数。根据二次规划最优解的KKT条件(Karush-Kuhn-Tucher conditions)与中位函数建立模型参数的嵌套迭代算法, 并通过算例与现有算法进行比较, 验证该算法的可行性。

    • 广泛使用的时间序列自回归AR模型为:

      $$\left\{ {\begin{array}{*{20}{c}} {{h_{m + 1}} = {a_1}{h_1} + {a_2}{h_2} + \cdots + {a_m}{h_m}{\rm{}}}\\ {{h_{m + 2}} = {a_1}{h_2} + {a_2}{h_3} + \cdots + {a_m}{h_{m + 1}}{\rm{}}}\\ \vdots \\ {{h_{m + n}} = {a_1}{h_n} + {a_2}{h_{n + 1}} + \cdots + {a_m}{h_{m + n - 1}}{\rm{}}} \end{array}} \right.$$ (1)

      式中, ${a_1}, {a_2} \cdots {a_m}$为m个模型参数; ${h_1}, {h_2} \cdots {h_{m + n}}$为m+n期观测值。

      AR模型表明, 同一观测元素在模型中多次出现。设观测值中含有的观测误差为$\left( {{v_1}{\rm{}}{v_2} \cdots {v_{m + n}}} \right)$, 将式(1)表示成G-M模型:

      $$\left\{ {\begin{array}{*{20}{c}} {\mathit{\boldsymbol{v}} = \mathit{\boldsymbol{Bx}} - \mathit{\boldsymbol{l}}}\\ {\mathit{\boldsymbol{v}} \sim \left( {0, \sigma _0^2{\mathit{\boldsymbol{Q}}_l}} \right)} \end{array}} \right.$$ (2)

      式中, v为观测向量l中含有的误差向量; B为系数矩阵; x为模型参数组成的参数向量; $\sigma _0^2$为单位权方差; ${\mathit{\boldsymbol{Q}}_l}$为协因数阵。

      vxlB的表达式如下:

      $$\mathit{\boldsymbol{v}} = {\left[ {\begin{array}{*{20}{c}} {{v_{m + 1}}}&{{v_{m + 2}}}& \cdots &{{v_{m + n}}} \end{array}} \right]^{\rm{T}}}$$
      $$\mathit{\boldsymbol{x}} = {\left[ {\begin{array}{*{20}{c}} {{a_1}}&{{a_2}}& \cdots &{{a_m}} \end{array}} \right]^{\rm{T}}}$$
      $$\mathit{\boldsymbol{l}} = {\left[ {\begin{array}{*{20}{c}} {{h_{m + 1}}}&{{h_{m + 2}}}& \cdots &{{h_{m + n}}} \end{array}} \right]^{\rm{T}}}$$
      $$\mathit{\boldsymbol{B}} = \left[ {\begin{array}{*{20}{c}} {{h_1}}&{{h_2}}& \cdots &{{h_m}}\\ {{h_2}}&{{h_3}}& \cdots &{{h_{m + 1}}}\\ \vdots & \vdots & \vdots & \vdots \\ {{h_n}}&{{h_{n + 1}}}& \cdots &{{h_{m + n - 1}}} \end{array}} \right]$$

      应用LS得到模型参数的估值$\mathit{\boldsymbol{\hat x}}$为:

      $$\left\{ {\begin{array}{*{20}{c}} {\mathit{\boldsymbol{vQ}}_l^{ - 1}\mathit{\boldsymbol{v}} = {\rm{min}}}\\ {\mathit{\boldsymbol{\hat x}} = {{(\mathit{\boldsymbol{BQ}}_l^{ - 1}\mathit{\boldsymbol{B}})}^{ - 1}}\mathit{\boldsymbol{BQ}}_l^{ - 1}\mathit{\boldsymbol{l}}} \end{array}} \right.$$ (3)

      考虑到时间序列的自回归G-M模型(式(2))中, 不仅模型的观测向量l中的元素由观测元素组成; 系数矩阵B中的元素也由观测元素组成, 其不可避免也会含有观测误差。此时, 建立G-M模型的前提假设并不成立, 因此, 由式(3)得到的参数估值不再具有最优、无偏的性质。为同时顾及系数矩阵与观测向量中的元素含有观测误差对参数求解的影响, 将G-M模型拓展为EIV模型:

      $$\left\{ {\begin{array}{*{20}{c}} {\mathit{\boldsymbol{v}} = \left( {\mathit{\boldsymbol{B}} + {\mathit{\boldsymbol{E}}_B}} \right)\mathit{\boldsymbol{x}} - \mathit{\boldsymbol{l}}}\\ {\left[ {\begin{array}{*{20}{c}} \mathit{\boldsymbol{v}}\\ \mathit{\boldsymbol{b}} \end{array}} \right] \sim \left( {\left[ {\begin{array}{*{20}{c}} 0\\ 0 \end{array}} \right], \sigma _0^2\left[ {\begin{array}{*{20}{c}} {{\mathit{\boldsymbol{Q}}_l}}&0\\ 0&{{\mathit{\boldsymbol{Q}}_B}} \end{array}} \right]} \right)} \end{array}} \right.$$ (4)

      式中, ${\mathit{\boldsymbol{E}}_B}$为系数矩阵B中含有的n×m阶误差矩阵; $\mathit{\boldsymbol{b}} = {\mathop{\rm vec}\nolimits} \left( {{\mathit{\boldsymbol{E}}_B}} \right)({\mathop{\rm vec}\nolimits} ( \cdot )$为矩阵的拉直符号, 表示将矩阵按列拉直得到列向量); QB为误差向量b对应的协因数阵。

      同时顾及模型系数矩阵与观测向量中观测误差的参数估计算法称为总体最小二乘算法(total least squares, TLS)。EIV模型是一种非线性参数估计模型, 目前广泛应用的参数估计算法是基于Lagrange函数的极值函数算法[11]

      需要指出的是, 区别于一般EIV模型, 时间序列的自回归EIV模型的系数矩阵与观测向量中多次出现同一观测值, 系数矩阵与观测向量中同一元素应具有相同的改正数, 但是传统算法无法保障。

    • 根据地表沉降的一般性规律, 同一观测点后期观测得到的高程值应小于前期观测得到的高程值。结合时间序列的自回归模型(式(1)), 建立不等式约束方程为:

      $$\left( {\begin{array}{*{20}{c}} {{h_2} - {h_1}{\rm{}}}&{{h_3} - {h_2}{\rm{}}}& \cdots &{{h_{m + 1}} - {h_m}{\rm{}}}\\ {{h_3} - {h_2}{\rm{}}}&{{h_4} - {h_3}{\rm{}}}& \cdots &{{h_{m + 2}} - {h_{m + 1}}{\rm{}}}\\ \vdots & \vdots &{}& \vdots \\ {{h_n} - {h_{n - 1}}}&{{h_{n + 1}} - {h_n}}& \cdots &{{h_{m + n - 1}} - {h_{m + n - 2}}} \end{array}} \right)\left( {\begin{array}{*{20}{c}} {{a_1}}\\ {{a_2}}\\ \vdots \\ {{a_m}} \end{array}} \right) < \left( {\begin{array}{*{20}{c}} 0\\ 0\\ \vdots \\ 0 \end{array}} \right)$$ (5)

      式(5)中, 定义G为约束方程的系数矩阵, w为约束方程的观测向量。与式(4)组成附有不等式约束的EIV函数模型:

      $$\left\{ {\begin{array}{*{20}{l}} {\mathit{\boldsymbol{v}} = \left( {\mathit{\boldsymbol{B}} + {\mathit{\boldsymbol{E}}_B}} \right)\mathit{\boldsymbol{x}} - \mathit{\boldsymbol{l}}}\\ {\mathit{\boldsymbol{Gx}} < \mathit{\boldsymbol{w}}} \end{array}} \right.$$ (6)

      根据AR模型(式(1)), 观测向量中的观测元素($\begin{array}{*{20}{c}} {{h_{m + 1}}}&{{h_{m + 2}}}& \cdots &{{h_{m + n}}} \end{array}$)在模型中出现的次数为(n, n-1…1), 模型中观测元素($h_{1}, h_{2} \cdots h_{m}$)在模型中出现的次数为(1, 2…m)次。在AR模型中, 同一观测元素在模型中多次出现, 在进行数据处理时, 同一观测元素应该得到相同的改正数。为达到这一目的, 文献[12]应用虚拟观测值方法建立的数据处理模型理论严密, 但是需要重新改写原模型, 算法较为复杂, 不便于程序化设计。本文根据同一观测元素的多个改正值, 应用中位函数确定其改正数, 并据此建立附不等式约束EIV模型的迭代算法。

      设观测元素($h_{1}, h_{2} \cdots h_{m+n}$)的改正数为($v_{1}, v_{2}\cdots v_{m+n}$), 由同一观测元素的改正值组成的向量为:

      $$\left\{ {\begin{array}{*{20}{l}} {{\mathit{\boldsymbol{ \boldsymbol{\varDelta} }}_1} = \left[ {{\mathit{\Delta }_1}} \right]}\\ {{\rm{}} \vdots }\\ {{\mathit{\boldsymbol{ \boldsymbol{\varDelta} }}_m} = \left[ {\begin{array}{*{20}{c}} {\mathit{\Delta }_m^1, {\rm{}}\mathit{\Delta }_m^2}& \cdots &{\mathit{\Delta }_m^m} \end{array}} \right]}\\ {{\mathit{\boldsymbol{ \boldsymbol{\varDelta} }}_{m + 1}} = \left[ {\begin{array}{*{20}{c}} {\mathit{\Delta }_{m + 1}^1, {\rm{}}\mathit{\Delta }_{m + 1}^2}& \cdots &{\mathit{\Delta }_{m + 1}^n} \end{array}} \right]}\\ {{\rm{}} \vdots }\\ {{\mathit{\boldsymbol{ \boldsymbol{\varDelta} }}_{m + n}} = \left[ {{\mathit{\Delta }_{m + n}}} \right]} \end{array}} \right.$$

      据此, 应用中位函数(med[·])得到观测元素的改正数为:

      $$\left\{ {\begin{array}{*{20}{l}} {{v_1} = {\rm{med}}\left[ {{\mathit{\Delta }_1}} \right]}\\ {{\rm{}} \vdots }\\ {{v_m} = {\rm{med}}\left[ {\begin{array}{*{20}{c}} {\mathit{\Delta }_m^1, {\rm{}}\mathit{\Delta }_m^2}& \cdots &{\mathit{\Delta }_m^m} \end{array}} \right]}\\ {{v_{m + 1}} = {\rm{med}}\left[ {\begin{array}{*{20}{c}} {\mathit{\Delta }_{m + 1}^1, \mathit{\Delta }_{m + 1}^2}& \cdots &{\mathit{\Delta }_{m + 1}^n} \end{array}} \right]}\\ {{\rm{}} \vdots }\\ {{v_{m + n}} = {\rm{med}}\left[ {{\mathit{\Delta }_{m + n}}} \right]} \end{array}} \right.$$

      应用中位函数使得模型中同一观测元素得到相同的改正数。同时, 中位函数可以降低含有粗差的观测值对参数估值的影响, 相关研究表明, 其崩溃污染率可以达到50%[13]

      文献[8]应用穷举法建立附有不等式约束的EIV模型, 使得此类模型的参数估计算法与应用被广泛关注。穷举法是鉴别不等式约束方程中有效约束与无效约束的可行方法, 适用于算法的程序化设计, 但是其也具有计算量随约束方程数量的增多呈指数级增长的缺陷。

      根据附不等式约束的EIV模型(式(6)), 将其转化为不等式约束的二次规划问题:

      $$ \min \boldsymbol{e}^{\mathrm{T}} \boldsymbol{Q}^{-1} \boldsymbol{e}, \text { 使得 } \boldsymbol{G x} \leqslant \boldsymbol{w} $$ (7)

      其中, $\mathit{\boldsymbol{e}} = \left[ {\begin{array}{*{20}{l}} \mathit{\boldsymbol{v}}\\ \mathit{\boldsymbol{b}} \end{array}} \right];\mathit{\boldsymbol{Q}} = \left[ {\begin{array}{*{20}{l}} {{\mathit{\boldsymbol{Q}}_l}}&0\\ 0&{{\mathit{\boldsymbol{Q}}_B}} \end{array}} \right]$。

      为建立不等式约束的二次规划KKT条件, 应用拉格朗日函数构建附不等式约束模型的目标方程如下:

      $$\begin{array}{*{20}{r}} {f\left( {\mathit{\boldsymbol{e}}, \mathit{\boldsymbol{x}}, \mathit{\boldsymbol{\mu }}, \mathit{\boldsymbol{\lambda }}} \right) = {\mathit{\boldsymbol{e}}^{\rm{T}}}{\mathit{\boldsymbol{Q}}^{ - 1}}\mathit{\boldsymbol{e}} + }\\ {2{\mathit{\boldsymbol{\mu }}^{\rm{T}}}\left( {\mathit{\boldsymbol{l}} - \mathit{\boldsymbol{Bx}} - {\mathit{\boldsymbol{E}}_B}\mathit{\boldsymbol{x}} + \mathit{\boldsymbol{v}}} \right) + }\\ {2{\mathit{\boldsymbol{\lambda }}^{\rm{T}}}\left( {\mathit{\boldsymbol{Gx}} - \mathit{\boldsymbol{w}} + \mathit{\boldsymbol{\eta }}} \right)} \end{array}$$ (8)

      式中, μλ分别为1×n维与1×(n-1)维拉格朗日乘数; η为非负松驰变量。对目标方程中的各变量求偏导数, 并令其等于零, 得到目标方程的法方程为[14]:

      $${(\mathit{\boldsymbol{B}} + {\mathit{\boldsymbol{E}}_B})^{\rm{T}}}{(\mathit{\boldsymbol{AQ}}{\mathit{\boldsymbol{A}}^{\rm{T}}})^{ - 1}}\left( {\mathit{\boldsymbol{l}} - \mathit{\boldsymbol{Bx}}} \right) - {\mathit{\boldsymbol{G}}^{\rm{T}}}\mathit{\boldsymbol{\lambda }} = 0$$ (9)

      式(9)中, $\mathit{\boldsymbol{A}} = \left[ {{{\mathit{\boldsymbol{\hat x}}}^{\rm{T}}} \otimes {\mathit{\boldsymbol{I}}_n}, {\rm{}}{\mathit{\boldsymbol{I}}_n}} \right]$。

      满足附不等式约束EIV模型参数最优解的KKT条件为:

      $$\left\{ {\begin{array}{*{20}{l}} {{{(\mathit{\boldsymbol{B}} + {\mathit{\boldsymbol{E}}_B})}^{\rm{T}}}{{({\bf{AQ}}{\mathit{\boldsymbol{A}}^{\rm{T}}})}^{ - 1}}\left( {\mathit{\boldsymbol{l}} - {\bf{Bx}}} \right) - {\mathit{\boldsymbol{G}}^{\rm{T}}}\mathit{\boldsymbol{\lambda }} = 0}\\ {\begin{array}{*{20}{l}} {{\mathit{\boldsymbol{G}}_i}\mathit{\boldsymbol{x}} - {w_i} \le 0}\\ {{\lambda _i}\left( {{\mathit{\boldsymbol{G}}_i}\mathit{\boldsymbol{x}} - {w_i}} \right) = 0, {\lambda _i} \ge 0} \end{array}} \end{array}} \right.$$ (10)

      进一步得到:

      $${(\mathit{\boldsymbol{B}} + {\mathit{\boldsymbol{E}}_B})^{\rm{T}}}{(\mathit{\boldsymbol{AQ}}{\mathit{\boldsymbol{A}}^{\rm{T}}})^{ - 1}}\left( {\mathit{\boldsymbol{l}} - \mathit{\boldsymbol{Bx}}} \right) - {\mathit{\boldsymbol{G}}^{\rm{T}}}\mathit{\boldsymbol{\lambda }} + \mathop \sum \limits_{i = 1}^{n - 1} \frac{{\partial \left[ {{\lambda _i}\left( {{\mathit{\boldsymbol{G}}_i}\mathit{\boldsymbol{x}} - {w_i}} \right)} \right]}}{{\partial x}} = 0$$ (11)
    • 附有不等式约束的地表沉降时间序列EIV模型参数估计迭代算法步骤如下:

      1) 计算模型参数初值$\mathit{\boldsymbol{\hat x}} = {\left( {{\mathit{\boldsymbol{B}}^{\rm{T}}}\mathit{\boldsymbol{Q}}_l^{ - 1}\mathit{\boldsymbol{B}}} \right)^{ - 1}} \bullet \left( {{\mathit{\boldsymbol{B}}^{\rm{T}}}\mathit{\boldsymbol{Q}}_l^{ - 1}\mathit{\boldsymbol{l}}} \right)$, 令拉格朗日向量λ的初值$\mathit{\boldsymbol{\hat \lambda }} = 0$;

      2) 由约束方程${\mathit{\boldsymbol{e}}^{\rm{T}}}{\mathit{\boldsymbol{Q}}^{ - 1}}\mathit{\boldsymbol{e}}$求解其Hessian矩阵H与梯度向量g;

      3) 根据参数估值$\mathit{\boldsymbol{\hat x}}$与式(11)判别约束条件是否为有效约束, 计算向量λ。若${\lambda _i} \ge 0$, 则为有效约束, 保留计算值; 若λi<0, 则为无效约束, 取λi=0, 剔除相应约束方程;

      4) 根据判别后的约束方程, 计算参数xλ估值的改正数[15]

      $$\left[ {\begin{array}{*{20}{c}} {{\rm{d}}\mathit{\boldsymbol{x}}}\\ {{\rm{d}}\mathit{\boldsymbol{\lambda }}} \end{array}} \right] = - {\left[ {\begin{array}{*{20}{c}} {\mathit{\boldsymbol{H}} + \mathit{\boldsymbol{\hat \lambda }}{\mathit{\boldsymbol{H}}_t}}&{{\mathit{\boldsymbol{G}}^{\rm{T}}}}\\ \mathit{\boldsymbol{G}}&0 \end{array}} \right]^{ - 1}}\left[ {\begin{array}{*{20}{c}} {\mathit{\boldsymbol{g}} + {\mathit{\boldsymbol{G}}^{\rm{T}}}\mathit{\boldsymbol{\hat \lambda }}}\\ {\mathit{\boldsymbol{G\hat x}} - \mathit{\boldsymbol{w}}} \end{array}} \right]$$

      式中, Ht为约束方程的Hessian矩阵。

      5) 根据模型参数估值计算误差向量e的估值$\mathit{\boldsymbol{\hat e}} = \mathit{\boldsymbol{Q}}{\mathit{\boldsymbol{A}}^{\rm{T}}}{(\mathit{\boldsymbol{AQ}}{\mathit{\boldsymbol{A}}^{\rm{T}}})^{ - 1}}\left( {\mathit{\boldsymbol{l}} - \mathit{\boldsymbol{B\hat x}}} \right)$ [14], 根据上文建立的中位数法计算观测元素的改正值, 并改正观测向量与系数矩阵中相应观测元素;

      6) 重复步骤2)~5), 迭代求解模型参数估值直至收敛。

    • 应用算例对本文建立的算法进行验证。首先, 为检验应用基于中位函数计算观测值改正数的EIV模型参数估计算法的可行性, 在无约束条件下, 分别应用AR模型参数估计的LS、奇异值分解算法(singular value decomposition, SVD)、文献[12]建立的基于虚拟观测值的总体最小二乘算法(TLS-virtul, TLS-V)与基于中位数的总体最小二乘算法(TLS-media, TLS-M)进行比较, 通过模型参数的求解精度验证算法的可行性。

      算例1采用文献[12]的模拟数据, 其模型参数模拟值分别为a1=-0.3、a2=-0.5、a3=0.8, 模拟数据列于表 1。算例2根据实测数据, 应用最小二乘算法(LS)、附不等式约束的最小二乘算法[16](inequality constraints LS, ICLS)与本文建立的附不等式约束EIV模型参数估计算法(inequality constraints TLS, ICTLS)进行比较, 并通过检核点的精度验证本文所建立算法的可靠性。

      表 1  算例1应用的模拟数据

      Table 1.  Simulation Data of Example 1

      序号 模拟值
      1 15.9189
      2 14.7974
      3 23.2462
      4 7.3572
      5 -10.8052
      6 -18.0692
      7 -10.5677
      8 2.7923
      9 14.6792
      10 13.1268
      11 1.1775
      12 -10.7078
      13 -13.8582
      14 -6.0010
      15 5.4601
      16 9.6419
      17 9.1258
      18 -1.9343
      19 -6.5519
      20 -9.4017
      21 -2.8574
      22 4.8899
      23 7.6215
      24 5.0390
      25 0.6078
      26 -6.3572
      27 -5.7925
      28 -1.3012
      29 2.1084
      30 5.8308

      应用上述LS算法、SVD算法、TLS-V算法与本文建立的TLS-M算法得到的参数估值与单位权中误差σ0列于表 2。算例结果表明, 由TLS-M算法得到的参数估值精度低于TLS-V算法, 高于LS、SVD算法, 但与TLS-V算法相比, TLS-M算法无须改写原AR模型, 便于程序化设计。

      表 2  4种算法得到的参数估值与其精度

      Table 2.  Model Parameters Estimation and Its Accuracy Based on the Four Methods

      算法 a1 a2 a3 精度
      LS -0.3246 -0.4667 0.7709 1.4160
      SVD -0.2725 -0.5461 0.8363 1.4455
      TLS-V -0.3068 -0.5084 0.8043 0.9042
      TLS-M -0.2903 -0.5077 0.8136 0.9203

      算例2采用某地区滑坡体形变监测中某一监测点观测数据, 该监测点共有28期形变监测数据, 其高程变化值列于表 3

      表 3  监测点形变观测数据

      Table 3.  Observation Data for Settlement of Monitoring Point

      周期 高差/mm
      1 -1.5
      2 9.6
      3 -2.4
      4 -1
      5 -5.6
      6 -1.5
      7 -8.5
      15 -3.6
      16 -2.9
      17 -10.4
      18 -14.3
      19 -3.1
      20 -9.4
      21 -0.6
      8 -3.3
      9 -4.1
      10 -5.3
      11 -1.2
      12 -4.4
      13 -1.6
      14 2
      22 -4.6
      23 -2.3
      24 -3.1
      25 -2.3
      26 -3.9
      27 -2.6
      28 -2.7

      监测结果表明, 受观测误差等因素的影响, 监测周期2、14观测到的高程变化值大于零, 不符合于滑坡体的形变垂直于坡度线向下、高程变化值应小于零的固有属性。本文分别应用LS、ICLS、本文建立的ICTLS时间序列模型对监测数据进行处理, 约束条件为监测点前一期观测值应大于后一期观测值, 即监测点高程变化值小于零。模型参数个数为3个, 即$a_{1}, a_{2}, a_{3}$, 应用前20期监测数据拟合模型参数。应用后8期监测数据计算均方根误差(root mean squared error, RMSE), 进行精度检核。3种算法求解的模型参数估值与检核精度列于表 4。从表 4可以看出来本文建立的附不等式约束的参数估计模型推估精度优于现有算法。

      表 4  3种算法求解的模型参数估值与检核精度

      Table 4.  Model Parameters Estimation and Accuracy by the Three Different Algorithms

      算法 a1 a2 a3 RMSE/mm
      LS 0.0277 0.3714 0.4470 0.8235
      ICLS 0.3234 0.4131 0.2401 0.7773
      ICTLS 0.2218 0.3747 0.4147 0.7173
    • 地表沉降可能对道路、房屋建筑、地下管线等造成破坏, 并可能引发次生地质灾害。由于地表沉降是随时间演变的地质演化过程, 应用时间序列的回归模型能够描述其演化特征。本文应用EIV模型讨论地表沉降时间序列模型, 根据时间序列模型与地表沉降先验信息建立具有不等式约束的地表沉降自回归EIV模型。应用非线性最优化理论, 构建地表沉降时间序列的附不等式约束模型的KKT条件, 并结合中位函数建立模型参数的嵌套迭代算法。通过算例对所建立的算法进行验证, 结果表明本文所建立的算法可行。应用地表沉降固有属性建立具有约束条件的参数估计模型能够克服观测误差等因素对地表沉降信息的扭曲, 客观反映监测对象的变形特征, 这一研究思路在以后的研究工作中需要继续关注。

参考文献 (16)

目录

    /

    返回文章
    返回