徐志敏，王利雙

(1.長(zhǎng)江勘測(cè)規(guī)劃設(shè)計(jì)研究有限責(zé)任公司，武漢 430010；2.長(zhǎng)江空間信息技術(shù)工程有限公司(武漢)，武漢 430010；3.湖北省水利信息感知與大數(shù)據(jù)工程技術(shù)研究中心，武漢 430010；4.國(guó)家知識(shí)產(chǎn)權(quán)局專利局專利審查協(xié)作湖北中心，武漢 430205)

0 引言

數(shù)字高程模型(digital elevation model，DEM)是國(guó)家基礎(chǔ)地理數(shù)據(jù)的重要組成部分，也是進(jìn)行三維空間分析和應(yīng)用的基礎(chǔ)。對(duì)DEM進(jìn)行全面的精度分析和質(zhì)量評(píng)價(jià)是保證其正確應(yīng)用的基本前提[1]。目前主要的DEM精度評(píng)價(jià)模型很多，其中高程中誤差模型[2]已經(jīng)成為國(guó)家測(cè)繪行業(yè)的采納標(biāo)準(zhǔn)(CH/T 1015.2—2007、CH/T 9009.2—2010)。然而研究表明，高程中誤差只能評(píng)價(jià)DEM采樣點(diǎn)的高程誤差，但無(wú)法評(píng)價(jià)DEM模擬表面對(duì)真實(shí)地表的逼近誤差，即使高程中誤差為0，DEM模擬表面也不一定能真實(shí)表達(dá)地表起伏。也正因此，學(xué)界提出了多種DEM誤差的評(píng)價(jià)模型，旨在對(duì)DEM誤差進(jìn)行更加科學(xué)的評(píng)價(jià)，包括重構(gòu)等高線模型[3]、地形信息熵模型[4]、Strahler積分模型[5]、填挖方模型[6]與地形描述精度模型[7]等。

其中，DEM地形描述精度模型由于其科學(xué)性和實(shí)用性，引起了學(xué)界的關(guān)注和研究。湯國(guó)安等[8]最早提出了地形描述誤差的概念和提取方法，并對(duì)地形描述誤差的空間結(jié)構(gòu)和定量計(jì)算模型做了研究。張勇等[9]以黃土丘陵溝壑區(qū)為例，采用統(tǒng)計(jì)和比較分析的方法，研究了不同空間尺度下的DEM地形描述誤差的數(shù)學(xué)轉(zhuǎn)換模型。劉春等[10]基于地形復(fù)雜因子概念，利用線性回歸分析方法探討了柵格DEM地形描述精度與地形復(fù)雜因子和平均高程之間的線性關(guān)系。賈敦新等[11]以模擬實(shí)驗(yàn)獲得的高精度DEM為數(shù)據(jù)源，分析坡度計(jì)算精度與DEM數(shù)據(jù)誤差和地形描述誤差的耦合關(guān)系。齊曉飛等[12]提出了一種基于散點(diǎn)圖矩陣建立多元回歸模型的方法，并依此構(gòu)建了DEM地形描述誤差與分辨率、相對(duì)高差、凹凸系數(shù)的關(guān)系公式。王春等[13-14]借鑒誤差橢圓和ε誤差帶建模的思想，提出了地形描述誤差的中誤差場(chǎng)和極值場(chǎng)，并研究了這些場(chǎng)的基本特性?？梢?jiàn)，DEM地形描述精度的定量模型一直是學(xué)界關(guān)注的熱點(diǎn)。

本文對(duì)地形描述誤差(以下簡(jiǎn)稱Et)與DEM空間分辨率和平均剖面曲率之間的定量關(guān)系進(jìn)行了重新研究。從地形描述精度的概念出發(fā)，分析了現(xiàn)有的地形描述精度定量模型的缺陷，提出了優(yōu)化的冪函數(shù)定量模型，并通過(guò)實(shí)驗(yàn)對(duì)比了本文所提出定量模型的擬合精度和泛化精度的優(yōu)勢(shì)。

1 現(xiàn)有Et定量模型的問(wèn)題

湯國(guó)安等首次提出了地形描述誤差的概念和計(jì)算方法，并以地形描述誤差概念和計(jì)算方法為基礎(chǔ)，選取了六種不同地貌類(lèi)型(平原、低丘、丘陵、中山、高山以及混合類(lèi)型)的DEM作為實(shí)驗(yàn)樣區(qū)，分別計(jì)算了各樣區(qū)在不同DEM空間分辨率下的地形描述誤差的均方差。首先，對(duì)各樣區(qū)的DEM地形描述誤差均方差和DEM空間分辨率之間的關(guān)系采用y=ax+b線性模型進(jìn)行了線性回歸；然后，對(duì)方程系數(shù)b與對(duì)應(yīng)樣區(qū)的平均剖面曲率V分別再次進(jìn)行線性回歸；最后，得到DEM地形描述誤差均方差與空間分辨率R和平均剖面曲率V的定量計(jì)算公式，如式(1)所示。

RMSEt=(0.006 3V+0.006 6)R-0.022V+0.241 5

(1)

對(duì)此，王光霞等[15]指出：上述研究中系數(shù)b與平均剖面曲率的線性相關(guān)性并不好，利用直線擬合勢(shì)必產(chǎn)生較大的誤差。因此，對(duì)系數(shù)b與平均剖面曲率采用了二次函數(shù)模型進(jìn)行了擬合，將參數(shù)a與平均剖面曲率的相關(guān)系數(shù)從0.989 9提高到0.998 9，將系數(shù)b與平均剖面曲率的相關(guān)系數(shù)從0.842 5提高到0.988 2，最終得到擬合精度更高的定量模型，計(jì)算方法如式(2)所示。該方法得到的定量模型與式(1)相比較具有更好的擬合效果，將在地勢(shì)較平緩的地貌類(lèi)型(如平原、低丘和丘陵地區(qū))上的擬合中誤差提高了一個(gè)數(shù)量級(jí)。

RMSEt=(0.006 1V+0.002 7R)+0.001V2-0.064 9V+0.569 5

(2)

事實(shí)上，根據(jù)湯國(guó)安等對(duì)地形描述誤差的定義，DEM地形描述誤差Et是假定DEM高程采樣誤差為0的條件下，模擬地面與實(shí)際地面之間的差異。

如圖1所示，弧面AB代表實(shí)際地形表面，面AB代表離散采樣點(diǎn)所表達(dá)的模擬表面，AB點(diǎn)在水平面上的投影點(diǎn)之間的距離為DEM的分辨率。在點(diǎn)C處的DEM模擬地面高程與實(shí)際地面高程差為C點(diǎn)的DEM地形描述誤差Etc。實(shí)際計(jì)算中，實(shí)際地面高程值一般采用4鄰域的DEM柵格點(diǎn)高程的平均值進(jìn)行替代。

圖1 地形描述誤差示意圖

那么，根據(jù)DEM地形描述誤差的定義和計(jì)算方法，在理想情況下，當(dāng)DEM平均剖面曲率趨向于0時(shí)，DEM模擬地面為近似平面，任意柵格點(diǎn)四鄰域的高程平均值與其本身高程差為0，地形描述誤差為0；而當(dāng)DEM分辨率趨向于0時(shí)，DEM模擬地面能夠完全表達(dá)實(shí)際地面，地形描述誤差也為0。因此，對(duì)于地形描述誤差Et與DEM空間分辨率R和平均剖面曲率V之間的定量模型RMSEt=f(R，V)應(yīng)該滿足式(3)。

(3)

然而，式(1)和式(2)均為帶常數(shù)項(xiàng)的多項(xiàng)式模型，在DEM分辨率R和平均剖面曲率V分別趨向于0時(shí)的極限均為帶常數(shù)的不定式，因而均不滿足式(3)。

圖2顯示了當(dāng)平均剖面曲率為0，即地形絕對(duì)平坦的情況下，式(1)、式(2)對(duì)應(yīng)模型的地形描述精度與DEM分辨率之間的函數(shù)關(guān)系。根據(jù)式(3)，當(dāng)?shù)匦谓^對(duì)平坦時(shí)，DEM地形描述誤差與DEM分辨率不相關(guān)，且恒為0。但式(1)、式(2)在平均剖面曲率為0時(shí)，地形描述誤差與DEM分辨率為線性關(guān)系，隨著地形分辨率的提升，地形描述誤差逐漸下降，且恒為某一大于0的值。

圖2 V為0時(shí)式(1)、式(2)中RMS Et與R的關(guān)系

這顯然于地形描述誤差的固有性質(zhì)不符。因此，式(1)、式(2)對(duì)應(yīng)的模型均欠缺一定的科學(xué)性。事實(shí)上，文獻(xiàn)[9]和文獻(xiàn)[11]研究所得的定量模型也存在類(lèi)似問(wèn)題。這些模型在對(duì)DEM誤差進(jìn)行定量建模的時(shí)候都沒(méi)有顧及到DEM誤差與DEM地形特征之間的內(nèi)在關(guān)系。因此，地形描述誤差與DEM空間分辨率和平均剖面曲率之間的定量模型需要基于式(3)這一前提進(jìn)行重新設(shè)計(jì)和研究。

2 基于冪函數(shù)的Et定量模型與應(yīng)用

考慮到上述DEM地形描述誤差的重要性質(zhì)，本文基于式(3)的前提，設(shè)計(jì)了基于冪函數(shù)的地形描述誤差Et與DEM空間分辨率R和平均剖面曲率V的定量關(guān)系模型，如式(4)所示。

RMSEt=aRbVc

(4)

式中：a、b、c為模型系數(shù)。

將(R，V，RMSEt)看作是三維空間中的離散點(diǎn)，將式(4)看作空間曲面方程。模型求解過(guò)程即可看作三維離散點(diǎn)的曲面擬合。由于式(4)與自變量的關(guān)系是非線性的，不能用線性最小二乘法求多元函數(shù)極值的方法獲得參數(shù)估計(jì)值，因此需要采取復(fù)雜的優(yōu)化算法進(jìn)行求解。常用的優(yōu)化求解算法包括最速下降算法、高斯-牛頓法和列文伯格-馬夸爾特算法。本文采用基于高斯-牛頓迭代的非線性最小二乘法，對(duì)離散三維點(diǎn)進(jìn)行曲面擬合。將三維離散點(diǎn)和三維曲面的擬合效果在MATLAB中顯示，結(jié)果如圖3所示。

圖3 曲面擬合結(jié)果

擬合得到的模型系數(shù)a、b、c的值分別為：0.005 713、1.056和0.943 1。因此，基于冪函數(shù)的地形描述誤差與DEM空間分辨率和平均剖面曲率之間的定量模型為式(5)所示。

RMSEt=0.005 713R1.056V0.943 1

(5)

進(jìn)一步地，計(jì)算并比較式(1)、式(2)和式(5)對(duì)樣本數(shù)據(jù)的擬合優(yōu)度相關(guān)指標(biāo)，包括誤差平方和SSE、確定系數(shù)R-square以及標(biāo)準(zhǔn)差RMSE，結(jié)果如表1所示?？梢?jiàn)，本文所提出定量模型的擬合效果很好，擬合精度與式(2)基本相當(dāng)，且各項(xiàng)指標(biāo)明顯優(yōu)于式(1)。

表1 擬合優(yōu)度

對(duì)式(5)進(jìn)行改寫(xiě)，可得DEM分辨率是地形描述誤差和平均剖面曲率的函數(shù)，即式(6)，以及平均剖面曲率是地形描述誤差和DEM分辨率的函數(shù)，即式(7)。DEM地形描述誤差定量模型在生產(chǎn)實(shí)踐中的實(shí)用價(jià)值在于：基于給定的DEM的限差，式(5)可以直接根據(jù)分辨率和剖面曲率檢驗(yàn)DEM是否符合精度要求；式(6)可以推算一定地貌類(lèi)型范圍內(nèi)符合精度要求的最低DEM分辨率；式(7)可以推算一定分辨率的DEM產(chǎn)品中符合精度要求的地貌類(lèi)型區(qū)域。

(6)

(7)

3 泛化精度對(duì)比實(shí)驗(yàn)與分析

擬合精度評(píng)價(jià)的是定量模型公式對(duì)已知的實(shí)驗(yàn)樣本數(shù)據(jù)的逼近程度，但是模型的擬合精度并不能保證模型在生產(chǎn)實(shí)踐中的實(shí)用性，還需要評(píng)價(jià)模型對(duì)未知樣本數(shù)據(jù)的估算效果，即泛化精度。

本文分別選取了平原、丘陵和中山地帶的0.5 m高分辨率DEM進(jìn)行了實(shí)驗(yàn)。數(shù)據(jù)來(lái)源于美國(guó)圣地亞哥超級(jí)計(jì)算機(jī)中心的開(kāi)放地形數(shù)據(jù)OpenTopography網(wǎng)站，這個(gè)網(wǎng)站開(kāi)放了高分辨率地形數(shù)據(jù)及處理工具(http：//opentopo.sdsc.edu/lidar?format=sd)。數(shù)據(jù)區(qū)域?yàn)槊绹?guó)加尼福尼亞州Napa Valley，地理范圍為1×104m2，高程分布在海拔-1～634 m之間，從西南到東北方向依次分布著平原丘陵和山地，包含了比較豐富的地形種類(lèi)。本文選取了該DEM范圍內(nèi)的平原、丘陵和中山的三塊地形作為實(shí)驗(yàn)樣區(qū)，源數(shù)據(jù)和實(shí)驗(yàn)樣區(qū)的分布如圖4所示。

圖4 源數(shù)據(jù)與實(shí)驗(yàn)樣區(qū)

采用文獻(xiàn)[7]關(guān)于地形描述誤差的計(jì)算方法，基于0.5 m分辨率的原始DEM，將地形描述精度計(jì)算窗口依次按照3×3、5×5、7×7擴(kuò)大，分別計(jì)算得到3種地形在分辨率為1 m、2 m、4 m三種情況下的地形描述誤差，總計(jì)九組數(shù)據(jù)。同時(shí)，采用式(1)、式(2)、式(5)分別計(jì)算九組不同DEM分辨率和平均剖面曲率情況下的地形描述誤差估算值，結(jié)果如表2所示。將地形描述精度的計(jì)算值與式(1)、式(2)和式(5)的估算值進(jìn)行統(tǒng)計(jì)比較?？梢钥闯?，式(5)的定量模型對(duì)于高分辨率DEM的地形描述誤差估算精度都在83.3%以上，同時(shí)隨著DEM分辨率的提高，模型的估算精度也呈現(xiàn)出增加的趨勢(shì)，可見(jiàn)本文的定量模型對(duì)于高分辨率DEM具有很好的泛化性能。而對(duì)于式(1)和式(2)，其定量模型的估算精度不太穩(wěn)定，式(1)的估算精度在45.5%和84.7%之間波動(dòng)，式(2)的估算精度在49.4%和98.9%之間波動(dòng)，甚至隨著DEM分辨率的提高，一些樣區(qū)地形描述誤差的均方差估算值產(chǎn)生了不符合邏輯的負(fù)值。

表2 泛化精度對(duì)比結(jié)果

因此，實(shí)驗(yàn)結(jié)果印證了本文指出的式(1)和式(2)的問(wèn)題：該定量模型在設(shè)計(jì)的過(guò)程中只追求對(duì)現(xiàn)有樣本數(shù)據(jù)的擬合精度，沒(méi)有考慮到地形描述精度的固有性質(zhì)并根據(jù)對(duì)應(yīng)自變量和因變量的固有數(shù)學(xué)關(guān)系來(lái)選擇合適的擬合模型，導(dǎo)致模型對(duì)原有樣本數(shù)據(jù)的擬合精度較好，但是對(duì)于更高分辨率DEM樣本的泛化性能不理想。

4 結(jié)束語(yǔ)

本文對(duì)地形描述精度與DEM分辨率和地形平均剖面曲率之間的定量關(guān)系模型進(jìn)行了研究。首先，從DEM地形描述精度的概念出發(fā)，分析了現(xiàn)有的帶常數(shù)項(xiàng)的多項(xiàng)式模型的缺陷；然后，提出了基于冪函數(shù)的定量模型，計(jì)算并對(duì)比了模型對(duì)樣本數(shù)據(jù)擬合精度；最后，基于 4 m、2 m和1 m的高分辨率DEM，比較了模型在新的樣區(qū)的泛化精度。分析結(jié)果表明，本文的冪函數(shù)模型在保證良好的擬合精度的同時(shí)，對(duì)高分辨率的DEM具有較好的泛化性能。

隨著測(cè)繪遙感技術(shù)的發(fā)展和計(jì)算機(jī)軟硬件性能的提升，高分辨率DEM生產(chǎn)和應(yīng)用由于能夠滿足用戶對(duì)三維地形精度和體驗(yàn)的要求，得到越來(lái)越廣泛的應(yīng)用。在這種趨勢(shì)下，本文的定量模型對(duì)于在生產(chǎn)實(shí)踐中，根據(jù)限差要求推算適宜的DEM分辨率具有更好的指導(dǎo)意義和實(shí)用價(jià)值。

然而，還需注意到，雖然本文基于冪函數(shù)的定量模型具有較好的擬合精度與泛化性能，但是滿足“當(dāng)DEM空間分辨率和平均剖面曲率各自趨向于0，地形描述誤差亦趨向于0這一客觀事實(shí)”的冪函數(shù)模型形式卻不限于本文這一種。事實(shí)上，將本文的冪函數(shù)定量模型進(jìn)行線性組合，所得到新模型都同樣滿足上述條件。因此，要找到更加精確、實(shí)用性更好的模型仍然需要進(jìn)一步的探索。