DEM 分辨率對提取東北地區坡度的影響研究1

劉佳軒 高 瑞 師黎靜

1）中國地震局工程力學研究所，地震工程與工程振動重點實驗室，哈爾濱 154200 2）黑龍江科技大學，建筑工程學院，哈爾濱 154200

引言

相比鉆孔數據受經費、環境及鉆孔技術等的限制，基于地形坡度的大區域場地分類方法得到大量應用。目前國內外研究多基于Wald 等（2006）建立的坡度?VS30模型進行場地類別劃分（史大成，2009；Thompson 等，2014；吳效勇，2019），其中坡度來源主要是基于30″分辨率數字高程數據（Digital Elevation Model，DEM）的提取。然而，提取的坡度受DEM 分辨率、算法及數據類型（網格點類型）等因素影響，且分辨率是最主要的影響參數。國內外大量研究表明（Walker 等，1999；Mukherjee 等，2013；Grohmann，2015；劉利峰，2018；師動等，2020），提取的坡度隨著DEM 分辨率的降低而減小，即坡度均值與標準差明顯降低，且不同地區分辨率的影響不同，這對建立坡度-剪切波速、坡度-覆蓋層厚度等關系模型具有重要影響，因此應選擇合適分辨率的DEM 進行數字地形分析（Thornley，1976）。李昕蕾（2020）通過對新疆地區不同分辨率DEM 坡度的提取，發現利用精度高或大區域小比例尺數據進行場地分類的效果較差，而30″分辨率效果較好；喬之軒（2019）研究發現，高分辨率DEM 提取的坡度具有更精細的局部變化，但有時會掩蓋區域坡度的整體變化趨勢和特征。受計算機能力與硬件顯示等條件限制，高分辨率DEM往往不能得到良好應用，有時需采用與研究比例尺相適應的分辨率。部分研究者在提取不同分辨率下的坡度時考慮地貌條件，發現不同地貌下的坡度尺度具有明顯差異性（Wang 等，2012；劉曉等，2017；土祥等，2018）。目前，對我國東北地區考慮地貌條件下的坡度研究相對較少，且不同研究區域可能搜集到不同分辨率DEM，需對提取的坡度進行統一，王英等（2019）、劉飛等（2019）以分形理論為基礎，建立不同分辨率下DEM 坡度轉換模型，分別以涇河區域和四川丘陵地區某流域為研究對象，研究成果不適用于我國東北地區。

基于上述問題，本文以美國地質勘探局（United States Geological Survey，USGS）發布的我國東北地區1″、3″、30″分辨率DEM 數據為基礎，分別提取研究區域內公里網格點的百分比坡度，分析不同分辨率下坡度差異原因。區分平原、丘陵和山地后進行坡度-頻率分布曲線統計，分析不同地貌單元下的坡度隨分辨率變化特征。通過線性函數模型、多項式函數模型及冪函數模型，擬合不同分辨率坡度轉換關系，并對比擬合優度R2和坡度分級平均相對差值Δmean，選取不同DEM 分辨率得到坡度最佳擬合公式。

1 研究區域及數據

1.1 研究區域

研究區域主要包括黑龍江省、吉林省、遼寧省及內蒙古自治區部分地區，地理坐標為116°E～135°E，41°N～53°N。東西向最大橫距約1 441 km，南北向最大縱距約1 656 km，面積約1 450 000 km2。東北地區海拔高度為?275～2 692 m，整個地勢大致分為3 環，外圍由黑龍江、烏蘇里江、鴨綠江等流域構成，整體地勢較低；中部為由大興安嶺、小興安嶺及長白山構成的山地與丘陵，整體地勢較高；內部為由松花江、嫩江等流域構成的東北平原，整體地勢較低。

1.2 DEM 高程數據

研究數據來源于USGS 官網發布的DEM，原始數據采樣間隔分別為地球等角坐標系的1″（約30 m）、3″（約90 m）和30″（約1 km），圖1 所示為3 種分辨率下高程數據模型，不同分辨率下DEM 表現出的整體地勢結構基本相同，高分辨率DEM 對研究區域的描述更詳細，具有更精確的高程數據范圍及更精細的地形刻畫，隨著分辨率的降低，區域內高程數據產生一定程度的過濾與平滑，使區域整體高程數據范圍減小，但整體地勢結構更突出。

圖1 東北地區不同分辨率下高程數據模型圖Fig. 1 Elevation maps of 1 arc second, 3 arc second, and 30 arc second resolution in the Northeast, China

2 坡度提取

2.1 提取方法

目前最常用的坡度提取方法為擬合曲面法（Olaya，2009），計算式為：

式中，slope 表示坡度，本文使用廣泛應用于大區域場地分類的百分比坡度；dz/dx，dz/dy分別表示計算點處2 個正交水平方向高程梯度變化率，無量綱，考慮復雜地貌坡度影響時，dz/dx，dz/dy優選算法建議使用三階反距離權算法，即Horn 算法（Horn，1981）。

圖2 所示為計算高程梯度變化率dz/dx，dz/dy的3×3 像元窗口，Horn 算法中高程梯度變化率dz/dx，dz/dy計算如下：

式中，w為像元寬度，即分辨率。

以1″和30″分辨率DEM 坡度計算過程為例，說明不同分辨率下坡度提取差異。當分辨率為1″時，式（2）和式（3）是對圖2 代表的約90 m×90 m 范圍內高程數據點的計算，最終通過代表中間柵格數據的高程梯度變化率得到坡度。當分辨率為30″時，式（2）和式（3）是對圖2 代表的約3 km×3 km 范圍內高程數據點的計算，進而計算得到坡度。雖然1″、30″分辨率下計算中心點可以相同，但坡度體現計算面積平均特性，隨著分辨率的增大，柵格計算面積變小，因此高程數據有所差異，導致不同分辨率下的坡度有時存在較大差距，且區域附近坡度變化越大，可能導致不同分辨率下坡度差異越大。

圖2 計算高程梯度變化率的像元窗口Fig. 2 3×3 unit schematic diagram of slope calculation point

2.2 提取結果

圖3 所示為通過式（1）及DEM 數據計算得到的不同分辨率下百分比坡度圖，由圖3 可知，高分辨率下的坡度圖對地表地形的描述更精細，坡度特征表現層次更豐富；隨著DEM 分辨率的降低，坡度值域逐漸減小，部分地區坡度趨于0，整體地形結構更突出，但缺少地形細節的表達。通過GIS 軟件建立與不同分辨率下坡度圖相同的漁網式公里網格點，并提取網格點對應的坡度。

圖3 東北地區不同分辨率下坡度圖Fig. 3 Slope maps of 1 arc second, 3 arc second, and 30 arc second resolution in the Northeast, China

圖4 所示為東北地區不同DEM 分辨率下坡度統計曲線，由圖4（a）可知，高分辨率下坡度累計頻率曲線較低分辨率更平緩，即高分辨率下坡度整體相對偏高；對于坡度值域而言，高分辨率下坡度值域范圍較低分辨率廣，1″、3″、30″分辨率下坡度分別約為0.6、0.4、0.2 m/m。由圖4（b）可知，隨著分辨率的降低，坡度均值和標準差明顯減小，且在高分辨率下減小速率較快，在低分辨率下減小速率較慢，標準差的降低說明了坡度值域整體向均值附近靠近。圖4 基本表現出東北地區整體坡度隨著DEM 分辨率的降低而減小的現象。

圖4 不同DEM 分辨率下的坡度統計Fig. 4 Slope statistics of different DEM resolutions

3 不同地貌坡度統計分析

區分不同地貌條件，對研究區域網格采樣點坡度進行統計分析，研究不同地貌下坡度隨分辨率變化特征，為擬合公式按照不同地貌條件提供理論基礎。

圖5 所示為我國東北地區地貌分布，平原、丘陵和山地地貌占比分別約為33.3%、12.2%、43.1%，其他地貌占比約為11.4%。通常，平原地貌較平坦或起伏較小；丘陵地貌由多個低矮山丘組成，具有一定起伏；山地地貌由眾多山脈區域組成，地面起伏較大。對于海拔高度而言，通常山地海拔較高，丘陵次之，平原較低。對于坡度而言，通常山地坡度較大、較陡，丘陵次之，平原坡度較小、較緩。

圖5 地貌分布圖Fig. 5 Landform distribution map

圖6 所示為不同分辨率下不同地貌單元坡度-頻率分布曲線。由圖6 可知，平原在坡度較小時頻率分布較大，山地在坡度較大時頻率分布相對較大；隨著坡度的增大，平原坡度-頻率分布曲線下降速率較山地快，而丘陵坡度-頻率分布曲線變化特征處于兩者中間；隨著分辨率的減小，平原坡度-頻率分布曲線峰值位置迅速向低坡度段移動，峰值區域變窄，坡度減小至0.05 m/m 以下，而山地坡度-頻率分布曲線峰值位置向低坡度移動速度相對較緩，峰值區域變窄程度相對較小，坡度減小至0.15 m/m 以下，丘陵坡度-頻率分布曲線變化特征處于平原與山地之間，坡度減小至0.1 m/m 以下。各地貌單元坡度減小程度不同，因此需考慮地貌條件，建立不同DEM 分辨率的坡度轉換關系。

圖6 不同分辨率下不同地貌單元坡度-頻率分布曲線Fig. 6 Slope-frequency distribution curves of different topography at different resolutions

4 轉換關系擬合回歸

4.1 擬合模型

本文采用常用函數模型對不同DEM 分辨率提取的坡度進行擬合，各函數模型如下：

線性函數模型：

式中，Si表示分辨率為i（其中i=1″、3″、30″）時計算得到的百分比坡度（m/m）；aj、bj（j=1、2、3）、c1、c2、d1分別表示不同函數模型的回歸參數。為使擬合公式計算結果更合理，需考慮邊界效應，即當1″、3″分辨率下坡度趨于0 時，30″下坡度也應趨于0。因此，當Si=0 時，S30=0，即a1、d1、c2=0。

4.2 擬合回歸

基于上述函數模型，分別對不同地貌單元1″、3″分辨率下坡度與30″分辨率下坡度進行擬合，坡度轉換關系參數如表1、2 所示。由表1、2 可知，平原回歸關系參數b最小，而山地回歸關系參數b最大，丘陵回歸關系參數b處于二者之間，基本符合隨著分辨率減小，平原坡度減小程度較大，而山地坡度減小程度相對較小的變化規律；多項式函數及冪函數模型擬合參數具有類似規律。

表1 不同地貌單元回歸模型擬合參數及評價參數（1″與30″分辨率下）Table 1 Fitting and evaluation parameters of regression models for different landforms (1 arc second and 30 arc second)

5 回歸方程評價方法

5.1 擬合優度

表2 不同地貌單元回歸模型擬合參數及評價參數（3″與30″分辨率下）Table 2 Fitting and evaluation parameters of regression models for different landforms (3 arc second and 30 arc second)

擬合優度R2計算如下：式中，SSR為回歸平方和，SST為總平方和，ESS為誤差平方和，TSS為總離差平方和。

由表1、2 可知，冪函數模型轉換的擬合優度R2最小，即冪函數模型擬合效果最差；線性函數模型擬合優度R2略小于多項式函數模型；整體來說，多項式函數模型擬合結果最優。

5.2 坡度分級平均相對差值

高分辨率DEM 可獲得高分辨率坡度，進而獲得高分辨率等效剪切波速或覆蓋層厚度，但場地分類主要使用30″分辨率DEM 提取的坡度，分類方法基本為建立坡度與等效剪切波速的分級對應關系，因此還需考慮轉換前后30″坡度分級結果。將坡度按0.01 m/m 寬度范圍進行分級統計，以30″分辨率DEM 提取的坡度分級作為參考值，分別計算不同擬合模型對高分辨率坡度轉換后的坡度分級與參考值的平均相對差值Δmean，平均相對差值越接近于0，說明該擬合公式向30″坡度轉化越好，平均相對差值越接近于1，說明該擬合公式向30″坡度轉化越差，平均相對差值計算如下：

式中，Pi,k表示在i（i=1″、3″）分辨率下轉換后第k組分段坡度，n表示坡度分級級數。

不同地貌單元不同分辨率下坡度轉換結果如圖7 所示，由圖7 可知，平原1″、3″與30″真值接近，基本將分布在0～0.2 m/m 的坡度轉換為0～0.03 m/m；丘陵1″向30″轉換時，線性函數模型轉換與多項式函數模型轉換結果較好，而3″向30″轉換時多項式函數模型轉換結果較好；山地線性函數模型轉換效果較好，而多項式函數模型與冪函數模型轉換效果較差。由表1、2 可知，線性函數模型分級平均相對差值Δmean最小，即僅考慮坡度分級時，線性函數擬合為最佳擬合。

圖7 不同地貌單元不同分辨率下坡度轉換結果Fig. 7 Conversion results of DEM slopes with different resolutions under different landforms

5.3 最佳擬合模型

由于分別考慮擬合優度R2和分級平均相對差值Δmean，各模型在坡度轉換上均有優劣，因此，以比例參數K綜合考慮R2、Δmeans，計算公式如下：

K值越大，說明擬合模型擬合效果越好。

根據表1、2 中的比例參數K，最終確定不同地貌單元不同分辨率下坡度轉換公式，如表3 所示。平原和山地最佳坡度擬合模型均為線性函數模型；丘陵1″與30″分辨率最佳坡度擬合模型為線性函數模型，3″與30″分辨率最佳坡度擬合模型為多項式函數模型。

表3 不同地貌單元不同分辨率下坡度轉換公式Table 3 The slope conversion formula of different resolution DEM for different landforms

6 結論

本文以我國東北地區1″、3″、30″分辨率DEM 數據為基礎，分別提取研究區域公里網格點百分比坡度，分析不同分辨率下坡度提取差異原因。區分平原、丘陵和山地進行坡度-頻率分布曲線統計，分析不同地貌單元坡度變化特征。通過線性函數模型、多項式函數模型及冪函數模型擬合不同分辨率坡度轉換關系，并計算擬合公式擬合優度R2及坡度分級平均相對差值Δmean，得出以下結論：

（1）不同分辨率提取的坡度不同，主要原因是坡度計算方法體現計算面積平均特性，隨著分辨率的增大，柵格計算面積變小，高程數據有所差異，導致不同分辨率下坡度存在較大差異。由于坡度計算平均特性，區域附近坡度變化越大，可能導致不同分辨率下獲得的坡度差異越大。

（2）我國東北地區坡度隨DEM 分辨率的降低而減小，高分辨率下坡度圖對地表地形的描述更精細，坡度特征表現層次更豐富，而低分辨率下整體地形結構更突出，但缺少地形細節的表達。隨著分辨率的降低，我國東北地區坡度均值、方差及變化范圍減小，坡度-頻率分布曲線向低坡度段移動。

（3）隨著分辨率的降低，不同地貌單元坡度-頻率分布曲線變化趨勢相同，變化程度不同。各地貌單元坡度-頻率分布曲線均向低坡度方向移動，且峰值區域變窄，其中平原地貌變化趨勢性顯著，丘陵次之，山地最小。因此，建立坡度轉換關系時需區分不同地貌條件。

（4）給出不同地貌單元下東北地區不同分辨率坡度轉換模型，為利用30″分辨率的坡度進行大區域場地分類提供轉換關系。平原和山地最佳坡度擬合模型均為線性函數模型；丘陵1″與30″分辨率最佳坡度擬合模型為線性函數模型，3″與30″分辨率最佳坡度擬合模型為多項式函數模型。

（5）本文給出的最佳坡度轉換公式主要應用于工程抗震領域中的場地分類。