CN110766308B - 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
CN110766308B
CN110766308B CN201910989918.6A CN201910989918A CN110766308B CN 110766308 B CN110766308 B CN 110766308B CN 201910989918 A CN201910989918 A CN 201910989918A CN 110766308 B CN110766308 B CN 110766308B
Authority
CN
China
Prior art keywords
assimilation
matrix
weight
year
normalized
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.)
Active
Application number
CN201910989918.6A
Other languages
Chinese (zh)
Other versions
CN110766308A (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—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

  • Business, Economics & Management (AREA)
  • Engineering & Computer Science (AREA)
  • Strategic Management (AREA)
  • Human Resources & Organizations (AREA)
  • Economics (AREA)
  • Theoretical Computer Science (AREA)
  • General Business, Economics & Management (AREA)
  • Entrepreneurship & Innovation (AREA)
  • General Physics & Mathematics (AREA)
  • Physics & Mathematics (AREA)
  • Tourism & Hospitality (AREA)
  • Marketing (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Primary Health Care (AREA)
  • Animal Husbandry (AREA)
  • General Health & Medical Sciences (AREA)
  • Mining & Mineral Resources (AREA)
  • Marine Sciences & Fisheries (AREA)
  • Agronomy & Crop Science (AREA)
  • Development Economics (AREA)
  • Game Theory and Decision Science (AREA)
  • Health & Medical Sciences (AREA)
  • Educational Administration (AREA)
  • Operations Research (AREA)
  • Quality & Reliability (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

The invention discloses a regional crop yield estimation method based on a set assimilation strategy, which comprises the steps of obtaining remote sensing L of predicted year Y and previous M years in a research regionAI data; 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‑HupAnd performing simulated estimation on the forecast year Y.

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 L AI data of a forecast year Y and a previous M year in the research area, wherein M is a natural number more than 3;
s2, carrying out spatial difference assimilation operation by using a crop model and the remote sensing L AI data according to the Y-1 year, wherein the spatial 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 L AI matrix and a normalized remote sensing inversion L AI matrix of an assimilation time point k, wherein H is an assimilation weight and has a value range of 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；
S6Repeating the steps S2-S5 for the 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.
According to an embodiment of the present invention, M in step S1 is 4, although other values may be selected as appropriate.
According to one embodiment of the invention, the crop model is a model of the MCW L a series, or alternatively other crop models.
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 L AI data may be selected from Copernicus L AI, G L ASS L AI and G L OBMAP L AI, or other remotely sensed L AI 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 simulation L AI matrix and the normalized remote sensing inversion L AI 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 BDA0002237926060000051
Figure BDA0002237926060000052
wherein,
Figure BDA0002237926060000053
representing the simulation L AI matrix at the point of assimilation k, DkRepresenting a corresponding remote sensing L AI matrix at the same time point k;
Figure BDA0002237926060000054
and
Figure BDA0002237926060000055
respectively represent by
Figure BDA0002237926060000056
And DkThe converted normalized simulation L AI matrix and normalized remote sensing L AI matrix 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
Reverse normalized transformation to obtain analytical L AI 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 shows the effect of the method according to the invention when different telemetric L AI products are used.
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, remote sensing L AI data of a forecast year Y and a previous year M in a research area are obtained, M is a natural number larger than 3, M is shown as 4, namely remote sensing L AI data of the Y-1, Y-2, Y-3 and Y-4 are obtained, wherein the remote sensing L AI data can be selected from Copernicus L AI, G L ASS L AI and G L OBMAP L AI, or other appropriate data.
And S2, carrying out differential assimilation by using a crop model (shown as MCW L A-Wheat model in the figure) and the remote sensing L AI data for the Y-1 year, wherein the differential assimilation comprises the following spatial differential assimilation operation equations:
Figure BDA0002237926060000071
wherein
Figure BDA0002237926060000072
Is a normalized analysis matrix of the assimilation time point k,
Figure BDA0002237926060000073
and
Figure BDA0002237926060000074
the method comprises the steps of respectively representing a normalized simulation L AI matrix and a normalized remote sensing inversion L AI matrix of an assimilation time point k, wherein H is an assimilation weight and ranges from 0.01 to 1.00, it should be understood that superscripts in the formula are used for distinguishing the matrixes, a is an analysis abbreviation and represents an analysis matrix, f is a forecast abbreviation and represents a forecast matrix, namely a model simulation matrix, and n is a normali abbreviationAnd (6) an abbreviation, indicating 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 L AI matrix and the corresponding remote sensing L AI matrix are normalized at a time point (called a assimilation time point) where remote sensing data exist, so that the original matrix is normalized into a normalized matrix (formulas (3) and (4) below) with a value range of (0, 1), and the normalization operation method is more, and for example, the normalization operation method can be obtained by dividing the L AI matrix by the maximum value of matrix elements:
Figure BDA0002237926060000081
Figure BDA0002237926060000082
wherein,
Figure BDA0002237926060000083
representing the simulation L AI matrix at the point of assimilation k, DkRepresenting a corresponding remote sensing L AI matrix at the same time point k;
Figure BDA0002237926060000084
and
Figure BDA0002237926060000085
respectively represent by
Figure BDA0002237926060000086
And DkThe converted normalized simulation L AI matrix and normalized remote sensing L AI matrix SfAnd SDIs a normalization function used to normalize the original matrix.
2) Assimilating the normalized simulated L AI matrix and the normalized remote sensing L AI matrix using assimilation equation (1):
Figure BDA0002237926060000087
wherein
Figure BDA0002237926060000088
The 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
Reverse normalized transformation to obtain analytical L AI matrix
Figure BDA00022379260600000811
And is composed of
Figure BDA00022379260600000812
The 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 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-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, expansion assimilation weight value-taking area based on S7Room Hlow-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 central part of North China plain was selected as a research area, the research time was 2008-Wheat, MCW L A-Wheat model was used, and the used remote sensing L AI data included Copernicius L0 AI (spatial resolution 1km × 1km, time frequency every 10 days), G L1 ASS L AI (spatial resolution 1km × 1km, time frequency every 8 days), and G L OBMAP L AI (spatial resolution 0.08 DEG × 0.08.08 DEG, time frequency every 8 days), each remote sensing L AI product was respectively assimilated with the MCW L A-Wheat model to estimate the yield of winter Wheat, the research area is shown in FIG. 2.
Based on the set assimilation strategy provided by the invention, the assimilation estimation implementation process is as follows (taking Copernicus L AI as an example, the implementation processes of G L ASS L AI and G L OBMAP L AI are consistent with Copernicus L AI):
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, the MCW L A-Wheat model and remote sensing L AI 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, so that 100 times of assimilation estimation operation can be driven, and 100 groups of estimation results are obtained, wherein the complete process of each assimilation estimation operation is as follows:
the method comprises the steps of simulating the growth of winter Wheat day by using MCW L A-Wheat from the sowing period of the winter Wheat, assimilating from the green turning period of the winter Wheat to the mature period of the winter Wheat, assimilating data of a model simulation L AI and remote sensing L AI by using a spatial difference assimilation algorithm at each time point (called assimilation time point) with remote sensing observation between the green turning period and the mature period, ending the assimilation and model simulation process at the mature period, and outputting the simulated yield of the winter Wheat grid by grid.
After one-time assimilation estimation operation is completed, the estimation results of the grid scale are aggregated to a county scale and compared with the statistical yield, and estimation accuracy is represented by a correlation coefficient r and a root mean square error rmse. 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.
For different remote sensing L AI data, the estimation mode is consistent with the process, and the difference is that when different remote sensing L AI data are used, the spatial resolution of the grid scale of the model operation needs to be consistent with the remote sensing L AI data, and the time point interval of assimilation needs to be consistent with the time frequency of the remote sensing L AI data.
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.
For example, for Copernicious L AI, the estimation accuracy can reach 98.1% of the estimation accuracy when the optimal assimilation weight is used, and the estimation accuracy respectively reaches 90.8% and 85.0% of the estimation accuracy when the optimal assimilation weight is used when G L ASS L AI and G L OBMAP L AI are used, so that relatively better accuracy can be obtained when the remote sensing product with higher resolution is used, and the estimation effect is closer to the estimation effect when the optimal assimilation weight is used.
In addition, the evaluation accuracy of the set assimilation strategy provided by the invention is stably superior to other evaluation methods in the assimilation evaluation performed by applying different remote sensing L AI data, and the method has universality.
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 L AI data of a forecast year y and the previous M years in the research area, wherein M is a natural number more than 3;
s2, performing differential assimilation aiming at the y-1 year, wherein the differential assimilation comprises the step of performing spatial differential assimilation operation by using a crop model and the remote sensing L AI data, and the spatial differential assimilation operation equation is as follows:
Figure FDA0002490997030000011
wherein
Figure FDA0002490997030000012
Is a normalized analysis matrix of the assimilation time point k,
Figure FDA0002490997030000013
and
Figure FDA0002490997030000014
respectively representing a normalized simulation L AI matrix and a normalized remote sensing inversion L AI matrix of an assimilation time point k, wherein H is an assimilation weight and has a value range of 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 FDA0002490997030000015
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 expanded assimilation weight value range and a lower limit of the expanded assimilation weight value rangeIs HlowThe upper limit is referred to as Hup，HlowObtained by a 10% expansion of the minimum, HupIs obtained by expanding the maximum value by 10 percent, and the value range of the expansion assimilation weight is not more than 0.99 and not less than 0.01, namely,
Hlow＝max(0.9*Hopmin,0.01)
Hup＝min(1.1*Hopmax,0.99)；
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 MCW L A 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 estimation of claim 4 based on collective assimilation strategies for regional crops, characterized by remote sensing L AI data selected from Copernicus L AI, G L ASS L AI and G L OBMAP L AI.
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 method based on collective assimilation strategy of claim 1, characterized in that in step S2, the normalized analysis matrix, the normalized simulation L AI matrix and the normalized remote sensing inversion L AI 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 FDA0002490997030000031
Figure FDA0002490997030000032
wherein,
Figure FDA0002490997030000033
representing the simulation L AI matrix at the point of assimilation k, DkRepresenting a corresponding remote sensing L AI matrix at the same time point k;
Figure FDA0002490997030000034
and
Figure FDA0002490997030000035
respectively represent by
Figure FDA0002490997030000036
And DkThe converted normalized simulation L AI matrix and normalized remote sensing L AI matrix 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 FDA0002490997030000037
by normalizing function SfIs inverse function of
Figure FDA0002490997030000038
Reverse normalized transformation to obtain analytical L AI matrix
Figure FDA0002490997030000039
And is composed of
Figure FDA00024909970300000310
The driving model runs to the next assimilation time point k + 1:
Figure FDA00024909970300000311
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 CN110766308A (en) 2020-02-07
CN110766308B true 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)

Families Citing this family (1)

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

Family Cites Families (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7702597B2 (en) * 2004-04-20 2010-04-20 George Mason Intellectual Properties, Inc. Crop yield prediction using piecewise linear regression with a break point and weather and agricultural parameters
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
CN104134095B (en) * 2014-04-17 2017-02-15 中国农业大学 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
CN109800921B (en) * 2019-01-30 2019-12-20 北京师范大学 Regional winter wheat yield estimation method based on remote sensing phenological assimilation and particle swarm optimization
CN110222870B (en) * 2019-05-05 2021-02-23 中国农业大学 Regional winter wheat yield estimation method by assimilating satellite fluorescence data and crop growth model

Also Published As

Publication number Publication date
CN110766308A (en) 2020-02-07

Similar Documents

Publication Publication Date Title
CN111027752B (en) Crop yield estimation method based on deep spatiotemporal feature joint learning
CN102651096B (en) Method for estimating yield of winter wheat by assimilating characteristics of leaf area index time-sequence curve
CN104134095B (en) Crop yield estimation method based on scale transformation and data assimilation
CN109800921B (en) Regional winter wheat yield estimation method based on remote sensing phenological assimilation and particle swarm optimization
CN109711102B (en) Method for rapidly evaluating crop disaster loss
CN110751094A (en) Crop yield estimation technology based on GEE comprehensive remote sensing image and deep learning method
CN110705182B (en) Crop breeding adaptive time prediction method coupling crop model and machine learning
CN111582551A (en) Method and system for predicting short-term wind speed of wind power plant and electronic equipment
CN113011372B (en) Automatic monitoring and identifying method for saline-alkali soil
CN110633841B (en) Provincial range plot scale data assimilation yield prediction method based on set sampling
CN111062526B (en) Winter wheat yield per unit prediction method and system
CN110766308B (en) Regional crop yield estimation method based on set assimilation strategy
CN109766871A (en) A kind of region Crop Estimation Method based on spatial diversity assimilation remotely-sensed data and crop modeling
CN115952421A (en) High-precision time-space simulation method for coupling ecological process model and machine learning algorithm
CN111915096B (en) Crop yield early-stage forecasting technology based on crop model, remote sensing data and climate forecasting information
CN117520783A (en) Crop yield prediction method
CN115907178B (en) Clean ecosystem CO 2 Exchange amount prediction method
CN113762768B (en) Agricultural drought dynamic risk assessment method based on natural gas generator and crop model
CN109359862A (en) A kind of real-time yield estimation method of cereal crops and system
CN114091344A (en) Power transmission line risk assessment model training method and device based on data coupling
JP6570929B2 (en) Characteristic estimation model generation apparatus and method, analysis target characteristic estimation apparatus and method
CN114692971A (en) Crop yield prediction method and device based on yield difference
CN112131954A (en) China northern winter wheat remote sensing estimated yield based on convolutional neural network
CN116701371B (en) Method and device for interpolating missing values of atmospheric temperature data under covariance analysis
CN111667167B (en) Agricultural grain yield estimation method and system

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