高志遠, 謝元禮, 劉 狀, 蔣廣鑫
(1.西北大學, 西安 710127; 2.陜西省地表系統與環境承載力重點實驗室, 西安 710127)
SRTM DEM(Shuttle Radar Topography Mission)是全球最常見的免費中分辨數字高程模型之一,在土壤侵蝕、水土流失、流域分析、地貌分析等地理學各個領域都有著十分廣泛的應用。SRTM DEM的精度相較于另一種常見的中分辨率數字高程模型ASTER GDEM,不但在絕大部分區域內精度優于ASTER GDEM[1],并且SRTM DEM在地形較為復雜的地區其誤差與地形具有較強的相關性[2],即其誤差隨地形變化的敏感性較ASTER GDEM更強,因此利用SRTM DEM作為數據源進行精度與地形因子研究較為合適。
現階段,SRTM DEM精度評價研究已取得了不少進展,Rodriguez等[3]利用USGS發布的NED(National Elevation Dataset)數據集作為高程標準數據,對全美SRTM DEM進行了精度評價。Manas等[4]利用動態全球定位系統,對全球6個大陸的SRTM DEM產品的絕對高程誤差和相對高程誤差進行了評估。而后Mukul等[5]又提出利用全球定位系統網絡IGS(International Global Navigation Satellite System Service Network)分別對X波段和C波段的SRTM DEM數據進行了全球性的精度比較和評估。Zhao等[6]人利用ICESat/GLAS數據對中國中部地區的ASTER GDEM與SRTM DEM數據進行精度比較分析。詹蕾等[7]以我國1∶5萬比例尺DEM數據作為參考數據,利用空間插值的方法對全陜西省的SRTM DEM數據進行了精度評價。在地形因子的研究中,張泉等[8]提出SRTM DEM誤差與不同地形和地表覆蓋密切相關,坡度增大,誤差由正變負,誤差絕對值增大。張朝忙等[9]在研究中國地區SRTM DEM精度時提出DEM數據高程精度受地形影響并存在一定的空間分布性。杜小平等[10]提出DEM的高程誤差隨地形起伏增加而增加。南希等[11]提出DEM數據垂直精度對坡度有較大的依賴性,且中誤差隨坡度有近似指數的曲線的增長趨勢。可以說不管是全球大尺度還是地區小尺度,SRTM DEM的研究工作已有不小進展,但針對影響SRTM DEM精度的地形因子的研究工作卻比較少,在以往的研究工作中,研究者也傾向于將研究精力投入在DEM精度的評估工作,而涉及到地形因子的影響性時也偏向于做定性的試驗性結論陳述,而未做詳細的定量研究。本文利用ICESat/GLAS GLA14陸地測高數據作為高程標準數據,結合不同的地形因子綜合分析,對1弧秒SRTM DEM精度隨地形因子的變化進行定量化研究。
研究區以陜北典型黃土溝壑地貌為主要研究地貌類型,其大致范圍為東經108.7°—110.5°,北緯37.0°—38.1°。包括了陜西省延安市的子長縣和安塞縣的一部分,以及榆林市的子洲縣、米脂縣全境和橫山縣、靖邊縣、佳縣與綏德縣的部分地區,總面積約為20 000 km2。黃土溝壑地區溝壑縱橫、地形破碎、地形起伏度較大,各種地形因子變化較為明顯,因此把黃土溝壑地區作為DEM地形因子評價的試驗樣區是比較合適的。
本文中研究數據分為兩種,一種是SRTM DEM數據,將其作為精度評估原始數據,具體為SRTM1 Arc-Second Global,其空間分辨率為30 m,下載自(https:∥earthexplorer.usgs.gov/)。第2種是GLA14測高數據,GLA14(Global Land Surface Altimetry Data)數據隸屬于ICESat/GLAS數據集中的二級產品,由GLA 05和GLA 06數據再生產,其傳感器每秒發射40次脈沖亮斑,星下光斑點的直徑為70 m,沿衛星軌跡相鄰兩光斑間隔為170 m,GLA14水平精度為±20 cm,垂直精度為±18 cm[8-9]。綜合比較兩種數據可以看出,GLA 14數據較SRTM DEM數據具有較高的精度,其垂直精度和水平精度都達到了cm級,因此可以作為評測SRTM DEM精度的標準對比數據,GLA 14數據來源于美國國家冰雪數據中心,下載自(http:∥nsidc.org/)。
在本文研究中,引用4個指標量作為精度評價及地形因子分析的參考指標,分別為高程誤差d、平均誤差Mean、標準偏差SD和中誤差RSME。高程誤差d代表了每一個GLA14數據點與對應SRTM DEM柵格點的差值(HSRTM-HGLAH14),平均誤差Mean是高程誤差的平均值,是數據集中誤差的平均體現。標準偏差SD表示高程誤差之間數據離散的程度。中誤差RSME是指測量值與真值之間的離散程度,是對測量值精度最直接的評價。4種指標公式如式(1—4)[10-13]所示。
d=HSRTM-HGLA14
(1)
(2)
(3)
(4)
本研究中,引入高程、坡度、坡向、總曲率、平面曲率以及剖面曲率6種地形因子對SRTM DEM精度進行相關性分析。除了常規的幾種地形因子高程、坡度和坡向外,還涉及到了3種曲率因子,其中平面曲率反映了地形在等高線方向的變化率,表達了坡向的變化;剖面曲率是一種對地形在最大坡度方向上高程的變化的度量,是地面坡度的變化率,高程變化的二階導,是坡度的坡度;總曲率是平面曲率、剖面曲率與正切曲率的和,是地表復雜度和破碎度的一種度量方式。假設z=f(x,y)是地形曲面函數,并且函數連續、二階可導,則剖面曲率KV和平面曲率KC的計算公式如式(5—6)[14-16]所示。
(5)
(6)
式中:p是指x方向上的高程變化率;q指y方向上的高程變化率;s表示x方向高程變化在y方向的變化率;t是指y方向高程變化的變化率;r是指x方向上高程變化的變化率。
為了探求上述6種地形因子與SRTM DEM精度之間的關系,本文引入地理探測器的手段進行因子的相關性分析。地理探測器是一種新型的空間統計方法,它的基本原理繼承了亞里士多德以來人類認識自然環境的主要方法。即空間分異表達,是一種能夠揭示背后驅動因子的統計學方法,不同于常規的相關性分析,地理探測器具有更高的空間敏感性,能夠更加直接地揭示自變量與因變量之間的關系。此方法不但能夠同時參與多因子的相關性分析,并且可以揭示自變量因子間的交互影響對因變量的影響。本文中應用王勁峰等[17]開發的代碼進行計算,代碼下載自(http:∥www.geodetectors.org)。
本文中所采用的是表面積—尺度的方法來計算地表的分形維數,表面積—尺度法也稱為投影覆蓋法,由Clarke[18]等人首次提出,其主要過程如下:首先對原始SRTM DEM進行重采樣,得到不同空間分辨率的數據。然后記錄每一分辨率DEM的尺度Ri并且計算對應尺度下的柵格表面積Si,再利用線性擬合原理擬合lgRi與lgSi,并在雙對數坐標軸上選取一段擬合度最好的一段無標度區[19],記錄下此時的擬合線性斜率K,最后由D=2-K計算得到分形維數。
在進行分維數的計算前,首先需要選擇樣區,在本文中共選取78個分形維數計算樣區,計算樣區面積大小為5 km×5 km,且樣區中包含數量較多的GLA14數據點。
對SRTM DEM數據進行精度評價,做出GLA 14與SRTM DEM散點圖,可以發現兩者之間存在明顯的線性遞增關系,做出線性擬合方程,發現其斜率為1.000 9,決定系數R2為0.997 2(圖1),其擬合性極強,對其做ANOVA檢驗,發現其p值約為0,顯著性明顯。可以看出,SRTM DEM精度總體上與GLA14數據存在較為明顯的相關性和一致性。圖2為高程誤差d的頻數直方圖,可以看出高程誤差正值的個數明顯多于負值的個數,總數上前者比后者多18.45%,為8 732個,這表明在試驗樣區內,多數SRTM DEM測量值較GLA14值偏大,其數據系統誤差為(0.470±9.520) m,數據精度為9.531 m。

圖1 GLA14與SRTM DEM擬合

圖2 高程誤差頻數
按照地理探測器的算法要求,將高程、坡度、坡向、總曲率、剖面曲率和平面曲率按照數量級進行分類。因為總曲率、剖面曲率以及平面曲率值的正負只代表該柵格表面開口的凹凸性,因此這里對三者進行絕對值處理,保證其值的大小代表了不同維度上地表的破碎程度。
設高程誤差絕對值為因變量Y,高程因子為自變量X1,坡度因子為自變量X2,坡向因子為自變量X3,總曲率因子為自變量X4,剖面曲率因子為自變量X5,水平曲率因子為自變量X6,具體詳情見表1。
將各地形因子分級后利用地理探測器進行解算,高程誤差絕對值Y的地形因子分析結果見表2。

表1 高程誤差絕對值Y各因子探測

表2 高程誤差絕對值Y各因子探測
在顯著性水平0.05的情況下(表1),6種地形因子皆具有顯著性的統計學差異,其中坡向因子X3的q統計值最小,僅有0.003 5,這里可以認為坡向因子對高程誤差絕對值的影響極小。而坡度因子X2,總曲率因子X4,剖面曲率因子X5與平面曲率因子X6的q統計值較其他因子明顯較大,四者之中,坡度因子X2的q值最大,為0.063 8,說明了有6.38%的數據點高程誤差絕對值Y受到坡度的影響,總曲率因子X4,剖面曲率因子X5和平面曲率因子X6的q值相當,分別是0.033 7,0 330,0.030 3,即3種因子分別解釋了3.37%,3.30%和3.03%數據點的高程誤差絕對值。而由表2可以看出,坡度因子與曲率因子之間的疊加以及各類曲率因子之間的疊加對高程誤差絕對值的影響較大,其中總曲率因子X4與坡度因子X2的疊加和高程誤差絕對值的相關性最強,q統計值為0.085 1,說明了上述兩者因子的疊加解釋了8.51%的數據點的高程誤差絕對值Y。

圖3 坡度因子q統計值與總曲率分級關系

圖4 3種曲率q統計值與坡度分級關系
研究發現,曲率和坡度是影響高程誤差較為顯著的兩種地形因子。由圖3可以看出,在總曲率絕對值小于0.8時,坡度因子的q統計值與總曲率絕對值呈現負相關關系,在總曲率絕對值為0~0.2時,坡度因子q統計值最大,為0.152,即坡度因子能夠解釋15.2%數據點的高程誤差絕對值,當總曲率絕對值大于0.8時,坡度因子統計值q大致維持在0.018左右,坡度因子對高程誤差絕對值的解釋性較差。3種曲率因子q統計值在不同坡度分級呈現不同的分異特征(圖4),其中剖面曲率因子和平面曲率因子q統計值隨坡度分級的變化,分異性不大,而總曲率因子q統計值在坡度小于10°時較大,當坡度為0°~5°和5°~10°時,總曲率因子統計值q分別為0.063,0.065,總體上,當坡度較小時,總曲率因子q統計值較大,當總曲率因子較小時,坡度因子q統計值較大。
計算出每個分形維數樣區內的平均誤差Mean和中誤差RSME,利用地理探測器進行計算,令分形維數D為自變量X,平均誤差為因變量Y1,中誤差為因變量Y2,結果見表3。

表3 分形維數因子探測表
在顯著性水平0.05的情況下,兩者的p檢驗值都小于0.05,具有顯著性的統計學差異,分形維數與平均誤差的q統計值為0.370 4,說明分維數可以解釋37.04%分維數樣區的平均誤差,而分形維數與中誤差的q統計值為0.672 6,說明分維數可以解釋67.26%分維數樣區的中誤差。可見相比較高程、坡度、坡向、總曲率、剖面曲率與水平曲率這些普通的地形因子,SRTM DEM誤差與分形維數的相關性更高。
為了進一步了解分形維數與平均誤差及中誤差的關系,對其進行擬合分析,結果與地理探測器的結果相一致,見圖5—6,兩種精度指標量均與分形維數存在明顯的相關關系,相比平均誤差,中誤差與分形維數的相關性更強,決定系數R2達到了0.683 3。平均誤差和中誤差均與分形維數存在二次多項式的擬合關系,在研究區內,平均誤差隨著分形維數的增大先減少后增大,其平均誤差最小值出現在分形維數為2.034左右,最小值約為0.260 m。中誤差和分形維數的擬合趨勢隨著分形維數的增大,中誤差呈現增大的趨勢,其最大值時的分形維數約為2.046,最大值為13.42 m,在研究區內其中誤差隨分形維數變化處于一直增加的狀態,增長逐漸趨于平穩,直到峰值。

圖5 平均誤差與分形維數擬合

圖6 中誤差與分形維數擬合
(1) 在黃土溝壑研究區域內,SRTM DEM高程平均誤差為0.470 m,高程誤差標準偏差為9.520 m,數據精度為9.531 m,在研究區內符合SRTM DEM理論精度。
(2) 單一地形因子對SRTM DEM精度的影響非常小,具體來講,高程因子、坡向因子對數據精度的影響最小,坡度因子對精度影響相對最大;地形因子兩兩疊加后對精度的影響較單因子變大,其中坡度因子與曲率因子綜合后對精度的影響相對最大,并且當坡度較小時,總曲率因子q統計值較大,當總曲率因子較小時,坡度因子q統計值較大。
(3) 地表粗糙度(分形維數)對數據精度的影響較為明顯,其作為地表自變量因子能夠解釋37.04%的樣區的平均誤差,對中誤差的解釋程度達到了總樣區數目的67.26%。分形維數與平均誤差和中誤差具有比較明顯的二次多項式擬合關系。