The Optimal Lid Method for the objective definition of cross-section limits in ephemeral gullies
划定沟道横截面边缘的最佳盖板法
切沟是在降雨时集中的地表径流流经易侵蚀的农田土壤而形成的小水沟。尽管利用常规农业耕作措施就能够填平,但是后续会在同一地点再次出现侵蚀沟。因此,实际工作中需要能够准确测量切沟横截面的可行性方法。
现阶段沟蚀测量通常包括三个主要阶段:(1)确定沟谷地形特征;(2)确定侵蚀沟边界;(3)计算侵蚀沟边界内的体积。首先要定义沟谷纵轴任意位置的两个横截面边缘(XSLL和XSLR,分别表示左边缘和右边缘)。然后就可以得出沟谷深度、顶宽和面积。因此沟道尺寸取决于沟道横截面限制的定义。对沟蚀量的估算主要采用横截面外推法,通常会沿沟道选择一些横截面,将相邻横截面的平均面积值乘以它们之间的距离,就能估算出侵蚀沟的体积。沟道体积的总和可估算出沟蚀造成的土壤流失总量。
尽管横截面是评估沟壑侵蚀率最常用的信息来源,基于横截面的沟道划分很少受到研究人员的关注。当前有两种主要的定义横截面界限的方法:坡度变化法(change-in-slope)和侵蚀前表面法(pre-erosion surface)。
坡度变化法的基础是由专业操作员目测识别出坡度变化剧烈的横截面位置(左右边缘),这可以在实地直接进行测量(例如使用卷尺测),也可以通过分析地形数据获得。这种方法的重要优势是操作人员可以借助额外的视觉信息,如水流迹象或土壤表明粗糙度的变化,来区分原有地表和侵蚀后的地表。然而坡度变化标准的一个主要缺点是在复杂的横截面中经常会出现很多次坡度剧烈变化。因此选择最突出的变化受到主观因素影响。
侵蚀前表面法包括估算一个理论上存在的侵蚀发生前的地表,从中减去侵蚀后的剖面,以量化特定时间跨度内的侵蚀体积。该方法可以使用栅格图层中的DEM来实现。这种方法在相关的科学文献中并不常见,但最近也是一种能够有效划分横截面边缘的方法。
所有基于人类对横截面边缘划分的方法都存在一定的主观性,且无法在建模程序中实施,而自动方法通常需要只适用于当地条件的阈值。由于学界缺乏统一规则,这些问题导致不同地点之间侵蚀估算值没有可比性。因此,开发一种可靠、可重复、客观的沟道横截面自动划定程序将有助于最大限度地减少不确定性,提高沟蚀估算的可比性。在这项工作中,主要解决了两个研究问题:(1)基于专业操作员的沟道横截面边缘估计涉及多少不确定性;(2)能否设计出一种自动划分沟道横截面边缘的方法,其准确性和精确性如何?
为了解决这些问题,本研究旨在设计一种新的方法(OLid法),通过在三种情况下进行校准测试,降低不确定性:(1)在已经确定的沟道上合成横截面;(2)在41个实地切沟上,应用三维摄影测量和实地目视识别来确定横截面边界;(3)针对学者先前仅使用二维地形数据对60个实地切沟进行评估。
图1 本研究使用的方法和数据集
沟谷横截面的左右边缘(XSLL和XSLR)确定了用于评估其尺寸的剖面部分(图2a)。从形态上看,沟谷横截面可分为两个区域:边缘和沟道,边缘通常坡度较缓(未受径流冲刷的原地表)。定义沟谷横截面的常用标准是宽度(W,Width)、深度(D,Depth)和面积(A,Area)。沟谷宽度定义为沟谷边界之间的水平距离(图2a)。深度被定义为横截面最低点与每个横截面边缘(XSLL和XSLR)之间的平均高差(DL和DR),横截面面积是连接两个界限的斜线与下方沟道之间的封闭表面。
图2
将沟谷横截面边缘确定为符合一系列地形特征的河谷两侧的点:(1)剖面上突出较高的位置;(2)与沟谷边缘平台保持一致;(3)确定侵蚀前的光滑表面。
OLid法基于一个迭代过程,在此过程中,这三个分量的总和(S)达到最大值,同时考虑到可能的左右XS边缘候选值的所有可能组合(公式 1)。突出位置分量由具有四个顶点(如图 2b 中的XSLL-ULL-LL-Th,图中阴影部分面积)的梯形面积表示(变量AT,在更突出的位置显示更大的值),它是该方法初步版本的主要变量;边缘对齐分量由边缘与连接极值和极限候选点的直线之间的平均高程差(e)表示(图2c,坡沟系统中坡面的高差),e越小,说明对齐程度越好;侵蚀前表面(pre-erosion surface),即连接极点与边缘候选点和沟盖的直线与水平面的夹角之和(β)(图2d,坡面坡度),β越小,说明侵蚀前表面越平滑。横截面边缘是突出部(梯形 "AT "面积最大)、边缘对齐(平均差 e 值最小)和侵蚀前表面光滑度(角度总和最小)的最佳组合。
对于每个组成部分,通过减去最小值并除以极差(最大值减去最小值),以及按项的数量(2或4)求平均值,对变量进行归一化处理(数据标准化)。在侵蚀前地表部分,沟盖角的系数(2)更大,以赋予沟道和沟边同等权重。OLid 公式(公式 1)中的加号表示每个分量的XS极限值偏高(正值,仅梯形区域)或偏低(负值,其余区域)的趋势。每项之前都有一个需要校准的权重参数("k")(第2.4.2节)。确定每点的S值的OLid公式为:
XS中的最大S值是通过迭代计算得出的,也就是说,在连续试验中探索所有可能的配对组合。迭代过程包括两个阶段:(1)对于左侧横截面边上的每个点,在右侧横截面边上寻找S变量最大的点(图3b);(2)在所有候选中找到S最大的一对(最终横截面边缘,图3c和d)。在每次迭代中,只有那些定义横截面边缘的线对在任何一点上都不与横截面剖面外轮廓相交,才会被认为对目标变量的计算有效。
图3
OLid方法是完全自动化的,因此很容易在代码脚本中实现。本研究使用Matlab代码脚本(Mathworks,Natick,MA,USA)进行分析,输入为包含横截面剖面坐标csv文件。输出的是OLid脚本生成一幅图像(图4a)和一个 csv 文件(图 4b),其中包含数据集中所有横截面的尺寸(宽、长、高、极限坐标等)。
图4
本研究使用了四个不同的XS数据集(表1)。前两个数据集(XSynOne、XSynTwo)是为模型筛选(公式1中参数的初步筛选)而合成的,第三个数据集(XSCor)用于校准和验证模型,最后一个数据集(XSNav)用于比较OLid和传统方法(专业操作人员实地测量)。
每个数据集的特点如下:
XSyn1:模拟1个侵蚀周期的沟道横截面(图5a)。这些是由直线段组成的简单横截面,呈现三角形、梯形或矩形,由两个边缘、两个堤岸和一个平坦的沟底组成。在Matlab脚本中使用均匀分布的随机数生成了总共100个XS,每个设计参数为:
沟顶宽度W为0.2-2m,宽深比WDR为1-7,沟底宽度BW为0.01-0.99Wm,左边缘坡度SL为1%-15%,右边缘坡度SR为1%-15%,且左右深度(图2a)与平均深度的差异范围在20%以内,两个边缘的斜率也是如此。
XSyn2:模拟2个侵蚀周期的沟道横截面(图5b)。这些横截面也是由直线段组成,其形状比前一种情况更复杂,包括中央通道上方的一个中间台阶。同样,在Matlab中使用随机数方法共生成了100个横截面XS。最低侵蚀周期的参数范围为:沟顶宽度W为0.2-2m,宽深比WDR为1-7,沟底宽度BW为0.01-0.99Wm。深度 DL 的最小值为0.3m。最顶端侵蚀周期的尺寸定义为最低水道宽度(WL)和落差高度(DL)的函数:顶端宽度WT为1.2-1.8倍的WL),左侧落差d为0.2-0.8倍DL,最小值为0.15 m。为了与水流侵蚀得出的水平几何形状趋势保持一致,并为了简化起见,假定中间台阶的左侧和右侧高原是平坦的。根据其他设计参数的定义,最顶端和最低端河岸的坡度被认为是相等的。一旦确定了左岸坡度d和其他参数,右岸坡度也就自动确定了(图5b,XS 23)。最后,左边缘和右边缘斜率在1%-15%之间浮动。
XSCor:在Campiña地区(距离西班牙科尔多瓦市20公里)种植了一年生作物的农田中,利用运动结构(Structure-from-Motion,SfM)技术重建了一条长约50米的切沟,以生成一个基准数据集,探索OLid在实际田地条件下的准确性。使用摄像机(分辨率为1200万像素)在延时模式下(间隔2秒)手动收集了总共300幅图像。为了给三维模型提供比例和方向,在土中挖出了一根已知尺寸的垂直水平钢筋(离地面0.5米)。以钢筋的长度和方向为基准对点云进行转换。通过使用CloudCompare中的"section"命令并绘制横向于沟谷纵向边缘的横截面,从三维云中获得了总共41个横截面(平均距离约为1.2米)。最后,在沟谷的三维彩色点云上目测确定横截面边缘,并将得出的尺寸作为实测值。横截面尺寸的范围包括:宽度W为0.29-1.26m、深度D为0.12-0.36m,而面积的变化系数VCA=0.48(图5c)。通过在所有横截面中随机选择,将该数据集进一步划分为校准子集(XSCorCal,n=27)和验证子集(XSCorVal,n=14)。
XSNav:通过使用激光轮廓仪或针式轮廓仪对Navarre地区的农田进行实地调查,共获得了60个切沟的横截面面积。调查发现,横截面的尺寸变化很大,沟宽0.2-2.0m不等,深度0.04-0.90m不等,面积变化系数也很大(VCA=1.06,图5d)。切沟剖面显示出截然不同的XS形状和形态,在某些情况下,边缘宽度很低,这使得确定横截面边缘对专业操作人员和OLid来说都是一项具有挑战性的任务。
图5
为了对横截面的形状多样性进行量化评估,本文使用了两个形态因子。宽度深度比(WDR)和形状因子(SF):
SF为1 表示横截面形状类似矩形(“[”),1-2之间表示是梯形或抛物线(“U”),2表示三角形(“V”),大于2则是典型的上宽下窄(“{”)。
图6显示了四个数据集的WDR和SF值以进行比较。合成的横截面数据集包括从矩形到三角形的所有剖面形态(SF值为1-2)和广泛的浅度范围(1-8),XSynTwo组中有一些括号型沟谷(SF>2)。最后一种类型在野外中最为常见,而野外又不包括矩形沟壑(最小SF值约为1.3)。四个数据集的形态复杂度呈递增趋势(XSynOne<XSynTwo<XSCor<XSNav)(图6b)。
图6
为了找到OLid参数的最佳方案,本文进行了逐步分析,包括筛选、校准和验证阶段。
4.1. 筛选
考虑了三个系数(ktrap、kfit和kang)的所有可能的排列组合,其值范围为0-1,间隔为0.1(11个值),即113=1331个案例。为避免冗余,只考虑了定义参数之间相同比例的排列组合(例如,当ktrap=1、kfit=1和kang=1时,去掉了ktrap=0.2、kfit=0.2和kang=0.2等组合),总共剩下 1064个案例。
为了在不同的几何条件下对OLid进行初步测试,并减少校准阶段的案例数量和处理时间,对两个合成数据集(XSynOne和XSynTwo)进行了分析。通过将实现横截面尺寸100% 精确度(即模型预测的尺寸与脚本中给出的设计值完全一致)定义为前提条件,筛选阶段可作为模型校准前的初步筛选。
4.2.模型校准:准确度和精确度
横截面边缘预测模型的两个验证参数是准确度(与误差相对)和精确度(与不确定性相对)。准确性是指模型对横截面量级的估计值接近参考技术的估计值。精确度是在定义估计值时具有较高的辨别能力,即预测的不确定性较低。
XSCorCal子集(n=27)用于校准OLid模型,只考虑成功通过筛选测试的参数组合。通过误差分析对模型的准确性进行了评估,以比较横截面尺寸的预测值和观测值(由现场 XSCorCal横截面的三维图像确定)。为了避免误差评估中的规模效应,选择了每个参数集的平均绝对相对误差(ErW)作为差异度量。在本分析中,选择沟宽(W)作为对误差最敏感的维度:
其中ErW代表宽度相对误差,Wpred和Wobs分别代表沟谷宽度的预测值和观测值。
除精度测试外,还对模型进行了精度分析。与使用误差评估精度的方法类似,我们使用 XS极限位置预测的不确定性来评估精度。OLid以最大化原则为基础,即横截面边缘由几何配置导致剖面中S值最大的一对点来定义。最有利于精确估算的情况是,有其他空间上接近XS极限(不确定性跨度小)的S幅值与最大值相似的数据。在这种情况下,地形或测量误差的微小变化对XS量级的影响很小(误差传播的可能性低)。为了评估 XS 极限值计算的精确度:(1)将S最大值的0.95分位数定义为阈值;(2)在剖面图中找到S超过该阈值的对,它们定义了不确定性的左右总跨度(LS和RS,图7);(3)将不确定性U定义为不确定性总跨度在XS剖面图总长度中所占的百分比:
图7
对模型的校准是通过确定哪种调整方式能使精度和准确度达到最佳综合效果来进行的。为此,根据误差和不确定性的大小从大到小对所有参数集进行排序,ErW和U平均值最小的参数集排在第一位(最佳性能),两个排序之和为总分,总分最小的参数集对应最佳调整(ktrap、kfit和kang的最佳组合)。
4.3.模型验证
利用校准模型,在14个实地XSCorVal断面上对模型进行了测试。同样,利用所有横截面的平均相对误差来评估模型的准确性。
5.比较OLid法和基于专业操作人员测量的方法:
在XSNav数据集中,OLid模型的性能还与专家小组方法(二维数据集中用于划定XS边缘的传统方法)进行了比较。对于专家小组(EP)方法,一组沟壑侵蚀研究人员利用收集到的实地数据进行了桌面划界测试。
专家们设计了Excel表格,用于记录操作员在图上目测到的横截面左右界限,并据此自动计算横截面尺寸。分析结果(每位专家60个Excel文件的宽度、深度和面积)由一位未参与划界工作的作者进行处理和分析。为了处理大量文件,编写了一个Matlab例程,用于打开和读取 Excel 文件中的结果,最后自动保存一个包含每位专家测量的所有XS尺寸的csv文件。此外,还利用变异系数(VC)及其与专家中位数观测值的相关趋势,评估了专家们对XS维度估计值的变异性。
1.定义和测试OLid:筛选、校准和验证:
在1331个案例中,去除冗余后剩下1064个案例(80%,图8b)。其中,820组成功预测了所有XSynOne数据集中横截面的宽度尺寸(62%,图8c),只有18%通过了XSynTwo数据集精确度测试(图8d)。
在筛选阶段的第一步(使用单周期侵蚀横截面),主要舍弃了高kang值和低ktrap值。然而,第二阶段的筛选最为严格(XSynTwo数据集),只将范围缩小到总数的五分之一(图8d)。大多数有效案例都集中在较宽的kang值范围内的高kfit和中等ktrap量级上。
图8
在大幅减少参数集的情况下,对模型校准进行了精度和准确度测试。图9描述了XSCorCal数据集中平均ErW和平均U之间的关系。发现准确度和精确度之间总体上呈正相关。虽然某些参数集的平均ErW接近100%(估计宽度值大于实际宽度的两倍),但表现最好的参数集有7%的宽度误差和2%的不确定性。表2显示了校准阶段表现最佳的10个参数集。总分最高(最低)的是集合474(ktrap=0.6、kfit=0.5和kang=0.9),它定义了最终校准的Olid模型(平均ErW=6.71%,U=2.30%)。ErA和ErD分别为6.04%和3.21%。
图9
通过优化组合,观测值和预测值之间的判定系数RW2为0.916(图 10)。同样,面积和深度的判定系数分别为RA2=0.977和RD2=0.968。研究还发现校准数据集有的高估趋势(回归线斜率为1.02)。图10b展示了误差最大的一个示例(XSCorCal 14),其中OLid将沟壑宽度高估了35%。
利用校准模型,在XSCorVal数据集(n=14)上进行了验证。验证结果显示与校准数据集的性能相似,平均宽度、面积和深度的相对误差分别为ErW=7.05%、ErA=4.66%和ErD=2.95%(RW2=0.921、RA2=0.986和RD2=0.993)。在这种情况下,发现有轻微的低估趋势(回归线斜率为0.981)。
图10
2. 对专家方法的评估:
2.1.专家小组在XS估算的差异
XSNav剖面图显示了专家对沟壑面积观测值(此处以XS尺寸为基准)的变异系数(图11a)。W、D和A估值的平均变异系数分别为16.7%、5.5%和12.2%。虽然在少数几个XS中,所有专家都提供了相同的估计值,但在其中一些XS中却有非常大的差异。图11b和c显示了受EP变异影响最大的XS(XS 20和35)。例如,对于XS 20,专家1在坡度变化最突然的位置确定了界限,而专家3则通过绘制准水平盖来选择更宽的XS。
为了说明专家组内部的差异和趋势,图 12 描述了高估和低估的趋势。由于XSNav数据集中没有 XS 维度的基本真实值,因此下文将这些趋势称为相对于组内中位数(中位数专家,表示对极端值的敏感度低于平均值)。一些专家倾向于报告高于专家中位值的观测值(专家1、6和7,其估算值在高估算值类别中所占比例较高),而另一些专家则提出了大量远低于平均值的观测值(如专家4和5)。专家3的情况最为平衡(与中位数的吻合度最高,为63%,"偏高"类别的数量最少,为27%),而专家5的情况恰恰相反,低估的情况占40%,与专家中位数吻合的情况占42%。
图11
图12
2.2.比较OLid和专家平均估计值
图13显示了XSNav数据集三个XS维度的OLid与专家估计中值之间的相关性。总体而言,宽度是差异最大的变量(R2=0.77),而深度和面积的一致性则非常高(RD2=0.98和RA2=0.96)。在所有情况下,最佳拟合线的斜率都小于1,与专家中值相比,呈现出数值较低的趋势,这主要是由于三个XS(XSNav 5、9和10,图13a)的差异较大。在回归中不考虑XSNav 9和10时,所有沟谷维度的最佳拟合线斜率都为1。图13d、e和f展示了这些情况。
图13
EP分析表明,在处理复杂的几何操作(如XSNav数据集中)时,基于人工的解释容易在XS维度上产生更大的不确定性。本研究中专家们的面积估计值(侵蚀评估的关键)在很多情况下相差超过20%,在某些特殊情况下甚至高达60%。这些不确定性的大小按比例影响着对沟壑侵蚀的估算,导致体积误差远远超出误差规格。
因此,不同的操作人员在界定横截面边缘时存在很大的主观性,因为他们在一些方面采用了不同的不明确标准,例如沟谷边缘是什么样子的、是划定第一个侵蚀周期还是以后的侵蚀周期,以及坡度的变化应该有多大才能界定沟谷的边缘。因此,根据这些内部规则,一些操作人员倾向于做出更保守的估算,从而得出更大的宽度和面积值,而另一些人则更严格地定义横截面边缘,并在剖面中寻找明显的区别。操作员之间的这些差异可能源于他们的专业背景、熟悉的侵蚀过程类型或个人偏好。
OLid是定义横截面边缘的自动替代方法,有助于避免专家小组方法的主观性。该方法的基础是优化横截面边缘的关键通用属性,这些属性可能适用于不同的地貌沟道,因此无需定义阈值。它的应用有助于为研究人员提供一个共同的基础,以便比较侵蚀沟体积。有了一个简单的OLid脚本,就可以批量模式划定沟谷横截面,从而简化操作流程。
通过对具有截然不同截面特征(图6)的复杂剖面进行筛选测试,证明OLid在各种情况下的精度都较高。此外,实地测试允许对模型进行最后调整,以产生准确的估计值(W、A和D的平均相对误差低于7%),同时保持这些预测的高精度。OLid精确度优于其他自动方法。
验证阶段的结果与校准阶段类似。将校准和验证结果结合在一起时,预测结果没有明显的低估或高估偏差。此外,OLid与专家评估中位数相近,这一策略用于消除专家操作员方法的隐含变异性。
尽管其性能良好,但很明显,仅靠地形信息最终无法成功确定横截面边缘。在沟壑的情况下,耕作作业产生的微地形可能与沟壑尺寸处于同一数量级(几厘米甚至几十厘米),这使得EG大小很难确定。因此,在实地收集边缘-沟谷界面的视觉信息(如土壤粗糙度和颜色)可能会为确保准确划定横截面边缘提供重要帮助。由于这些原因,我们认为这种方法更适用于较大的沟道,因为沟道和边缘微地形尺度之间的比例更大,不过今后还需要开展进一步的研究来探讨这个问题。
我们在专家小组的基础上进行了一项实验,测试了人为对沟谷横截面边缘的视觉判读的变异性,在这项实验中没有给出先前的划界标准。发现操作者之间存在显著差异:宽度是受影响最大的变量,其次是面积。该分析表明,采用这一程序沟壑面积的不确定性约为20-40%,这可能导致对侵蚀率的估算不准确。一些专家似乎更容易根据主观划分标准做出偏高或偏低的估计。
在此提出的最优上盖法(OLid)是一种不依赖阈值的自动方法,它为划界提供了一种可靠的替代方法,并避免了与人类主观性的影响。通过使用以下S变量的最大值,能够有效调整该模型,使其在合成和实地切沟的横截面中具有最佳的准确度和精确度性能:
尽管事实证明,与人工方法相比,OLid是一种可靠的方法,但并不能说二维沟壑划定问题已接近最终解决。在许多情况下,要可靠地界定沟谷的边缘,可能需要的不仅仅是几何信息:例如,实地收集现场证据,如对地形、径流侵蚀痕迹或植被存在的视觉解读。此外,由于世界各地的沟壑侵蚀形式多种多样,只有未来的研究才能确定OLid能在多大程度上对不同的沟壑情况做出适当的反应。还需要采取进一步措施来阐明这个复杂的问题,例如设计特定的实验室或实地实验装置,对侵蚀前和侵蚀后的情况进行评估,以验证这一新方法。
作者 | 锐鹏
编辑 | 回毅滢