基于分层抽样的水稻种植面积回归估计

申克建　王飞　裴志远

摘要：随着农业统计需求和遥感技术的发展，中国农业部即将使用中分辨率遥感数据开展中国主要农作物的全覆盖遥感调查，如何修正全覆盖调查结果使其更接近农作物面积真值是需要研究的问题。以水稻为例，在区域尺度上利用基于样本无人机数据解译得到的水稻种植面积修正基于中分辨率数据解译得到的水稻种植面积，分别比较分层分别回归估计和分层合并回归估计，利用基于全覆盖无人机数据解译得到的水稻种植面积验证回归估计值。在给定5%误差和95%置信水平下，抽样比为036%（最小样本量公式换算）的抽样结果表明：（1）2种方法的估计精度都大于95%；（2）合并比估计的精度好于分别比估计。

关键词：水稻种植面积；中分辨率数据；无人机数据；分层抽样；回归估计

中图分类号： S127文献标志码：

文章编号：1002-1302（2017）16-0206-05

收稿日期：2016-10-11

基金项目：国家自然科学基金（编号：41301506）；国家重点研发计划（编号：2016YFB0501505）。

作者简介：申克建（1982—），男，河北廊坊人，博士，工程师，研究方向为农业遥感。E-mail：ashenkejian@126com。

农作物面积是各级管理部门进行生产指导、农业补贴、农业保险核查等决策的重要依据，也是农作物产量预测和环境管理等研究的重要依据。这些管理和研究需要解决的重要问题是提高农作物面积估算精度。遥感技术结合地理信息系统（GIS）和全球定位系统（GPS）等现代地球信息技术能提供及时、客观的农作物面积信息。目前，30 m格网尺度、20～30个类别的土地覆盖类型监测总体精度大约在80%～85%之间，高分辨数据的分类精度超过95%[2-3]。随着地球资源卫星对地观测能力不断提高，中分辨率数据（10～60 m）正逐渐可以保证大尺度全覆盖，但是由于混合像元等问题不能满足农业统计等行业应用的精度要求[4]。高分辨率数据（<5 m）可以得到符合要求的行业应用精度，但它受幅宽和重访周期的局限，无法满足调查区域的影像全覆盖，甚至不能提供满足抽样要求的影像[5]。在3S技术（指遥感、地理信息系统、全球定位系统）的支持下，遥感和抽样调查相结合的空间抽样调查方法得到了广泛应用，地面调查是获取抽样样本的重要手段，但地面调查存在强度大、费用高等不足。近些年来，无人机被广泛应用于民用领域，它具有机动、灵活的特点，有降低地面调查强度和保障高分辨数据获取的优势。Laliberte 等应用无人机遥感数据进行牧场监测，证明无人机遥感数据可以补充或者替代部分地面调查数据[8]。Breckenridge研究证明，无人机可以在较短的时间内采集大范围影像，比地面调查的效率高很多，他预见无人机对牧草地监测的影响，将像全球定位系统影响导航和野外数据采集一样[9]。

目前遥感与抽样相结合的空间抽样方案设计被广泛应用[10-13]。例如，美国农业部现行农作物面积估算方法是将遥感全覆盖结果和6月调查结果进行回归分析，回归模型是Y=a+b×X。其中Y是估计面积；X是全覆盖作物遥感分类面积（cropland data layer classified acres，簡称CDL），6月调查数据是全美11 000个地块的地面调查信息（约占国土面积 250%，耕地地块约为161 km2，城市区域地块约为 016 km2，未利用地块为 644～1288 km2）；a和b是基于6月调查数据和CDL样本数据使用最小二乘法得到的[14]。孙佩军等使用分层分别回归估计估算河南省冬小麦面积，使用无人机解译的冬小麦面积样本进行事后分层，辅助信息来自全覆盖中分辨率数据（Landsat TM5和HJ）解译的冬小麦面积[15]。

随着农业统计需求和遥感技术的发展，中国农业部即将使用中分辨率遥感数据开展中国主要农作物的全覆盖遥感调查，如何修正全覆盖调查结果使其更接近农作物面积的真实值是需要研究的问题。本研究以水稻作物为例，研究在区域尺度上利用基于样本无人机数据解译得到的水稻种植面积修正基于中分辨率数据解译得到的水稻种植面积，旨在为修正中分辨率全覆盖遥感调查结果提供参考抽样方法。

1材料与方法

11研究数据

研究数据位于江苏省盐城市，选取理由：（1）该区域位于水稻主产区；（2）在水稻生长季节内获取了无云覆盖的光学遥感数据SPOT5；（3）该区域不是无人机禁飞区。

研究数据主要包括无人机数据和SPOT5多光谱数据。无人机数据拍摄于2014年9月24日至2014年9月27日，经过处理后覆盖面积12 16600 hm2，一共5条样带，每条样带长度约300 km，宽度07～10 km，通过人工目视解译得到水稻种植面积5 23774 hm2，基于样本无人机数据解译得到的水稻解译结果简称UAV-Rice，SPOT5多光谱数据获取于2014年8月3日，该数据经过裁剪后覆盖和无人机相同的范围，通过人工目视解译得到水稻种植面积6 00504 hm2，基于SPOT5数据解译得到的水稻解译结果简称SPOT5-Rice，数据示意见图1，数据用途说明见表1。

12研究方法

121总体思路

为了用样本UAV-Rice修正SPOT5-Rice，本研究首先基于SPOT5-Rice构建抽样框，并按水稻种植面积规模分为6层；然后按回归估计最小样本量公式计算最小样本量，按等比例分配将样本分配到各层；随后用分别比估计和联合比估计进行回归估计，最后用UAV-Rice进行回归估计精度评估（图2）。

122抽样框准备

基于SPOT5-Rice构建100 m×100 m格网，将格网和SPOT5-Rice进行空间分析，统计每个格网内的SPOT5-Rice种植面积。

具体是通过ArcGIS Identity工具， 用SPOT5-Rice矢量

数据叠合成100 m×100 m矢量格网，这样对100 m×100 m内的每个地块都给定了格网编号；最后用Dissolve工具以格网编号为依据，使每个100 m×100 m中的水稻合并成一个整体；最后统计每个100 m×100 m中水稻的种植面积，完成抽样框构建，N=8 627，这样每个抽样单元都有来自SPOT5-Rice的种植面积。同理将UAV-Rice矢量数据叠合抽样框后的数据再融合，这样每个抽样单元都有来自UAV-Rice的种植面积。

最后，抽样框中的每个抽样单元属性表包含抽样单元编码、来自SPOT5-Rice的种植面积和来自UAV-Rice的种植面积，如图1所示。

123抽样框分层标志与层数的确定

本研究选择以抽样单元中水稻种植面积为分层标志。一般分层数越多，抽样统计越准确，但抽样方差的降低是与分层数的平方成反比的，分层抽样理论要求每层至少必须抽取2个样本单元，层数不超过总样本量的一半，但在层数大于6层时方差的减少幅度将大为减缓，理论与实践研究表明，层数以不超过6层为宜[16]。本试验将层数定为6层，分层界限确定采用累计频数直方图法。

124抽样样本量设计

在给定研究变量均值允许误差（本研究设定为5%）和置信水平1-α（本研究设定为α=005）下，按最小样本量[17]计算公式：

[JZ（]n=[SX（]n01+[SX（]n0N[SX）][SX）]；n0=[SX（]u2α/2·S2y·（1-ρ2）Δ2[SX）]。[JZ）][JY]（1）

式中：n为无放回抽样的样本量；n0为有放回抽样得到的样本量；N为总体单元数；uα/2为正态分布的上侧分位数；Δ为研究变量均值允许误差；S2y=[SX（]1N-1[SX）]·∑[DD（]Ni=1[DD）]（Yi-Y[TX-]）2，为研究变量（无人机水稻种植面积）总体方差；ρ=[SX（]SxySx·Sy[SX）]，为总体相关系数，其中辅助变量（中分辨率SPOT5水稻种植面积）总体标准差为Sx=[KF（][SX（]1N-1[SX）]·∑[DD（]Ni=1[DD）]（Xi-X[TX-]）2[KF）]，总体协方差为Sxy=[SX（]1N-1[SX）]·∑[DD（]Ni=1[DD）]（Xi-X[TX-]）·（Yi-Y[TX-]），其中X[TX-]为辅助变量总体均值，Y[TX-]为研究变量总体均值。

125分层随机抽样下的回归估计与精度评估

分层回归估计按分别回归估计和合并回归估计2种方法分别计算[17]，2种方法各层样本量都按等比例分配，2种方法都按相同精度评估方法评估。

1251分别回归估计

设第i层的目标变量和辅助变量分别为yi和xi，其均值分别为y[TX-5]i和x[TX-5]i，层权重为Wi=Ni/N，辅助变量层总体为Ni，每一层的回归系数为bi，辅助变量层样本方差为s2xi=[SX（]1n-1[SX）]·∑[DD（]ni=1[DD）]（xi-x[TX-5]i）2，sxiyi=[SX（]1n-1[SX）]·∑[DD（]ni=1[DD）]（xi-x[TX-5]i）·（yi-y[TX-5]i）为层样本协方差。

各層均值的回归估计量：

式中：r为估计误差；Y[DD（-1][HT6]^[DD）]为分别比估计或合并比估计的估计值；Y为真值（UAV-Rice）。

2结果与分析

在给定均值允许误差5%和95%置信水平下，根据公式（1）计算得到n为31，将31个样本按等比例分配到6层中，按2种回归估计的公式计算得到表2所示结果，详细计算过程见表3。

从表2可以看出，|025%|<|-103%|，|-333%|<|-457%|，说明合并比估计结果都好于分别比估计。同时2种方法估计精度都大于95%。

SPOT5-Rice解译水稻存在漏分，通过无人机数据对比

分析，漏分面积为18709 hm2，漏分面积占实际面积 5 23774 hm2 的357%。357%与合并比估计的误差 -333% 很接近，说明中分辨率数据的漏分精度直接影响合并比回归估计的精度，那么在实际中中分辨率数据解译水稻时要特别减少漏分。

3结论与讨论

本研究使用中分辨率分类数据构建抽样框和分层，目的是在实际应用中可指导无人机样带的布设。本研究在以下几个方面值得探讨：抽样框覆盖误差，本研究方法采用中分辨率数据（SPOT5）解译的水稻空间分布设计抽样框，为了达到95%的农业统计精度要求，根据上一段合并比回归估计误差和中分辨率数据（SPOT5）漏分误差接近的分析，其漏分误差保守估计要≤5%。

样本量的确定：本研究采用基于无人机分类数据的总体方差S2y，计算得到最小样本量是31个，转为最小抽样比是036%。基于无人机分类数据的总体方差S2y在实际抽样调查中常是难以获取的，由于中分辨率结果和无人机结果具有很强的相关性，研究用中分辨率（SPOT5）的总体方差S2x替代无人机的总体方差时，计算得到的最小样本量是34个，比原来的31个增加了3个，增涨了968%。

样本量的分配：本研究采用等比例法分配样本，将来试验可以考虑奈曼分配和调查费用的最优分配。

分别比估计和合并比估计比较：本研究中分别比估计精度不如合并比估计精度的结论，部分验证了倪加勋的结论[17]：分层分别回归的效率高于分层合并回归的估计效率，但在样本量较小的情况下，分层分别回归估计可能使偏差增大。

合并比估计回归系数的确定：本研究中可以采用基于总体数据得到的准确回归系数，但考虑到本研究方法实际用于推广，采用了基于样本数据的回归系数估计值。

本研究设计了基于中分辨率数据和无人机数据水稻种植面积分层回归估计方法，该方法在抽样框设计、样本量确定、样本量分配、分别比估计和合并比估计的选用、合并比回归估计系数确定方面都给出了参考，可以为农作物面积遥感抽样调查方案设计提供参考。

参考文献：

Chauhan H J，Arora M K，Agarwal A Estimating land cover class area from remote sensing classification[J] Journal of Applied Remote Sensing，2008，2（1）：183-198

张磊，吴炳方 关于土地覆被遥感监测的几点思考[J] 国土资源遥感，2011，3（1）：15-20

[3]刘慧平，朱启疆 应用高分辨率遥感数据进行土地利用与覆盖变化监测的方法及其研究进展[J] 资源科学，1999，21（3）：25-29

[4]Gallego F J Remote sensing and land cover area estimation[J] International Journal of Remote Sensing，2004，25（15）：3019-3047

[5]吴炳方，李强子 基于两个独立抽样框架的农作物种植面积遥感估算方法[J] 遥感学报，2004，8（6）：551-569

[6]吴炳方，蒙继华，李强子 国外农情遥感监测系统现状与启示[J] 地球科学进展，2010，25（10）：1003-1012

[7]范承啸，韩俊，熊志军，等 无人机遥感技术现状与应用[J] 测绘科学，2009，34（5）：214-215

[8]Laliberte A S，Herrick J E，Rango A，et al Acquisition，orthorectification，and object-based classification of unmanned aerial vehicle （UAV） imagery for rangeland monitoring[J] Photogrammetric Engineering & Remote Sensing，2010，76（6）：661-672

[9]Breckenridge R P Improving rangeland monitoring and assessment： integrating remote sensing，GIS，and unmanned aerial vehicle systems[D] Idaho Falls： University of Idaho，2007

[10]张焕雪，李强子，文宁，等 农作物种植面积遥感抽样调查的误差影响因素分析[J] 农业工程学报，2014（13）：176-184

[11]刘国栋，邬明权，牛铮，等 基于GF-1卫星数据的农作物种植面积遥感抽样调查方法[J] 农业工程学报，2015，31（5）：160-166[HJ17mm]

[12]Gallego F J，Kussul N，Skakun S，et al Efficiency assessment of using satellite data for crop area estimation in Ukraine[J] International Journal of Applied Earth Observation and Geoinformation，2014，29：22-30

[13]王迪，周清波，陳仲新，等 玉米种植面积空间抽样调查方案优化设计[J] 农业工程学报，2014，30（8）：117-125

[14]Claire B，Yang Z W，Rick M，et al Monitoring US agriculture： the US department of agriculture，national agricultural statistics service，cropland data layer program[J] Geocarto International，2011，26（5）：341-358[HJ]

[15]孙佩军，张锦水，潘耀忠，等 基于无人机样方事后分层的作物面积估算[J] 中国农业资源与区划，2016，37（2）：1-10

[16]杜子芳 抽样技术及其应用[M] 北京：清华大学出版社，2005：187-191

[17]倪加勋 抽样调查[M] 桂林：广西师范大学出版社，2003：61-84

[18]刘建红，朱文泉 耕地变化空间抽样调查方案的精度与效率分析[J] 农业工程学报，2010，26（10）：331-336