常用穩(wěn)健估計(jì)方法應(yīng)用于概率積分法參數(shù)求取時(shí)的有效性研究

王曉輝，姜佃高

(太原理工大學(xué)測(cè)繪科學(xué)與技術(shù)系，山西太原030024)

一、引 言

開采沉陷預(yù)計(jì)可以定量地研究開采影響的地表移動(dòng)變形在時(shí)間和空間上的分布規(guī)律，對(duì)“三下”開采也具有重要指導(dǎo)作用。而開采沉陷預(yù)計(jì)方法中應(yīng)用最為廣泛的是概率積分法，其預(yù)計(jì)結(jié)果的精度主要取決于預(yù)計(jì)參數(shù)的可靠性［1］。因此，準(zhǔn)確確定概率積分法參數(shù)對(duì)開采沉陷理論和生產(chǎn)實(shí)踐都具有重要意義。

概率積分法的預(yù)測(cè)參數(shù)一般根據(jù)地表移動(dòng)觀測(cè)站實(shí)測(cè)資料通過(guò)最小二乘曲線擬合確定，最小二乘法有極好的配賦誤差的能力，但它易受粗差的影響，使參數(shù)估值失真［2-4］。當(dāng)觀測(cè)值中不可避免地存在粗差時(shí)［5-6］，穩(wěn)健估計(jì)求參技術(shù)具有有效抵御粗差或異值對(duì)參數(shù)干擾的特性，能夠克服最小二乘法擬合求參時(shí)常出現(xiàn)的結(jié)果發(fā)散問(wèn)題，保證求參結(jié)果的可靠性［7-8］。然而，不同穩(wěn)健估計(jì)方法的穩(wěn)健特性不同。Mitra等［9］研究表明，考慮到重尾噪聲分布，Huber法和Andrews法的穩(wěn)健性要優(yōu)于L1法。Li等［10］指出Fair法的穩(wěn)健性要優(yōu)于 L1法。Pennacchi［11］通過(guò)算例表明，在迭代 100次的情況下，Cauchy法比German-McClure法、Welsch法和Tukey法效果更好。方俊濤等［12］研究認(rèn)為，Welsch法的穩(wěn)健性最好，其次是Tukey法，Huber法的穩(wěn)健性最差。李惠芬等［13］分析指出，Andrews法、Tukey法和IGG方案均優(yōu)于Fair法。

本文在嘉樂(lè)泉煤礦實(shí)測(cè)資料基礎(chǔ)上進(jìn)行人工異值干擾求參試驗(yàn)，以觀測(cè)值中含有不同粗差值為例，比較了13種常用穩(wěn)健估計(jì)方法求取概率積分法參數(shù)時(shí)的穩(wěn)健特性。

二、概率積分法參數(shù)的穩(wěn)健估計(jì)

1.概率積分法參數(shù)的穩(wěn)健估計(jì)模型

概率積分法參數(shù)的穩(wěn)健估計(jì)模型詳細(xì)見文獻(xiàn)［14］。

2.改進(jìn)算法麥夸爾特法

用曲線擬合求取實(shí)測(cè)參數(shù)的基本算法是泰勒級(jí)數(shù)展開法，但其存在一個(gè)缺點(diǎn)，即當(dāng)各參數(shù)初始值選的偏離其真值過(guò)大時(shí)，各次求得的各參數(shù)值可能不逐漸趨近于真值(稱為迭代收斂)，而是越來(lái)越偏離真值(稱為迭代發(fā)散)。為了克服這個(gè)缺點(diǎn)，放寬對(duì)初始值的限制，采用其改進(jìn)算法，即麥夸爾特法。

首先選取參數(shù) bi(i=1，2，…，m)的初始值可從經(jīng)驗(yàn)值中選取，也可根據(jù)本礦區(qū)已有的實(shí)參數(shù)據(jù)求出，則有(Δi為第i個(gè)參數(shù)的修正值)

將f(Xk;B)附近按泰勒級(jí)數(shù)展開，并假設(shè)與bi足夠接近，使Δi足夠小，并可在展開式中略去Δi的二次及二次以上的項(xiàng)，則有

式中

根據(jù)最小二乘原理VΤV=min，則有

即

寫出m個(gè)未知數(shù)Δi的線性方程組，并進(jìn)行整理、移項(xiàng)，則上式可寫作

式中

上述中的d稱為阻尼因子，是一個(gè)≥0的可調(diào)整的常數(shù)，其作用是適當(dāng)選擇d可使每次迭代后求得的Q'值均較迭代前的值小。d值的選擇原則是:在迭代收斂時(shí)，d選取較小的值，以減少迭代次數(shù);在迭代后的Q'值大于迭代前的值時(shí)，將d值放大，重新進(jìn)行迭代，直至迭代后的Q'值比迭代前的值小為止，這樣能保證經(jīng)過(guò)一次迭代Q'就減少一次，bi的值就向其真值接近一步，從而大大放寬了對(duì)接近bi的要求，最終可求得各參數(shù)值bi。對(duì)應(yīng)于概率積分求參下的參數(shù)有下沉系數(shù)q，主要影響角的正切 tanβ，主要影響傳播角 θ0和點(diǎn)偏移距 s1、s2、s3、s4等參數(shù)。

三、常用穩(wěn)健估計(jì)方法和穩(wěn)健性比較

1.常用穩(wěn)健估計(jì)方法

13種常用的穩(wěn)健估計(jì)方法及它們的調(diào)和系數(shù)見文獻(xiàn)［15］。

2.兩種參數(shù)估計(jì)方法的比較

為了比較任意兩種參數(shù)估計(jì)方法在人工異值干擾試驗(yàn)下哪種方法更有效，給出參數(shù)的相對(duì)減益的概念，簡(jiǎn)述如下:

定義:在概率積分求參中，當(dāng)觀測(cè)值中不存在粗差時(shí)，通過(guò)最小二乘曲線擬合法求得的參數(shù)理論上是最優(yōu)值，而當(dāng)觀測(cè)值中存在粗差時(shí)通過(guò)各種穩(wěn)健估計(jì)方法求得同一參數(shù)的值與無(wú)粗差的最小二乘曲線擬合法求得的參數(shù)值會(huì)發(fā)生變化，二者之間的差值與無(wú)粗差的最小二乘曲線擬合法求得的參數(shù)比值的絕對(duì)值的百分?jǐn)?shù)就稱作參數(shù)的相對(duì)減益。當(dāng)用穩(wěn)健估計(jì)方法求得參數(shù)的相對(duì)減益越大，說(shuō)明這種方法越不穩(wěn)健;反之，其值越小則越穩(wěn)健。

在概率積分求參中選取3個(gè)主要參數(shù):下沉系數(shù)q、主要影響角的正切tanβ和拐點(diǎn)偏移距D=(H0為平均采深)。令下沉系數(shù)q的相對(duì)減益為

主要影響角的正切tanβ的相對(duì)減益為

拐點(diǎn)偏移距D的相對(duì)減益為

上式中，當(dāng)觀測(cè)值中不存在粗差時(shí)，通過(guò)最小二乘曲線擬合法求得的3個(gè)主要參數(shù)值為下沉系數(shù)q0、主要影響角的正切tanβ0和拐點(diǎn)偏移距D0;當(dāng)觀測(cè)值中存在粗差時(shí)，通過(guò)各種穩(wěn)健估計(jì)方法求得的3個(gè)主要參數(shù)值為下沉系數(shù)qi、主要影響角的正切tan βi和拐點(diǎn)偏移距 Di(i=1，2，3，…，13)，i為所選的第i種穩(wěn)健估計(jì)方法。

四、結(jié)果與討論

基于嘉樂(lè)泉煤礦某工作面地表移動(dòng)走向觀測(cè)線實(shí)測(cè)資料，通過(guò)人為地增加一些異值點(diǎn)的辦法，進(jìn)行人工干擾穩(wěn)健求參試驗(yàn)。該工作面采用走向長(zhǎng)壁采煤法采煤，全部陷落法管理頂板，上覆巖層巖性綜合評(píng)價(jià)為中硬，地表沉陷規(guī)律基本符合概率積分模型，地表線A是半無(wú)限開采，其實(shí)測(cè)下沉值見表1。本文采用將走向觀測(cè)線的拐點(diǎn)和最大下沉點(diǎn)處第一次加上100 mm粗差，第二次加上200 mm的粗差，經(jīng)計(jì)算得出了13種常用穩(wěn)健估計(jì)方法與無(wú)粗差下最小二乘法3個(gè)主要參數(shù)的相對(duì)減益，其結(jié)果如圖1—圖3所示。

表1 走向觀測(cè)線地表下沉值

圖1 不同穩(wěn)健估計(jì)方法相對(duì)于無(wú)粗差的最小二乘法下沉系數(shù)的相對(duì)減益

圖2 不同穩(wěn)健估計(jì)方法相對(duì)于無(wú)粗差的最小二乘法拐點(diǎn)偏移距的相對(duì)減益

圖3 不同穩(wěn)健估計(jì)方法相對(duì)于無(wú)粗差的最小二乘法的主要影響角正切的相對(duì)減益

由圖1可看出，當(dāng)粗差取100 mm時(shí)，L1法和German-McClure法下沉系數(shù)的相對(duì)減益均為0%，其次是IGGIII方案，其相對(duì)減益為1%，其余方法的結(jié)果都≥2%;當(dāng)粗差取200 mm時(shí)，L1法和IGGIII方案相對(duì)減益均為6%，而German-McClure法下沉系數(shù)的相對(duì)減益為7%，其余方法的結(jié)果都≥8%。

由圖2可看出，當(dāng)粗差取100 mm時(shí)，L1法和German-McClure法拐點(diǎn)偏移距的相對(duì)減益均為0%，而IGGIII方案的相對(duì)減益為1%，其余方法的結(jié)果都≥2%;當(dāng)粗差取200 mm時(shí)，L1法和IGGIII方案相對(duì)減益均為4%，而German-McClure法拐點(diǎn)偏移距的相對(duì)減益為5%，其余方法的結(jié)果都≥6%。總體上看各個(gè)穩(wěn)健估計(jì)方法的拐點(diǎn)偏距的相對(duì)減益相差不是很大，都在9%以內(nèi)。

由圖3可看出，當(dāng)粗差取100 mm時(shí)，L1法、L1-L2法和German-McClure法主要影響角正切的相對(duì)減益均為 0%，而是 Huber法、Andrews法、Welsch法、Tukey法、Danish法、IGG 方案、IGGIII方案和Cauchy法主要影響正切的相對(duì)減益為1%，各方法的結(jié)果都特別接近，都不大于3%;當(dāng)粗差取200 mm時(shí)，L1法、German-McClure法和IGGIII方案相對(duì)減益均為4%，其余方法的結(jié)果都≥6%。

五、結(jié)束語(yǔ)

綜合圖1、圖2和圖3的試驗(yàn)結(jié)果，可知當(dāng)給走向觀測(cè)線的拐點(diǎn)和最大下沉點(diǎn)處同時(shí)加上100 mm的粗差時(shí)，對(duì)于各穩(wěn)健估計(jì)方法3個(gè)主要參數(shù)的相對(duì)減益相差都不是很大，其中 L1法、German-McClure法和IGGIII方案相對(duì)更穩(wěn)健;當(dāng)給走向觀測(cè)線的拐點(diǎn)和最大下沉點(diǎn)處同時(shí)加上200 mm的粗差時(shí)，各穩(wěn)健估計(jì)主要參數(shù)相對(duì)減益值差異逐漸增大，而 L1法、German-McClure法和IGGIII方案的3個(gè)主要參數(shù)相對(duì)減益值比其他穩(wěn)健估計(jì)方法的都小。對(duì)于概率積分法參數(shù)的穩(wěn)健估計(jì)方法，L1、German-McClure方法和IGGIII方案較其他穩(wěn)健估計(jì)方法相對(duì)更穩(wěn)健。

［1］何國(guó)清，楊倫，凌賡娣，等.礦山開采沉陷學(xué)［M］.徐州:中國(guó)礦業(yè)大學(xué)出版社，1994.

［2］王素立，劉永.基于波動(dòng)相關(guān)性及主分量變換的多元線性回歸模型研究［J］.統(tǒng)計(jì)與決策，2012(22):18-21.

［3］高山，李孝軍.關(guān)于最小二乘法的穩(wěn)健性分析［J］.統(tǒng)計(jì)與決策，2006(8):125-126.

［4］何霞，劉衛(wèi)鋒.基于全最小一乘準(zhǔn)則的灰色GM(1，1)模型參數(shù)估計(jì)［J］.統(tǒng)計(jì)與決策，2006(8):30-33.

［5］GRIEP M I，WAKELING I N，VANKEERBERGHEN P，et al.Comparison of Semi-robust and Robust Partial Least Squares Procedures［J］.Chemometrics and Intelligent Laboratory Systems，1995，29(1):37-50.

［6］GUO JF，OU J K，WANG H T.Robust Estimation for Correlated Observations:Two Local Sensitivity-based Down Weighting Strategies ［J］.Journal of Geodesy，2010，84(4):243-250.

［7］吳健平，張立.地理數(shù)據(jù)線性回歸中的穩(wěn)健估計(jì)方法［J］.干旱區(qū)地理，1994，17(1):83-88.

［8］馮守平，石澤，鄒瑾.一元線性回歸模型中參數(shù)估計(jì)的幾種方法比較［J］.統(tǒng)計(jì)與決策，2008(24):152-153.

［9］MITRA S，MITRA A，KUNDU D.Genetic Algorithm and M-Estimator Based Robust Sequential Estimation of Parameters of Nonlinear Sinusoidal Signals［J］.Communications in Nonlinear Science and Numerical Simulation，2010，16(7):2796-2809.

［10］LI H F，JIANG X Q，LI Z.Robust Estimation in Gaussian Filtering for Engineering Surface Characterization［J］.Precision Engineering，2004，28(2):186-193.

［11］PENNACCHI P.Robust Estimate of Excitations in Mechanical Systems Using M-Estimators-Theoretical Background and Numerical Applications［J］.Journal of Sound and Vibration，2008，310(4-5):923-946.

［12］方俊濤，何楨，宋琳曦.響應(yīng)曲面建模的穩(wěn)健M-回歸方法［J］.工業(yè)工程，2012，15(3):98-103.

［13］李惠芬，蔣向前，李柱.高斯濾波穩(wěn)健性能的研究與改進(jìn)［J］.儀器儀表學(xué)報(bào)，2004，25(5):633-637.

［14］郭廣云，汪云甲.概率積分法參數(shù)的穩(wěn)健估計(jì)模型及其應(yīng)用研究［J］.測(cè)繪學(xué)報(bào)，2000，29(2):162-165.

［15］GE Y H，YUAN Y，JIA N N.More Efficient Methods among Commonly Used Robust Estimation Methods for GPS Coordinate Transformation ［J］.Survey Review，2013，45(330):229-234.