留言板

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

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

病态方程基于Liu估计的一种迭代估计新方法

姜兆英 刘国林 于胜文

姜兆英, 刘国林, 于胜文. 病态方程基于Liu估计的一种迭代估计新方法[J]. 武汉大学学报 ● 信息科学版, 2017, 42(8): 1172-1178. doi: 10.13203/j.whugis20150218
引用本文: 姜兆英, 刘国林, 于胜文. 病态方程基于Liu估计的一种迭代估计新方法[J]. 武汉大学学报 ● 信息科学版, 2017, 42(8): 1172-1178. doi: 10.13203/j.whugis20150218
JIANG Zhaoying, LIU Guolin, YU Shengwen. A New Iterative Estimator Method Based on Liu Estimator for Ill-Posed Equations[J]. Geomatics and Information Science of Wuhan University, 2017, 42(8): 1172-1178. doi: 10.13203/j.whugis20150218
Citation: JIANG Zhaoying, LIU Guolin, YU Shengwen. A New Iterative Estimator Method Based on Liu Estimator for Ill-Posed Equations[J]. Geomatics and Information Science of Wuhan University, 2017, 42(8): 1172-1178. doi: 10.13203/j.whugis20150218

病态方程基于Liu估计的一种迭代估计新方法

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

国家自然科学基金 41274007

国家自然科学基金 41404003

山东省自然科学基金 ZR2012DM001

高等学校博士学科点专项科研基金 20123718110001

详细信息
    作者简介:

    姜兆英, 博士, 主要从事InSAR数据处理和模型算法研究。15964287012@163.com

    通讯作者: 刘国林, 博士, 教授。gliu@sdust.edu.cn
  • 中图分类号: P207

A New Iterative Estimator Method Based on Liu Estimator for Ill-Posed Equations

Funds: 

The National Natural Science Foundation of China 41274007

The National Natural Science Foundation of China 41404003

Shandong Province Natural Science Foundation of China ZR2012DM001

Specialized Research Fund for the Doctoral Program of Higher Education 20123718110001

More Information
    Author Bio:

    JIANG Zhaoying, PhD, specializes in InSAR data processing and model algorithm. E-mail:15964287012@163.com

    Corresponding author: LIU Guolin, PhD, professor. E-mail: gliu@sdust.edu.cn
  • 摘要: 当线性回归模型的设计矩阵病态时,最小二乘(least square,LS)估值方差大且不稳定,已不是一种优良估计。为了减弱病态性,许多有偏估计法如岭估计、主成分估计、Liu估计等被提出。基于Liu估计,引入迭代的思想,提出了一种新的有偏估计法—迭代估计法。借助对称正定矩阵的谱分解,将迭代公式转化为便于解算的解析表达式,并证明迭代公式在修正因子d∈[-1,1]是收敛的。基于Liu估计中修正因子d的确定方法,在均方误差最小的情况下给出最优修正因子d的确定公式。最后,分别利用LS估计、岭估计、Liu估计和提出的迭代估计对两个算例进行计算并给出实验结果。在第一个算例中,对观测向量添加不同的扰动,结果表明迭代估计法具有更强的抗干扰能力;第二个算例的结果表明,迭代估计法所得结果更接近于真值,即迭代估计法在均方误差意义下优于LS估计、岭估计和Liu估计。
  • 表  1  ε=0.1时不同方法所得结果

    Table  1.   Results of Different Methods when ε=0.1

    较值 真值 LS估计 岭估计 Liu估计 本文方法
    1 1.606 5 1.489 1 1.050 6 1.050 6
    1 -16.160 0 -5.304 3 1.053 7 1.053 9
    1 69.522 0 12.243 6 1.050 9 1.050 2
    1 294.290 0 12.176 4 1.049 1 1.046 1
    $\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$ 1 -2 082.400 0 -6.846 8 1.021 0 1.042 3
    1 3 765.700 0 -19.525 6 1.077 7 1.039 1
    1 -1 717.100 0 -16.402 1 1.018 6 1.036 2
    1 -1 428.500 0 2.268 8 1.019 1 1.033 8
    1 1 128.000 0 32.849 2 1.043 2 1.031 6
    ‖Δ$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$‖ 0 4.986 9×103 45.749 3 0.139 6 0.129 9
    MSE($\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$) 0 2.486 9×107 20.930 0.019 48 0.016 8
    下载: 导出CSV

    表  2  ε=0.01时不同方法所得结果

    Table  2.   Results of Different Methods when ε=0.01

    较值 真值 LS估计 岭估计 Liu估计 本文方法
    1 1.060 7 1.048 9 1.005 1 1.005 1
    1 -7.160 3 -3.695 7 1.005 4 1.005 4
    1 7.852 2 2.124 4 1.005 1 1.005 0
    1 30.329 0 2.117 6 1.004 9 1.004 6
    $\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$ 1 -207.340 0 2.153 2 1.002 1 1.004 2
    1 377.470 0 -1.052 6 1.007 8 1.003 9
    1 -170.810 0 -7.402 1 1.001 9 1.003 6
    1 -141.950 0 1.126 9 1.001 9 1.003 4
    1 113.700 0 4.184 9 1.004 3 1.003 2
    ‖Δ$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$‖ 0 4.986 9×102 4.574 9 0.013 97 0.012 99
    MSE($\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$) 0 2.486 9×105 2.093 0 1.951 1×10-4 1.686 8×10-4
    下载: 导出CSV

    表  3  ε=0.001时不同方法所得结果

    Table  3.   Results of Different Methods when ε=0.001

    较值 真值 LS估计 岭估计 Liu估计 本文方法
    1 1.006 1 1.004 9 1.000 5 1.000 5
    1 -0.828 4 0.936 9 1.000 5 1.000 5
    1 1.685 2 1.002 4 1.000 5 1.000 5
    1 3.932 9 1.111 8 1.000 5 1.000 5
    $\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$ 1 -1.983 4 0.921 5 1.000 2 1.000 4
    1 3.864 7 0.794 7 1.000 8 1.000 4
    1 -1.618 1 0.825 9 1.000 2 1.000 4
    1 -1.329 5 1.012 7 1.000 2 1.000 3
    1 1.227 0 1.318 5 1.000 4 1.000 3
    ‖Δ$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$‖ 0 4.986 9 0.457 5 1.406 5×10-3 1.298 8×10-3
    MSE($\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$) 0 2.486 9×103 0.209 30 1.978 2×10-6 1.686 8×10-6
    下载: 导出CSV

    表  4  控制点的坐标和观测距离

    Table  4.   Coordinate and Observation Distance of Control Points

    点号 坐标/m 观测距离/m
    x y z di, 10 di, 11
    P1 23.000 10.000 0.010 25.078 69 16.765 17
    P2 10.000 9.990 0.000 14.134 51 17.719 65
    P3 35.000 10.010 -0.010 36.415 88 28.442 94
    P4 100.000 19.990 0.005 101.479 43 93.168 39
    P5 -36.000 10.005 0.000 37.364 22 43.299 05
    P6 0.000 10.010 -0.005 10.010 04 8.600 60
    P7 56.000 9.995 0.010 56.996 06 49.256 18
    P8 -15.000 10.105 -0.010 18.035 90 22.559 66
    P9 -1.7000 10.008 -0.015 10.150 63 10.043 82
    下载: 导出CSV

    表  5  LS估计、岭估计、Liu估计以及本文方法的结果比较

    Table  5.   Comparison of the Results Among LS Estimator、Ridge Estimator、Liu Estimator and Our Method

    较值 真值 LS估计 岭估计 Liu估计 本文方法
    0 -0.036 8 -0.041 1 -0.038 0 -0.038 4
    0 0.051 8 0.027 5 0.022 7 0.022 0
    0 9.363 1 0.086 0 0.375 9 0.055 9
    $\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$ 7 7.048 8 6.974 6 6.975 9 6.973 3
    10 5.596 3 10.036 3 9.866 2 10.019 6
    -5 -4.648 2 -4.873 5 -4.904 4 -4.911 7
    参数 0 0 0.482 7 0.137 5 0.137 5
    ‖Δ$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$‖ 0 10.353 3 0.166 7 0.214 2 0.118 3
    MSE($\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$) 0 18.281 3 0.219 6 0.170 9 0.014 0
    下载: 导出CSV
  • [1] Stein C. Inadmissibility of the Usual Estimator for Mean of Multivariate Normal Distribution[C]. 3th Berkley Symposium on Mathematical and Statistics Probability, Berkley, 1956
    [2] Malthouse E C. Shrinkage Estimation and Direct Marketing Scoring Model[J]. Journal of Interactive Marketing, 1999, 13(4):10-23 doi:  10.1002/(SICI)1520-6653(199923)13:4<10::AID-DIR2>3.0.CO;2-3
    [3] Massy W F. Principal Components Regression in Exploratory Statistical Research[J]. Journal of the American Statistical Association, 1965, 60(309):234-266 doi:  10.1080/01621459.1965.10480787
    [4] Hoerl A E, Kennard R W. Ridge Regression:Biased Estimation for Non-Orthogonal Problems[J]. Technometrics, 1970, 12(1):55-67 doi:  10.1080/00401706.1970.10488634
    [5] Hoerl A E, Kennard R W. Ridge Regression:Application for Non-Orthogonal Problems[J]. Technometrics, 1970, 12(1):69-82 doi:  10.1080/00401706.1970.10488635
    [6] 郭杰, 归庆明, 郭淑妹, 等.利用复共线性诊断确定偏差矫正项的截断型岭估计[J].武汉大学学报·信息科学版, 2015, 40(6):785-789 http://ch.whu.edu.cn/CN/abstract/abstract3277.shtml

    Guo Jie, Gui Qingming, Guo Shumei, et al. Truncation Ridge Estimation Based on the Biased Corrected Theory by Multicollinearity Diagnosis[J].Geomatics and Information Science of Wuhan University, 2015, 40(6):785-789 http://ch.whu.edu.cn/CN/abstract/abstract3277.shtml
    [7] Liu Kejian. A New Class of Biased Estimate in Linear Regression[J].Communications in Statistics-Theory and Methods, 1993, 22(2):393-402 doi:  10.1080/03610929308831027
    [8] Liu Kejian. Using Liu-Type Estimator to Combat Collinearity[J].Communications in Statistics-Theory and Methods, 2003, 32(5):1009-1020 doi:  10.1081/STA-120019959
    [9] Mansson K. Developing a Liu Estimator for the Negative Binomial Regression Model:Method and Application[J].Journal of Statistical Computation and Simulation, 2013, 83(9):1773-1780 doi:  10.1080/00949655.2012.673127
    [10] 王新洲, 刘丁酉, 张前勇, 等.谱修正迭代法及其在测量数据处理中的应用[J].黑龙江工程学院学报, 2001, 15(2):3-6 http://www.cnki.com.cn/Article/CJFDTOTAL-JTGZ200102000.htm

    Wang Xinzhou, Liu Dingyou, Zhang Qianyong, et al. The Iteration by Correcting Characteristic Value and Its Application in Surveying Data Processing[J]. Journal of Heilongjiang Institute of Technology, 2001, 15(2):3-6 http://www.cnki.com.cn/Article/CJFDTOTAL-JTGZ200102000.htm
    [11] 刘斌, 龚健雅, 江万寿, 等.基于岭参数的谱修正迭代法及其在有理多项式参数求解中的应用[J].武汉大学学报·信息科学版, 2012, 37(4):399-402 http://ch.whu.edu.cn/CN/abstract/abstract176.shtml

    Liu Bin, Gong Jianya, Jiang Wanshou, et al. Improvement of the Iteration by Correcting Characteristic Value Based on Ridge Estimation and Its Application in RPC Calculating[J]. Geomatics and Information Science of Wuhan University, 2012, 37(4):399-402 http://ch.whu.edu.cn/CN/abstract/abstract176.shtml
    [12] 郭秋英, 胡振琪.遗传算法在GPS快速定位病态方程解算中的应用[J].武汉大学学报·信息科学版, 2009, 34(2):240-243 http://ch.whu.edu.cn/CN/abstract/abstract1181.shtml

    Guo Qiuying, Hu Zhenqi. Application of Genetic Algorithm to Solve Ill-Conditioned Equations for GPS Rapid Positioning[J].Geomatics and Information Science of Wuhan University, 2009, 34(2):240-243 http://ch.whu.edu.cn/CN/abstract/abstract1181.shtml
    [13] 陈希孺, 王松桂.近代回归分析:原理方法及应用[M].合肥:安徽教育出版社, 1987

    Chen Xiru, Wang Guisong. Modern Regression Analysis:Theory and Application[M]. Hefei:Anhui Education Press, 1987
    [14] Hansen P C, O'Leary D P. The Use of the L-Curve in the Regularization of Discrete Ill-Posed Problems[J]. SIAM Journal on Scientific Computing, 1993, 14(6):1487-1503 doi:  10.1137/0914086
    [15] 王振杰, 欧吉坤.用L-曲线法确定岭估计中的岭参数[J].武汉大学学报·信息科学版, 2004, 29(3):235-238 http://ch.whu.edu.cn/CN/abstract/abstract4657.shtml

    Wang Zhenjie, Ou Jikun. Determining the Ridge Parameter in a Ridge Estimation Using L-Curve Method[J]. Geomatics and Information Science of Wuhan University, 2004, 29(3):235-238 http://ch.whu.edu.cn/CN/abstract/abstract4657.shtml
    [16] 王家映.地球物理反演理论[M].北京:高等教育出版社, 2002

    Wang Jiaying. Inverse Theory in Geophysics[M]. Beijing:Higher Education Press, 2002
    [17] 黄维彬.近代平差理论及其应用[M].北京:解放军出版社, 1992

    Huang Weibin. Modern Adjustment Theory and Its Application[M]. Beijing:People's Liberation Army Press, 1992
    [18] Sakallioglu S, Kaciranlar S. A New Biased Estimator Basedon Ridge Estimation[J].Statistics Papers, 2008, 49:669-689 doi:  10.1007/s00362-006-0037-0
    [19] Sakallioglu S, Kaciranlar S, Akdeniz F. Mean Squared Error Comparisons of Some Biased Regression Estimators[J].Communications in Statistics:Theory and Methods, 2001, 30(2):347-361 doi:  10.1081/STA-100002036
    [20] 刘丁酉.矩阵分析[M].武汉:武汉大学出版社, 2006

    Liu Dingyou. Matrix Analysis[M]. Wuhan:Wuhan University Press, 2006
    [21] Golam K B M. Applications of Some Improved Estimators in Linear Regression[J].Journal of Modern Applied Statistical Methods, 2006, 5(2):367-380 http://digitalcommons.wayne.edu/cgi/viewcontent.cgi?article=1198&context=jmasm
    [22] 郭建峰. 测量平差系统病态性的诊断与处理[D]. 郑州: 信息工程大学, 2002

    Guo Jianfeng. Diagnostics and Processing of Ill-Posed in Surveying Adjustment System[D].Zhengzhou:Information Engineering University, 2002
  • [1] 韩松辉, 宫轶松, 李建文, 马朝忠, 李新娜, 郭淑妹.  基于ARMA模型的BDS卫星钟差异常值探测及其短期预报 . 武汉大学学报 ● 信息科学版, 2021, 46(2): 244-251. doi: 10.13203/j.whugis20190232
    [2] 李豪, 顾勇为, 韩松辉.  基于信噪比检验的双截断奇异值估计 . 武汉大学学报 ● 信息科学版, 2019, 44(2): 228-232, 239. doi: 10.13203/j.whugis20170051
    [3] 周跃寅, 潘国荣, 吴廷, 汪大超.  数据融合模型选取对工业测量整体平差结果的影响分析 . 武汉大学学报 ● 信息科学版, 2017, 42(12): 1840-1846. doi: 10.13203/j.whugis20150450
    [4] 林东方, 朱建军.  附有奇异值修正限制的改进的岭估计方法 . 武汉大学学报 ● 信息科学版, 2017, 42(12): 1834-1839. doi: 10.13203/j.whugis20150581
    [5] 李永明, 归庆明, 顾勇为, 韩松辉.  有偏卡尔曼滤波及其算法 . 武汉大学学报 ● 信息科学版, 2016, 41(7): 946-951. doi: 10.13203/j.whugis20140072
    [6] 郭杰, 归庆明, 郭淑妹, 张宁.  利用复共线性诊断确定偏差矫正项的截断型岭估计 . 武汉大学学报 ● 信息科学版, 2015, 40(6): 785-789. doi: 10.13203/j.whugis20130465
    [7] 任锴, 宋小勇, 贾小林, 欧阳桂崇.  一种区域卫星导航系统导航策略 . 武汉大学学报 ● 信息科学版, 2012, 37(11): 1356-1359.
    [8] 刘斌, 龚健雅, 江万寿, 祝小勇.  基于岭参数的谱修正迭代法及其在有理多项式参数求解中的应用 . 武汉大学学报 ● 信息科学版, 2012, 37(4): 399-402.
    [9] 陈刘成, 李静, 马瑞, 陈娉娉.  工程化广播星历参数拟合算法与接口设计 . 武汉大学学报 ● 信息科学版, 2011, 36(1): 18-23.
    [10] 张永军, 吴磊, 林立文, 赵家平.  摄影测量中病态问题的条件数指标分析 . 武汉大学学报 ● 信息科学版, 2010, 35(3): 308-312.
    [11] 袁修孝, 林先勇.  基于岭估计的有理多项式参数求解方法 . 武汉大学学报 ● 信息科学版, 2008, 33(11): 1130-1133.
    [12] 王振杰, 欧吉坤, 柳林涛.  一种解算病态问题的方法——两步解法 . 武汉大学学报 ● 信息科学版, 2005, 30(9): 821-824.
    [13] 王振杰, 欧吉坤.  L-曲线法确定岭估计中的岭参数 . 武汉大学学报 ● 信息科学版, 2004, 29(3): 235-238.
    [14] 徐天河, 杨元喜.  均方误差意义下正则化解优于最小二乘解的条件 . 武汉大学学报 ● 信息科学版, 2004, 29(3): 223-226.
    [15] 王新洲.  在无偏估计类中改进最小二乘估计的方法 . 武汉大学学报 ● 信息科学版, 1995, 20(1): 46-50.
    [16] 吴子安.  大坝变形监测数据回归分析中的因子选择 . 武汉大学学报 ● 信息科学版, 1993, 18(1): 20-26.
    [17] 张方仁, 汪晓庆.  平差参数的岭估计和压缩估计 . 武汉大学学报 ● 信息科学版, 1989, 14(3): 46-58.
    [18] 刘丁酉.  秩亏网伪逆平差值的修正估计 . 武汉大学学报 ● 信息科学版, 1988, 13(3): 90-100.
    [19] 黄幼才.  岭估计及其应用 . 武汉大学学报 ● 信息科学版, 1987, 12(4): 64-73.
    [20] 李德仁.  自检校平差中过度参数化的克服 . 武汉大学学报 ● 信息科学版, 1986, 11(3): 95-104.
  • 加载中
计量
  • 文章访问数:  1599
  • HTML全文浏览量:  54
  • PDF下载量:  342
  • 被引次数: 0
出版历程
  • 收稿日期:  2015-09-29
  • 刊出日期:  2017-08-05

病态方程基于Liu估计的一种迭代估计新方法

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

    国家自然科学基金 41274007

    国家自然科学基金 41404003

    山东省自然科学基金 ZR2012DM001

    高等学校博士学科点专项科研基金 20123718110001

    作者简介:

    姜兆英, 博士, 主要从事InSAR数据处理和模型算法研究。15964287012@163.com

    通讯作者: 刘国林, 博士, 教授。gliu@sdust.edu.cn
  • 中图分类号: P207

摘要: 当线性回归模型的设计矩阵病态时,最小二乘(least square,LS)估值方差大且不稳定,已不是一种优良估计。为了减弱病态性,许多有偏估计法如岭估计、主成分估计、Liu估计等被提出。基于Liu估计,引入迭代的思想,提出了一种新的有偏估计法—迭代估计法。借助对称正定矩阵的谱分解,将迭代公式转化为便于解算的解析表达式,并证明迭代公式在修正因子d∈[-1,1]是收敛的。基于Liu估计中修正因子d的确定方法,在均方误差最小的情况下给出最优修正因子d的确定公式。最后,分别利用LS估计、岭估计、Liu估计和提出的迭代估计对两个算例进行计算并给出实验结果。在第一个算例中,对观测向量添加不同的扰动,结果表明迭代估计法具有更强的抗干扰能力;第二个算例的结果表明,迭代估计法所得结果更接近于真值,即迭代估计法在均方误差意义下优于LS估计、岭估计和Liu估计。

English Abstract

姜兆英, 刘国林, 于胜文. 病态方程基于Liu估计的一种迭代估计新方法[J]. 武汉大学学报 ● 信息科学版, 2017, 42(8): 1172-1178. doi: 10.13203/j.whugis20150218
引用本文: 姜兆英, 刘国林, 于胜文. 病态方程基于Liu估计的一种迭代估计新方法[J]. 武汉大学学报 ● 信息科学版, 2017, 42(8): 1172-1178. doi: 10.13203/j.whugis20150218
JIANG Zhaoying, LIU Guolin, YU Shengwen. A New Iterative Estimator Method Based on Liu Estimator for Ill-Posed Equations[J]. Geomatics and Information Science of Wuhan University, 2017, 42(8): 1172-1178. doi: 10.13203/j.whugis20150218
Citation: JIANG Zhaoying, LIU Guolin, YU Shengwen. A New Iterative Estimator Method Based on Liu Estimator for Ill-Posed Equations[J]. Geomatics and Information Science of Wuhan University, 2017, 42(8): 1172-1178. doi: 10.13203/j.whugis20150218
  • 最小二乘估计自Gauss-Markov定理建立后,便作为一个良好的估计而广泛应用。然而,在应用中人们逐渐发现,法方程系数阵呈病态时,模型设计阵的列向量间存在复共线性关系(multicollinearity),这往往会引起LS估计性质显著变坏,模型的解不稳定。为了解决法矩阵病态的问题,研究学者提出了有偏估计这一理论和方法。目前常用的有偏估计有Stein压缩估计[1, 2]、主成分估计[3]、岭估计[4-6]、Liu估计[7-9]等,它们在均方误差(mean squared error, MSE)准则下均局部或一致地改进了参数的LS估计。但当存在较严重的复共线性时,岭估计、Liu估计不能完全克服病态问题。除此之外,谱修正迭代法[10, 11]、遗传算法[12]等也是常用的解决病态问题的方法。

    考虑线性模型:

    $$ \mathit{\boldsymbol{L}} = \mathit{\boldsymbol{AX + \boldsymbol{\varDelta} }} $$ (1)

    式中,设计阵A为已知的m×n(m>n)阶列满秩阵;Xn×1阶未知参数向量;Lm×1阶观测向量;Δm×1阶误差向量;且E(Δ)=0, Cov(Δ)=σ2I。其经典最小二乘估值为:

    $$ {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_{{\rm{LS}}}} = {\left( {\mathit{\boldsymbol{A'A}}} \right)^{ - 1}}\mathit{\boldsymbol{A'L}} $$ (2)

    AA病态,即它的特征根中至少有一个近似于零,(AA)-1中的元素变得很大,一方面导致最小二乘估值不稳定且偏离真值;另一方面会放大模型观测中干扰的影响。所以,若模型(1) 病态,LS估计已不再适应。判断矩阵AA是否“病态”或“病态”的程度,可以采用条件数来衡量。AA (因A是列满秩的,所以AA为对称正定矩阵)的条件数为[13]

    $$ {\rm{cond}}\left( {\mathit{\boldsymbol{A'A}}} \right) = \frac{{{r_{\max }}}}{{{r_{\min }}}} $$ (3)

    式中,rmaxrmin表示矩阵AA的最大、最小特征值。条件数越大,矩阵AA病态性越严重。为了解决AA病态问题,Hoerl和Kennard[4]在1970年提出了岭估计,则模型(1) 的岭估计为:

    $$ {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_\lambda } = {\left( {\mathit{\boldsymbol{A'A}} + \lambda \mathit{\boldsymbol{I}}} \right)^{ - 1}}\mathit{\boldsymbol{A'L}} $$ (4)

    式中,λ>0为岭参数,其确定方法有L-曲线法[14, 15]、岭迹法[16]、广义交叉效验(generalized cross-validation,GCV)法[17]等。这里AA +λI的条件数为:

    $$ \frac{{{r_{\max }} + \lambda }}{{{r_{\min }} + \lambda }} \buildrel \Delta \over = f\left( \lambda \right) $$ (5)

    显然,λ>0时有如下关系:

    $$ {\rm{cond}}\left( {\mathit{\boldsymbol{A'A}} + \lambda \mathit{\boldsymbol{I}}} \right) < {\rm{cond}}\left( {\mathit{\boldsymbol{A'A}}} \right) $$ (6)

    $$ f'\left( \lambda \right) = \frac{{{r_{\min }} - {r_{\max }}}}{{{{\left( {{r_{\min }} + \lambda } \right)}^2}}} < 0 $$ (7)

    故cond(AA +λI)是λ的递减函数。所以,当AA严重病态时,可选取足够大的λ使AA +λI的条件数降到预期值。但参数λ太大,虽然可以稳定地求解,但该解与模型(1) 的实际解相差甚远,是一个相当糟糕的逼近。但若λ太小,则对问题病态性的改善不起作用,解的不稳定性依然存在。1993年,Liu结合岭估计和Stein估计提出了一种新的有偏估计:

    $$ {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_d} = {\left( {\mathit{\boldsymbol{A'A}} + \mathit{\boldsymbol{I}}} \right)^{ - 1}}\left( {\mathit{\boldsymbol{A'L}} + \mathit{\boldsymbol{d}}{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_{{\rm{LS}}}}} \right) $$ (8)

    式中,0 < d < 1为一修正因子。特别当d=1时,Liu估计为经典的最小二乘估计。同时Liu还给出了其推广形式:

    $$ {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_d} = {\left( {\mathit{\boldsymbol{A'A}} + \mathit{\boldsymbol{I}}} \right)^{ - 1}}\left( {\mathit{\boldsymbol{A'L}} + \mathit{\boldsymbol{D}}{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_{{\rm{LS}}}}} \right) $$ (9)

    式中,D =diag(di),为修正因子di的对角矩阵,且0 < di < 1(i=1, 2, …, n)。文献[7]指出Liu估计可得到比LS估计更小的MSE,但相较于岭估计,Liu估计得到的MSE不一定优于岭估计。由于${\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} _d}$是参数d的线性函数,故d的确定比岭参数的确定要容易些。2008年,Sakallioglu[18]基于岭估计又提出了另外一种有偏估计:

    $$ {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_{\lambda ,d}} = {\left( {\mathit{\boldsymbol{A'A}} + \mathit{\boldsymbol{I}}} \right)^{ - 1}}\left( {\mathit{\boldsymbol{A'L}} + d{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_\lambda }} \right) $$ (10)

    式中,-∞ < d < +∞。该估计就是将Liu估计中的LS估计用岭估计代替,并将修正因子d推广到可取任意实数。在文献[18]中,Sakallioglu证明了该方法在MSE意义下优于最小二乘、岭估计和Liu估计。本文在Liu估计的基础上引入迭代求解的思想,提出了一种迭代估计方法。通过实验验证了该方法在均方误差衡量准则下优于LS估计、岭估计和Liu估计。

    • 虽然Liu估计能得到比LS更小的MSE,但相较于岭估计,只有在满足一定条件下Liu估计得到的MSE才会小于岭估计得到的[19]。为此,本文基于Liu估计,引入迭代的思想,提出一种有偏迭代估计法。

      对于初值的选取可以采用两种方式,一是LS估计,另一种为岭估计。下面以LS估计作为初值为例给出迭代公式。

      取$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$(0)= $\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$LS,为了方便书写,记q =(AA + I)-1,则:

      $$ \begin{array}{l} {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 1 \right)}} = {\left( {\mathit{\boldsymbol{A'A}} + \mathit{\boldsymbol{I}}} \right)^{ - 1}}\left( {\mathit{\boldsymbol{A'L}} + d{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_{{\rm{LS}}}}} \right) = \\ \;\;\;\;\;\;\;\;\;\;\mathit{\boldsymbol{qA'L}} + d\mathit{\boldsymbol{q}}{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 0 \right)}} \end{array} $$ (11)

      将$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$(1)代替式(11) 中的$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$(0),得第二次估值:

      $$ \begin{array}{l} {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 2 \right)}} = \mathit{\boldsymbol{qA'L}} + \mathit{\boldsymbol{q}}d{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 1 \right)}} = \\ \;\;\;\;\;\;\;\;\;\mathit{\boldsymbol{qA'L}} + \mathit{\boldsymbol{q}}d\left( {\mathit{\boldsymbol{qA'L}} + \mathit{\boldsymbol{q}}d{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 0 \right)}}} \right) = \\ \;\;\;\;\;\;\;\;\;\mathit{\boldsymbol{qA'L}} + d{\mathit{\boldsymbol{q}}^2}\mathit{\boldsymbol{A'L}} + {d^2}{\mathit{\boldsymbol{q}}^2}{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 0 \right)}} \end{array} $$ (12)

      同样将$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$(2)代替式(12) 中的$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$(1),得第三次估值:

      $$ \begin{array}{l} {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 3 \right)}} = \mathit{\boldsymbol{qA'L}} + \mathit{\boldsymbol{q}}d{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 2 \right)}} = \\ \;\;\;\;\;\;\;\;\;\;\mathit{\boldsymbol{qA'L}} + d\mathit{\boldsymbol{q}}\left( {\mathit{\boldsymbol{qA'L}} + d{\mathit{\boldsymbol{q}}^2}\mathit{\boldsymbol{A'L}} + } \right.\\ \;\;\;\;\;\;\;\;\;\;\left. {{d^2}{\mathit{\boldsymbol{q}}^2}{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 0 \right)}}} \right) = \mathit{\boldsymbol{qA'L}} + d{\mathit{\boldsymbol{q}}^2}\mathit{\boldsymbol{A'L}} + \\ \;\;\;\;\;\;\;\;\;\;{d^2}{\mathit{\boldsymbol{q}}^3}\mathit{\boldsymbol{A'L}} + {d^3}{\mathit{\boldsymbol{q}}^3}{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 0 \right)}} \end{array} $$ (13)

      依次类推,即得迭代公式为:

      $$ \begin{array}{l} {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( k \right)}} = \mathit{\boldsymbol{qA'L}} + d\mathit{\boldsymbol{q}}{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( {k - 1} \right)}} = \left( {\mathit{\boldsymbol{q}} + d{\mathit{\boldsymbol{q}}^2} + {d^2}{\mathit{\boldsymbol{q}}^3} + \cdots + {d^{k - 1}}{\mathit{\boldsymbol{q}}^k}} \right)\mathit{\boldsymbol{A'L}} + {d^k}{\mathit{\boldsymbol{q}}^k}{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 0 \right)}} = \\ \;\;\;\;\;\;\;\;\;\;\sum\limits_{j = 0}^{k - 1} {{d^j}{{\left[ {{{\left( {\mathit{\boldsymbol{A'A}} + \mathit{\boldsymbol{I}}} \right)}^{ - 1}}} \right]}^{j + 1}}\mathit{\boldsymbol{A'L}}} + {d^k}{\left[ {{{\left( {\mathit{\boldsymbol{A'A}} + \mathit{\boldsymbol{I}}} \right)}^{ - 1}}} \right]^k}{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_{{\rm{LS}}}} = \\ \;\;\;\;\;\;\;\;\;\;\left\{ {\sum\limits_{j = 0}^{k - 1} {{d^j}{{\left[ {{{\left( {\mathit{\boldsymbol{A'A}} + \mathit{\boldsymbol{I}}} \right)}^{ - 1}}} \right]}^{j + 1}}\mathit{\boldsymbol{A'A}}} + {d^k}{{\left[ {{{\left( {\mathit{\boldsymbol{A'A}} + \mathit{\boldsymbol{I}}} \right)}^{ - 1}}} \right]}^k}} \right\}{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_{{\rm{LS}}}} \end{array} $$ (14)

      式中,k=1, 2, …。

    • 矩阵AA为实对称正定阵,故存在正交矩阵Q (即QQ = I)使得AA = QΛQ′,这里Λ = diag{s1, s2, …, sn}是对角矩阵,且s1s2≥…≥sn>0为矩阵AAn个特征值。AA + I = QΛQ′ + QIQ′ = Q (Λ + I) Q′,则:

      $$ \mathit{\boldsymbol{q}} = {\left( {\mathit{\boldsymbol{A'A}} + \mathit{\boldsymbol{I}}} \right)^{ - 1}} = \mathit{\boldsymbol{Q}}{\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} + I}}} \right)^{ - 1}}\mathit{\boldsymbol{Q'}} = \mathit{\boldsymbol{QDQ'}} $$

      式中,

      $$ \mathit{\boldsymbol{D}} = {\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} + I}}} \right)^{ - 1}} = {\rm{diag}}\left\{ {1/\left( {{s_1} + 1} \right), \cdots ,1/\left( {{s_n} + 1} \right)} \right\} $$

      进一步有:

      $$ \begin{array}{*{20}{c}} {{d^j}{{\left[ {{{\left( {\mathit{\boldsymbol{A'A}} + \mathit{\boldsymbol{I}}} \right)}^{ - 1}}} \right]}^{j + 1}}\mathit{\boldsymbol{A'A}} = }\\ {{d^j}\left( {\mathit{\boldsymbol{QDQ'}} \cdot \mathit{\boldsymbol{QDQ'}} \cdot \cdots \cdot \mathit{\boldsymbol{QDQ'}}} \right)\mathit{\boldsymbol{Q \boldsymbol{\varLambda} Q' = }}}\\ {\mathit{\boldsymbol{Q}}\left( {{d^j}{\mathit{\boldsymbol{D}}^{j + 1}}\mathit{\boldsymbol{ \boldsymbol{\varLambda} }}} \right)\mathit{\boldsymbol{Q'}}} \end{array} $$ (15)

      将式(15) 代入式(14),并整理得:

      $$ \begin{array}{*{20}{c}} {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( k \right)}} = \left[ {\mathit{\boldsymbol{Q}}\left( {\sum\limits_{j = 0}^{k - 1} {{d^j}{\mathit{\boldsymbol{D}}^{j + 1}}\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + {d^k}{\mathit{\boldsymbol{D}}^k}} } \right)\mathit{\boldsymbol{Q'}}} \right]{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_{{\rm{LS}}}} = }\\ {\mathit{\boldsymbol{Q}}{\mathit{\boldsymbol{D}}_{k,d}}\mathit{\boldsymbol{Q'}}{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_{{\rm{LS}}}}} \end{array} $$ (16)

      式中,${\mathit{\boldsymbol{D}}_{k, d}} = \sum\limits_{j = 0}^{k-1} {{d^j}{\mathit{\boldsymbol{D}}^{j + 1}}\mathit{\boldsymbol{ \boldsymbol{\varLambda} }}} + {d^k}{\mathit{\boldsymbol{D}}^k}$,为一个n阶对角矩阵,第i个对角元素为:

      $$ {d_i} = \sum\limits_{j = 0}^{k - 1} {{{\left( {\frac{d}{{{s_i} + 1}}} \right)}^j} \cdot \frac{{{s_i}}}{{{s_i} + 1}} + {{\left( {\frac{d}{{{s_i} + 1}}} \right)}^k}} $$ (17)

      式中,i=1, …, nkN+N+表示正整数。

    • 将式(12) 与(11) 相减,得:

      $$ {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 2 \right)}} - {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 1 \right)}} = d\mathit{\boldsymbol{q}}\left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 1 \right)}} - {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 0 \right)}}} \right) $$ (18)

      将式(13) 与(12) 相减,得:

      $$ {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 3 \right)}} - {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 2 \right)}} = d\mathit{\boldsymbol{q}}\left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 2 \right)}} - {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 1 \right)}}} \right) $$ (19)

      将式(18) 代入式(19),得:

      $$ {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 3 \right)}} - {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 2 \right)}} = {d^2}{\mathit{\boldsymbol{q}}^2}\left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 1 \right)}} - {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 0 \right)}}} \right) $$ (20)

      以此类推,得:

      $$ \begin{array}{l} {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( {k + 1} \right)}} - {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( k \right)}} = {d^k}{\mathit{\boldsymbol{q}}^k}\left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 1 \right)}} - {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 0 \right)}}} \right) = \\ {d^k}{\left[ {{{\left( {\mathit{\boldsymbol{A'A}} + \mathit{\boldsymbol{I}}} \right)}^{ - 1}}} \right]^k}\left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 1 \right)}} - {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}^{\left( 0 \right)}}} \right) \end{array} $$ (21)

      定理1 [20]  设n阶任意方阵M,则方阵M幂收敛于零矩阵的充要条件为M的所有特征值的模都小于1。

      s(s>0) 为方阵AA的某一特征值,则矩阵d(AA + I)-1相应的特征值为$\frac{d}{{s + 1}}$。所以,若-1≤d≤1,则有$\left| {\frac{d}{{s + 1}}} \right| < 1$。根据定理1可知,当k→+∞,dk[(AA + I)-1]k→0,从而$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$(k+1)-$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$(k)→0,即迭代公式(14) 在-1≤d≤1是收敛的。

    • 在迭代公式(14) 中,若k=1则为Liu估计。为此,式(14) 中d的确定方法依据Liu估计中d的确定思路。

      将模型(1) 改写为典则形式:

      $$ \mathit{\boldsymbol{L}} = \mathit{\boldsymbol{BY + \boldsymbol{\varDelta} }} $$ (22)

      式中,B = AQY = QX。则有:

      $$ \mathit{\boldsymbol{B'B = Q'A'AQ = Q'Q \boldsymbol{\varLambda} Q'Q}} = \mathit{\boldsymbol{ \boldsymbol{\varLambda} }} $$ (23)

      则模型(22) 的最小二乘估值为$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{Y}}}$LS= Λ-1BL,岭估值为$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{Y}}}$λ=(Λ +λI)-1 BL,Liu估值为:

      $$ \begin{array}{*{20}{c}} {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over Y} }}}_d} = {{\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + \mathit{\boldsymbol{I}}} \right)}^{ - 1}}\left( {\mathit{\boldsymbol{B'L + }}d{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over Y} }}}_{{\rm{LS}}}}} \right) = }\\ {{{\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + \mathit{\boldsymbol{I}}} \right)}^{ - 1}}\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + d\mathit{\boldsymbol{I}}} \right){{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over Y} }}}_{{\rm{LS}}}}} \end{array} $$ (24)

      则有如下的线性关系[21]

      $$ {{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_d} = \mathit{\boldsymbol{Q}}{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over Y} }}}_d},{\rm{MSE}}\left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}_d}} \right) = {\rm{MSE}}\left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over Y} }}}_d}} \right) $$

      由式(24) 知:

      $$ \begin{array}{*{20}{c}} {E\left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over Y} }}}_d}} \right) = {{\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + \mathit{\boldsymbol{I}}} \right)}^{ - 1}}\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + d\mathit{\boldsymbol{I}}} \right)E\left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over Y} }}}_{{\rm{LS}}}}} \right) = }\\ {{{\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + \mathit{\boldsymbol{I}}} \right)}^{ - 1}}\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + d\mathit{\boldsymbol{I}}} \right)\mathit{\boldsymbol{Y}}} \end{array} $$ (25)
      $$ \begin{array}{*{20}{c}} {{\mathop{\rm cov}} \left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over Y} }}}_d}} \right) = {\sigma ^2}{{\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + \mathit{\boldsymbol{I}}} \right)}^{ - 1}}\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + d\mathit{\boldsymbol{I}}} \right){\mathit{\boldsymbol{ \boldsymbol{\varLambda} }}^{ - 1}} \cdot }\\ {\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + d\mathit{\boldsymbol{I}}} \right){{\left( {\mathit{\boldsymbol{ \boldsymbol{\varLambda} }} + \mathit{\boldsymbol{I}}} \right)}^{ - 1}}} \end{array} $$ (26)

      所以,Liu估值$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{Y}}}$d的均方误差为:

      $$ \begin{array}{*{20}{c}} {{\rm{MSE}}\left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over Y} }}}_d}} \right) = {{\left\| {E\left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over Y} }}}_d}} \right) - \mathit{\boldsymbol{Y}}} \right\|}^2} + {\rm{trcov}}\left( {{{\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over Y} }}}_d}} \right) = }\\ {\sum\limits_{i = 1}^n {\frac{{{{\left( {d - 1} \right)}^2}y_i^2}}{{{{\left( {{s_i} + 1} \right)}^2}}}} + {\sigma ^2}\sum\limits_{i = 1}^n {\frac{{{{\left( {{s_i} + d} \right)}^2}}}{{{s_i}{{\left( {{s_i} + 1} \right)}^2}}} = f\left( d \right)} } \end{array} $$ (27)

      所以,

      $$ \begin{array}{*{20}{c}} {f'\left( d \right) = 2\left( {d - 1} \right)\sum\limits_{i = 1}^n {\frac{{y_i^2}}{{{{\left( {{s_i} + 1} \right)}^2}}} + } }\\ {2{\sigma ^2}\sum\limits_{i = 1}^n {\frac{{{s_i} + d}}{{{s_i}{{\left( {{s_i} + 1} \right)}^2}}}} } \end{array} $$

      $$ f'\left( 1 \right) = 2{\sigma ^2}\sum\limits_{i = 1}^n {\frac{1}{{{s_i}{{\left( {{s_i} + 1} \right)}^2}}}} > 0 $$ (28)

      因此一定存在0 < d < 1,满足MSE($\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{Y}}}$d) < MSE($\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{Y}}}$LS)。使MSE达到最小的是使f′(d)=0的解,解得使MSE最小的修正因子为:

      $$ {d_{{\rm{opt}}}} = \frac{{\sum\limits_{i = 1}^n {\left( {y_i^2 - {\sigma ^2}} \right)/{{\left( {{s_i} + 1} \right)}^2}} }}{{\sum\limits_{i = 1}^n {\left( {{s_i}y_i^2 + {\sigma ^2}} \right)/{s_i}{{\left( {{s_i} + 1} \right)}^2}} }} $$ (29)

      式中,σ2yi(i=1, …, n)未知,可以用LS估计法或者岭估计法得到的$\overset\frown{{{\sigma }^{2}}}$、$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{Y}}}$替代。

    • 算例1  将11×9阶的Hilbert矩阵作为线性模型AX = L的观测矩阵A,其元素aij=$\frac{1}{i+j-1}$。设未知参数的真值均为1,则没有误差的常数项为:

      L =[2.829 0  1.929 0  1.519 9  1.269 9  1.096 8  0.968 2  0.868 2  0.787 9  0.721 7  0.666 1  0.618 8]′。

      根据式(3),求出法矩阵的条件数为8.985 9×1016,属于严重病态问题。在具体的处理过程中,对L的每个分量分别添加不同的扰动误差ε(ε=0.1, ε=0.01, ε=0.001)。分别用LS估计、岭估计(参数的确定选用L-曲线法)、Liu估计和本文的迭代估计公式进行计算,比较其估值$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$、$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$的均方误差以及$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$与X差的范数‖Δ$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$‖=‖ $\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$ -X ‖的大小,其计算结果见表 1~3

      表 1  ε=0.1时不同方法所得结果

      Table 1.  Results of Different Methods when ε=0.1

      较值 真值 LS估计 岭估计 Liu估计 本文方法
      1 1.606 5 1.489 1 1.050 6 1.050 6
      1 -16.160 0 -5.304 3 1.053 7 1.053 9
      1 69.522 0 12.243 6 1.050 9 1.050 2
      1 294.290 0 12.176 4 1.049 1 1.046 1
      $\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$ 1 -2 082.400 0 -6.846 8 1.021 0 1.042 3
      1 3 765.700 0 -19.525 6 1.077 7 1.039 1
      1 -1 717.100 0 -16.402 1 1.018 6 1.036 2
      1 -1 428.500 0 2.268 8 1.019 1 1.033 8
      1 1 128.000 0 32.849 2 1.043 2 1.031 6
      ‖Δ$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$‖ 0 4.986 9×103 45.749 3 0.139 6 0.129 9
      MSE($\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$) 0 2.486 9×107 20.930 0.019 48 0.016 8

      表 2  ε=0.01时不同方法所得结果

      Table 2.  Results of Different Methods when ε=0.01

      较值 真值 LS估计 岭估计 Liu估计 本文方法
      1 1.060 7 1.048 9 1.005 1 1.005 1
      1 -7.160 3 -3.695 7 1.005 4 1.005 4
      1 7.852 2 2.124 4 1.005 1 1.005 0
      1 30.329 0 2.117 6 1.004 9 1.004 6
      $\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$ 1 -207.340 0 2.153 2 1.002 1 1.004 2
      1 377.470 0 -1.052 6 1.007 8 1.003 9
      1 -170.810 0 -7.402 1 1.001 9 1.003 6
      1 -141.950 0 1.126 9 1.001 9 1.003 4
      1 113.700 0 4.184 9 1.004 3 1.003 2
      ‖Δ$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$‖ 0 4.986 9×102 4.574 9 0.013 97 0.012 99
      MSE($\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$) 0 2.486 9×105 2.093 0 1.951 1×10-4 1.686 8×10-4

      表 3  ε=0.001时不同方法所得结果

      Table 3.  Results of Different Methods when ε=0.001

      较值 真值 LS估计 岭估计 Liu估计 本文方法
      1 1.006 1 1.004 9 1.000 5 1.000 5
      1 -0.828 4 0.936 9 1.000 5 1.000 5
      1 1.685 2 1.002 4 1.000 5 1.000 5
      1 3.932 9 1.111 8 1.000 5 1.000 5
      $\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$ 1 -1.983 4 0.921 5 1.000 2 1.000 4
      1 3.864 7 0.794 7 1.000 8 1.000 4
      1 -1.618 1 0.825 9 1.000 2 1.000 4
      1 -1.329 5 1.012 7 1.000 2 1.000 3
      1 1.227 0 1.318 5 1.000 4 1.000 3
      ‖Δ$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$‖ 0 4.986 9 0.457 5 1.406 5×10-3 1.298 8×10-3
      MSE($\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$) 0 2.486 9×103 0.209 30 1.978 2×10-6 1.686 8×10-6

      表 1~3可以得出如下结论:

      1) 在法矩阵严重病态的情况下,LS的解与真值偏差很大,该估计法已不适应,且偏差随着观测量扰动误差大小的变化而变化,扰动越大,偏差越大。

      2) 在不同扰动误差的情况下,岭估计有效改善了LS估计,但在添加扰动ε=0.1时,岭估值与真值的偏差依然很大。所以岭估计虽然能改善法方程病态时的解算效果,但在扰动误差比较大时其求解效果仍不理想。

      3) 在加入不同程度的干扰后,Liu估计法和迭代估计法所得结果均与真值偏差不大,都是可靠的。干扰越小,偏差就越小,估值越接近于真值,但本文所提的迭代法在精度上略高于Liu估计法。

      4) 在该算例中,3种情况的迭代次数均为2次,这是因为Liu估计就是它的第一次迭代,精度已经很高,所以迭代完第二次后,MSE几乎不发生变化,故在这里只取了迭代两次的结果进行对比。

      算例2  取自文献[22]中第47页的例子。设P1, P2, …, P9为9个已知点,其坐标列于表 4。9个已知点到两个未知点P10, P11(假设模拟真值分别为(0, 0, 0) 和(7, 10, -5))的观测距离也列于表 4。两个未知点之间的观测距离为d10, 11=13.107 85 m。要求根据19个观测距离确定两个未知点的坐标,其中各距离为等精度观测,中误差为±0.000 1 m。

      表 4  控制点的坐标和观测距离

      Table 4.  Coordinate and Observation Distance of Control Points

      点号 坐标/m 观测距离/m
      x y z di, 10 di, 11
      P1 23.000 10.000 0.010 25.078 69 16.765 17
      P2 10.000 9.990 0.000 14.134 51 17.719 65
      P3 35.000 10.010 -0.010 36.415 88 28.442 94
      P4 100.000 19.990 0.005 101.479 43 93.168 39
      P5 -36.000 10.005 0.000 37.364 22 43.299 05
      P6 0.000 10.010 -0.005 10.010 04 8.600 60
      P7 56.000 9.995 0.010 56.996 06 49.256 18
      P8 -15.000 10.105 -0.010 18.035 90 22.559 66
      P9 -1.7000 10.008 -0.015 10.150 63 10.043 82

      该算例中法方程系数阵的条件数为cond(A′·A)=4.601 5×103,呈病态。分别用LS估计、岭估计、Liu估计和本文的迭代估计进行计算,比较计算的估值$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$、$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$的均方误差以及$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$与X差范数‖Δ$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$‖=‖ $\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$ -X ‖的大小,计算结果见表 5

      表 5  LS估计、岭估计、Liu估计以及本文方法的结果比较

      Table 5.  Comparison of the Results Among LS Estimator、Ridge Estimator、Liu Estimator and Our Method

      较值 真值 LS估计 岭估计 Liu估计 本文方法
      0 -0.036 8 -0.041 1 -0.038 0 -0.038 4
      0 0.051 8 0.027 5 0.022 7 0.022 0
      0 9.363 1 0.086 0 0.375 9 0.055 9
      $\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$ 7 7.048 8 6.974 6 6.975 9 6.973 3
      10 5.596 3 10.036 3 9.866 2 10.019 6
      -5 -4.648 2 -4.873 5 -4.904 4 -4.911 7
      参数 0 0 0.482 7 0.137 5 0.137 5
      ‖Δ$\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$‖ 0 10.353 3 0.166 7 0.214 2 0.118 3
      MSE($\mathit{\boldsymbol{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{X}}}$) 0 18.281 3 0.219 6 0.170 9 0.014 0

      表 5中的“参数”是指在不同估计方法中所需要确定的参数,如岭估计法中指的是岭参数λ,Liu估计法和本文所提的方法中指的是修正因子d

      表 5可知,岭估计、Liu估计和本文所提出的迭代方法所得结果与真值比较接近,均提高了LS解的精度。进一步对比估值的第三个分量,迭代估计法得到的值比Liu估计法得到的值更接近于真值,说明迭代估计法在病态问题求解中比Liu估计更可靠准确。此次实验中,在采用Liu估计的修正因子d后,迭代法所得的MSE和偏差均小于LS估计、岭估计和Liu估计,表明新方法对病态方程的求解具有一定的有效性,所得结果具有可信度。

    • 在GPS数据处理、工程控制网平差、形变观测分析、大地测量反演等测量数据处理领域,系统的病态性问题是常见的,并且病态性的危害作用非常严重。在病态问题的求解中,统计学家们提出了许多有偏估计方法。本文在Liu估计或Sakallioglu提出的有偏估计法的基础上,引入迭代的思想提出了一种新的有偏迭代估计法。若初值选为LS估值,则第一次迭代式就是Liu估计;若初值选为岭估计,则第一次迭代就是Sakallioglu提出的估计法。针对迭代公式,证明了修正因子在[-1, 1]范围内的收敛性。借助矩阵的谱分解给出迭代公式的解析形式,并给出修正因子的选取公式。通过两个不同情况例子的实验结果来看,迭代估计法优于LS估计和岭估计,求解精度略高于Liu估计。还需要说明的是,本文的两个例子均采用初值是LS估值。事实上不管初值选为LS估值还是岭估计值,最后迭代公式总是达到同一个精度。

      本文所提的有偏迭代估计法受修正因子d∈[-1, 1]的限制,不能使该方法所得到的估值达到更高的精度。故如何拓宽修正因子的取值范围以及每次迭代步骤中是否可以采用不同的最优的修正因子都是下一步需要研究的问题。

参考文献 (22)

目录

    /

    返回文章
    返回