基于旋翼无人机低空遥感的高原鼠兔危害等级划分技术研究

花蕊，周睿，2，包达尔罕，董克池，唐庄生，花立民*

（1.甘肃农业大学草业学院，草业生态系统教育部重点实验室，国家林业和草原局高寒草地鼠害防控工程技术研究中心，甘肃 兰州 730070；2.青海大学畜牧兽医科学院，青海 西宁 810016）

青藏高原具有丰富的草地资源，不仅为我国畜牧业发展提供了保障，在发挥生态屏障功能及维持生物多样性等方面也具有至关重要的地位［1］。近年来，随着全球气候不断变化和过度放牧等因素影响，青藏高原草地逐渐退化［2］。同时在景观层面上出现了破碎化，部分严重地区出现大面积黑土滩，对草地生态安全和区域经济的发展产生了重大影响［3］。草地退化提高了啮齿动物生境适合度，使其繁殖率与存活率同步提升［4-5］，导致种群密度增大。当种群密度超过一定阈值后，发生草地鼠害，进而使得区域内草地破坏及退化加重［6-7］。

高原鼠兔（Ochotona curzoniae）是我国青藏高原主要害鼠之一。对草地的破坏表现在频繁的挖掘活动使得草地景观破碎化、土地斑块化［8］，导致土壤养分减少、加剧水土流失［9］；大量采食活动使得植被生物量及盖度变化明显［10］。但高原鼠兔的存在对维系草地生态系统良性循环也具有一定的积极作用［11］。高原鼠兔作为组成草地生态系统食物链中的一部分，它们的种群变化对天敌的种群变化有重要影响［12］。同时，翻土掘地行为造成微地形改变，形成植被斑块，可以提高植被及土壤的异质性［13-14］。鉴于高原鼠兔在草地生态系统中的多重功能，在盲目灭鼠或提出保护前，都应该综合分析其对草地生态系统的干扰，科学划分其危害程度，进而根据不同危害程度提出差异化的防治措施。

我国草地鼠害调查主要以人工实地调查为主，最常用的有堵洞开洞法、夹捕法及洞口系数法等［15］。这些调查方法操作简单直观，对鼠害调查人员技术要求低［16］。但是高原鼠兔分布范围广，且栖息地具有强选择性。在大范围调查中，传统方法的人力资源和经济投入需求较高且调查结果的准确性较差。而且我国对于高原鼠兔危害程度的不同等级界定比较模糊。常采用的高原鼠兔危害确定标准为洞口数是否达到150个·hm-2（NY/T 1240-2006）。已有研究认为这种危害评价方法指标单一，难以评估不同草地类型、不同程度的高原鼠兔危害实际情况［17］。对高原鼠兔危害评价的指标应从单一的洞口数向多指标转变，同草地植被群落相结合［17］。可见评价指标的选择对危害评价结果的准确性影响较大。

近年来，随着无人机技术的发展，无人机遥感在鼠害监测中已有初步应用［18］。无人机技术为鼠害监测调查提供了一种新型手段，在中小尺度监测中具有十分显著的优势，可以弥补传统调查中的不足［19］。本研究以无人机为载体结合地面调查，开展基于无人机低空遥感的高原鼠兔危害监测及危害等级划分技术研究，形成一套适用于高原鼠兔危害的低空遥感监测方法。重点解决人工监测费时费力、精度低、范围小等问题。使鼠害及草地调查更加高效，以期为高原鼠兔危害状况监测提供技术支持。建立的高原鼠兔危害等级划分体系可作为参考指标，对研究区高原鼠兔危害程度进行明确的定位，使得该区域鼠害防治做到“分级而治”。

1 材料与方法

1.1 研究区概况

本研究调查区域选择在青藏高原东部甘肃省甘南藏族自治州玛曲县（100°45′45″E-102°29′00″E，33°06′30″N-34°30′15″N）。地势西部高东部低，总体东南偏低，总面积达100190.8 km2，平均海拔为3600 m。县内约有91.06万hm2草地，占全县总面积的89.4%［20］。植被类型包括山地针叶林、高寒常绿灌木林、高寒草甸、高寒草原等，其中高寒草甸分布面积最广，约占天然草地总面积的68.35%。天然草地植被共有47科413种植物，垂穗披碱草（Elymus nutans）、早熟禾（Poa annua）、鹅绒委陵菜（Potentilla ansrina）和紫花针茅（Stipa purpurea）等为优势物种（图1）［21］。

1.2 数据采集

1.2.1高原鼠兔分布点调查 于2019年5月进行高原鼠兔的野外调查。选取玛曲县高原鼠兔主要分布乡镇，将其网格化，每个网格大小为3 km×3 km。在每个乡镇内随机抽取鼠兔危害面积大于10%的网格作为调查样地，调查取样点见图1。使用大疆精灵3无人机在20 m高度进行航拍，每个样地拍摄照片数量约为230张，覆盖面积为100 m×100 m，用GPS记录相应空间位置。

图1 高原鼠兔调查取样点Fig.1 Sampling sites of plateau pika in Maqu county

1.2.2地面数据采集 每块调查样地中，用步测法调查样地内植被总盖度［22］，在每个样地内设置3个0.25 m2的样方，用样方框刈割法测定地上生物量，并在65℃下烘干称重。依据《中国草地饲用植物资源》［23］辨识样地可食牧草种类，采用样圆法［24］测定可食牧草频度；统计样地鼠洞个数。在调查时记录每个样地经纬度（WGS84坐标系）和海拔等信息［25］。同时在样地内选择植被、鼠洞、裸地等典型地物，用相机拍摄实际照片，以便后期将无人机图片中的地物与实际地物进行对比观察。

1.3 数据分析

1.3.1生物量估算模型的建立 本研究通过植被指数来估算植物生长状态，基于无人机航拍得到的可见光影像，选取3种不同植被指数进行计算。分别为红、绿、蓝3个波段参与的可见光波段差值植被指数（visible-band difference vegetation index，VDVI）［26-27］；红、绿 波 段 参 与 的 归 一 化 绿 红 差 值 指 数NGRDI（normalized green-red difference index）［28-29］；绿、蓝波段参与的归一化绿蓝差值指数（normalized green-blue difference index，NGBDI）［30］。不同植被指数计算方式如下：

式中：VDVI、NGRDI、NGBDI分别为不同植被指数，ρRed、ρGreen、ρBlue分别为红、绿、蓝3个波段的反射率或像元值。

在ENVI 5.1软件中对植被指数进行计算后，利用ArcGIS提取每个调查点的植被指数值，将实测生物量分别与不同植被指数进行相关性分析，选取相关性最强的植被指数；然后，将该植被指数作为自变量，生物量作为因变量，建立一元、二次、指数、对数模型，进而模拟分析其之间的关系，选取最优模型对研究区地上生物量进行反演计算。在模型建立过程中，将70%的数据用于建模，其余30%的数据用于精度验证。计算与实测值的平均相对误差（mean relative error，MRE）、均方根差（mean square error，RMSE）及相关系数R，综合评价模型精度：

式中：MRE为平均相对误差；y为研究区生物量地面调查实测值；yi为生物量预测值；n为样本数量；RMSE为均方根差。

1.3.2可食牧草比例计算 对样地利用样圆法进行物种调查，分别记录每一物种出现的频度。统计后计算可食牧草比例（edible forage proportion，EFP），公式如下：

式中：r为可食牧草总数；n为物种总数；A i为群落中某一种可食牧草频度；B i为群落中某一物种的频度。

1.3.3高原鼠兔危害等级划分 在计算前需要对原始各指标进行标准化处理。本研究采用的数据归一化方法为最小-最大标准化法（min-max normalization），转换函数为：

式中：X i为指标进行标准化后的值；x i为初始指标；xmax为整体样本实测值的最大值；xmin为整体样本实测值的最小值。

为了降低评价者的主观误差，采用主成分分析法（principal components analysis，PCA）来确定各指标在退化评价中所占的权重［31］。基于主成分分析法确定权重，首先需要检验数据是否适合于主成分分析。对数据进行KMO（Kaiser-Meyer-Olkin）和Bartlett检验，KMO越接近1，说明变量之间相关性越强，原有变量适合做主成分分析；Bartlett的球度检验值小于显著水平0.05，说明变量之间存在相关关系。得到各指标的权重后，进行加权，得到草地危害指数DRI（damage rangeland index）［32］。DRI值越大，说明草地状况越好，鼠害程度越轻；DRI值越小，代表草地状况越差，受危害程度越重，计算公式为：

式中：DRI为草地危害指数；X i为各指标实测数据的标准化值；W i为各变量的权重；n为变量的个数。

1.3.4洞口、植被盖度信息解译 使用草地鼠洞数字化识别软件v1.0［19］对鼠兔洞口数及样地植被盖度进行提取。

2 结果与分析

2.1 生物量估算模型建立

利用1.3.1中公式（1）～（3）分别计算样地植被指数。以其中3个样地为例，得到植被指数空间分布（图2）。VDVI及NGBDI指数的植被与非植被区亮度区别较明显，白色较亮区域为草地，裸地区域则比较暗。3种植被指数NGRDI、NGBDI和VDVI范围为［-1，1］。

图2 各植被指数计算结果Fig.2 Results of each vegetation index

选取NGRDI、NGBDI和VDVI为自变量，地上生物量为因变量，建立草地生物量的估算模型方程（图3）。通过对比，VDVI指数与实测地上生物量具有较好的正相关关系（R=0.87），其次为NGRDI（R=0.74），三者中最差为NGBDI（R=0.66）。据此，本研究选取可见光波段差异植被指数VDVI来反映研究区生物量状况。

图3 生物量与植被指数相关性关系Fig.3 Corr elation between biomass and vegetation index

进一步分析表明，基于VDVI指数构建的二次多项式与指数模型生物量估测模型较其他几个植被指数建立的模型相关性高（图4）。考虑到构建模型的简单化，在本试验中，二次多项式作为估测草地生物量的较优模型，建立的回归方程为：

图4 生物量与可见光波段差异植被指数的拟合Fig.4 Corr elation between biomass and VDVI

式中：y为草地地上生物量（g·m-2），x为可见光波段差异植被指数VDVI。

2.2 模型精度检验

模型构建后对选取的函数模型进行检验，采用未参与建模的13个样地进行验证，将实测值与估算模型估测得到的地上生物量进行了对比分析，并计算了其绝对误差以及精度。估产精度大于80%的样方数据有10个，平均精度为86.54%。平均相对误差（MRE）为0.183，均方根差（RMSE）为0.109，这说明利用VDVI指数建立的二次项模型来估算研究区地上生物量是可行的（表1）。

表1 预测值与实测值误差及估产精度Table 1 The absolute error and yield estimation accuracy between simulate data and measuring data

2.3 可食牧草比例计算

研究区可食牧草比最大值为98.5%，最小值为63.6%，有11块样地可食牧草比例为78%～81%，所有样地平均数值为81%（图5）。

图5 可食牧草比例频度Fig.5 Frequency of edible forage proportion

2.4 洞口及植被盖度提取

使用草地鼠洞数字化识别软件v1.0，对所有研究样地洞口数量进行提取，随机选取10块样地将实测洞口数与解译洞口数作精度验证（图6），平均相对误差（MRE）为0.05，相关系数（R）为0.95。同理，对研究样地植被盖度进行提取，将实测盖度与解译盖度作精度验证。通过验证发现实测值与预测值之间MRE为0.07，R为0.81。

图6 影像解译结果与实测值相关性分析Fig.6 Correlation analysis between image interpretation results and measured values

2.5 高原鼠兔危害等级划分

通过计算得：KMO（Kaiser-Meyer-Olkin）检验=0.612、P=0.001，符合分析条件，可以进一步进行主成分分析。从总方差解释表（表2）能够得出，最大的主成分对应的特征根值为1.868＞1，前两个主成分的累计方差贡献率达到70.315%，超过70%，因此，鼠害信息基本上可由鼠洞个数与VDVI来反映，原来的鼠洞个数、VDVI、植被盖度、可食牧草比例4个指标可由这两个主成分替换。然后计算两个主成分载荷数及各指标在两个主成分线性组合中的系数（表3）。

表2 主成分特征根值及方差贡献率Table 2 The pr incipal component eigenvalue and var iance contribution rate

第一主成分与第二主成分对原来4个指标（α1，α2，α3，α4）的载荷数见表3。根据各成分载荷数和特征根值，进而得到两个主成分组合如下：

表3 各指标载荷数和系数Table 3 Load number and coefficient of each indicator

对指标在不同主成分线性组合中的系数及方差贡献率进行计算，指标系数可以看成是以鼠洞个数与VDVI的方差贡献率为权重，对指标在鼠洞个数与VDVI两个成分线性组合中的系数作加权平均得到的方差贡献率。最后在综合得分系数的基础上进行归一化即得到各指标权重（表4）。

表4 各指标贡献率和权重Table 4 Contribution rate and weight of each indicator

将综合指数定义为草地生物量、洞口个数、盖度和可食牧草比例的加权，利用公式（8）计算危害指数DRI。采用聚类分析方法依据DRI对研究区高原鼠兔危害状况划分等级，共分为5类。I类为无明显鼠害地，Ⅱ类为轻度鼠害地，Ⅲ类为中度鼠害地，Ⅳ类为重度鼠害地，V类为极度鼠害地（表5～6）。不同等级中指标变化为，洞口数：I＜Ⅱ＜V＜Ⅲ＜Ⅳ；地上生物量：V＜Ⅳ＜Ⅲ＜Ⅱ＜I；植被盖度：V＜Ⅳ＜Ⅲ＜Ⅱ＜I；可食牧草比例：V＜Ⅳ＜Ⅲ＜Ⅱ＜I。

表5 高原鼠兔危害级别划分标准Table 5 Demarcation standard of damage gr ade for plateau pika

3 讨论

在高原鼠兔危害评估中，通常以鼠洞个数为评估的主要指标。但除了洞口个数，植被状况、土壤养分等多种因子均可反映鼠兔破坏对草地产生的影响。本研究选择高原鼠兔洞口个数、草地总盖度、地上生物量及可食牧草比例这4个因子作为高原鼠兔危害评价的指标，指标选取考虑的首要因素是能否从无人机影像中解译得到，其次是考虑鼠害对草地造成危害的重要表现。通过各指标权重计算，得出地上生物量所占权重最大，为56%；其次为植被盖度，所占比重为55%；再次为可食牧草比例，为29%；最小的为鼠洞个数，为5%。可见草地退化的最直观表现为植被盖度和生物量的显著降低。随着危害程度的上升，地上生物量、植被盖度、可食牧草比例等指标均呈降低趋势，主要因为高原鼠兔频繁的挖掘采食活动能直接干扰植被的生长，使得植被破碎化严重，形成斑块景观，导致裸地面积增加，活动干扰大的地区植物地上生物量也随之降低。

洞口个数在危害评价中权重占比最小。同时，随着危害程度的递增，研究区高原鼠兔洞口数表现出起初数量低、后有一定升高、再降低的变化规律。可能由于随着危害程度的增加，区域范围内高原鼠兔密度不断增加，洞口数量同步增多。当高原鼠兔密度发展到一定量时，栖息地各种条件不断恶化，鼠兔死亡率及迁移率增高。当鼠兔离开后，由于降水、家畜踩踏等因素使得很多废弃洞口塌陷被填埋，因此极度危害时其洞口数量又呈下降趋势。通过分析可见，以往用单一的洞口数来评价其危害程度是不准确的。

表6 危害级别中各指标的极值Table 6 Ranges of damage indexes in all damage levels

鼠害的发生并非必然，通常导致高原鼠兔种群数量增多的因素为种群自身繁殖率上升和种群大量迁入两种，但造成这两种情况的共同因素均为某一地区的生境条件适合其生存繁殖［33］。很多研究均证实因为草地植被的退化，产生了适合高原鼠兔生存的环境，进而使其种群数量急剧增长［34］。只有客观的评价鼠害，了解啮齿动物在草地生态系统中的积极作用，才能合理科学地进行防治。以往将有鼠兔分布的区域就定义为鼠害地区的做法并不客观，在处于轻度危害时，生态系统内部可以通过生物机制来自我控制较高种群数量的增加，无须人为治理干预。只有当其发生到一定程度，危害到生态系统稳定性或草地生产力时才能算作危害。因此，一定要根据实际危害程度进行治理，在危害较轻地区采用监测防控措施，在危害严重区域采用科学的治理措施并对未来种群变化加以预测，使得鼠害防治真正的有效化和科学化。

最后，对比无人机航拍监测法与传统野外人工调查法。以1 hm2样地为例，传统人工调查洞口数需要2～3人，且容易重复计数，不适用于大面积多重复的调查。对于生物量监测，1 hm2样地平均需要3人，耗时在1 h左右，且传统方法需要对草地造成一定破坏。而无人机监测仅需1人，耗时约15 min，在不破坏草场的情况下极大地提高了工作效率。同时，通常生物量调查主要在盛草期进行，青藏高原紫外线极强，严酷的工作环境会对人工监测产生一定影响，无人机监测仅操纵飞行器进行航拍，无需大量体力劳动，后续工作均可室内完成，很大程度的解放了劳动力。

本研究也存在一定的不足，通过无人机遥感技术可以快速提取计算植被盖度、生物量、鼠洞等指标，但却无法对草地植物群落中的物种进行很好的识别。随着计算机智能学习技术的不断发展，后期可以在植物物种识别方面展开研究，在识别过程中对不同物种形态及不同光谱反射率的分析及综合特征库的建立是关键。或许可通过无人机搭载高光谱成像仪，针对高光谱图像的不同光谱特征来设计不同物种自动分类算法，让无人机监测代替人工进行分类，使得单张影像可获取信息增多，完善危害评价指标。

4 结论

本研究以高原鼠兔为例，分别选择处于不同危害程度的典型样地，利用无人机低空遥感监测，与同一时期的地面调查数据相结合，对样地高原鼠兔洞口数、地上生物量、草地盖度及可食牧草比例进行调查。通过自主设计软件对无人机部分指标进行解译，最后综合多种指标对研究区不同危害等级进行划分，通过研究得到以下结论：

1）本研究探讨了基于无人机低空遥测技术的草地啮齿动物危害监测新方法。通过试验证实无人机监测的各项指标与实地调查结果基本一致，精度较高。说明在草地啮齿动物危害监测中，无人机可以取代部分传统人工监测指标。试验从外业监测到内业处理整体的高效性符合现代用机器智能代替人力的发展新形势，使得鼠害及草场调查更加轻简化与高精度化。

2）通过归一化绿红差值指数NGRDI、归一化绿蓝差值指数NGBDI和可见光波段差值植被指数VDVI 3种植被指数对比，可见光波段差值植被指数VDVI与实测地上生物量具有较好的正相关关系。上述研究结果均表明基于可见光波段的植被指数可以准确地提取无人机影像中植被信息，可以较好地估算地上生物量。

3）试验证明单靠鼠洞数量来评估鼠害程度的方法具有一定的弊端，多项指标的综合测定才是鼠害程度的最准确评估。本研究基于综合指数将研究区鼠害状况划分为5个等级，后续在鼠害防治中要避免一刀切式灭鼠，应做到因级而治。