作物多元数轮作周期轮作补贴的遥感估算方法

技术领域  本发明涉及农业及相关领域，用于估算在大地理区域中农作物轮作补贴。

背景技术  通常将在同一地块上，一种作物种植一定的时间后，再种植另一种作物一定的时间，并且一直交替下去的过程称为作物轮作，这两种时间分别称为从一种作物到另一种作物的轮作时间，也称为轮作周期。轮作是作物主要的耕作栽培模式之一，具有明显的生态效益和经济效益，不仅可以大大改善土壤理化性状和生态环境，而且还可以有效地破坏病虫、杂草的生存环境，显著降低农药的使用量，减少农药污染，提高农作物的产量。我国不仅在南方地区水旱以及旱旱轮作广泛应用，而且在北方地区轮作也广泛应用，如玉米与大豆的轮作就具有明显的生态和经济效益，国家每年投资近亿元用于补贴实施大豆与玉米轮作的农户。在世界范围内，轮作同样也正在成为提高农业生产的生态和经济效益的关键技术，特别是研究对实施轮作的农户如何发放轮作补贴，已成为现代农业生产及管理中必须要解决的关键问题。根据WTO的《农业协议》，对休耕、休渔、土地轮作、培肥地力的农业生产单位和农户，可以在农业生产结构调整补贴项目下发放生产补贴。因此，轮作补贴应属于“绿箱”政策下的对结构调整提供的投资补贴，它是一种直接支付制度，符合WTO《农业协议》的精神。我国已经加入WTO，各项补贴都应在WTO框架下实施，应参照各发达国家以及其他发展中国家的经验，建立起一套国内支农体系，通过加大“绿箱”政策的支持力度，提高我国农业的整体素质和国际竞争力，轮作补贴完全符合这一指导思想。事实上，如何科学的发放轮作补贴已经成为影响农业生产的生态和经济效益的重要因素。因此，在大地理范围内研究如何科学的发放轮作补贴具有十分重要的意义。

通常，将农作物不进行轮作或轮作周期太长所产生的病虫害称为农作物轮作病虫害，由于农作物不进行轮作或轮作周期太长都会使轮作病虫害加剧，造成使用农药的数量增加，并且使农作物产量下降，而科学的发放轮作补贴是鼓励农户积极进行轮作的有效措施，但是如何建立一种科学发放轮作补贴的模型，用于在大地理区域范围内为农户发放轮作补贴一直是农业和其它相关领域中面临的一个难题，国内外的学者对此已进行了大量的探索，但没有找到解决此问题的有效方法。

导致对建立一种科学发放轮作补贴的模型研究进展缓慢的主要原因之一是：并不知道如何准确的表示和描述农作物轮作周期，而农作物轮作自身的特点已说明，传统的用单一数值描述农作物之间的轮作周期方法是不科学的，首先农作物之间的轮作存在明显的对称性，其次这种农作物之间的轮作是相互的，第三在农作物轮作周期之间存在明显的不对称性，即它们数值的大小并不一定相等，而所有这三个关键特征根本就无法用单一数值来描述，采用单一数值很难对轮作补贴进行客观估算，从而制约了用科学的方法对建立一种科学发放轮作补贴的模型进行深入研究，所以必须寻找新的方法。

本发明的目的是采用一种新的方法运用对遥感影像的解译结果，并通过建立在统计学以及随机过程基础上的数学模型，来分析估算大地理区域范围内的农作物的实际多元数轮作周期，然后再将农作物的实际多元数轮作周期与最小必须的多元数轮作周期进行比较，来估算该农作物轮作应该发放的补贴数量，该方法具有高效、简单、客观、易于应用等特点，可广泛用于在区域或全球范围内对农作物的轮作补贴进行估算等领域。

发明内容  本发明将从有限多个在年份上两两连续的卫星遥感影像上解译出来的在大地理范围内分布的轮作区域中的所有同一地块上不同年份的农作物轮作的基本数据，根据地理范围的大小不同，分别按市、县、乡或村为单位来计算农作物的轮作周期，再根据从一种作物到另一种作物轮作顺序不同，他们的轮作周期也不尽相同的结论，运用在轮作为平稳随机过程基础上建立的公式来估算不同顺序农作物的实际轮作周期，最终估算该农作物的实际多元数轮作周期。然后，再将农作物的实际多元数轮作周期与最小必须的多元数轮作周期进行比较，来估算该农作物轮作应该发放的补贴数量。因此，将有限多个在年份上两两连续的卫星影像的解译结果用于高精度估算大地理区域范围内农作物实际多元数轮作周期以及将农作物的实际多元数轮作周期与最小必须的多元数轮作周期进行比较，来估算该农作物轮作应该发放的补贴数量的方法成为本发明的重要特征。

本发明作物多元数轮作周期轮作补贴的遥感估算方法的技术方案是：

首先给出农作物多元数轮作周期的基本概念，然后对获取的覆盖要研究轮作区域的农作物的遥感影像进行解译，获得每个乡或镇(市或县)的农作物轮作数据；再运用农作物轮作周期估算公式对该轮作数据进行分析，获得农作物的实际轮作周期和实际多元数轮作周期；再将农作物的实际多元数轮作周期与最小必须的多元数轮作周期进行比较，根据相关公式即可计算出该农作物轮作应该发放的补贴数量。

本发明的农作物遥感影像解译主要包括四个步骤，首先对原始影像进行几何校正，线性拉伸增强等预处理；其次根据野外实地调查确定轮作农作物以及其他作物的解译标志；然后再通过人工目视解译，从影像中提取轮作农作物的信息并生成含轮作农作物图斑的解译图件，之后，对该图件实地验证检查，对有问题的图斑进行修正，使图件的解译准确可靠；最后生成包含轮作农作物的矢量化解译图件，并通过地理信息系统对解译图件进行空间分析处理，获得轮作农作物的轮作数据。

研究农作物轮作的多元数轮作周期，通常必须对轮作区域进行多年观察，并且对所获得的数据进行统计分析，才能获得令人信服的结果，因此，本发明利用从有限多个在年份上两两连续的遥感卫星影像中获得的数据来估算农作物的多元数轮作周期。有限多个在年份上两两连续的遥感卫星影像是指选取的多年遥感卫星影像在年份上每两幅影像是连续的，但连续的每两幅影像之间在年份上可以是间断的。

本发明设计的农作物多元数轮作周期估算公式，除了利用有限多个在年份上两两连续的遥感卫星影像监测外，采用了对轮作区域内所有乡镇的农作物轮作水平同时观察的方法，来估算一个理想化的典型乡镇的农作物轮作水平(对市、县可同样的讨论)，相当于将乡一级的随机试验样本数量扩大到多个，以更加准确地推算出该区域的农作物的多元数轮作周期，其依据为：各个乡镇以及各个农户之间在轮作及耕作栽培管理上具有明显的相互独立性，基本符合统计学上对样本独立性的要求；农作物轮作的周期性反映在用随机过程描述的作物轮作过程时，该过程具有平稳性，在统计学上就相当于每年农作物轮作的数量应大致相同，可根据每年农作物轮作的数量估算出完成整个农作物轮作所需要的时间或周期；有限多个在年份上两两连续的多个乡镇的独立样本在数量上对于研究具有平稳性的乡一级的农作物轮作过程的统计特征或周期来说，应该是足够大。

本发明设计的作物多元数轮作周期轮作补贴的遥感估算方法适用于所有农作物是本发明的重要特征。

现在，首先以棉花与水稻轮作为例来说明农作物多元数轮作周期估算公式，但本公式适用于所有农作物多元数轮作周期的估算是本发明的重要特征。

在有限多个在年份上两两连续的遥感卫星影像监测的相邻两年中，用“棉/稻”表示上一年种棉花的地块在下一年种了水稻的面积。对于给定的地块，如果“棉_i/稻_i+1”表示在此地块上，第i年种棉花而第i+1年种水稻的面积，则下一年种水稻占上一年种棉花面积的比例即“水稻的面积_i+1/棉花的面积_i”被定义为此地块i年到i+1年的棉稻轮作因子或简称为此地块i年的棉稻轮作因子，并且用CRRF_i(Cotton-Rice Rotation Factor)来表示，

$>>>CRRF>i>>=>>>NRA>i>>>CCA>i>>>;>$>

其中CCA_i为对于给定的地块在第i年种棉花的面积；NRA_i为CCA_i中，第i+1年种水稻的面积。

假设仅有棉花与水稻在某一给定的地块上进行轮作，对此地块连续监测N年，这一轮作过程可以用具有明显周期性的平稳随机过程来描述，而且在此地块上，棉稻轮作因子为CRRF_i(i＝1，......，N-1)，则棉花与水稻的轮作周期CRRP(Cotton-Rice Rotation Period)为

$>>CRRP>=>>lim>>N>\to >\infty >\; >>1>>N>->1>>over>>\Sigma >>i>=>1>>>N>->1>>>>1>>CRRF>i>>>;>$>

由于实际上观察的年数N是有限的，因此当N足够大时，棉稻轮作周期CRRP的近似值：

$>ver>>CRRP>‾>>=>>1>>N>->1>>over>>\Sigma >>i>=>1>>>N>->1>>>>1>>CRRF>i>>>.>$>

由于棉花与水稻轮作过程可以用具有明显周期性的平稳随机过程来描述，因此每年棉花轮作的数量在统计学上应该大致相同，而且棉花与水稻轮作周期是存在并可计算的。

因为仅有棉花与水稻在给定的地块上进行轮作，所以根据CRRF_i的定义，在给定“上一年棉花面积”CCA_i的情况下，如果每年有“棉/稻”面积NRA_i改种水稻，则年可将“上一年棉花面积”改种完，故CRRP_i＝1/CRRF_i成立。

假设所要求的周期CRRP是随机变量ξ的数学期望E(ξ)，那么近似确定CRRP的方法是对ξ进行N次重复抽样，产生相互独立的ξ值的序列ξ₁，......，ξ_N(ξ₁＝1/CRRF₁，........，ξ_N＝1/CRRF_N)，则其算术平均值

$>ver>>>\xi >N>>‾>>=>>1>N>over>>\Sigma >>i>=>1>>N>>>\xi >i>>;>$>

根据强大数定理，

$>>P>>(>>lim>>N>\to >\infty >\; ver>>>\xi >N>>‾>>=>CRRP>)>>=>1>;>$>

因此，当N充分大时，

$>ver>>>\xi >N>>‾>>\approx >E>>(>\xi >)>>=>CRRP>;>$>

成立的概率等于1，因此可用作为所求量CRRP的估算值。

根据上述定理，可以利用对轮作区域的M个乡镇连续N年卫星监测数据，计算该轮作区域的棉稻轮作周期，

$>ver>>CRRP>‾>>=>>1>>N>->1>>over>>\Sigma >>i>=>1>>>N>->1>>>>1>>CRRF>i>>>=>>1>>N>->1>>over>>\Sigma >>i>=>1>>>N>->1>>>>(over>>\Sigma >>j>=>1>>M>>>>>CCA>ij>>2>>>>NRA>ij>over>>\Sigma >>j>=>1>>M>>>CCA>ij>>>>)>>;>$>

其中M为乡镇总数；NRA_ij为第j个乡镇在第i年种棉花的面积中，在第i+1年种水稻的面积；CCA_ij为第j个乡镇在第i年种棉花的面积。

同理，对于水稻与棉花轮作的情况，可推出计算轮作区域的稻棉轮作周期公式

$>ver>>RCRP>‾>>=>>1>>N>->1>>over>>\Sigma >>i>=>2>>N>>>1>>RCRF>i>>>=>>1>>N>->1>>over>>\Sigma >>i>=>2>>N>>>(over>>\Sigma >>j>=>1>>M>>>>>CCA>ij>>2>>>>NCA>ij>over>>\Sigma >>j>=>1>>M>>>CCA>ij>>>>)>>;>$>

其中M为乡镇总数；RCRF_i为i年种棉花的面积中有多少在i-1年种水稻；NCA_ij为第j个乡镇在第i-1年种水稻的面积中，在第i年种棉花的面积；CCA_ij为第j个乡镇在第i年种棉花的面积。

推广到一般情况，称用上述方法求得的农作物轮作周期为农作物的统计学轮作周期，简称为农作物的统计轮作周期。因此，在大地理范围内，作物x与y的轮作可以用它们之间的统计轮作周期来描述，根据上述分析，在一般情况下，作物x与作物y的轮作周期与作物y与作物x的轮作周期是不同的，这种在农作物统计轮作周期之间的不对称性与通常农作物轮作周期之间的不对称性的物理含义是相同的。因此，农作物x和y之间的轮作周期可以用下述m阶多元数来描述：

(x₁，...，x_m，y₁，...，y_m，f₁(x₁，...，x_m，y₁，...，y_m)，...，f_p(x₁，...，x_m，y₁，...，y_m))

其中：x_i(i＝1，...，m)为用m种不同的方法所求得的农作物x与农作物y的轮作周期；y_i(i＝1，...，m)为用m种不同的方法所求得的农作物y与农作物x的轮作周期；f_j(x₁，...，x_m，y₁，...，y_m)(j＝1，...，p)为用p种不同的方法所求的反映在轮作周期x_i和y_i(i＝1，...，m)之间的差异函数，因此可以将多元数轮作周期中的元分为三类，第一类与农作物x与y的轮作有关，第二类与农作物y与x的轮作有关，第三类与第一类和第二类元有关。

现在假定仅用上述稻棉轮作周期与棉稻轮作周期公式求得农作物x与y之间的轮作周期分别为x₁和y₁，描述x₁和y₁之间的周期差异函数为f₁(x₁，y₁)＝y₁-x₁，则描述农作物x与y轮作的一阶多元数轮作周期为(x₁，y₁，y₁-x₁)。通过与传统的用单一数值描述的农作物x与y之间的轮作周期相比较，农作物多元数轮作周期提供了更加全面的信息，是一种更为科学的描述方法，并且可以得到许多重要的性质，以一阶多元数轮作周期(x₁，y₁，y₁-x₁)为例：如果y₁-x₁＞0，则说明轮作使农作物x的收益比作物y的收益更大，这是因为如果对于农作物x与y的轮作来说，轮作是使作物x具有更高的生态和经济效益，那么农作物x与y轮作周期的长短通常与生态和经济效益有密切的关系，农作物x与y的轮作周期短些以及农作物y与x的轮作周期长些，通常对农作物x更为有利，即对于两个多元数轮作周期的比较，在x₁相同的情况下，更要考虑差异函数值y₁-x₁的大小，即对y₁值进行比较，y₁值越大则轮作的效果越好。

综上所述，以棉花与水稻轮作为例来说明农作物的一阶多元数轮作周期估算公式如下，但本公式适用于所有农作物一阶多元数轮作周期的估算是本发明的重要特征，对于农作物的m阶多元数轮作周期可作类似分析。农作物x与y轮作的一阶多元数轮作周期：

(x₁，y₁，y₁-x₁)

其中： $>>>x>1>>=ver>>CRRP>‾>>;>>y>1>>=ver>>RCRP>‾>>;>>y>1>>->>x>1>>=ver>>CRRP>‾>>-ver>>RCRP>‾>>.>$>

在以下分析中，如果不加以特别说明，所有多元数轮作周期均是指一阶多元数轮作周期。

本发明设计的作物多元数轮作周期轮作补贴的遥感估算方法如下：

假定农作物i与农作物j进行轮作，用于提高农作物生产的生态和经济效益，则可以用下述公式来估算在农作物i与农作物j进行轮作时，相关的农户应该得到的轮作补贴δ_ijs：

δ_ijs＝f(x_a，x_b，x_c)；

其中：f(x_a，x_b，x_c)为求δ_ijs值的函数；变量x_a为农作物i的实际多元数轮作周期：x_a＝(x₁₁，y₁₁，y₁₁-x₁₁)；变量x_b为农作物i的最小必须多元数轮作周期：x_b＝(x₁₂，y₁₂，y₁₂-x₁₂)；变量x_c为与农作物i、农作物j、农作物轮作病虫害、耕作栽培方式、轮作区域的地理条件以及其它生态和经济因素有关的量。

f(x_a，x_b，x_c)可以有多种形式，例如可以用多元回归方法求δ_ijs值，此时有：

δ_ijs＝f(x_a，x_b，x_c)＝a₁₁x₁₁+b₁₁y₁₁+c₁₁(y₁₁-x₁₁)+a₁₂x₁₂+b₁₂y₁₂+c₁₂(y₁₂-x₁₂)；

其中：回归系数a₁₁，...，c₁₂与农作物i、农作物j、农作物轮作病虫害、耕作栽培方式、轮作区域的地理条件以及其它生态和经济因素有关。可以通过在农作物i的轮作区域，根据实际多元数轮作周期的不同，设置若干个样方，再在每个样方中设定农户的实际补贴δ_ijs，并根据相对应的x₁₁，...，(y₁₂-x₁₂)的值，用回归的方法确定系数a₁₁，...，c₁₂的值。

通过上述多元回归方程对样方中的实际数据的处理，并且令x_a＝(CRTR_ij，RCTP_ij，RCTP_ij-CRTP_ij)，x_b＝(CRTN_ij，RCTN_ij，RCTN_ij-CRTN_ij)，回归系数a₁₁，...，c₁₂用系数S_ij、C_ij和D_ij来取代，可以得到求δ_ijs值的更简单的公式如下：

δ_ijs＝S_ij-{C_ij×(CRTP_ij-CRTN_ij)+D_ij×(RCTP_ij-RCTN_ij)}；

其中：(1)S_ij为农作物i和农作物j轮作，应该发放的最大补贴数，即规定δ_ijs的最大值为S_ij，0≤δ_ijs≤S_ij，实际上δ_ijs也可为负值用于对轮作水平太低进行惩罚。

(2)C_ij和D_ij为与农作物i、农作物j、农作物轮作病虫害、耕作栽培方式、轮作区域的地理条件以及其它因素有关的系数，用于对δ_ijs进行规范化，使δ_ijs更加合理。

(3)最小多元数轮作周期是指能够最大限度的发挥农作物i和农作物j轮作的经济和生态效益同时又能防止农作物轮作产量损失所必须的最小多元数轮作周期，可以通过科学试验、测定，或根据长期积累的经验来确定。

(4)1≤CRTP_ij＜∞，1≤CRTN_ij＜∞，1≤RCTP_ij＜∞，1≤RCTN_ij＜∞。

根据估算的农作物i与农作物j实施轮作的实际多元数轮作周期，应用上式即可估算出应该发放的轮作补贴。