王吉?jiǎng)P

（中國(guó)人民解放軍61175 部隊(duì)，江蘇 南京 210000）

0 引言

數(shù)字高程模型是在二維空間上對(duì)三維地形表面的客觀表達(dá)，它在二維平面上對(duì)研究區(qū)域進(jìn)行網(wǎng)格劃分，形成格網(wǎng)空間架構(gòu)，以覆蓋整個(gè)區(qū)域。地理學(xué)分析中的數(shù)字高程模型作為等高線地形圖的替代產(chǎn)品，其是科學(xué)研究、經(jīng)濟(jì)建設(shè)和國(guó)防建設(shè)的基礎(chǔ)數(shù)據(jù)和有力工具，也是地理學(xué)分析、生物學(xué)等區(qū)域性科學(xué)基本參數(shù)的提供者和科學(xué)研究成果的表現(xiàn)形式。同時(shí)，各種不同分辨率的數(shù)字高程模型是數(shù)字繪圖和綜合必需數(shù)據(jù)[1]。在數(shù)字高程模型建立過(guò)程中的重要內(nèi)容是格網(wǎng)點(diǎn)的高程內(nèi)插計(jì)算，內(nèi)插計(jì)算是構(gòu)建數(shù)字高程模型的核心內(nèi)容，貫穿了數(shù)字高距離模型的生產(chǎn)、品質(zhì)控制、精度評(píng)價(jià)以及分析應(yīng)用的所有環(huán)節(jié)[2]。在DEM 數(shù)據(jù)內(nèi)插研究過(guò)程中，研究人員不僅對(duì)各種傳統(tǒng)方法在地形表面建模中的可行性、特點(diǎn)進(jìn)行了充分而深入的研究，而且也引進(jìn)相關(guān)學(xué)科的內(nèi)插方法（例如克里金、有限元等內(nèi)插技術(shù)），同時(shí)根據(jù)地形表面的特點(diǎn)發(fā)展地形數(shù)據(jù)的各種內(nèi)插方法。

該文利用實(shí)際測(cè)量的GPS 數(shù)據(jù)，運(yùn)用廣義回歸神經(jīng)網(wǎng)絡(luò)算法、反距離權(quán)法、普通克里格法以及多項(xiàng)式插值法對(duì)地形高程插值進(jìn)行研究。

1 地表高程插值

1.1 基本概念

廣義回歸神經(jīng)網(wǎng)絡(luò)是基于徑向基網(wǎng)絡(luò)的一種演變算法。該算法基于非參數(shù)回歸算法，使用樣本集充當(dāng)驗(yàn)證數(shù)據(jù)，運(yùn)用Parzen 非參數(shù)估計(jì)，以最大概率原則輸出廣義回歸神經(jīng)網(wǎng)絡(luò)的運(yùn)算結(jié)果。非線性逼近能力十分優(yōu)異，與徑向基相比，訓(xùn)練更方便，廣義回歸神經(jīng)網(wǎng)絡(luò)非常適合解決曲線擬合的問(wèn)題。

x、y為隨機(jī)變量，其聯(lián)合概率密度為f（x，y），假設(shè)x的觀測(cè)值為x0，對(duì)y進(jìn)行求導(dǎo)，那么y對(duì)x的回歸E（y|x0）如公式（1）所示。

那么在輸入為x0的情況下，y的輸出預(yù)測(cè)即為y（x0）。執(zhí)行Parzen 非參數(shù)估計(jì)，得到樣本數(shù)據(jù)集{xi，yi}，根據(jù)公式（2）估計(jì)密度函數(shù)f（x0，y）。

式中：n為樣本容量；p為隨機(jī)變量x的維數(shù)；σ為光滑因子，也就是高斯函數(shù)的標(biāo)準(zhǔn)差。

密度函數(shù)f（x0，y）公式中導(dǎo)數(shù)d（x0，xi）和d（y，yi）分別為x數(shù)據(jù)集xi的均方根誤差和y數(shù)據(jù)集yi的方差，如公式（3）所示。

將公式（3）的d（x0，xi）代入密度函數(shù)，并交換求和與積分順序，如公式（4）所示。

式中：分子為所有訓(xùn)練樣本算得的yi值的加權(quán)和；e-d（x0，xi）為權(quán)值。

廣義回歸神經(jīng)網(wǎng)絡(luò)不需要訓(xùn)練，但是平滑因子的值對(duì)網(wǎng)絡(luò)性能有很大的影響，并且需要優(yōu)化值。如果光滑因子取值非常大，那么d（x0，xi）趨近于0，y（x0）近似于所有樣本因變量的平均值。如果光滑因子趨近于0，那么y（x0）與訓(xùn)練樣本的值非常接近，當(dāng)需要預(yù)測(cè)的點(diǎn)在訓(xùn)練樣本時(shí)，算得的預(yù)測(cè)值與樣本中的期望輸出非常接近，但是一旦給定新的輸入，預(yù)測(cè)的結(jié)果就會(huì)急劇變差，使網(wǎng)絡(luò)失去推廣能力，這種現(xiàn)象就是過(guò)學(xué)習(xí)[3]。

廣義回歸神經(jīng)網(wǎng)絡(luò)由4 層構(gòu)成，即輸入層、模式層、求和層和輸出層。設(shè)輸入層為X=[x1，x2，...，xn]T，其輸出層為Y=[y1，y2，...，yn]T。

1.1.1 輸入層

輸入層神經(jīng)元的數(shù)量等于學(xué)習(xí)樣本中輸入向量的維數(shù)，每個(gè)神經(jīng)元是簡(jiǎn)單的分布單位，并且將輸入變量直接傳遞到模式層[4]。

1.1.2 模式層

模式層神經(jīng)元的數(shù)目等于學(xué)習(xí)樣本的數(shù)目n，每個(gè)神經(jīng)元對(duì)應(yīng)不同的樣本，而模式層神經(jīng)元傳遞函數(shù)如公式（6）所示。

式中：σ為樣本標(biāo)準(zhǔn)差；X為網(wǎng)絡(luò)輸入變量；Xi為第i個(gè)神經(jīng)元對(duì)應(yīng)的學(xué)習(xí)樣本[4]。

神經(jīng)元i的輸出為輸入變量與其對(duì)應(yīng)的樣本X之間Eucd距離平方的指數(shù)平方Di2=（X-Xi）T（X-Xi）的指數(shù)形式。

1.1.3 求和層

在求和層中，使用2 種神經(jīng)元進(jìn)行求和。

第一種求和如公式（7）所示。

它對(duì)所有模式層神經(jīng)元的輸出進(jìn)行算術(shù)求和，其模式層與各神經(jīng)元的連接權(quán)值為1，傳遞函數(shù)如公式（8）所示[4]。

式中：SD為傳遞函數(shù)。

第二種求和如公式（9）所示。

它對(duì)所有模式層的神經(jīng)元進(jìn)行加權(quán)求和，模式層中第i個(gè)神經(jīng)元與求和層中第j個(gè)分子求和神經(jīng)元之間的連接權(quán)值為第i個(gè)輸出樣本Yi中的第j個(gè)元素，傳遞函數(shù)如公式（10）所示[4]。

式中：k為無(wú)限值（因?yàn)槭街写嬖趎，所以用k代指）。

1.1.4 輸出層

輸出層中的神經(jīng)元數(shù)目等于學(xué)習(xí)樣本中輸出向量的維數(shù)k，各神經(jīng)元與求和層的輸出相除，神經(jīng)元j的輸出對(duì)應(yīng)估計(jì)結(jié)果（X）的第j個(gè)元素，如公式（11）所示[4]。

1.2 數(shù)據(jù)處理

使用MATLAB 編程實(shí)現(xiàn)廣義回歸神經(jīng)網(wǎng)絡(luò)對(duì)高程進(jìn)行預(yù)測(cè)的功能，輸入歸一化后的信號(hào)X、Y通過(guò)中間節(jié)點(diǎn)（隱層點(diǎn)）作用于輸出節(jié)點(diǎn)，經(jīng)過(guò)非線形變換產(chǎn)生輸出信號(hào)Z，并反歸一化X、Y值[5]。廣義回歸神經(jīng)網(wǎng)絡(luò)的所有訓(xùn)練樣本為輸入值X、輸入值Y、期望值Z0以及輸出值Z與期望值Z0的偏差。調(diào)節(jié)隱層與輸入節(jié)點(diǎn)的連接強(qiáng)度、隱層與輸出節(jié)點(diǎn)的連接強(qiáng)度以及閾值，使偏差沿梯度下降，通過(guò)不斷訓(xùn)練學(xué)習(xí)，最終確定與最小誤差相適應(yīng)的權(quán)值和閾值。經(jīng)過(guò)學(xué)習(xí)訓(xùn)練的神經(jīng)網(wǎng)絡(luò)算法可以自動(dòng)對(duì)相似的樣本信息進(jìn)行計(jì)算，輸出擬合程度最高的結(jié)果。

該程序設(shè)置學(xué)習(xí)率為0.01，訓(xùn)練精度要求為0.000 1，最大訓(xùn)練次數(shù)為10 000 次，網(wǎng)絡(luò)輸入向量的最大值和最小值為歸一化后X、Y的最大值和最小值，設(shè)置每層神經(jīng)元個(gè)數(shù)，共設(shè)3 層神經(jīng)元，每層個(gè)數(shù)分別為15、40 和1。誤差反向傳播神經(jīng)網(wǎng)絡(luò)算法的隱含層的傳遞函數(shù)是“S”形函數(shù)，輸出層是線性函數(shù)。前兩層使用雙正切“S”形函數(shù)傳輸函數(shù)，“S”形函數(shù)（Sigmoid function）是誤差反向傳播神經(jīng)網(wǎng)絡(luò)算法中常用的非線性作用函數(shù)，由于誤差反向傳播神經(jīng)網(wǎng)絡(luò)算法的傳遞函數(shù)必須可微，因此一般使用“S”形函數(shù)或者線性函數(shù)作為傳遞函數(shù)。第三層為線性傳輸函數(shù)。

選定Levenberg-Marquardt 訓(xùn)練函數(shù)，萊文拜格－馬夸特方法（Levenberg–Marquardt algorithm）能計(jì)算數(shù)值非線性最小化的解。該方法可以通過(guò)在計(jì)算時(shí)修改參數(shù)來(lái)實(shí)現(xiàn)梯度下降法和高斯算法的優(yōu)點(diǎn)，并優(yōu)化這2 種方法的缺點(diǎn)。核心思想就是用雅可比矩陣代替H 矩陣的計(jì)算，從而提高優(yōu)化效率。萊文拜格－馬夸特方法的優(yōu)點(diǎn)是可以優(yōu)化調(diào)節(jié)，一旦下降過(guò)快，就可以改變參數(shù)，使用小的阻尼系數(shù)λ，使該算法與高斯法類似，一旦下降過(guò)慢，就可以改變參數(shù)，使用更大的λ，使算法與梯度下降法類似。

使用該神經(jīng)網(wǎng)絡(luò)輸出預(yù)期值Z與實(shí)際值Z0進(jìn)行比較，得出預(yù)測(cè)精度，并使用神經(jīng)網(wǎng)絡(luò)得出該地區(qū)所有高程值并可視化，得到的平面地形圖和三維地形圖如圖1 和圖2 所示。

圖1 帶點(diǎn)地形圖

圖2 三維地形圖

2 預(yù)測(cè)精度評(píng)價(jià)

該文采用廣義回歸神經(jīng)網(wǎng)絡(luò)法、反距離權(quán)重法、克里格法以及多項(xiàng)式法的地表高程插值得到預(yù)測(cè)結(jié)果，并進(jìn)行精度驗(yàn)證和評(píng)價(jià)分析，得出預(yù)測(cè)值與驗(yàn)證數(shù)據(jù)高程值對(duì)比，如圖3 所示。

由圖3 可知，數(shù)據(jù)點(diǎn)由低到高的離散程度逐漸變大，廣義回歸網(wǎng)絡(luò)方法在整個(gè)高程范圍內(nèi)的數(shù)據(jù)點(diǎn)距擬合線平均距離最近，近似呈線性分布，擬合程度最高，其他方法隨著高程的增加逐漸變高，數(shù)據(jù)點(diǎn)離散程度提高，擬合程度逐漸降低。綜上所述，GRNN 方法擬合程度最高，其插值數(shù)據(jù)逐漸接近驗(yàn)證數(shù)據(jù)。

圖3 離散回歸點(diǎn)線圖

該文用均方根誤差（RMSE）、方差（SSE）以及確定系數(shù)（R2）來(lái)評(píng)價(jià)預(yù)測(cè)精度，見(jiàn)表1。

表1 精度指數(shù)

對(duì)比發(fā)現(xiàn)，在4 種方法中，廣義回歸神經(jīng)網(wǎng)絡(luò)算法的均方根誤差、方差最小，確定系數(shù)最接近1。說(shuō)明在樣本數(shù)據(jù)的數(shù)目和分布狀況均相同的條件下，廣義回歸神經(jīng)網(wǎng)絡(luò)算法的插值精度最高，且使用廣義回歸神經(jīng)網(wǎng)絡(luò)算法得到的均方根誤差的值達(dá)到比較理想的值（5）附近。綜上所述，廣義回歸神經(jīng)網(wǎng)絡(luò)算法是4 種方法中最適合該地區(qū)進(jìn)行地形擬合的插值算法。

3 結(jié)語(yǔ)

該文以陜西省延安市安塞地區(qū)GNSS 實(shí)測(cè)高程數(shù)據(jù)為基礎(chǔ)，分別采用GRNN、O-Kriging、IDW 以及Polynome 方法對(duì)樣本點(diǎn)進(jìn)行高程插值及精度驗(yàn)證。結(jié)果表明，在樣本數(shù)據(jù)一定的情況下，GRNN 方法的插值精度明顯高于其他方法，是最適合該地形區(qū)域的插值算法。試驗(yàn)證明該高程插值方法能夠較準(zhǔn)確地插值未測(cè)量地區(qū)的高程，但是不能滿足對(duì)工程施工的精度要求，目前只能作為非精密測(cè)量區(qū)域的地形參考。且各種地表高程插值方法適用的地形特征不同，在實(shí)際使用過(guò)程中可能需要嘗試多種方法，以比較預(yù)測(cè)精度。