CN110766308A - Regional crop yield estimation method based on set assimilation strategy - Google Patents

Regional crop yield estimation method based on set assimilation strategy Download PDF

Info

Publication number
CN110766308A
CN110766308A CN201910989918.6A CN201910989918A CN110766308A CN 110766308 A CN110766308 A CN 110766308A CN 201910989918 A CN201910989918 A CN 201910989918A CN 110766308 A CN110766308 A CN 110766308A
Authority
CN
China
Prior art keywords
assimilation
lai
weight
matrix
year
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
Application number
CN201910989918.6A
Other languages
Chinese (zh)
Other versions
CN110766308B (en
Inventor
陶福禄
陈一
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Institute of Geographic Sciences and Natural Resources of CAS
Original Assignee
Institute of Geographic Sciences and Natural Resources of CAS
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Institute of Geographic Sciences and Natural Resources of CAS filed Critical Institute of Geographic Sciences and Natural Resources of CAS
Priority to CN201910989918.6A priority Critical patent/CN110766308B/en
Publication of CN110766308A publication Critical patent/CN110766308A/en
Application granted granted Critical
Publication of CN110766308B publication Critical patent/CN110766308B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • 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/067—Enterprise or organisation modelling
    • 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—Systems or methods specially adapted for specific business sectors, e.g. utilities or tourism
    • G06Q50/02—Agriculture; Fishing; Mining

Abstract

The invention discloses a regional crop production estimation method based on a set assimilation strategy, which comprises the following steps: acquiring remote sensing LAI data of a forecast year Y and a previous year M in a research area; carrying out differential assimilation aiming at the Y-1 year to obtain a plurality of grid scale simulation yields, a plurality of corresponding groups of correlation coefficients r and a root mean square error rmse; calculating a fitting index F according to r and rmse to obtain a plurality of fitting indexes; selecting the maximum fitting index as the optimal result, and recording the corresponding assimilation weight as the optimal assimilation weight H of the (y-1) th yearop_y‑1(ii) a Repeating the steps S2-S5 for the y-2, y-3op_y‑2，Hop_y‑3......Hop_y‑M(ii) a Expanded assimilation weight value-taking interval Hlow‑HupThen, the assimilation weight value-taking interval H is utilizedlow‑HupSimulating for the predicted year YAnd (6) estimating yield.

Description

Regional crop yield estimation method based on set assimilation strategy
Technical Field
The invention relates to the technical field of agricultural remote sensing, in particular to a regional crop production estimation method based on a set assimilation strategy.
Background
Crop models often face the problem of insufficient input data when performing region-scale crop yield estimation. There are often significant spatial differences in surface features, near-surface environment, and crop management measures, and some necessary model input data, such as crop phenology data, meteorological data, crop management information, etc., are recorded at site scale, which makes it difficult for a crop model to obtain enough data to represent spatial heterogeneity of key factors such as initial conditions, crop parameters, growth process, etc., when applied to regional scale. The satellite remote sensing data provides continuous monitoring data of large-range earth surface information and can reflect the space continuity and time sequence change characteristics of the earth surface information. The advantage effectively supplements the weakness of the application of the crop model on the regional scale, so that the crop model and remote sensing data assimilation technology becomes an important way for improving the regional simulation precision of the crop model at present and is one of the hot directions in the field of crop assessment research in recent years.
In the process of assimilating the crop model and the remote sensing data, one factor influencing the key related to the assimilation result is uncertainty of the remote sensing data. The main effect of the uncertainty is to quantify the assimilation weight of the model analog value and the remote sensing data. Whether the assimilation weight is reasonable or not directly influences the final assimilation result precision, once the assimilation weight is set to have deviation, the assimilation result can also generate obvious deviation, and the application effect of the assimilation technology is seriously influenced. The best method for quantifying the uncertainty of the remote sensing data is to compare and analyze the remote sensing data with the earth surface observation result, but the method requires earth surface observation data with space-time density, so that the labor cost and the time cost are very high, and the method is difficult to popularize. The existing crop model and remote sensing data assimilation technology usually adopt a simplified method, for example, directly setting an assimilation weight value to reduce the requirement on earth surface observation. This approach has achieved some positive findings, but has limited assistance for practical assimilation assessment applications. The reason for this is that in practical applications, the assimilation weights need to exist as prior parameters to drive the assimilation process, and the existing research does not relate to how to accurately determine the prior of the assimilation weights. In fact, since the assimilation weights have many influence factors and are varied in the year, the optimal assimilation weights cannot be directly derived through historical experience, which greatly reduces the practicability of the assimilation technology. At present, the fact that the optimal assimilation weight cannot be known a priori is a main limitation of assimilation estimation technology towards practical application, how to break through the limitation is that the assimilation estimation technology develops the assimilation estimation and obtains an estimation result with sufficient precision on the premise that the optimal assimilation weight is unknown, the technology is still blank, and breakthrough is needed urgently.
Accordingly, new techniques are needed to at least partially address the above-described limitations of the prior art.
Disclosure of Invention
The assimilation estimation scheme based on the set assimilation strategy is provided, the estimation precision of the assimilation estimation scheme is close to the estimation precision based on the optimal assimilation weight, namely, a sufficiently accurate estimation result can be obtained under the condition that the optimal assimilation weight is unknown, and the practicability of the assimilation estimation technology is greatly improved.
According to one aspect of the invention, a regional crop estimation method based on a set assimilation strategy is provided, and comprises the following steps:
s1, acquiring remote sensing LAI data of a forecast year Y and a previous M years in a research area, wherein M is a natural number more than 3;
s2, aiming at the Y-1 year, carrying out space difference assimilation operation by using a crop model and the remote sensing LAI data, wherein the space difference assimilation operation equation is as follows:
Figure BDA0002237926060000031
wherein
Figure BDA0002237926060000032
Is a normalized analysis matrix of the assimilation time point k,
Figure BDA0002237926060000033
and
Figure BDA0002237926060000034
respectively representing a normalized simulation LAI matrix and a normalized remote sensing inversion LAI matrix of an assimilation time point k, wherein H is an assimilation weight and the value range is 0.01-1.00;
in the value range of 0.01-1.00, H takes values by step length S respectively, and assimilation estimation operation is carried out under the condition of 1/S different assimilation weight values, so that 1/S corresponding grid scale simulation yield is obtained;
s3, based on corresponding yield recorded data in the research area, performing precision evaluation on each grid scale simulation yield, representing by using a correlation coefficient r and a root mean square error rmse, and obtaining a 1/S group r and rmse in total;
s4, calculating a fitting index F according to r and rmse obtained in S3 to obtain 1/S fitting indexes;
the calculation formula of the fitting index F is as follows:
Figure BDA0002237926060000035
where Fi is the fitting index of the ith assimilation estimate, riAnd rmseiCorrelation coefficient and root mean square error, r, representing the ith assimilation estimatemax,rmin,rmsemaxAnd rmseminRespectively represent the maximum value and the minimum value of r and the maximum value and the minimum value of rmse in the 1/S group r and rmse;
s5, selecting the maximum fitting index as the optimal result, and recording the corresponding assimilation weight as the optimal assimilation weight H of the (y-1) th yearop_y-1；
S6, repeating the steps S2-S5 for y-2, y-3op_y-2，Hop_y-3......Hop_y-M；
S7 based on Hop_y-1，Hop_y-2，Hop_y-3......Hop_y-MMaximum value of (1)opmaxAnd a minimum value HopminDefining an expansion assimilation weight value range, wherein the lower limit of the expansion assimilation weight value range is called HlowThe upper limit is referred to as Hup，HlowObtained by a 10% expansion of the minimum, HupThe maximum value is obtained by expanding 10% outwards, and the value range of the expansion assimilation weight is not more than 0.99 and not less than 0.01;
s8, expanding assimilation weight value section H based on S7low-HupAnd performing simulated assessment on the predicted year Y, wherein the simulated assessment comprises the following steps: value range H is taken by expanding assimilation weightlow-HupUnder the condition of replacing the value range of 0.01-1.00, the step S2 is repeated, and the space difference assimilation operation is expandedComputing, proceeding together (H)up-Hlow) (S + 1) assimilation operations, from which the corresponding (H) is obtainedup-Hlow) And simulating the yield by the/S +1 grid scales, and taking the average value as the simulated yield of the predicted year Y.
According to an embodiment of the present invention, M in step S1 is 4, although other values may be selected as appropriate.
According to an embodiment of the invention, the crop model is a MCWLA series model, or other crop models may be selected.
According to one embodiment of the invention, the step size S is 0.01, or other values are chosen as the case may be, for example 0.02, etc.
According to one embodiment of the invention, the remotely sensed LAI data may be selected from the group consisting of copernius LAI, GLASS LAI and GLOBMAP LAI, or other remotely sensed LAI products.
According to one embodiment of the invention, the crop is wheat, or another crop such as corn.
According to an embodiment of the present invention, in step S2, the normalized analysis matrix, the normalized simulated LAI matrix and the normalized remote sensing inversion LAI matrix are obtained by a normalization method, so that the original matrix is normalized to a normalized matrix with a value range of (0, 1), which is represented by the following formulas (3) and (4):
Figure BDA0002237926060000051
wherein the content of the first and second substances,
Figure BDA0002237926060000053
representing the simulated LAI matrix at the point of assimilation k, DkRepresenting the corresponding remote sensing LAI matrix at the same time point k;
Figure BDA0002237926060000054
andrespectively represent by
Figure BDA0002237926060000056
And DkThe converted normalized simulated LAI matrix and the normalized remote sensing LAI matrix are obtained; sfAnd SDIs a normalization function used to normalize the original matrix.
According to an embodiment of the invention, in step S2, the differential assimilation further comprises:
Figure BDA0002237926060000057
by normalizing function SfIs inverse function of
Figure BDA0002237926060000058
Inverse normalized transformation to obtain analytical LAI matrix
Figure BDA0002237926060000059
And is composed of
Figure BDA00022379260600000510
The driving model runs to the next assimilation time point k + 1:
Figure BDA00022379260600000511
compared with the prior art, the method integrates the idea of collective operation on the basis of the assimilation technology of remote sensing data and crop models. A solution strategy is provided for the practical problem of how to carry out assimilation estimation and obtain an estimation result with sufficient precision under the condition that the optimal assimilation weight cannot be known a priori. By applying the assimilation estimation scheme based on the set assimilation strategy, the estimation precision is close to the estimation precision based on the optimal assimilation weight, namely, a sufficiently accurate estimation result can be obtained under the condition that the optimal assimilation weight is unknown. The invention makes up the deficiency of the assimilation assessment technology when oriented to business application, and greatly improves the practical application capability of the assimilation assessment technology.
Drawings
The same reference numbers in the drawings identify the same or similar elements or components. The objects and features of the present invention will become more apparent in view of the following description taken in conjunction with the accompanying drawings, in which:
FIG. 1 is a schematic flow diagram of a regional crop estimation method based on collective assimilation strategies according to one embodiment of the invention.
FIG. 2 is a schematic illustration of an investigation region for carrying out the method of the invention according to one embodiment of the invention;
fig. 3 is a diagram of the effect of the method according to the invention when applied to different remote sensing LAI products.
Detailed Description
For a clear description of the solution according to the invention, preferred embodiments are given below and are described in detail with reference to the accompanying drawings. The following description is merely exemplary in nature and is not intended to limit the present disclosure, application, or uses
It should be understood that the crop models and telemetry data referenced in the present invention are known per se, e.g., individual sub-modules of the models, various parameters, operating mechanisms, etc., and the present invention focuses on the collective assimilation process between the crop models and the telemetry data.
FIG. 1 is a schematic illustration according to an embodiment of the present invention. Referring to fig. 1, the regional crop estimation method based on collective assimilation strategy according to the present invention includes the following steps:
s1, acquiring remote sensing LAI data of a forecast year Y and a previous M years in a research area, wherein M is a natural number more than 3; the diagram shows M being 4, i.e., remotely sensed LAI data of the (y-1, y-2, y-3, y-4) th is acquired, wherein the remotely sensed LAI data may be selected from Copernicius LAI, GLASS LAI and GLOBMAP LAI, or other suitable data. That is, the method can be applied to various remote sensing LAI data, and the spatial resolution and the time frequency of assimilation estimation operation are consistent with the spatial-temporal resolution of the remote sensing LAI data.
And S2, carrying out differential assimilation by using a crop model (shown as MCWLA-Wheat model in the figure) and the remote sensing LAI data for the Y-1 year, wherein the differential assimilation comprises a spatial differential assimilation operation equation:
Figure BDA0002237926060000071
wherein
Figure BDA0002237926060000072
Is a normalized analysis matrix of the assimilation time point k,
Figure BDA0002237926060000073
and
Figure BDA0002237926060000074
respectively representing a normalized simulation LAI matrix and a normalized remote sensing inversion LAI matrix of an assimilation time point k, wherein H is an assimilation weight and the value range is 0.01-1.00; it should be understood that the superscripts in the formula are for the purpose of distinguishing matrices, a is an analysis abbreviation and denotes an analysis matrix; f is a forecast abbreviation and represents a forecast matrix, namely a model simulation matrix; n is a normalization abbreviation and represents normalization.
And in the value range of 0.01-1.00, H takes values respectively according to step length S (such as 0.01), and under the condition of 1/S (such as 100) different assimilation weight values, assimilation estimation operation is respectively carried out, so that 1/S (100) grid scale simulation yield is obtained correspondingly.
More specifically, the primary assimilation process of differential assimilation includes:
1) the simulated LAI matrix and the corresponding remote sensing LAI matrix are normalized at a time point (called a assimilation time point) where remote sensing data exists, so that the original matrix is normalized to a normalized matrix (following formulas (3) and (4)) with a value range of (0, 1), and the normalization operation method is many, and can be obtained by dividing the LAI matrix by the maximum value of matrix elements:
Figure BDA0002237926060000081
Figure BDA0002237926060000082
wherein the content of the first and second substances,
Figure BDA0002237926060000083
representing the simulated LAI matrix at the point of assimilation k, DkRepresenting the corresponding remote sensing LAI matrix at the same time point k;and
Figure BDA0002237926060000085
respectively represent by
Figure BDA0002237926060000086
And DkThe converted normalized simulated LAI matrix and the normalized remote sensing LAI matrix are obtained; sfAnd SDIs a normalization function used to normalize the original matrix.
2) Assimilating the normalized simulated LAI matrix and the normalized remote sensing LAI matrix using assimilation equation (1):
Figure BDA0002237926060000087
whereinThe method is characterized in that the method is a normalized analysis matrix of an assimilation time point k, H is an assimilation weight, the value range is 0.01-1.00, and the values of H are respectively obtained by step length 0.01, so that 100 different assimilation weight values are obtained.
3)
Figure BDA0002237926060000089
By normalizing function SfIs inverse function of
Figure BDA00022379260600000810
Inverse normalized transformation to obtain analytical LAI matrixAnd is composed ofThe driving model runs to the next assimilation time point k + 1:
Figure BDA00022379260600000813
wherein, the one-time complete assimilation estimation operation process in the step comprises the following steps: the model starts day by day simulation from the seeding stage, starts assimilation with remote sensing data from the assimilation starting point, the assimilation process continues to the assimilation end point, and the model continues to run to the maturation stage to output grid scale simulation yield.
S3, based on corresponding yield recorded data in the research area, performing precision evaluation on each grid scale simulation yield, representing by using a correlation coefficient r and a root mean square error rmse, and obtaining a 1/S group r and rmse in total; the estimation precision evaluation is based on the spatial scale of actual statistical yield, and the estimation result obtained by assimilation operation is a grid scale and needs to be aggregated to the scale of the statistical yield, such as county level.
S4, calculating a fitting index F according to r and rmse obtained in S3 to obtain 1/S fitting indexes;
the calculation formula of the fitting index F is as follows:
Figure BDA0002237926060000091
wherein, FiIs the fitting index, r, of the ith assimilation estimateiAnd rmseiCorrelation coefficient and root mean square error, r, representing the ith assimilation estimatemax,rmin,rmsemaxAnd rmseminRespectively represent the maximum value and the minimum value of r and the maximum value and the minimum value of rmse in the 1/S group r and rmse; a larger F value represents a higher accuracy of the estimated result.
S5, selecting the maximum fitting index as the optimal result, and recording the corresponding assimilation weight as the y-1 yearOptimal assimilation weight Hop_y-1；
S6, repeating the steps S2-S5 for y-2, y-3op_y-2，Hop_y-3......Hop_y-MIn the figure, M is shown as 4, i.e. Hop_y-4；
S7 based on Hop_y-1，Hop_y-2，Hop_y-3......Hop_y-MMaximum value of (1)opmaxAnd a minimum value HopminDefining an expansion assimilation weight value range, wherein the lower limit of the expansion assimilation weight value range is called HlowThe upper limit is referred to as Hup，HlowObtained by a 10% expansion of the minimum, HupThe maximum value is obtained by expanding 10% outwards, and the value range of the expansion assimilation weight is not more than 0.99 and not less than 0.01; that is to say that the first and second electrodes,
Hlow＝max(0.9*Hopmin,0.01) (6)
Hup＝min(1.1*Hopmax,0.99) (7)
s8, expanding assimilation weight value section H based on S7low-HupAnd performing simulated assessment on the predicted year Y, wherein the simulated assessment comprises the following steps: value range H is taken by expanding assimilation weightlow-HupUnder the condition of replacing the value range of 0.01-1.00, repeating the step S2, and carrying out space difference assimilation operation (H)up-Hlow) (S + 1) (101) assimilation operations whereby the corresponding (H) is obtainedup-Hlow) (ii) simulating the yield (i.e. n sets of results in the figure) by using/S +1 grid scales, and taking the average value as the simulated yield of the predicted year Y.
The process according to the invention is further illustrated below with reference to specific examples:
the specific application of the method of the invention is exemplified below by taking the assessment of winter wheat in the plain area of North China as an example.
The middle part of North China plain was selected as the study area, and the study time was 2008-2015. By adopting an MCWLA-Wheat model, the used remote sensing LAI data comprises Copernicius LAI (spatial resolution is 1km multiplied by 1km, time frequency is one period every 10 days), GLASS LAI (spatial resolution is 1km multiplied by 1km, time frequency is one period every 8 days) and GLOBMAP LAI (spatial resolution is 0.08 degrees multiplied by 0.08 degrees, time frequency is one period every 8 days), and each remote sensing LAI product is respectively assimilated with the MCWLA-Wheat model to estimate the yield of winter Wheat. The study area is shown in FIG. 2.
Based on the set assimilation strategy provided by the invention, the assimilation and estimation implementation process is as follows (taking Copernics LAI as an example, the implementation processes of GLASS LAI and GLOBMAP LAI are consistent with Copernics LAI):
in this example, 2012-2015 is taken as the estimated production period. For each year of assessment, assimilation assessments are made with their first 4 years as historical periods to provide a priori knowledge about optimal assimilation weights. That is, assimilation assessment for year y, it is necessary to perform assimilation assessment first for its first 4 years (i.e., y-1, y-2, y-3, y-4) to provide prior knowledge about optimal assimilation weights for year y. For example, for the assimilation assessment in 2012, the assimilation assessment operation needs to be performed in advance in 2008-2011 to obtain the optimal assimilation weight year by year. Assimilation is performed using a spatial differential assimilation algorithm.
The operation steps of the present invention will be described in detail with reference to FIG. 1.
As shown in the figure, in the year y-1, a MCWLA-wheel model and remote sensing LAI data are used for carrying out spatial difference assimilation operation, the assimilation weight H is 0.01-1.00, the step length is 0.01, namely, the H has 100 values, and therefore 100 times of assimilation estimated production operation can be driven, and 100 groups of estimated production results are obtained. The complete process of each assimilation estimation operation is as follows:
the growth of the winter Wheat is simulated day by using MCWLA-Wheat from the sowing period of the winter Wheat. The assimilation starting point is the green turning period of the winter wheat, and the assimilation end point is the mature period of the winter wheat. Between the green turning period and the mature period, the model simulated LAI and the remote sensing LAI data are assimilated by using a spatial difference assimilation algorithm at each time point (referred to as an assimilation time point) having the remote sensing observation. And (5) ending the assimilation and model simulation process in the mature period, and outputting the winter wheat simulation yield grid by grid.
After one assimilation estimation operation is completed, aggregating estimation results of grid scales to county scale and comparing with statistical yield to obtain a resultThe correlation coefficient r and the root mean square error rmse represent the estimated production accuracy. 100 sets of r and rmse can be obtained by repeating the assimilation operation 100 times by using different assimilation weights H, and each assimilation weight H corresponds to one set of r and rmse values. On the basis, calculating a fitting index F of each assimilation estimation result, wherein the calculation formula is the formula (2): wherein r ismax,rmin,rmsemaxAnd rmseminRepresents the maximum and minimum values of r and rmse in 100 groups of r and rmse, respectively.
The 100 fitting indexes F corresponding to the 100 assimilation weights H are obtained by calculation of formula (2). Selecting the maximum fitting index as the optimal fitting index, and taking the corresponding assimilation weight H as the optimal assimilation weight of year y-1, named as Hop_y-1。
The process of calculating the Hop _ y-1 in y-1 year needs to be repeated in y-2, y-3 and y-4 years respectively to obtain the annual optimal assimilation weight Hop_y-2，Hop_y-3，Hop_y-4。
op_y-1，Hop_y-2，Hop_y-3，Hop_y-4Maximum value of (1)opmaxAnd a minimum value HopminAnd (3) defining an assimilation weight value range, wherein the value range is formed by outwards expanding the maximum value range and the minimum value range by 10%, and the value range is not more than 0.99 and not less than 0.01. The lower limit of the interval is called HlowThe upper limit is referred to as HupSee, in particular, the above equations (6) and (7).
After the calculation is finished, carrying out collective assimilation estimation operation on the year y, wherein the assimilation weight H adopted by the operation is based on the obtained HlowAnd HupThe value range [ H ] enclosed by the twolow,Hup]The values are incremented by 0.01 steps. Co-operation (H)up-Hlow) 0.01+1 assimilation operations. The mean of the estimates obtained from the set operation is used as the final estimate, and the simulation results are shown in fig. 3.
Aiming at different remote sensing LAI data, the estimation mode is consistent with the process, and the difference is that when different remote sensing LAI data are used, the spatial resolution of the grid scale of the model operation needs to be consistent with the remote sensing LAI data, and the assimilation time point interval needs to be consistent with the time frequency of the remote sensing LAI data. In addition, when different remote sensing LAI data are used, the optimal assimilation weight obtained each year changes.
In order to verify the technical effects of the present invention, the following methods were used to evaluate the effects of the present invention:
firstly, the optimal assimilation weight calculation method applied to the years y-1 to y-4 is also applied to the years y, and the following parameters of the years y are obtained: comprising rmax，rmin，rmsemax，rmseminAnd maximum fitting index F, referred to as Fop. Then, the collective assimilation estimation result of year y is compared with the actual statistical yield, r and rmse are calculated, and r is combined with the abovemax，rmin，rmsemax，rmseminCalculating the F value of the collective assimilation estimate, called Fs. According to FopAnd FsCalculating the accuracy index pc:
Figure BDA0002237926060000131
the accuracy index pc is used for representing the percentage of the estimated production accuracy degree of the set assimilation compared with the estimated production accuracy degree of the optimal assimilation weight, and the larger the pc is, the higher the accuracy is.
In addition, in order to compare the estimation accuracy, an accuracy index pc of some typical estimation methods can be calculated, and the methods include:
a. using the optimal assimilation weight H of the year y-1 to drive the assimilation estimation operation of the year y;
b. using the average value of the optimal assimilation weight H of the year y-1 to year y-2 to drive the assimilation estimation operation of the year y;
c. using the average value of the optimal assimilation weight H of the year y-1 to year y-3 to drive the assimilation estimation operation of the year y;
d. using the average value of the optimal assimilation weight H of the year y-1 to year y-4 to drive the assimilation estimation operation of the year y;
e. using the optimal assimilation weight H of the year y-1 to year y-4 to drive the assimilation estimation operation of the year y and taking the result mean value as an estimation result;
f. the maximum and minimum value interval of the optimal assimilation weight H of the year y-1 to the year y-4, namely [ Hopmin,Hopmax]Step length is 0.01, assimilation estimation operation of year y is driven, and the result mean value is used as an estimation result;
g. the method provided by the invention is used for estimating the yield of the year y.
The methods a-g are used for carrying out assimilation estimation operation in 2012-2015 respectively, and the average value of the accuracy index pc of the four-year estimation result is used for measuring the comprehensive implementation effect of the method. The results are shown in fig. 3, which shows that the collective assimilation strategy proposed by the present invention shows excellent implementation effects when different remote sensing LAI products are applied. For example, for CoperniciusLAI, the estimation accuracy may be up to 98.1% of the estimation accuracy using the optimal assimilation weights, and for GLASS LAI and GLOBMAP LAI, the estimation accuracy may be up to 90.8% and 85.0% of the estimation accuracy using the optimal assimilation weights, respectively. The use of remote sensing products with higher resolution can achieve relatively better accuracy, closer to the effect of estimation when optimal assimilation weights are used.
In addition, the evaluation accuracy of the set assimilation strategy provided by the invention is stably superior to other evaluation methods in assimilation evaluation performed by applying different remote sensing LAI data, and the method has universality. The results show that the set assimilation strategy provided by the invention can carry out assimilation estimation under the condition of unknown optimal assimilation weight and obtain an estimation result with sufficient precision. The invention makes up the deficiency of the assimilation assessment technology when oriented to business application, and greatly improves the practical application capability of the assimilation assessment technology.
The principles and embodiments of the present invention have been described herein using specific examples, which are presented solely to aid in the understanding of the apparatus and its core concepts; meanwhile, for a person skilled in the art, according to the idea of the present invention, there may be variations in the specific embodiments and the application scope, and in summary, the content of the present specification should not be construed as a limitation to the present invention.

Claims (8)

1. A regional crop production estimation method based on a set assimilation strategy comprises the following steps:
s1, acquiring remote sensing LAI data of a forecast year Y and a previous M years in a research area, wherein M is a natural number more than 3;
s2, aiming at the Y-1 year, carrying out differential assimilation, including carrying out spatial differential assimilation operation by using a crop model and the remote sensing LAI data, wherein the spatial differential assimilation operation equation is as follows:
Figure FDA0002237926050000011
wherein
Figure FDA0002237926050000012
Is a normalized analysis matrix of the assimilation time point k,
Figure FDA0002237926050000013
andrespectively representing a normalized simulation LAI matrix and a normalized remote sensing inversion LAI matrix of an assimilation time point k, wherein H is an assimilation weight and the value range is 0.01-1.00;
in the value range of 0.01-1.00, H takes values by step length S respectively, and assimilation estimation operation is carried out under the condition of 1/S different assimilation weight values, so that 1/S corresponding grid scale simulation yield is obtained;
s3, based on corresponding yield recorded data in the research area, performing precision evaluation on each grid scale simulation yield, representing by using a correlation coefficient r and a root mean square error rmse, and obtaining a 1/S group r and rmse in total;
s4, calculating a fitting index F according to r and rmse obtained in S3 to obtain 1/S fitting indexes;
the calculation formula of the fitting index F is as follows:
Figure FDA0002237926050000015
where Fi is the fitting index of the ith assimilation estimate, riAnd rmseiCorrelation coefficient and root mean square error, r, representing the ith assimilation estimatemax,rmin,rmsemaxAnd rmseminRespectively represent the maximum value and the minimum value of r and the maximum value and the minimum value of rmse in the 1/S group r and rmse;
s5, selecting the maximum fitting index as the optimal result, and recording the corresponding assimilation weight as the optimal assimilation weight H of the (y-1) th yearop_y-1；
S6, repeating the steps S2-S5 for y-2, y-3op_y-2，Hop_y-3......Hop_y-M；
S7 based on Hop_y-1，Hop_y-2，Hop_y-3......Hop_y-MMaximum value of (1)opmaxAnd a minimum value HopminDefining an expansion assimilation weight value range, wherein the lower limit of the expansion assimilation weight value range is called HlowThe upper limit is referred to as Hup，HlowObtained by a 10% expansion of the minimum, HupThe maximum value is obtained by expanding 10% outwards, and the value range of the expansion assimilation weight is not more than 0.99 and not less than 0.01;
s8, expanding assimilation weight value section H based on S7low-HupAnd performing simulated assessment on the predicted year Y, wherein the simulated assessment comprises the following steps: value range H is taken by expanding assimilation weightlow-HupUnder the condition of replacing the value range of 0.01-1.00, repeating the step S2, and carrying out space difference assimilation operation (H)up-Hlow) (S + 1) assimilation operations, from which the corresponding (H) is obtainedup-Hlow) And simulating the yield by the/S +1 grid scales, and taking the average value as the simulated yield of the predicted year Y.
2. The regional crop estimation method based on collective assimilation strategy of claim 1, wherein in step S1, M is 4.
3. The regional crop estimation method based on collective assimilation strategy of claim 1, characterized in that the crop model is the MCWLA series model.
4. The regional crop estimation method based on collective assimilation strategy of claim 1, characterized in that the step size S is 0.01.
5. The regional crop yield assessment method based on collective assimilation strategies of claim 4 wherein the remotely sensed LAI data is selected from Copernicus LAI, GLASS LAI and GLOBMAP LAI.
6. The regional crop estimation method based on collective assimilation strategy of claim 1, characterized in that the crop is wheat.
7. The regional crop estimation and production method based on collective assimilation strategy of claim 1, wherein in step S2, the normalized analysis matrix, the normalized simulated LAI matrix and the normalized remote sensing inversion LAI matrix are obtained by a normalization method, and the original matrix is normalized to a normalized matrix with a value range of (0, 1), which is represented by the following formulas (3) and (4):
Figure FDA0002237926050000031
Figure FDA0002237926050000032
wherein the content of the first and second substances,
Figure FDA0002237926050000033
representing the simulated LAI matrix at the point of assimilation k, DkRepresenting the corresponding remote sensing LAI matrix at the same time point k;
Figure FDA0002237926050000034
and
Figure FDA0002237926050000035
respectively represent by
Figure FDA0002237926050000036
And DkThe converted normalized simulated LAI matrix and the normalized remote sensing LAI matrix are obtained; sfAnd SDIs a normalization function used to normalize the original matrix.
8. The regional crop estimation method based on collective assimilation strategy of claim 1, characterized in that the differential assimilation further comprises:
Figure FDA0002237926050000037
by normalizing function SfIs inverse function of
Figure FDA0002237926050000038
Inverse normalized transformation to obtain analytical LAI matrix
Figure FDA0002237926050000039
And is composed of
Figure FDA00022379260500000310
The driving model runs to the next assimilation time point k + 1:
Figure FDA00022379260500000311
CN201910989918.6A 2019-10-17 2019-10-17 Regional crop yield estimation method based on set assimilation strategy Active CN110766308B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910989918.6A CN110766308B (en) 2019-10-17 2019-10-17 Regional crop yield estimation method based on set assimilation strategy

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910989918.6A CN110766308B (en) 2019-10-17 2019-10-17 Regional crop yield estimation method based on set assimilation strategy

Publications (2)

Publication Number Publication Date
CN110766308A true CN110766308A (en) 2020-02-07
CN110766308B CN110766308B (en) 2020-07-10

Family

ID=69332455

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910989918.6A Active CN110766308B (en) 2019-10-17 2019-10-17 Regional crop yield estimation method based on set assimilation strategy

Country Status (1)

Country Link
CN (1) CN110766308B (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111667167A (en) * 2020-06-03 2020-09-15 福建慧政通信息科技有限公司 Agricultural grain yield estimation method and system

Citations (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050234691A1 (en) * 2004-04-20 2005-10-20 Singh Ramesh P Crop yield prediction
CN101044823A (en) * 2007-04-17 2007-10-03 华南师范大学 Method for estimating crop energy utilization rate and predetermining the yield
CN102162850A (en) * 2010-04-12 2011-08-24 江苏省农业科学院 Wheat yield remote sensing monitoring and forecasting method based on model
CN103955860A (en) * 2014-04-17 2014-07-30 中国农业大学 Regional crop yield estimation method based on ensemble Kalman filter assimilation
CN104134095A (en) * 2014-04-17 2014-11-05 中国农业大学 Crop yield estimation method based on scale transformation and data assimilation
CN102651096B (en) * 2012-04-28 2015-06-03 中国农业大学 Method for estimating yield of winter wheat by assimilating characteristics of leaf area index time-sequence curve
CN107941713A (en) * 2017-10-17 2018-04-20 河海大学 A kind of rice yield estimation method based on coupling crop modeling assimilation spectral reflectivity
CN109766871A (en) * 2019-01-30 2019-05-17 中国科学院地理科学与资源研究所 A kind of region Crop Estimation Method based on spatial diversity assimilation remotely-sensed data and crop modeling
CN109800921A (en) * 2019-01-30 2019-05-24 北京师范大学 A kind of Regional Fall Wheat yield estimation method based on remote sensing phenology assimilation and particle swarm optimization algorithm
CN110222870A (en) * 2019-05-05 2019-09-10 中国农业大学 Assimilate the Regional Fall Wheat yield estimation method of satellite fluorescence data and crop growth model

Patent Citations (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050234691A1 (en) * 2004-04-20 2005-10-20 Singh Ramesh P Crop yield prediction
CN101044823A (en) * 2007-04-17 2007-10-03 华南师范大学 Method for estimating crop energy utilization rate and predetermining the yield
CN102162850A (en) * 2010-04-12 2011-08-24 江苏省农业科学院 Wheat yield remote sensing monitoring and forecasting method based on model
CN102651096B (en) * 2012-04-28 2015-06-03 中国农业大学 Method for estimating yield of winter wheat by assimilating characteristics of leaf area index time-sequence curve
CN103955860A (en) * 2014-04-17 2014-07-30 中国农业大学 Regional crop yield estimation method based on ensemble Kalman filter assimilation
CN104134095A (en) * 2014-04-17 2014-11-05 中国农业大学 Crop yield estimation method based on scale transformation and data assimilation
CN107941713A (en) * 2017-10-17 2018-04-20 河海大学 A kind of rice yield estimation method based on coupling crop modeling assimilation spectral reflectivity
CN109766871A (en) * 2019-01-30 2019-05-17 中国科学院地理科学与资源研究所 A kind of region Crop Estimation Method based on spatial diversity assimilation remotely-sensed data and crop modeling
CN109800921A (en) * 2019-01-30 2019-05-24 北京师范大学 A kind of Regional Fall Wheat yield estimation method based on remote sensing phenology assimilation and particle swarm optimization algorithm
CN110222870A (en) * 2019-05-05 2019-09-10 中国农业大学 Assimilate the Regional Fall Wheat yield estimation method of satellite fluorescence data and crop growth model

Non-Patent Citations (7)

* Cited by examiner, † Cited by third party
Title
CHEN YEPEI，ETC: "Cotton Growth Monitoring and Yield Estimation Based on Assimilation of Remote Sensing Data and Crop Growth Model", 《IEEE XPLORE》 *
JUANA PAUL MOIWO AND FULU TAO: "Contributions of precipitation, irrigation and soil water to evapotranspiration in (semi)-arid regions", 《ROYAL METEOROLOGICAL SOCIETY》 *
LI HE，ETC: "Improving Winter Wheat Yield Estimation from the CERES-Wheat Model to Assimilate Leaf Area Index with Different Assimilation Methods and Spatio-Temporal Scales", 《REMOTE SENSING》 *
TAO FULU，ETC: "Improving regional winter wheat yield estimation through assimilation of phenology and leaf area index from remote sensing data", 《EUROPEAN JOURNAL OF AGRONOMY》 *
李妍: "基于WOFOST-HYDRUS耦合模型的玉米遥感估产研究", 《CNKI中国博士学位论文全文数据库农业科技辑》 *
王品，等: "气候变暖背景下水稻低温冷害和高温热害的研究进展", 《资源科学》 *
黄健熙，等: "基于时间序列LAI 和ET 同化的冬小麦遥感估产方法比较", 《农业工程学报》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111667167A (en) * 2020-06-03 2020-09-15 福建慧政通信息科技有限公司 Agricultural grain yield estimation method and system
CN111667167B (en) * 2020-06-03 2022-12-06 福建慧政通信息科技有限公司 Agricultural grain yield estimation method and system

Also Published As

Publication number Publication date
CN110766308B (en) 2020-07-10

Similar Documents

Publication Publication Date Title
CN109800921B (en) Regional winter wheat yield estimation method based on remote sensing phenological assimilation and particle swarm optimization
CN110751094A (en) Crop yield estimation technology based on GEE comprehensive remote sensing image and deep learning method
CN111027752A (en) Crop yield estimation method based on deep spatiotemporal feature joint learning
CN109711102B (en) Method for rapidly evaluating crop disaster loss
CN106408132A (en) Method and device of crop yield prediction based on plantation device
CN110705182B (en) Crop breeding adaptive time prediction method coupling crop model and machine learning
CN106483147B (en) Long-time sequence passive microwave soil moisture precision improvement research method based on multi-source data
CN110633841B (en) Provincial range plot scale data assimilation yield prediction method based on set sampling
Terando et al. Observed and modeled twentieth-century spatial and temporal patterns of selected agro-climate indices in North America
CN111062526B (en) Winter wheat yield per unit prediction method and system
CN113011372B (en) Automatic monitoring and identifying method for saline-alkali soil
CN111160680A (en) Agricultural drought assessment method based on information assimilation and fusion
CN110766308B (en) Regional crop yield estimation method based on set assimilation strategy
CN109272144B (en) BPNN-based prediction method for NDVI (normalized difference of variance) in northern grassland area of China
CN111915096B (en) Crop yield early-stage forecasting technology based on crop model, remote sensing data and climate forecasting information
CN115759524B (en) Soil productivity grade identification method based on remote sensing image vegetation index
CN109359862A (en) A kind of real-time yield estimation method of cereal crops and system
CN109766871A (en) A kind of region Crop Estimation Method based on spatial diversity assimilation remotely-sensed data and crop modeling
JP6570929B2 (en) Characteristic estimation model generation apparatus and method, analysis target characteristic estimation apparatus and method
CN114091344A (en) Power transmission line risk assessment model training method and device based on data coupling
Higgins Limitations to seasonal weather prediction and crop forecasting due to nonlinearity and model inadequacy
CN113762768B (en) Agricultural drought dynamic risk assessment method based on natural gas generator and crop model
CN113378381B (en) Method for calculating winter wheat crop coefficient based on air temperature distribution and surface heat
CN109886497A (en) Surface air temperature interpolation method based on the improved inverse distance weight of latitude
CN116701371B (en) Method and device for interpolating missing values of atmospheric temperature data under covariance analysis

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