CN110008621A - 基于双重流依赖集合平方根滤波同化算法和基于该算法的作物模型遥感同化估产方法 - Google Patents
基于双重流依赖集合平方根滤波同化算法和基于该算法的作物模型遥感同化估产方法 Download PDFInfo
- Publication number
- CN110008621A CN110008621A CN201910300249.7A CN201910300249A CN110008621A CN 110008621 A CN110008621 A CN 110008621A CN 201910300249 A CN201910300249 A CN 201910300249A CN 110008621 A CN110008621 A CN 110008621A
- Authority
- CN
- China
- Prior art keywords
- algorithm
- assimilation
- ensrf
- yield
- crop
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000004422 calculation algorithm Methods 0.000 title claims abstract description 128
- 238000001914 filtration Methods 0.000 title claims abstract description 40
- 238000000034 method Methods 0.000 title claims abstract description 30
- 230000009977 dual effect Effects 0.000 title claims abstract description 23
- 230000012010 growth Effects 0.000 claims abstract description 43
- 238000004458 analytical method Methods 0.000 claims abstract description 28
- 238000004088 simulation Methods 0.000 claims abstract description 20
- 230000036962 time dependent Effects 0.000 claims abstract description 9
- 239000002689 soil Substances 0.000 claims abstract description 8
- 238000005457 optimization Methods 0.000 claims abstract description 4
- 239000011159 matrix material Substances 0.000 claims description 21
- 230000008859 change Effects 0.000 claims description 11
- 238000012545 processing Methods 0.000 claims description 7
- 230000002265 prevention Effects 0.000 claims description 4
- 230000004807 localization Effects 0.000 claims description 3
- 230000008635 plant growth Effects 0.000 claims description 3
- 241000209140 Triticum Species 0.000 description 65
- 235000021307 Triticum Nutrition 0.000 description 63
- 238000011160 research Methods 0.000 description 14
- 239000000523 sample Substances 0.000 description 12
- 238000012360 testing method Methods 0.000 description 12
- 230000000694 effects Effects 0.000 description 7
- 238000005259 measurement Methods 0.000 description 7
- 230000008569 process Effects 0.000 description 6
- 238000004519 manufacturing process Methods 0.000 description 4
- 235000013339 cereals Nutrition 0.000 description 3
- 230000007812 deficiency Effects 0.000 description 3
- 238000011161 development Methods 0.000 description 3
- 230000018109 developmental process Effects 0.000 description 3
- 238000009826 distribution Methods 0.000 description 3
- 239000000203 mixture Substances 0.000 description 3
- 238000012544 monitoring process Methods 0.000 description 3
- 241001269238 Data Species 0.000 description 2
- 240000008042 Zea mays Species 0.000 description 2
- 235000002017 Zea mays subsp mays Nutrition 0.000 description 2
- 230000008901 benefit Effects 0.000 description 2
- 238000004364 calculation method Methods 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000011156 evaluation Methods 0.000 description 2
- 238000003306 harvesting Methods 0.000 description 2
- 230000007246 mechanism Effects 0.000 description 2
- 238000010899 nucleation Methods 0.000 description 2
- 238000001556 precipitation Methods 0.000 description 2
- 230000005855 radiation Effects 0.000 description 2
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 description 1
- 229920000742 Cotton Polymers 0.000 description 1
- 241000196324 Embryophyta Species 0.000 description 1
- 241000607479 Yersinia pestis Species 0.000 description 1
- 235000005824 Zea mays ssp. parviglumis Nutrition 0.000 description 1
- 235000016383 Zea mays subsp huehuetenangensis Nutrition 0.000 description 1
- 230000032683 aging Effects 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 238000009395 breeding Methods 0.000 description 1
- 230000001488 breeding effect Effects 0.000 description 1
- 235000005822 corn Nutrition 0.000 description 1
- 201000010099 disease Diseases 0.000 description 1
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 1
- 235000013399 edible fruits Nutrition 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 235000009973 maize Nutrition 0.000 description 1
- 238000013507 mapping Methods 0.000 description 1
- 230000035800 maturation Effects 0.000 description 1
- 239000008267 milk Substances 0.000 description 1
- 210000004080 milk Anatomy 0.000 description 1
- 235000013336 milk Nutrition 0.000 description 1
- 239000002245 particle Substances 0.000 description 1
- 238000012827 research and development Methods 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 238000010206 sensitivity analysis Methods 0.000 description 1
- 239000007787 solid Substances 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
- 230000017260 vegetative to reproductive phase transition of meristem Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q10/00—Administration; Management
- G06Q10/06—Resources, workflows, human or project management; Enterprise or organisation planning; Enterprise or organisation modelling
- G06Q10/063—Operations research, analysis or management
- G06Q10/0639—Performance analysis of employees; Performance analysis of enterprise or organisation operations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q50/00—Information and communication technology [ICT] specially adapted for implementation of business processes of specific business sectors, e.g. utilities or tourism
- G06Q50/02—Agriculture; Fishing; Forestry; Mining
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Business, Economics & Management (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Human Resources & Organizations (AREA)
- Economics (AREA)
- Data Mining & Analysis (AREA)
- Strategic Management (AREA)
- Development Economics (AREA)
- Entrepreneurship & Innovation (AREA)
- Mathematical Analysis (AREA)
- Computational Mathematics (AREA)
- Pure & Applied Mathematics (AREA)
- Operations Research (AREA)
- General Engineering & Computer Science (AREA)
- Mathematical Optimization (AREA)
- Mathematical Physics (AREA)
- Marketing (AREA)
- Educational Administration (AREA)
- Life Sciences & Earth Sciences (AREA)
- Tourism & Hospitality (AREA)
- General Business, Economics & Management (AREA)
- Animal Husbandry (AREA)
- Primary Health Care (AREA)
- General Health & Medical Sciences (AREA)
- Game Theory and Decision Science (AREA)
- Health & Medical Sciences (AREA)
- Quality & Reliability (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Mining & Mineral Resources (AREA)
- Marine Sciences & Fisheries (AREA)
- Evolutionary Biology (AREA)
- Agronomy & Crop Science (AREA)
- Geometry (AREA)
- Probability & Statistics with Applications (AREA)
- Evolutionary Computation (AREA)
- Algebra (AREA)
- Computer Hardware Design (AREA)
Abstract
提出一种双重流依赖集合平方根滤波同化算法，包括：S1，将EnSRF算法中的控制变量投影于观测空间，形成4DEnSRF算法；S2，在4DEnSRF算法中引入流依赖时间窗口形成DFD‑EnSRF算法。本发明还提出基于该算法的作物模型遥感同化估产方法，包括：T1，将气象数据、作物参数以及土壤参数输入到WOFOST作物生长模型中，获得相应时刻不同模型参数运行下的模拟LAI集合作为背景场数据；T2，利用DFD‑EnSRF算法同化模型模拟LAI集合和外部观测LAI集合，得到分析LAI集合；T3，针对分析LAI集合，使用SCE‑UA优化算法得到最终同化产量。本发明为提高作物同化估产精度提供新的思路。
Description
技术领域
本发明涉及遥感数据同化技术领域，更具体地，涉及一种双重流依赖集合平方根滤波同化算法，以及基于该算法的同化估产方法。
背景技术
自20世纪60年代起，作物生长模型的研究随着农业科学、计算机技术的发展以及人类对作物生长机理认识的不断加深得到了快速的发展，目前已进入了应用研究阶段。但是，当作物模拟从单点研究发展到区域尺度时，空间尺度的增加会带来地表或近地表环境的非均匀性，使得作物生长模型中的参数获取以及区域化等问题的解决遇到了困难。卫星遥感具有快速、宏观、动态等优点，能够监测作物宏观状况，反映环境因子对作物影响的综合作用结果，但遥感数据的获取受卫星运行周期和天气等因素影响，在作物整个生育期内只能获得时间有限、离散的观测数据，这些数据不能有效支持对于作物生长发育和产量形成动态过程的研究和应用。因此，可将遥感信息引入作物生长模型进行数据同化，这不仅可以在单点尺度上揭示作物生长发育和产量形成的内在机理，同时解决了宏观数据的获取和参数区域化的困难。
数据同化技术作为能将作物生长过程模型与卫星遥感数据结合的最佳技术手段，为解决区域尺度高精度作物长势监测和产量预测提供了可靠、有效的方法和途径，目前已成为定量遥感和精准农业的研究重点和发展趋势，越来越受到研究人员的关注。但由于作物生长过程的复杂性和数学理论的局限性，数据同化系统只是对真实作物生长过程的近似描述，因此同化系统中必然存在着诸多误差和不确定性。尽管利用常规的各种观测数据可降低同化系统误差，提高同化系统输出的准确性和可预测性，但观测数据也存在着诸如观测不足、观测误差等问题。
为了解决以上问题，通常的做法是在同化系统中引入预报场作为分析场的初始场，即背景场，从而使数据同化转换为估计作物生长状态与背景场间的偏差，且估计偏差时需要考虑背景场以及观测场的误差及不确定性。目前在顺序同化、变分同化等大多数同化算法中，通过背景场误差协方差和观测场误差协方差分别表示背景场和观测场的误差和不确定性。一般情况下，不同观测点的观测误差是独立互不相关的，因此可由观测数据直接计算得出观测误差；而背景误差协方差在很大程度上决定了分析场的质量，且难以直接计算，因此在数据同化系统中就显得十分重要。在作物生长模型与遥感数据同化系统中，受气象、土壤以及田间管理等因素的影响，作物产量的空间相关往往是非均匀、各向异性的。而目前存在的大部分作物生长模型与遥感数据同化系统中，尤其是四维同化系统中，常将背景误差假设为均匀各向同性，并不能满足高精度作物长势监测和产量估算的同化应用需求。同时，在如集合卡尔曼滤波等顺序同化算法中，观测资料按时间顺序逐一进入同化系统，在一次分析中只能同化一个时次的观测资料，非分析时刻的观测资料不能被同化系统有效吸收，从而影响同化系统最终的同化精度。
发明内容
针对背景技术中，为了更好的提高同化算法效率和同化结果精度，可将顺序同化算法控制变量投影于观测空间，使顺序同化算法扩展为四维同化算法，扩展后的四维同化算法可以吸收原算法中被抛弃的非分析时刻观测资料，使新算法更优于原顺序同化算法。同时，对于作物生长模型与遥感数据同化系统来说，受卫星重访周期或遥感成像质量的制约，在作物关键生长期内仅能获得数量有限的高分辨率晴空遥感影像，而难以形成充足、稳定的观测数据。而在四维同化算法中需要同化一个时间窗口内所有的观测数据，设置合理的时间窗口大小对于同化系统输出精度具有十分重要的影响。在目前存在的大部分作物生长模型与遥感数据同化系统中常将时间窗口大小设为固定值或利用滤波算法对观测时序数据进行滤波处理，在观测数据不足或不稳定的情况下，同化系统难以高精度高效率的模拟作物参数或产量。
对此，本发明首先将集合平方根滤波(EnSRF)算法扩展为四维集合平方根滤波(4DEnSRF)算法，并在此基础上引入流依赖时间窗口，提出一种双重流依赖集合平方根滤波(DFD-EnSRF)算法，该算法中背景误差协方差和同化时间窗口均随时空发生变化，可提高集合平方根滤波同化系统的稳定性和输出精度；最后将所提新同化算法应用于作物生长模型与遥感数据同化系统中，以期为提高作物同化估产精度提供新的思路和技术支撑。
一方面，本发明提出一种双重流依赖集合平方根滤波同化算法，包括：
S1，将集合平方根滤波(EnSRF)算法中的控制变量投影于观测空间，形成四维集合平方根滤波(4DEnSRF)算法；
S2，在4DEnSRF算法中引入流依赖时间窗口形成双重流依赖集合平方根滤波(DFD-EnSRF)算法，令背景误差协方差和同化时间窗口均随时空发生变化。
另一方面，本发明提出一种基于双重流依赖集合平方根滤波的作物模型遥感同化估产方法，包括：
T1，将气象数据、作物参数以及土壤参数输入到WOFOST作物生长模型中，按照作物生长模型演进方程，获得相应时刻不同模型参数运行下的模拟LAI集合作为背景场数据；
T2，利用DFD-EnSRF算法同化模型模拟LAI集合和外部观测LAI集合，得到分析LAI集合；
T3，针对分析LAI集合，使用SCE-UA优化算法，通过迭代WOFOST作物生长模型寻找最优参数组合并得到最终同化产量。
本发明的有益效果为：本发明首先将集合平方根滤波算法扩展为四维集合平方根滤波算法，并在此基础上引入流依赖时间窗口，提出一种双重流依赖集合平方根滤波(DFD-EnSRF)算法，该算法中背景误差协方差和同化时间窗口均随时空发生变化，可提高集合平方根滤波同化系统的稳定性和输出精度；最后将所提新同化算法应用于作物生长模型与遥感数据同化系统中，以期为提高作物同化估产精度提供新的思路和技术支撑。
附图说明
为了更容易理解本发明，将通过参照附图中示出的具体实施方式更详细地描述本发明。这些附图只描绘了本发明的典型实施方式，不应认为对本发明保护范围的限制。
图1为本发明的方法的一个实施过程的流程图。
图2显示了本发明的方法的试验研究区。
图3显示了基于DFD-EnSRF和EnSRF算法的同化效果对比。
图4显示了基于DFD-EnSRF算法的衡水地区冬小麦遥感同化估产空间结果图。
图5显示了基于DFD-EnSRF算法的衡水地区冬小麦遥感同化估产精度验证结果。
具体实施方式
下面参照附图描述本发明的实施方式，其中相同的部件用相同的附图标记表示。在不冲突的情况下，下述的实施例及实施例中的技术特征可以相互组合。
双重流依赖集合平方根滤波同化算法
在顺序同化算法中，观测资料按时间顺序逐一进入同化系统，在一次分析中只能同化一个时次的观测资料。为了更好的提高同化算法效率和同化结果精度，可以将顺序同化算法控制变量投影于观测空间，使顺序同化算法扩展为四维同化算法。从理论角度看，扩展后的四维同化算法等价于原顺序同化算法；而从应用角度来看，扩展后的四维同化算法可以吸收原算法中被抛弃的非分析时刻观测资料，更优于原顺序同化算法。
本发明的双重流依赖集合平方根滤波(DFD-EnSRF)算法包括：
S1，将EnSRF算法中的控制变量投影于观测空间，形成4DEnSRF算法；
S11，将4DEnSRF算法的分析样本分为集合平均和样本扰动两个部分，并构造了集合平均的分析公式，具体包括步骤S111-S113：
S111，预防滤波发散考虑4DEnSRF算法的局地化，构造Kk矩阵；
S112，使用集合预报样本近似构造Gk矩阵；
S113，使用新构造的Kk矩阵和Gk矩阵替换EnSRF算法中的增益矩阵。
S12，定义了新的样本扰动分析公式(式(14)-式(16))。
S2，在4DEnSRF算法中引入流依赖时间窗口形成DFD-EnSRF算法，令背景误差协方差和同化时间窗口均随时空发生变化。其中，优选地，利用权值合理的调控时间窗口的大小。
下面详细描述本发明的双重流依赖集合平方根滤波算法。
控制变量投影于观测空间的目标函数为：
其中，K为总观测次数；观测算子被重新定义为G，其元素为(G1,G2,…,Gk,…)，k≤K，Gk是一个m×p维的矩阵，m是模式自由度，p是观测数据长度，它的作用是将初始模式状态投影到k时刻的观测空间，Gk＝HkMk；Hk为线性化的标准观测算子，Mk表示切线模式从初始时刻积分到k时刻；w是一个K×p维的控制变量，其元素为(w1 T,w2 T,…,wk T,…)，k≤K，且wk的长度与观测y长度一致。
该投影算法将分析过程限制于观测空间内，不仅减小了控制变量的维数，而且避免了对B矩阵的求逆。对该目标函数的极小化实际上等价于
wk就可以定义为：
wk可以被认为是某种观测空间的增量，它对分析结果的贡献可以表示为：
将式(3)代入式(4)，同时考虑多时刻观测就可以得到：
利用式(5)形式将EnSRF算法扩展为四维集合平方根滤波(4DEnSRF)算法，扩展后的4DEnSRF算法中将自带流依赖背景误差，
4DEnSRF集合平均分析公式为：
对应每个观测时刻，都存在一个卡尔曼增益矩阵，预防滤波发散考虑4DEnSRF的局地化，那么：
其中，分析时刻背景误差协方差可通过使用集合预报样本进行近似的方式投影到任意时刻的观测空间，从而得出：
在实际应用中，可由k时刻的预报样本代替，即可得到背景误差协方差的估计值：
在4DEnSRF算法中，通过式(6)对预报中的集合平均进行分析，样本扰动分析公式更新为：
其中，
增益矩阵(式(7)和式(15))中的背景误差项则由式(10)和式(11)进行估计。
在目前存在的大部分作物生长模型与遥感数据同化系统中常将时间窗口大小设为固定值或利用滤波算法对观测时序数据进行滤波处理，在观测数据不足或不稳定的情况下，同化系统难以高精度高效率的模拟作物参数或产量。因此本发明在4DEnSRF同化算法的基础上引入流依赖时间窗口，提出一种双重流依赖集合平方根滤波(DFD-EnSRF)算法，即背景误差协方差和同化时间窗口均随时空发生变化，以提高同化系统稳定性和同化系统输出精度。
在4DEnSRF算法中，时间窗口的大小决定了同化运算所吸收观测数据的数量和顺序。若时间窗口设置过小，同化运算所吸收的观测数据过少甚至只有一个，四维拓展后的4DEnSRF同化系统可能退化回原顺序同化系统；时间窗口设置过大，同化系统所吸收的观测数据过多，观测数据之间互相产生影响，增加同化系统的同化误差，尤其当所吸收的观测数据误差一致时，同化系统输出精度受到的影响越大。因此，本发明利用两次相邻观测时间差表示随时空发生变化的时间窗口以达到在同化算法中引入流依赖时间窗口的目的，并利用权值调控时间窗口大小，本发明中DFD-EnSRF算法的时间窗口定义为：
其中，wL为权值，可根据遥感卫星重访周期和作物生育周期经试验得出。为ti时刻观测数据在整个分析场时间序列的序号，为ti-1时刻观测数据在整个分析场时间序列的序号。
基于双重流依赖集合平方根滤波的作物模型遥感同化估产方法WOFOST作物生长模型是当前应用较多，影响力较大的作物生长模型之一。
作物LAI是作物生长模型中重要的状态变量，也是开展作物生长模型与遥感数据同化的前提数据基础。因此本发明选择冬小麦LAI作为观测算子基于WOFOST作物生长模型和双重流依赖集合平方根滤波算法对研究区作物LAI和产量进行模拟。如图1所示，基于双重流依赖集合平方根滤波的作物同化估产方法包括：
S1，将每个作物种植面积比例大于50％的像元格网的气象数据、作物参数以及土壤参数输入到WOFOST作物生长模型中，按照作物生长模型的演进方程，获得相应时刻不同模型参数运行下的模拟LAI集合作为背景场数据。
S2，利用双重流依赖集合平方根滤波(DFD-EnSRF)算法同化模型模拟LAI集合(S1中得到的背景场LAI)和外部观测LAI集合，得到分析LAI集合。本发明中观测误差由不同物候期LAI地面实测值或遥感反演值标准差计算而来。
S3，针对分析LAI集合，使用SCE-UA优化算法，通过迭代WOFOST作物生长模型寻找最优参数组合。当代价函数达到给定范围内的最小值时，基于该最优参数组合输出作物产量即为最终的同化产量结果。
本发明以冬小麦为例，阐述如何利用所提双重流依赖集合平方根滤波算法结合WOFOST作物生长模型进行作物模型遥感同化估产研究。研究区位于河北省衡水市，是我国典型的冬小麦种植区和重要的商品粮生产基地。研究区位于暖温带大陆性季风气候区，且地形较为平坦，主要作物种植制度为冬小麦-夏玉米一年两熟制。夏收作物主要为冬小麦，秋收作物主要为玉米和棉花等。其中，冬小麦从每年10月上旬播种至次年6月上中旬收获：11月下旬到次年2月中下旬为越冬期，次年2月底至3月上旬进入返青期，3月中下旬至4月中旬为起身拔节期，4月中下旬至5月上旬为孕穗抽穗期，5月中旬为开花灌浆期，6月上中旬冬小麦成熟收获。具体研究区分布如图2所示。
(1)单点冬小麦同化估产试验
在单点试验中，本发明选择冬小麦出苗期作为数据同化开始时间，单点冬小麦同化估产试验所需要数据主要包括：
1)典型试验区4个关键生育期(2014年3月23日(起身拔节期)，2014年4月13日(拔节孕穗期)，2014年5月12日(抽穗开花期)以及2014年5月29日(乳熟期))59个地面样方点所获取的冬小麦LAI实测值。
使用LAI-2200植被冠层分析仪测量采样点冬小麦LAI。LAI-2200植被冠层分析仪中的天空光设置为1，测量目标值设置为6，探头镜头盖张角设置为45°。测量过程中，观测者背向太阳以避免因太阳光直射所引起的LAI测量误差，同时LAI-2200植被冠层分析仪的探头需加盖镜头盖。测量时，首先沿第一垄冬小麦垄行45°角方向向前，逐渐靠向第二垄冬小麦，当移至第二垄冬小麦时同样沿垄行45°角方向向前再转向回第一垄冬小麦，每个采样点测量3次，3次测量结果平均后作为该采样点的冬小麦LAI测量结果。
2)典型试验区59个地面样方点获取的冬小麦单产实测值。
冬小麦完全成熟前的1～2天进行数据采集。首先在地面样方内使用采样框选取50cm×50cm大小区域的冬小麦植株，数出实际株数和总有效麦穗数并计算亩穗数。然后，对冬小麦麦穗进行脱粒处理并数出穗粒数，随机选取其中的1000粒干麦粒，称重3次取平均得到1000粒干麦粒的总重，即千粒重。最后，通过冬小麦单产计算公式，冬小麦单产＝亩穗数×穗粒数×千粒重/1000,000以及折算系数计算冬小麦单产。
3)2013～2014年衡水地区的逐日气象数据，其中包括日最高气温、日最低气温、太阳辐射强度、水汽压、平均风速以及降水量等。这些数据可从中国气象局国家气象信息中心(http://www.nmic.cn/)直接下载或通过计算获得。
4)土壤参数和田间管理参数等。其中一些参数对作物生长发育状态模拟结果影响不大，使用模型默认值或参考其他文献即可；另一些参数则对作物生长发育状态模拟结果影响较大，需要在合理的参数范围内根据实测数据，结合“试错法”进行调节，以获得适用于研究区作物生长发育状态模拟的参数组合。
数据同化优化参数的选择遵循以下原则：所选参数与作物生长发育和产量形成敏感程度较高且参数数值变化范围较大。由PCSE/WOFOST作物生长模型参数敏感性分析结果，初始根长(RDI)、初始总干物质(TDWI)以及叶片衰老系数(SPAN)3个参数对作物生长发育和产量形成的敏感程度较高且参数数值变化范围较大。同时，在理论上TDWI的变化能够很大程度地影响作物LAI的增长速度和最大值；SPAN不仅决定了叶片衰老速度和衰老时间，也可以显示作物受营养物质和病虫害的影响。考虑到上述2个作物参数的重要性和敏感性，本研究选取TDWI和SPAN作为待优化作物参数进行冬小麦同化估产研究。
为验证DFD-EnSRF算法在同化多时刻稀疏观测资料的可行性和有效性，本发明基于DFD-EnSRF算法和EnSRF算法进行了冬小麦LAI估算和产量模拟试验。图3显示了DFD-EnSRF算法和EnSRF算法同化效果的对比，其中EnSRF算法中粒子维数设置为150、扰动方差设置为0.15。由图3可以看出，DFD-EnSRF算法和EnSRF算法均可以获得较准确的冬小麦LAI估算和产量模拟结果。对比两种同化算法，基于DFD-EnSRF同化算法的冬小麦LAI估算和产量模拟精度均高于EnSRF同化算法。冬小麦LAI模拟值与实测值间的决定系数R2，调整确定系数Adjusted R2以及均方根误差RMSE分别为0.997，0.996以及0.138；冬小麦产量模拟值与实测值间的决定系数R2，调整确定系数Adjusted R2以及均方根误差RMSE分别为0.711，0.706以及173.6kg·ha-1。基于EnSRF同化算法的冬小麦LAI模拟值与实测值间的决定系数R2，调整确定系数Adjusted R2以及均方根误差RMSE分别为0.970，0.970以及0.394；冬小麦产量模拟值与实测值间的决定系数R2，调整确定系数Adjusted R2以及均方根误差RMSE分别为0.690，0.685以及224.1kg·ha-1。以上试验结果证明了DFD-EnSRF算法在冬小麦同化估产中的有效性和可行性。DFD-EnSRF同化算法不仅将顺序同化算法拓展至四维空间，有利于同化算法吸收原算法中被抛弃的非分析时刻观测资料，而且在4DEnSRF同化算法中引入了流依赖的时间窗口，使同化算法更好的适应作物生长模型与遥感数据同化系统中参数和变量的空间、时间流依赖特性，因此该算法具有较高的同化精度和较好的同化效果。
(2)衡水地区冬小麦同化估产试验
本发明选择冬小麦出苗期作为数据同化开始时间，衡水地区冬小麦同化估产试验所需要数据主要包括：
1)通过遥感影像反演获得的时序冬小麦LAI数据。本研究使用研究区冬小麦返青至成熟期(2014年3月至6月)间所获得的HJ-1、GF-1以及Landsat-8多光谱遥感数据，基于PROSAIL辐射传输模型反演得到研究区时序冬小麦LAI数据，且冬小麦LAI反演数据与地面实测数据间的决定系数R2均大于0.9，冬小麦LAI反演结果可满足数据同化研究输入数据的精度要求。
2)衡水地区县市级冬小麦产量官方统计数据。
3)2013～2014年衡水地区的逐日气象数据，其中包括日最高气温、日最低气温、太阳辐射强度、水汽压、平均风速以及降水量等。这些数据可从中国气象局国家气象信息中心(http://www.nmic.cn/)直接下载或通过计算获得。
4)衡水地区土壤参数和田间管理参数等。其中一些参数对作物生长发育状态模拟结果影响不大，使用模型默认值或参考其他文献即可；另一些参数则对作物生长发育状态模拟结果影响较大，需要在合理的参数范围内根据实测数据，结合“试错法”进行调节，以获得适用于研究区作物生长发育状态模拟的参数组合。
5)其他辅助数据。主要包括行政边界矢量图、2014年衡水地区冬小麦空间分布数据等，其中2014年衡水地区冬小麦空间分布数据由农业部农业遥感重点实验室提供。
使用DFD-EnSRF同化算法构建基于WOFOST作物生长模型的区域冬小麦遥感同化估产方案，将遥感数据反演的时序冬小麦LAI数据作为外部观测数据，在衡水地区进行基于DFD-EnSRF同化算法的区域冬小麦遥感同化估产研究，并利用衡水地区县市级冬小麦产量官方统计数据对区域冬小麦遥感同化估产结果进行精度评价。
图4显示了衡水地区冬小麦遥感同化估产空间结果。图5显示了衡水地区冬小麦遥感同化估产精度验证结果。
利用2014年衡水地区县市级冬小麦产量官方统计数据对基于DFD-EnSRF同化算法的冬小麦遥感同化估产结果进行精度评价。由图5精度验证结果可以看出，衡水地区冬小麦产量模拟值为6910kg·ha-1，与官方统计值的均方根误差RMSE为6.46％；分县来看，衡水地区11个县的冬小麦遥感同化产量模拟值与官方统计值的均方根误差RMSE均在13％以内，总体来看冬小麦遥感同化估产取得了较高的估产精度；冬小麦产量模拟值与官方统计值均方根误差最小的3个县依次为衡水市、景县以及武强县，均方根误差RMSE分别为1.74％、2.22％以及3.80％，说明以上3个县的冬小麦遥感同化估产效果最好；冬小麦产量模拟值与官方统计值均方根误差最大的3个县依次为冀州市、武邑县以及故城市，均方根误差RMSE分别为12.21％、11.05％以及9.55％，说明以上3个县的遥感同化估产效果最差。精度验证结果表明基于DFD-EnSRF同化算法的作物生长模型与遥感数据同化系统在衡水地区冬小麦同化估产应用中获得了较高估产精度和较好的制图效果，进一步证明了DFD-EnSRF同化方案在区域尺度上作物估产的可行性和有效性。
本发明虽然仅以冬小麦为例，阐述和实施了如何利用所提DFD-EnSRF算法进行作物模型遥感同化估产的思路和方法，但该方法对于一般作物估产具有普遍的适用性。
以上所述的实施例，只是本发明较优选的具体实施方式，本领域的技术人员在本发明技术方案范围内进行的通常变化和替换都应包含在本发明的保护范围内。
Claims (10)
1.一种双重流依赖集合平方根滤波同化算法，其特征在于，包括：
S1，将集合平方根滤波EnSRF算法中的控制变量投影于观测空间，形成四维集合平方根滤波4DEnSRF算法；
S2，在4DEnSRF算法中引入流依赖时间窗口形成双重流依赖集合平方根滤波DFD-EnSRF算法，令背景误差协方差和同化时间窗口均随时空发生变化。
2.根据权利要求1所述的双重流依赖集合平方根滤波算法，其特征在于，步骤S1包括：
S11，将4DEnSRF算法的分析样本分为集合平均和样本扰动两个部分，并构造了集合平均的分析公式；
S12，定义了新的样本扰动分析公式。
3.根据权利要求1所述的双重流依赖集合平方根滤波同化算法，其特征在于，步骤S11包括：
S111，预防滤波发散考虑4DEnSRF算法的局地化，构造Kk矩阵；
S112，使用集合预报样本近似构造Gk矩阵；
S113，使用新构造的Kk矩阵和Gk矩阵替换EnSRF算法中的增益矩阵。
4.根据权利要求3所述的双重流依赖集合平方根滤波同化算法，其特征在于，
在步骤S112中，分析时刻背景误差协方差可通过使用集合预报样本进行近似的方式投影到任意时刻的观测空间，从而构造Gk矩阵。
5.一种基于双重流依赖集合平方根滤波算法的作物模型遥感同化估产方法，其特征在于，包括：
T1，将气象数据、作物参数以及土壤参数输入到WOFOST作物生长模型中，按照作物生长模型演进方程，获得相应时刻不同模型参数运行下的模拟LAI集合作为背景场数据；
T2，利用DFD-EnSRF算法同化模型模拟LAI集合和外部观测LAI集合，得到分析LAI集合；
T3，针对分析LAI集合，使用SCE-UA优化算法，通过迭代WOFOST作物生长模型寻找最优参数组合并得到最终同化产量。
其中，DFD-EnSRF算法为：S1，将EnSRF算法中的控制变量投影于观测空间，形成4DEnSRF算法；S2，在4DEnSRF算法中引入流依赖时间窗口形成DFD-EnSRF同化算法，令背景误差协方差和同化时间窗口均随时空发生变化。
6.根据权利要求5所述的同化估产方法，其特征在于，还包括：
T4，当代价函数达到给定范围内的最小值时，基于该最优参数组合输出作物产量即为最终的同化产量结果。
7.根据权利要求5所述的同化估产方法，其特征在于，还包括：
在T1中，将每个作物种植面积比例大于50％的像元格网的气象数据、作物参数以及土壤参数输入到WOFOST作物生长模型中。
8.根据权利要求5所述的同化估产方法，其特征在于，步骤S1包括：
S11，将4DEnSRF算法的分析样本分为集合平均和样本扰动两个部分，并构造了集合平均的分析公式；
S12，定义了新的样本扰动分析公式。
9.根据权利要求8所述的同化估产方法，其特征在于，步骤S11包括：
S111，预防滤波发散考虑4DEnSRF算法的局地化，构造Kk矩阵；
S112，使用集合预报样本近似构造Gk矩阵；
S113，使用新构造的Kk矩阵和Gk矩阵替换EnSRF算法中的增益矩阵。
10.根据权利要求9所述的同化估产方法，其特征在于，
在步骤S112中，分析时刻背景误差协方差可通过使用集合预报样本进行近似的方式投影到任意时刻的观测空间，从而构造Gk矩阵。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910300249.7A CN110008621B (zh) | 2019-04-15 | 2019-04-15 | 基于双重流依赖集合平方根滤波同化算法的作物模型遥感同化估产方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910300249.7A CN110008621B (zh) | 2019-04-15 | 2019-04-15 | 基于双重流依赖集合平方根滤波同化算法的作物模型遥感同化估产方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110008621A true CN110008621A (zh) | 2019-07-12 |
| CN110008621B CN110008621B (zh) | 2023-01-10 |
Family
ID=67171940
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910300249.7A Active CN110008621B (zh) | 2019-04-15 | 2019-04-15 | 基于双重流依赖集合平方根滤波同化算法的作物模型遥感同化估产方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110008621B (zh) |
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111259324A (zh) * | 2020-01-07 | 2020-06-09 | 南京大学 | 卫星数据同化在垂直方向的适应性局地化方法及集合卡曼滤波天气同化预报方法 |
| CN111898922A (zh) * | 2020-08-10 | 2020-11-06 | 深圳市大德众和科技有限公司 | 一种多尺度农作物产量评估方法及系统 |
| CN113945997A (zh) * | 2021-09-26 | 2022-01-18 | 天津大学 | 一种基于解析四维集合变分的提高海洋预报精度的方法 |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103955860A (zh) * | 2014-04-17 | 2014-07-30 | 中国农业大学 | 一种基于集合卡尔曼滤波同化的区域作物估产方法 |
| CN104134095A (zh) * | 2014-04-17 | 2014-11-05 | 中国农业大学 | 一种基于尺度转换和数据同化的农作物产量估测方法 |
-
2019
- 2019-04-15 CN CN201910300249.7A patent/CN110008621B/zh active Active
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103955860A (zh) * | 2014-04-17 | 2014-07-30 | 中国农业大学 | 一种基于集合卡尔曼滤波同化的区域作物估产方法 |
| CN104134095A (zh) * | 2014-04-17 | 2014-11-05 | 中国农业大学 | 一种基于尺度转换和数据同化的农作物产量估测方法 |
Non-Patent Citations (3)
| Title |
|---|
| WANG SHIZHANG: "A four‐dimensional asynchronous ensemble square‐root filter (4DEnSRF) algorithm and tests with simulated radar data", 《QUARTERLY JOURNAL OF THE ROYAL METEOROLOGICAL SOCIETY》 * |
| 姜志伟等: "粒子滤波同化方法在CERES-Wheat作物模型估产中的应用", 《农业工程学报》 * |
| 王航等: "基于更新和同化策略相结合的遥感信息与水稻生长模型耦合技术的研究", 《生态学报》 * |
Cited By (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111259324A (zh) * | 2020-01-07 | 2020-06-09 | 南京大学 | 卫星数据同化在垂直方向的适应性局地化方法及集合卡曼滤波天气同化预报方法 |
| CN111259324B (zh) * | 2020-01-07 | 2022-07-12 | 南京大学 | 卫星数据同化在垂直方向的适应性局地化方法及集合卡曼滤波天气同化预报方法 |
| CN111898922A (zh) * | 2020-08-10 | 2020-11-06 | 深圳市大德众和科技有限公司 | 一种多尺度农作物产量评估方法及系统 |
| CN113945997A (zh) * | 2021-09-26 | 2022-01-18 | 天津大学 | 一种基于解析四维集合变分的提高海洋预报精度的方法 |
| CN113945997B (zh) * | 2021-09-26 | 2023-10-10 | 天津大学 | 一种基于解析四维集合变分的提高海洋预报精度的方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN110008621B (zh) | 2023-01-10 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109829234B (zh) | 一种基于高分辨率遥感数据和作物模型的跨尺度高精度动态作物长势监测和估产方法 | |
| Fang et al. | Corn‐yield estimation through assimilation of remotely sensed data into the CSM‐CERES‐Maize model | |
| Casanova et al. | Monitoring rice reflectance at field level for estimating biomass and LAI | |
| CN106372592B (zh) | 一种基于冬小麦面积指数的冬小麦种植面积计算方法 | |
| Huang et al. | BESS-Rice: A remote sensing derived and biophysical process-based rice productivity simulation model | |
| Dong et al. | Deriving maximum light use efficiency from crop growth model and satellite data to improve crop biomass estimation | |
| Ramarohetra et al. | How satellite rainfall estimate errors may impact rainfed cereal yield simulation in West Africa | |
| Zhou et al. | Detection of phenology using an improved shape model on time-series vegetation index in wheat | |
| CN110008621A (zh) | 基于双重流依赖集合平方根滤波同化算法和基于该算法的作物模型遥感同化估产方法 | |
| Srivastava et al. | Comparison of actual evapotranspiration of irrigated maize in a sub-humid region using four different canopy resistance based approaches | |
| McNider et al. | A real-time gridded crop model for assessing spatial drought stress on crops in the Southeastern United States | |
| Liu et al. | Winter wheat yield estimation based on assimilated Sentinel-2 images with the CERES-Wheat model | |
| CN116308866A (zh) | 一种基于冠层反射光谱的水稻穗生物量估算方法及系统 | |
| CN115439402A (zh) | 基于SAFY-Sugar模型的糖料蔗单产遥感估算方法和装置 | |
| Della Nave et al. | Calibrating APSIM for forage sorghum using remote sensing and field data under sub-optimal growth conditions | |
| Wu et al. | Adaptation of cotton production to climate change by sowing date optimization and precision resource management | |
| CN106203673A (zh) | 考虑水分胁迫的农作物生物量遥感估算方法 | |
| Allies et al. | Assimilation of multisensor optical and multiorbital SAR satellite data in a simplified agrometeorological model for rapeseed crops monitoring | |
| Jones | Remote sensing of plant stresses and its use in irrigation management | |
| He et al. | Characteristics of winter wheat evapotranspiration in Eastern China and comparative evaluation of applicability of different reference evapotranspiration models | |
| CN112667955B (zh) | 基于遥感估算区域尺度玉米潜在产量与产量差的方法与应用 | |
| Soleymani | Modeling the water requirement of wheat and safflower using the ET-HS model | |
| Ren et al. | Based on historical weather data to predict summer field-scale maize yield: Assimilation of remote sensing data to WOFOST model by ensemble Kalman filter algorithm | |
| Du et al. | Evaluating of simulated carbon flux phenology over a cropland ecosystem in a semiarid area of China with SiBcrop | |
| Samanta et al. | Deriving PAR use efficiency of wet season rice from bright sunshine hour data and canopy characteristics |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |