CN102156886A - 基于统计数据和遥感影像数据的区域化肥投入空间化方法 - Google Patents

Abstract

本发明公开了一种基于统计数据和遥感影像数据的区域化肥投入空间化方法，包括以下步骤：S1、获取归一化植被指数NDVI时间序列数据，并进行平滑处理；S2、利用所述数据提取分类特征，并采用决策树分类方法对种植模式进行分类；S3、在不同种植模式下进行NDVI时间序列主成分分析，并提取多个主成分；S4、对化肥投入数据和其中两个主成分进行尺度匹配，并利用尺度分配的结果分析不同种植模式下化肥投入数据与这两个主成分的关系，选取能够反映化肥投入所造成的生产力的差异的主成分进行化肥投入的空间化。本发明能够从本质上全面地反映化肥投入的空间分布，且能够获得很好的区分度。

Description

基于统计数据和遥感影像数据的区域化肥投入空间化方法

技术领域

本发明涉及农业非点源污染和遥感空间化技术领域。尤其涉及一种统计数据和遥感影像数据相结合的区域化肥投入空间化方法。

背景技术

化肥是农业生产中一项至关重要的生产性资源，在粮食增产中起着重要的作用。化肥投入的大幅度增加，在不断提高作物产量的同时，还带来了水体富营养化、土壤物理和化学属性恶化、地下水污染、生物多样性丧失等一系列严重的生态环境问题，是农业非点源污染的主要类型之一。加上人们越来越关注环境保护和食品安全问题，将来随着生产环境标准和食品检疫标准的逐渐严格，化肥的使用将受到严格限制。因此，研究区域化肥投入存在的空间差异性，不仅是对农田化肥投入进行优化配置，实现农田集约化、科学化管理的基础，更是农业非点源污染研究中建立分布式模型的需要。

基于统计数据的农业化肥投入属于社会经济数据，社会经济数据空间化是当前地理科学、社会科学共同关注的焦点之一。社会经济数据以区域为单位进行统计，而这些统计单元一般都是一些大小各异的行政单元，在这些行政单元内一个单元赋予同一个数值，生成区域内指标的均匀分布，不同区域间突变分布的特点，这显然不能代表统计数据的实际分布情况。在实际的空间分布中，统计单元内的数据常常不是均匀分布的，相邻统计单元的边界处分布也是连续的。这些数据妨碍了经济地理学与自然地理学的数据共享，抑制了人文地理学理论的发展及其对GIS等新技术的运用。要解决这些问题，迫切需要建立一定分辨率的基础地理单元，将社会经济数据转化到一个可以方便操作、分析的数据平台。社会经济数据数据空间化就是最好的解决方法之一，即将以行政区域为单元的统计数据展布到一定尺寸的地理网格上。

社会经济数据空间化方法的一个重要特征是单元内的容量不变性，即空间分配过程不影响原有属性值对整个单元总量或平均值的表达。空间化的方法可以归纳为以下四种：(1)密度分配：通过分布函数、控制区或控制函数的方法将统计指标的数量以密度分配的方式分配到空间上。(2)邻域关系：利用地理变量在空间上的相互依赖关系，通过邻域的值来估计所在位置的值。可通过简单的插值方法或空间自相关方法来估计。(3)地理变量之间的相关关系：通过一个变量的分布差异去估计另一个变量的空间差异。这种关系可能是非线性和不确定的。(4)分形方法：分形具有自相似和自仿射特征，其局部通过变换与整体相等或相似。根据分形分析检验地理实体分布的分形特征并提取参数，然后通过分形模拟的方法可能实现对密度分布的估计。社会经济数据根据属性变量性质不同，可以将变量分为广延量(extensive)和强度量(intensive)，当区域被分成一系列子区域，如果区域的值是子区域的加权平均值，这个变量就叫广延量，如果区域的值是子区域的和，这个变量就叫强度量。

社会经济数据空间化的研究目前主要集中在人口数据的空间化、GDP(Gross Domestic Product，国内生产总值)空间化以及空间化方法的改良上。处理这一类问题的方法主要是使用面域加权，实际上是多边形叠加(polygon overlay)，就是将目标区和源区叠加，分别计算各交叉区域的属性值，再按目标区计算，也有学者将其称为“比例分配”。由于该方法不能只考虑面积作为分配的权重，于是学者们开始使用面积以外辅助的信息，即“借力”，使用控制区的面域加权或者是使用回归关系的面域加权。目前世界上这种“借力”方法的研究最主要集中在人口数据的空间化上。社会经济数据具有区别于自然地理数据的特点，应用辅助数据研究面插值是未来发展的趋势，特别是遥感数据作为辅助数据的研究有着很大的发展潜力。

在对化肥投入的空间变异研究中，学者们研究的着重点在于大尺度上(国家尺度)的化肥施用对粮食产量的影响，有的是直接进行定性的分区，有的是在行政界线的基础上对投入进行分区。还有的是划分为高投入、中投入、低投入，每种投入都有具体的定义，可以认为是一种半定量或者定性的方式。在将化肥投入空间化到一定的网格尺度方面，目前研究几乎还处于空白。农田系统是一个非常复杂的系统，农业投入受到诸多因素的影响，它有别于人口数据，简单地应用辅助数据进行面域加权的数学方法还不能够从本质上全面地反映化肥投入的空间分布。

发明内容

(一)要解决的技术问题

本发明要解决的技术问题是：如何提供一种能够从本质上全面地反映化肥投入的空间分布，且获得较好的区分度的化肥投入空间化方法。

(二)技术方案

为解决上述技术问题，本发明提供了一种基于统计数据和遥感影像数据的区域化肥投入空间化方法，包括以下步骤：

S1、获取归一化植被指数NDVI时间序列数据，并对该数据进行平滑处理；

S2、利用所述数据提取分类特征，并根据该分类特征采用决策树分类方法对作物种植模式进行分类；

S3、在不同类型的种植模式下进行NDVI时间序列主成分分析，并提取多个主成分；

S4、对统计获得的化肥投入数据和其中前2个主成分进行尺度匹配，并利用尺度分配的结果分析不同种植模式下化肥投入数据与这两个主成分的关系，利用前2个主成分中能够反映化肥投入所造成的农用地生产力的差异的主成分进行化肥投入的空间化，所述化肥投入数据包括投入的化肥种类数据和数量数据，主成分按照与化肥投入存在的关系的显著性从大到小排列。

其中，步骤S3中的所述主成分分析针对各种植模式的NDVI时间序列分别进行，分析时首先针对每种种植模式提取6个主成分，这6个主成分从不同方面反映了农用地生产力的状况，然后对各种植模式下前2个主成分的特征向量进行分析来认识前2个主成分分别代表的物理意义。

其中，在步骤S4中使用各种植模式下化肥投入数据与该种植模式6个主成分中的前2个主成分进行相关分析，找到各种植模式下与化肥投入存在一定关系的主成分。

其中，步骤S4中进行尺度匹配的方法为：将化肥投入数据和前2个主成分都匹配到像元尺度或乡镇尺度上。

其中，步骤S1中，利用权重滑动平均滤波法Savizky-Golay对MODIS NDVI时间序列数据进行平滑滤波。

其中，利用中分辨率成像光谱仪获得NDVI时间序列数据。

(三)有益效果

本发明通过统计数据和遥感影像数据相结合进行区域化肥投入空间化，试验证明，该方法从本质上全面地反映化肥投入的空间分布，且能够获得很好的区分度。

附图说明

图1是本发明的方法流程图；

图2是NDVI时间序列平滑前后对比图；

图3是北京北部山区2000年作物种植模式分类图；

图4是对不同种植模式的主成分分析得到的前2个成分的特征向量曲线图；

图5是不同种植模式的第一主成分(也称为主成分一)与化肥投入散点图；

图6是不同种植模式的第二主成分(也称为主成分二)与化肥投入散点图；

图7是三种种植模式下化肥投入空间分配结果图。

具体实施方式

下面结合附图和实施例，对本发明的具体实施方式作进一步详细说明。以下实施例用于说明本发明，但不用来限制本发明的范围。

如图1所示，本发明包括以下步骤：

首先利用MODIS NDVI时间序列数据(为遥感影像数据)进行作物种植模式的分类。包括MODIS NDVI时间序列数据获取及处理、作物种植模式的分类和分类精度评价。

S1、MODIS NDVI时间序列数据的获取与处理

植被的遥感数据中包含了叶子及其生长状况等综合信息，这些信息直接指示作物干物质的积累以及生物量的多少。一些研究指出NDVI与作物产量、生物量有较好的相关性，并且用于表征农田生产力。因此，本发明采用中分辨率成像光谱仪(MODIS)的NDVI时间序列表达作物一年初级生产力的分布规律。

首先，获取2000年MODIS 250米分辨率的8d(8天)合成地表反射率产品，该数据来源于NASA的MODIS数据产品分发网站(http://modis.gsfc.nasa.gov/)，包括了第一、第二波段两个波段。利用MODIS产品专门处理工具MRT(Modis Reprojection Tool)软件对下载的产品进行文件格式、空间镶嵌和地图投影的转换，最后进行掩模得到研究区域的影像图。

其次，利用NDVI计算公式对该影像图进行计算得到MODIS NDVI时间序列数据。

NDVI是近红外与红色通道反射比值的一种变换形式。许多研究表明NDVI与LAI、绿色生物量、植被覆盖度、光合作用等植被参数有关，是植被生长状况及植被覆盖度的最佳指示因子。其计算公式为：

$\mathrm{NDVI}=\frac{{\mathrm{\ρ}}_{\mathrm{NIR}}-{\mathrm{\ρ}}_R}{{\mathrm{\ρ}}_{\mathrm{NIR}}+{\mathrm{\ρ}}_R}$

式中：ρ_NIR代表近红外波段的反射率，ρ_R代表红光波段的反射率；由公式可以看出NDVI的取值范围在(-11)开区间。负值代表地面覆盖为云、水、雪等，对可见光高反射；0代表有岩石或裸土等，ρ_NIR和ρ_R近似相等；正值表示有植被覆盖，且随覆盖度增大而增大。在ENVI(The Environment for Visualizing Images)软件中采用band math功能对所有时期的反射率影像进行计算，并进行波段叠加，可以得到研究区每个像元的NDVI时间序列影像数据集。

最后，对NDVI时间序列进行平滑处理。理论上NDVI时间序列的变化可以作为地表条件变化尤其是植被变化的指示，但是研究表明，利用MODIS、AVHRR(先进的高分辨率辐射计)等中低分辨率的传感器进行观测时，太阳光照角度、观测视角以及云的条件随时间变化都很大，因此得到的观测值是来自地表的双向反射率信息，包含很多不可预测的噪音，严重影响了土地覆被和陆地生态系统的变化监测结果。所以需要对NDVI时间序列进行重建，即通过一系列预处理方法降低遥感NDVI合成数据中的噪声水平，提供更为可靠的时序数据集。本发明利用现有技术Savitzky-Golay滤波法对MODIS-NDVI时间序列进行去噪和重构处理。该方法有两点假设：第一，NDVI数据的时间序列变化对应于植被的生长与衰落；第二，云和大气的影响总是降低NDVI的数值，即NDVI与植被变化变缓过程不一致的突降都应该作为噪声消除。基于以上假设，该方法可应用于不同时间尺度，对传感器也没有限制。

本发明在ENVI中的IDL(Interactive Data Language)平台中应用Savizky-Golay滤波方法，编写程序对北京北部山区2000年NDVI时序数据逐点进行重构处理，平滑前后时序数据效果见图2，图中是某单个像元NDVI值按照一定时间顺序排列，以时间为横坐标轴排，得到的该像元NDVI时序的动态曲线。从图2中平滑前后曲线对比可以看出：平滑前时序曲线有严重的锯齿状波动，如果直接用于时间序列分析很难比较准确地捕捉作物生长趋势，对建立分类模型也造成一定的难度；重构后的生长曲线有效消除了云和缺失数据的影响，但同时基本上保持了原有曲线的基本形状，较为真实地恢复了植被指数曲线。因此，重构后的生长曲线可以更清楚地反应作物生长曲线变化趋势和规律，有利于分类模型的建立。

S2、作物种植模式分类

本发明主要利用决策树分类器对北京市北部山区的作物种植模式进行了研究，详细步骤论述如下：

为可以快速、实时提供大中尺度土地覆盖信息，需要简单、运行速度快且有精度保证的分类方法。对于经典的监督分类算法，训练样本的选择是决定分类精度的关键，训练样本不仅要全面代表研究区各地类，而且样本的数量和纯度均要有所保证。对全球、区域尺度的土地覆盖遥感分类时，训练样本的选择则更为费时、费力，因此监督分类法不适于大中尺度的业务化土地覆盖遥感分类。与监督分类法相比，非监督分类方法更适用于大中尺度时间序列数据土地覆被业务化分类，它无需选择训练样本，只需定义类别数、循环迭代次数等几个有限参数。然而，区域或全球尺度的空间范围大，地物类型复杂，如果只用非监督分类方法对全区域分类，随机因素将加入统计特征中，定义的分类数目要得到保证，且需要不断的调整参数。这样不仅增加分类运行时间及分类后的处理工作量，而且影响最终的分类精度。除上述两种传统分类方法外，结合专家知识的遥感分类方法--决策树分类方法，将地物由主到次，由复杂到简单，逐级分解为均质的子集，该方法以其非参数和树结构特性在全球及区域土地覆盖/植被覆盖遥感制图中具有良好的稳健性和鲁棒性，得到广泛应用。

本发明根据平滑后2000年具有39个波段的NDVI数据进行分类，采用决策树分类方法，在分类过程中，参照2000年土地利用分类图和各地物NDVI时间特征曲线，提取各种NDVI指标作为分类特征，如NDVI最大值、NDVI总和，并比较不同时期的NDVI值，用少数阈值对分类过程加以控制，通过计算边分类边学习建立并改进分类决策。分类得到的最初结果要多于最终需要的类别个数，最后根据分类参考图进行类别的识别及归并，并赋予每一个类别特定的含义。最终将北京北部地区分成了5类：水域及建成区、林地、果园、一年一熟农用地和一年两熟农用地，提取果园、一年一熟农用地和一年两熟农用地三种种植模式(图3)。

接下来还可以进行分类精度评价：

精度评价是遥感分类中不可缺少的工作，理想的精度评价是在抽样理论的基础上，进行实地考察，从野外获取检验数据进行实地检验，但该方法比较昂贵，特别在严格的统计设计要求下进行实地检验更加困难，对于低空间分辨率的影像分类，最简单的精度评价方法是通过获取该地区高分辨率的分类图评价低分辨率影像的相对精度。本发明即是采用后者方法，抽取北京北部山区30m分辨率遥感影像TM对低分辨率影像分类结果作精度检验。

用已预处理过的(几何纠正、拼接、裁切等)的2000年的30m分辨率的TM遥感影像进行精度检验，影像包括昌平、密云、怀柔、延庆、平谷五个区县，通过目视解译选取土地利用类型的训练样区，对2000年基于MODIS的土地利用类型分类图进行精度检验。训练样区选取时量均匀选取，覆盖全区；不仅选取土地利用类型单一的训练区，并要选取多种土地利用类型混杂的训练区，当其中某土地利用类型的面积占此混杂训练区面积的1/2时即分为此类型；此外在TM上选取的训练样区的长宽需要大于9个象元，这样才可以使MODIS一个象元的真实地表面积与TM上的训练样区的真实地表面积相匹配。在ENVI软件中作精度检验，总体分类精度为90.01％，Kappa系数为0.86。

S3、在本发明的优选实施例中，针对不同种植模式，提取NDVI时间序列主成分，并进行分析。

采用主成分变换可将NDVI年序列中初级生产力的有用信号压缩到少量的前几个主分量中，在信息损失最小的前提下，用较少的分量替代原来的高维数据，生成的主成分影像既包含了时间序列信息，也包含了空间信息，每个主成分表达的主题需要利用专业领域知识解释。

本发明假设这些主成分从不同的方面反映了同一作物种植模式下初级生产力的变化特征，它包含了土地自然质量、化肥投入、灌水措施等管理措施的影响。经过主成分分析，这些影响因素所造成的区域生产力的差异保持相对独立。在同一种植模式下，通常情况下土地自然质量是造成初级生产力空间差异最大的主导因素，第二主导因素是化肥投入，因此，NDVI年序列前2个主成分应分别对应土地自然质量和化肥投入在时间和空间上对年初级生产力的响应特征和空间差异信息。

分别对三种种植模式(果园、一年一熟农用地和一年两熟农用地)的NDVI时间序列做主成分分析。经过主成分分析，每种种植模式的前6个主成分在总方差中的累计贡献率如表1所示。根据方差的大小排序，每种种植模式的前6个主成分的方差最大，包含了信息量的大部分，累计达到86％以上，占有绝对优势。不同作物种植模式下，第一主成分反映了初级生产力的30.93％～46.65％，果园最大，一年两熟农用地最小。果园模式下第二主成分的贡献最小，为11.78％，一年一熟农用地最大占到了25.14％。

表1不同种植模式下主成分变换前6个主成分累积贡献率

各种植模式的前2个主成分特征向量曲线如图4所示。纵观各类种植模式的第一主成分特征向量图，3种种植模式下的第一主成分特征向量在1年中大部分时间为正值，而在作物关键生长期到收获期期间表现为负的特征值，通常这些时期是作物经济产量形成的关键时期。果园(图4中的a)表现为第26～29个8d(9月中旬-10月中旬)，一年两熟表现为第20～25个8d(7月底-九月初)，第34～39个8d(10底-12月初)，一年一熟表现为第23～30个8d(7月初-9月底)。这些时期的NDVI以负权重表征在一年的NDVI合成上，使得关键生长期的累积NDVI得到对比突出，增强了作物的生长成熟期与其他时间的生物量对比情况，可见3种种植模式下的第一主成分增强了作物的生长关键期对最终产量(尤其经济产量)的影响。

3种种植模式下的第二主成分特征向量曲线与第一主成分不同(图4)，但基本表现出在第一主成分所刻画的作物关键生长期前期有一个高峰上升。通常假设施肥发生在作物关键期的初期，因此可以认为第二主成分增强了施肥前后所导致的生产力的变化。

S4、化肥投入数据作为一种统计数据，与NDVI时间序列主成分数据的尺度不同，为进一步分析不同种植模式下化肥投入数据与其中两个主成分的关系，需要得到不同种植模式的化肥投入数据，并实现化肥投入数据和遥感数据在尺度上一致。本发明实施例中提出利用两种方法可实现其尺度匹配，并在此基础上分析不同种植模式下化肥投入与NDVI前两个主成分的关系，选取能够反映化肥投入所造成的生产力的差异的主成分，利用其关系进行化肥的空间化。

第一种方法，将两种数据都匹配到像元尺度上。通过在不同种植模式主成分影像上选取有代表性的纯净像元，开展不同种植模式下像元内典型农户地块上施肥品种、数量等行为的调查与监测，获取不同种植模式像元尺度上的化肥投入数据，与影像数据的尺度相匹配，进而分析像元尺度上不同种植模式NDVI前两个主成分与化肥投入的关系。如果存在显著地相关关系，则利用这些样区建立化肥投入的回归模型，实现化肥投入的空间分配。

第二种方法，将两种数据都匹配到统计数据的行政单元尺度上。化肥投入数据可来源于统计年鉴中的以乡镇为最小单元的化肥投入量，如果是这种来源的统计数据，先计算乡镇尺度上不同种植模式的化肥投入数据，然后将不同种植模式NDVI前两个主成分遥感数据按照乡镇级进行整合，以便与化肥投入信息相匹配，分析乡镇尺度上不同种植模式下化肥投入与前两个主成分的关系，并利用两者关系，按乡镇界线汇总的不同种植模式化肥投入统计数据分配到辖区上基于像元分辨率的农业用地上。

优选地，鉴于第二种方法的数据可获得性，本发明实施例中利用第二种方法实现两种数据的尺度匹配，进而分析乡镇尺度上不同种植模式下化肥投入与前两个主成分的关系，最后进行化肥的空间化分配，步骤如下：

第一，整理北京北部五个山区县的统计年鉴(2000年)，得到以乡镇为单元的化肥投入数据。

第二，计算基于乡镇单元的不同种植模式化肥投入数据。所用资料是2000年北京市北部山区按乡镇统计的化肥折纯量，并没有区分不同种植模式下的投入在乡镇之间的比例，也就是应该区别不同种植模式化肥投入的差异。可通过几种方法得到：第一种方法，对各乡镇进行实地走访调研，确定不同种植模式化肥投入的比例，然后根据统计年鉴上以乡镇为单元的化肥投入统计数据计算各种植模式下化肥的投入量。第二种方法，如果没有实地调研数据，则可以先确定每种种植模式下种植的农作物类型，再按照统计资料中的单位面积农作物的化肥投入量或者化肥投入成本来计算并汇总不同种植模式下的化肥投入量。

优选地，由于缺乏实地调研数据，本发明实施例中采用第二种方法计算基于乡镇单元的不同种植模式化肥投入量。按照2000年《全国农产品成本收益资料汇编》有关数据：2000年小麦的化肥投入为1137.45元/hm²，玉米的化肥投入为763.35元/hm²，果园的化肥投入为1642.05元/hm²。采用同一比例模式来说明研究区不同种植模式化肥投入的差异，将不同种植模式的单位面积化肥投入比例简化为：一年一熟农用地单位面积化肥投入：果园单位面积化肥投入：一年两熟单位面积农用地化肥投入＝玉米/果园/(小麦+玉米)＝1/2.15/2.49＝k。按照乡镇界线对各类种植模式的面积进行统计，每个乡镇将得到一年一熟农用地面积：果园面积：一年两熟农用地面积＝m，k*m可得到各个乡镇的各类种植模式按面积分配的比例，按照这个比例，可以得到各个乡镇每类种植模式的化肥投入总值，这个总值也一定是按照跟面积成比例的，同时各个不同种植模式的化肥投入总值之间也等于比值k。

第三，将基于像元的遥感数据分析结果转换成基于乡镇单元的汇总数据。由于NDVI时间序列主成分分析结果是像元尺度上的数值，需要按照乡镇级进行整合，以便与化肥投入信息相匹配。根据上述主成分分析结果，三种种植模式的第一主成分的特征向量基本上均为正值，对第一主成分的值为正贡献。由于乡镇单元的化肥投入数据是整个乡镇的化肥总投入，直接使用研究区域的乡镇界限图对三种种植模式的每个像元的第一主成分按照乡镇加和汇总，得到乡镇尺度上不同种植模式的第一主成分汇总值。

对于第二主成分，一年一熟农用地和一年两熟农用地这两种模式的第二主成分的特征向量基本上均为负值，对第二主成分的值为负的贡献，第二主成分绝对值越大，应该是代表的生产力越高，因此将这两种种植模式的第二主成分的影像取绝对值。然后使用研究区域的乡镇界限图，对三种种植模式的每个像元的第二主成分(绝对值)按照乡镇加和汇总，得到乡镇尺度上不同种植模式的第二主成分(绝对值)汇总值。

第四，按照不同的种植模式可以对乡镇化肥投入与按乡镇汇总的第一和第二主成分进行相关分析(表2)，绘制散点图(图5和图6)，可看出各类种植模式下乡镇化肥投入与按乡镇汇总的第二主成分均有较高的相关性，并且为0.01水平下显著正相关，而与第一主成分的相关性很低。因此，各种种植模式下的第二主成分可以代表由化肥投入所造成的生产力的差异。

一年两熟种植模式下的化肥投入与由化肥投入所造成的生产力线性相关的关系最显著，这是因为在本实施例的研究区一年两熟种植模式涉及的作物类型和品种组合相对稳定，对化肥的需求也成稳定的相应关系；而一年一熟种植模式下的作物类型除了有春玉米以外，还有其他的类型，二者的化肥生产力响应显然不同，因此图6明显地出现2条线性关系，本发明中采用二者的平均线性关系作为一年一熟种植模式对化肥用量的响应；对于果树，本实施例研究区的果树类型存在差异，由化肥投入所造成的生产力响应也有一定的分散性。

表2不同种植模式NDVI的前两个主成分与化肥投入的相关系数

表2中，**代表显著水平0.01，*代表显著水平0.05

第五，利用乡镇尺度上各种植模式第二主成分与化肥投入之间的关系，进行化肥投入的空间分配，具体实施如下：通常化肥投入与作物生产力间关系是抛物线特征，而上述分析中表明，在本发明实施例研究区内各种植模式的第二主成分可以反映由于化肥投入所造成的生产力差异，且表现出线性关系。因此，在该实施例中假设化肥投入与由化肥投入所造成的生产力之间是一种线性相关的关系，按照第二主成分的大小对每个区县各类种植面积进行化肥分配。具体来说，各类种植模式都有一个第二主成分的影像，并且有每个乡镇不同种植模式的第二主成分(绝对值)汇总值，对于每一种种植模式，第二主成分影像的每个像元值与所属乡镇的第二主成分汇总值的比值就是用于分配的比例，乘以其所属乡镇该种植模式下的化肥投入，就得到各类种植模式下每个像元的化肥投入。这样既保证了每个乡镇的化肥投入分配汇总后与统计数据相等，又使得不同种植模式下的化肥投入之间保持应有的投入比例，还能反映相同投入方式下不同的像元之间的差异。最后将三种种植模式的化肥投入空间分配图合并得到三种种植模式的化肥投入空间分配图(图7)。

以上实施方式仅用于说明本发明，而并非对本发明的限制，有关技术领域的普通技术人员，在不脱离本发明的精神和范围的情况下，还可以做出各种变化和变型，因此所有等同的技术方案也属于本发明的范畴，本发明的专利保护范围应由权利要求限定。

Claims (6)

1.一种基于统计数据和遥感影像数据的区域化肥投入空间化方法，其特征在于，包括以下步骤：

S1、获取归一化植被指数NDVI时间序列数据，并对该数据进行平滑处理；

S2、利用所述数据提取分类特征，并根据该分类特征采用决策树分类方法对作物种植模式进行分类；

S3、在不同类型的种植模式下进行NDVI时间序列主成分分析，并提取多个主成分；

S4、对统计获得的化肥投入数据和其中前2个主成分进行尺度匹配，并利用尺度分配的结果分析不同种植模式下化肥投入数据与这两个主成分的关系，利用前2个主成分中能够反映化肥投入所造成的农用地生产力的差异的主成分进行化肥投入的空间化，所述化肥投入数据包括投入的化肥种类数据和数量数据，主成分按照与化肥投入存在的关系的显著性从大到小排列。

2.如权利要求1所述的基于统计数据和遥感影像数据的区域化肥投入空间化方法，其特征在于，步骤S3中的所述主成分分析针对各种植模式的NDVI时间序列分别进行，分析时首先针对每种种植模式提取6个主成分，这6个主成分从不同方面反映了农用地生产力的状况，然后对各种植模式下前2个主成分的特征向量进行分析来认识前2个主成分分别代表的物理意义。

3.如权利要求2所述的基于统计数据和遥感影像数据的区域化肥投入空间化方法，其特征在于，在步骤S4中使用各种植模式下化肥投入数据与该种植模式6个主成分中的前2个主成分进行相关分析，找到各种植模式下与化肥投入存在一定关系的主成分。

4.如权利要求2所述的基于统计数据和遥感影像数据的区域化肥投入空间化方法，其特征在于，步骤S4中进行尺度匹配的方法为：将化肥投入数据和前2个主成分都匹配到像元尺度或乡镇尺度上。

5.如权利要求1所述的基于统计数据和遥感影像数据的区域化肥投入空间化方法，其特征在于，步骤S1中，利用权重滑动平均滤波法Savizky-Golay对MODIS NDVI时间序列数据进行平滑滤波。

6.如权利要求1～5任一项所述的基于统计数据和遥感影像数据的区域化肥投入空间化方法，其特征在于，利用中分辨率成像光谱仪获得NDVI时间序列数据。

Images