DEM地形描述精度的幂函数定量模型

徐志敏，王利双

(1.长江勘测规划设计研究有限责任公司，武汉 430010；2.长江空间信息技术工程有限公司(武汉)，武汉 430010；3.湖北省水利信息感知与大数据工程技术研究中心，武汉 430010；4.国家知识产权局专利局专利审查协作湖北中心，武汉 430205)

0 引言

数字高程模型(digital elevation model，DEM)是国家基础地理数据的重要组成部分，也是进行三维空间分析和应用的基础。对DEM进行全面的精度分析和质量评价是保证其正确应用的基本前提[1]。目前主要的DEM精度评价模型很多，其中高程中误差模型[2]已经成为国家测绘行业的采纳标准(CH/T 1015.2—2007、CH/T 9009.2—2010)。然而研究表明，高程中误差只能评价DEM采样点的高程误差，但无法评价DEM模拟表面对真实地表的逼近误差，即使高程中误差为0，DEM模拟表面也不一定能真实表达地表起伏。也正因此，学界提出了多种DEM误差的评价模型，旨在对DEM误差进行更加科学的评价，包括重构等高线模型[3]、地形信息熵模型[4]、Strahler积分模型[5]、填挖方模型[6]与地形描述精度模型[7]等。

其中，DEM地形描述精度模型由于其科学性和实用性，引起了学界的关注和研究。汤国安等[8]最早提出了地形描述误差的概念和提取方法，并对地形描述误差的空间结构和定量计算模型做了研究。张勇等[9]以黄土丘陵沟壑区为例，采用统计和比较分析的方法，研究了不同空间尺度下的DEM地形描述误差的数学转换模型。刘春等[10]基于地形复杂因子概念，利用线性回归分析方法探讨了栅格DEM地形描述精度与地形复杂因子和平均高程之间的线性关系。贾敦新等[11]以模拟实验获得的高精度DEM为数据源，分析坡度计算精度与DEM数据误差和地形描述误差的耦合关系。齐晓飞等[12]提出了一种基于散点图矩阵建立多元回归模型的方法，并依此构建了DEM地形描述误差与分辨率、相对高差、凹凸系数的关系公式。王春等[13-14]借鉴误差椭圆和ε误差带建模的思想，提出了地形描述误差的中误差场和极值场，并研究了这些场的基本特性。可见，DEM地形描述精度的定量模型一直是学界关注的热点。

本文对地形描述误差(以下简称Et)与DEM空间分辨率和平均剖面曲率之间的定量关系进行了重新研究。从地形描述精度的概念出发，分析了现有的地形描述精度定量模型的缺陷，提出了优化的幂函数定量模型，并通过实验对比了本文所提出定量模型的拟合精度和泛化精度的优势。

1 现有Et定量模型的问题

汤国安等首次提出了地形描述误差的概念和计算方法，并以地形描述误差概念和计算方法为基础，选取了六种不同地貌类型(平原、低丘、丘陵、中山、高山以及混合类型)的DEM作为实验样区，分别计算了各样区在不同DEM空间分辨率下的地形描述误差的均方差。首先，对各样区的DEM地形描述误差均方差和DEM空间分辨率之间的关系采用y=ax+b线性模型进行了线性回归；然后，对方程系数b与对应样区的平均剖面曲率V分别再次进行线性回归；最后，得到DEM地形描述误差均方差与空间分辨率R和平均剖面曲率V的定量计算公式，如式(1)所示。

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

(1)

对此，王光霞等[15]指出：上述研究中系数b与平均剖面曲率的线性相关性并不好，利用直线拟合势必产生较大的误差。因此，对系数b与平均剖面曲率采用了二次函数模型进行了拟合，将参数a与平均剖面曲率的相关系数从0.989 9提高到0.998 9，将系数b与平均剖面曲率的相关系数从0.842 5提高到0.988 2，最终得到拟合精度更高的定量模型，计算方法如式(2)所示。该方法得到的定量模型与式(1)相比较具有更好的拟合效果，将在地势较平缓的地貌类型(如平原、低丘和丘陵地区)上的拟合中误差提高了一个数量级。

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

(2)

事实上，根据汤国安等对地形描述误差的定义，DEM地形描述误差Et是假定DEM高程采样误差为0的条件下，模拟地面与实际地面之间的差异。

如图1所示，弧面AB代表实际地形表面，面AB代表离散采样点所表达的模拟表面，AB点在水平面上的投影点之间的距离为DEM的分辨率。在点C处的DEM模拟地面高程与实际地面高程差为C点的DEM地形描述误差Etc。实际计算中，实际地面高程值一般采用4邻域的DEM栅格点高程的平均值进行替代。

图1 地形描述误差示意图

那么，根据DEM地形描述误差的定义和计算方法，在理想情况下，当DEM平均剖面曲率趋向于0时，DEM模拟地面为近似平面，任意栅格点四邻域的高程平均值与其本身高程差为0，地形描述误差为0；而当DEM分辨率趋向于0时，DEM模拟地面能够完全表达实际地面，地形描述误差也为0。因此，对于地形描述误差Et与DEM空间分辨率R和平均剖面曲率V之间的定量模型RMSEt=f(R，V)应该满足式(3)。

(3)

然而，式(1)和式(2)均为带常数项的多项式模型，在DEM分辨率R和平均剖面曲率V分别趋向于0时的极限均为带常数的不定式，因而均不满足式(3)。

图2显示了当平均剖面曲率为0，即地形绝对平坦的情况下，式(1)、式(2)对应模型的地形描述精度与DEM分辨率之间的函数关系。根据式(3)，当地形绝对平坦时，DEM地形描述误差与DEM分辨率不相关，且恒为0。但式(1)、式(2)在平均剖面曲率为0时，地形描述误差与DEM分辨率为线性关系，随着地形分辨率的提升，地形描述误差逐渐下降，且恒为某一大于0的值。

图2 V为0时式(1)、式(2)中RMS Et与R的关系

这显然于地形描述误差的固有性质不符。因此，式(1)、式(2)对应的模型均欠缺一定的科学性。事实上，文献[9]和文献[11]研究所得的定量模型也存在类似问题。这些模型在对DEM误差进行定量建模的时候都没有顾及到DEM误差与DEM地形特征之间的内在关系。因此，地形描述误差与DEM空间分辨率和平均剖面曲率之间的定量模型需要基于式(3)这一前提进行重新设计和研究。

2 基于幂函数的Et定量模型与应用

考虑到上述DEM地形描述误差的重要性质，本文基于式(3)的前提，设计了基于幂函数的地形描述误差Et与DEM空间分辨率R和平均剖面曲率V的定量关系模型，如式(4)所示。

RMSEt=aRbVc

(4)

式中：a、b、c为模型系数。

将(R，V，RMSEt)看作是三维空间中的离散点，将式(4)看作空间曲面方程。模型求解过程即可看作三维离散点的曲面拟合。由于式(4)与自变量的关系是非线性的，不能用线性最小二乘法求多元函数极值的方法获得参数估计值，因此需要采取复杂的优化算法进行求解。常用的优化求解算法包括最速下降算法、高斯-牛顿法和列文伯格-马夸尔特算法。本文采用基于高斯-牛顿迭代的非线性最小二乘法，对离散三维点进行曲面拟合。将三维离散点和三维曲面的拟合效果在MATLAB中显示，结果如图3所示。

图3 曲面拟合结果

拟合得到的模型系数a、b、c的值分别为：0.005 713、1.056和0.943 1。因此，基于幂函数的地形描述误差与DEM空间分辨率和平均剖面曲率之间的定量模型为式(5)所示。

RMSEt=0.005 713R1.056V0.943 1

(5)

进一步地，计算并比较式(1)、式(2)和式(5)对样本数据的拟合优度相关指标，包括误差平方和SSE、确定系数R-square以及标准差RMSE，结果如表1所示。可见，本文所提出定量模型的拟合效果很好，拟合精度与式(2)基本相当，且各项指标明显优于式(1)。

表1 拟合优度

对式(5)进行改写，可得DEM分辨率是地形描述误差和平均剖面曲率的函数，即式(6)，以及平均剖面曲率是地形描述误差和DEM分辨率的函数，即式(7)。DEM地形描述误差定量模型在生产实践中的实用价值在于：基于给定的DEM的限差，式(5)可以直接根据分辨率和剖面曲率检验DEM是否符合精度要求；式(6)可以推算一定地貌类型范围内符合精度要求的最低DEM分辨率；式(7)可以推算一定分辨率的DEM产品中符合精度要求的地貌类型区域。

(6)

(7)

3 泛化精度对比实验与分析

拟合精度评价的是定量模型公式对已知的实验样本数据的逼近程度，但是模型的拟合精度并不能保证模型在生产实践中的实用性，还需要评价模型对未知样本数据的估算效果，即泛化精度。

本文分别选取了平原、丘陵和中山地带的0.5 m高分辨率DEM进行了实验。数据来源于美国圣地亚哥超级计算机中心的开放地形数据OpenTopography网站，这个网站开放了高分辨率地形数据及处理工具(http：//opentopo.sdsc.edu/lidar?format=sd)。数据区域为美国加尼福尼亚州Napa Valley，地理范围为1×104m2，高程分布在海拔-1～634 m之间，从西南到东北方向依次分布着平原丘陵和山地，包含了比较丰富的地形种类。本文选取了该DEM范围内的平原、丘陵和中山的三块地形作为实验样区，源数据和实验样区的分布如图4所示。

图4 源数据与实验样区

采用文献[7]关于地形描述误差的计算方法，基于0.5 m分辨率的原始DEM，将地形描述精度计算窗口依次按照3×3、5×5、7×7扩大，分别计算得到3种地形在分辨率为1 m、2 m、4 m三种情况下的地形描述误差，总计九组数据。同时，采用式(1)、式(2)、式(5)分别计算九组不同DEM分辨率和平均剖面曲率情况下的地形描述误差估算值，结果如表2所示。将地形描述精度的计算值与式(1)、式(2)和式(5)的估算值进行统计比较。可以看出，式(5)的定量模型对于高分辨率DEM的地形描述误差估算精度都在83.3%以上，同时随着DEM分辨率的提高，模型的估算精度也呈现出增加的趋势，可见本文的定量模型对于高分辨率DEM具有很好的泛化性能。而对于式(1)和式(2)，其定量模型的估算精度不太稳定，式(1)的估算精度在45.5%和84.7%之间波动，式(2)的估算精度在49.4%和98.9%之间波动，甚至随着DEM分辨率的提高，一些样区地形描述误差的均方差估算值产生了不符合逻辑的负值。

表2 泛化精度对比结果

因此，实验结果印证了本文指出的式(1)和式(2)的问题：该定量模型在设计的过程中只追求对现有样本数据的拟合精度，没有考虑到地形描述精度的固有性质并根据对应自变量和因变量的固有数学关系来选择合适的拟合模型，导致模型对原有样本数据的拟合精度较好，但是对于更高分辨率DEM样本的泛化性能不理想。

4 结束语

本文对地形描述精度与DEM分辨率和地形平均剖面曲率之间的定量关系模型进行了研究。首先，从DEM地形描述精度的概念出发，分析了现有的带常数项的多项式模型的缺陷；然后，提出了基于幂函数的定量模型，计算并对比了模型对样本数据拟合精度；最后，基于 4 m、2 m和1 m的高分辨率DEM，比较了模型在新的样区的泛化精度。分析结果表明，本文的幂函数模型在保证良好的拟合精度的同时，对高分辨率的DEM具有较好的泛化性能。

随着测绘遥感技术的发展和计算机软硬件性能的提升，高分辨率DEM生产和应用由于能够满足用户对三维地形精度和体验的要求，得到越来越广泛的应用。在这种趋势下，本文的定量模型对于在生产实践中，根据限差要求推算适宜的DEM分辨率具有更好的指导意义和实用价值。

然而，还需注意到，虽然本文基于幂函数的定量模型具有较好的拟合精度与泛化性能，但是满足“当DEM空间分辨率和平均剖面曲率各自趋向于0，地形描述误差亦趋向于0这一客观事实”的幂函数模型形式却不限于本文这一种。事实上，将本文的幂函数定量模型进行线性组合，所得到新模型都同样满足上述条件。因此，要找到更加精确、实用性更好的模型仍然需要进一步的探索。