CN104134095B - Crop yield estimation method based on scale transformation and data assimilation - Google Patents
Crop yield estimation method based on scale transformation and data assimilation Download PDFInfo
- Publication number
- CN104134095B CN104134095B CN201410155913.0A CN201410155913A CN104134095B CN 104134095 B CN104134095 B CN 104134095B CN 201410155913 A CN201410155913 A CN 201410155913A CN 104134095 B CN104134095 B CN 104134095B
- Authority
- CN
- China
- Prior art keywords
- lai
- crop
- crops
- parameter
- remote sensing
- 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
Links
Landscapes
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
The invention belongs to the field of agricultural remote sensing, and relates to a crop yield estimation method based on scale transformation and data assimilation and application of the crop yield estimation method to the guiding of crop production. The method comprises the following concrete steps that: 1, parameters are collected for completing spatialization of a WOFOST crop model; 2, a crop type distribution map and a purity percentage map of the crops to be tested are obtained; 3, TMLAI in a range of 30m is obtained; 4, a secondary scale conversion model is built, and a time sequence scale regulating LAI is generated; 5, the crop model error and the remote sensing observation error in a key phonological period of the crops to be tested are obtained; 6, a four-dimensional variation cost function is built, and an optimized crop model parameter is obtained; and 7, the per unit yield of the crops to be tested in a county region is output. The method provided by the invention has the advantages that the problem of scale mismatching between a remote sensing observation picture element and a crop model simulation unit is solved; the precision of a data assimilation model is improved; and the method is suitable for the crop yield estimation in the county region scale, and is particularly suitable for the winter wheat yield estimation in the county region scale.
Description
Technical field
The invention belongs to agricultural remote sensing field and in particular to a kind of based on spatial scaling data assimilation crop yield
Estimating and measuring method.
Background technology
Using data assimilation remote-sensing inversion Parameter fusion to crop mechanism process model, it is that current improved confinement is made
Thing growth simulation precision and the important channel improving agricultural output assessment precision.On maturity in field, based on crop photosynthesis, breathing, steam
Rise, the crop growth model of the mechanism process such as nutrition rely in it physical process and dynamical mechanism, can be with accurate simulation
Crop object with " my god " growth and development state of crop and yield on continuous evolution as time step and single-point yardstick.And as
When thing model extension is applied to regional scale, due to earth's surface, near surface environment non-uniform, result in initial in crop modeling
The difficulty that condition, soil parameters, crop parameter, the uncertainty of spatial distribution of meteorological forcing factors and data obtain.Satellite
Remote sensing has that space is continuous and the advantage of time dynamic, can effectively solve the problem that area crops parameter acquiring this bottle difficult
Neck.But remote sensing earth observation is due to being restricted by factors such as satellite spatial and temporal resolutions, can't really disclose plant growth and send out
Educate and the intrinsic procedure mechanism of yield composition, individual growth development condition and its relation with Meteorological Elements in China, and this is exactly
The advantage of crop modeling is located.Data assimilation passes through to couple remote sensing observations and crop modeling, is capable of both advantages
Complementation, has huge potentiality to improving area crops Granule weight.
Can be flutterred due to high time resolution remotely-sensed data and catch crop growth change and phenology information, therefore to farming
Thing Growing state survey and the yield by estimation have very important significance, however, high time resolution remote sensing satellite number in orbit at present
All very low according to often spatial resolution, such as AVHRR, MODIS, MERIS and SPOT vegetation spatial resolution is
250m-1km yardstick.This low spatial resolution remotely-sensed data, further increases pixel heterogeneity, causes remote sensing observations chi
Mismatch between degree and crop modeling simulation yardstick, strongly limit the precision of data assimilation model.In addition, standard
MODIS LAI product is towards the exploitation of Global Scale all of vegetation pattern, is not to design and develop for crops.Research
Show, MODIS LAI exists in winter wheat region and serious underestimates phenomenon.Therefore, directly assimilation standard MODIS LAI can lead to
LAI after assimilation and yield are unpractical low.Based on geo-statistic（Variation function）Method and pixel analysis that yardstick is corrected
The method of NO emissions reduction is to process the heterogeneous two kinds of means with scale effect of pixel.However, the one of both approaches common
Weak point be a series of Crop growing stage of requirement high and medium between resolution data, and to heterogeneous pixel modeling be one
It is desirable to strict approximate and some prioris, this is difficult acquisition to individual extremely complex problem in actual applications.One kind has
The method of potentiality is based on ground surface sample LAI, middle high-resolution remotely-sensed data and low resolution remotely-sensed data, builds two grades of yardsticks
Conversion method, generates the revised LAI of yardstick.The assimilation revised LAI of yardstick, it will larger raising Regional Fall Wheat yield
Estimation.To overcome the error that direct assimilation MODIS LAI product brings to be urgently to solve in the yield by estimation practice of application data assimilation model
Key issue certainly.
Content of the invention
For solving following problem present in prior art：" how to build the scale-transformation method in data assimilation model,
The revised LAI of assimilation time sequence scale to crop growth model, to reduce between remote sensing observations and crop modeling yardstick not
Mate the data assimilation model error causing ", the present invention provides a kind of crop yield based on the assimilation of spatial scaling data
Estimating and measuring method.
The present invention provides a kind of crop yield estimating and measuring method based on the assimilation of spatial scaling data, and concrete steps are such as
Under：
Meteorologic parameter in S1 collection research area, crop parameter, soil parameters and management parameters simultaneously carry out WOFOST with it
The single-point yardstick of crop modeling is demarcated；Afterwards compartmentalization is carried out based on meteorologic parameter, complete the spatialization of WOFOST crop modeling；
The remote sensing image based on multidate for the S2, the agrotype sample point that combined ground investigation obtains, in conjunction with farming to be measured
The phenology feature of thing, builds classifying ruless, obtains agrotype scattergram and crops percent purity figure to be measured；
S3 is filtered to crops period of duration time serieses MODIS LAI curve to be measured, to eliminate shortage of data and cloud
The impact of pollution；Survey actual measurement LAI and the multi-temporal vegetation indices of sample point based on ground（TM VI）Build regression calculation mould
Type, obtains the TM LAI in 30 meters of region；
S4 merges ground actual measurement sample point LAI, multidate TM LAI and filtered MODIS LAI information, builds two grades
Spatial scaling model, generates time serieses rescaling LAI；
The LAI based on the ground sample point key phenological period for the S5, calculates remote sensing observations and the standard deviation of crop modeling simulation, obtains
Obtain remote sensing observations error and the crop modeling error in crops key phenological period to be measured；
Rescaling LAI and crop modeling simulation LAI that S6 is obtained with S4, and introduce remote sensing observations error and crop modeling
Error, builds four-dimensional variation cost function, minimizes cost function using optimization algorithm, by successive ignition, obtains and optimizes
Crop modeling parameter；
The optimization crop modeling parameter that S6 obtains is substituted into crop modeling by S7, selects crops purity to be measured to be more than 60%
Pixel, runs crop modeling simulation yield by pixel unit, is then aggregated into county domain administration cell, output county domain crops to be measured
Per unit area yield, instructs grain-production.
Wherein, meteorologic parameter described in S1 is that the highest temperature, the lowest temperature, total radiation, vapour pressure, wind speed, precipitation etc. are joined
Number.
Wherein, crop parameter described in S1 is crops distributed constant, phenology feature of crops etc..
Wherein, described in S1, compartmentalization is carried out based on meteorologic parameter, the spatialization completing crop modeling refers to reference to remote sensing shadow
As the parameter collected is carried out with the coupling of locus, to the accumulated temperature parameter needed for meteorologic parameter and crop modeling, adopt
Kriging interpolation algorithm, completes parameter regionization and demarcates.
Wherein, the remote sensing image of multidate described in S2 is the Landsat TM of multidate.
Wherein, obtain the decision-making that agrotype scattergram and crops percent purity figure to be measured are using C5.0 described in S2
Tree classification algorithm, the phenology feature in conjunction with crops and spectral signature, build classifying ruless, obtain agrotype scattergram, adopt
With 1km grid fit, generate the crops percent purity spatial distribution map to be measured of 1km, in order to reduce non-crops area to be measured
The impact to data assimilation result for the domain, the pixel selecting crops purity to be measured to be more than 60%, estimate as follow-up County Scale yield
The pixel surveyed.
Wherein, it is filtered described in S3, using coenvelope line（Savitzky-Golay, S-G）Filtering algorithm, its expression formula
See formula（1）：
In formula：Yj+iRepresent the value in one piece of window on original LAI curve, m is the radius of window, N is convolution number,
The width of window is 2m+1；Yj *Represent filtered LAI value；C represents the filter factor of i-th LAI value.
Wherein, it is filtered described in S3, on the basis of improved to S-G filtering, form coenvelope line filtering method, pass through
LAI value is divided into true value and falsity by LAI timing variations trend, with the mode of local circulation iteration, S-G filter value is replaced falsity,
The LAI smoothed curve new with true value synthesis, progressively matching forms the coenvelope line of LAI time series data.
Wherein, build two grades of spatial scaling models described in S4 to be made up of 2 committed steps：
1）The TM LAI of 3 key developmental stages is selected from 30 meters of crops TM LAI to be measured of S3 formation zone；
2）Using the ratio coefficient between the key developmental stages TM LAI and the S-G MODIS LAI in corresponding period of choosing,
Adjust other key developmental stages S-G MODIS LAI, in ascent stage and the decline stage of crops period of duration to be measured, profit respectively
It is fitted with Logistic curve,
Logistic curvilinear equation is shown in formula（2）：
Wherein, t is MODIS LAI seasonal effect in time series index, and y (t) is t time corresponding LAI value, a and b is matching ginseng
Number, c+d is maximum LAI value, and d is LAI initial value, i.e. first value in LAI time serieses.
Wherein, the remote sensing observations error in crops key phenological period to be measured described in S5 is based on actual measurement LAI and yardstick correction
Standard variance between LAI calculates；Crop modeling source of error described in S5 is in 2 model parameters optimizing, TDWI and SPAN.
Wherein, minimize cost function using optimization algorithm described in S6 and adopt SCE_UA algorithm to four-dimensional variation cost letter
Number is minimized, and four-dimensional variation cost function presses formula（3）Calculate：
Wherein, k represents the number of Optimized model parameter in cost function, xkThe WOFOST parameter of representing optimized is in interval model
The numerical value enclosing, xk0The empirical value of the WOFOST number of parameters of representing optimized, B represents model error, with 2 Optimized model parameters
Error represents, the transposition of T representing matrix, and N represents the number of times of time serieses remote sensing observations data, yiRepresent remote sensing observations LAI, Hi
X () represents WOFOST modeling LAI, Q0Represent the error of remote sensing LAI.
Wherein, S7 concretely comprises the following steps：Select the pixel more than 60% for the crops purity to be measured by continuous iteration, again first
2 crop modeling parameters of beginningization are so that the LAI of crop modeling output and yield constantly change, then in cost function
Contrast MODIS LAI and crop modeling export under the minimum difference of LAI, can terminate when the three below condition of convergence meets one
Assimilation.Obtain the numerical value of parameters optimization：
1. after continuous 5 circulations, parameter value to be optimized has been retracted to the codomain scope specified；
2. target function value cannot improve 0.0001% after 5 circulations；
3. the number of times of calculation cost function is more than 10000 times；
Model parameter after optimizing is substituted into crop modeling run, obtain yield, in conjunction with the crops distribution to be measured of pixel
Percentage ratio, is collected by county domain Administrative boundaries, the per unit area yield result on output area.
Wherein, described crops to be measured are preferably winter wheat.
The present invention is instructing farming based on the crop yield estimating and measuring method of spatial scaling data assimilation described in also providing
Application in thing production.
The present invention compared with prior art, has the beneficial effect that：
Instant invention overcomes the unmatched problem of yardstick between remote sensing observations pixel and crop modeling analogue unit, improve
The precision of data assimilation model, is suitable for the winter wheat yield monitoring of County Scale.
Specific embodiment
With reference to embodiment, the specific embodiment of the present invention is described in further detail.Following examples are used for
The present invention is described, but is not limited to the scope of the present invention.
Embodiment 1
Technical scheme is expanded on further taking the winter wheat yield monitoring of Baoding and Hengshui Prefecture as a example.This
The flow process of the Regional Fall Wheat Granule weight method based on the assimilation of spatial scaling data of embodiment includes：
Meteorologic parameter in S1 collection research area, crop parameter, soil parameters and management parameters simultaneously carry out WOFOST with it
The single-point yardstick of crop modeling is demarcated；Afterwards compartmentalization is carried out based on meteorologic parameter, complete the spatialization of crop modeling.Specifically such as
Under：
Baoding and Hengshui Prefecture is selected to be survey region, this area is the important winter wheat producing region in Hebei province.Obtain with
Lower data：6 gas needed for the highest/lowest temperature of 6 meteorological site of selection, global radiation, vapour pressure, wind speed, dewatering model
As key element；The soil parameters of agricultural weather experimental station collection and crop parameter etc. in research area；The Hebei province Fen Xian winter in 2009 is little
Wheat output statistics data；MODIS15A2 data product, the coordinate projection of unified all data is Albers Conical Equal
Area, georeferencing is WGS1984.
For insensitive crop parameter in WOFOST model, but or the less ginseng of the higher excursion of sensitivity
Number, by consulting literatures or is directly determined using the default value of WOFOST model（As developed the long DLO of the most suitable light, emerge
Low temperature lower limit TBASE, the maintenance of leaf breathes with RML etc.）；For the crop parameter that sensitivity is high, excursion is big, energy
Enough by consulting the pertinent literature in this research area or utilizing measured data to calculate acquisition（As accumulated temperature TSUM1, TSUM2, Ye Tong
Compound conversion ratio CVL, specific leaf area SLATB etc.）.The priori in area etc. can also be studied and determine the possible span of parameter
And numerical value（As the vital stage SPAN when 35 DEG C for the blade, maximum CO2 Assimilation rate AMAXTB, single blade assimilates the luminous energy profit of CO2
With rate EFFTB etc.）.
The observation data of 6 meteorological site in research on utilization area, in conjunction with gram in golden spatial interpolation methods to 2008～2009 years
Day by day the highest temperature, the lowest temperature, total radiation, vapour pressure, wind speed, six key elements of precipitation are entered row interpolation, and are converted
For the reference format of WOFOST model, obtain studying the corresponding meteorological data of each lattice point of area.For easy in spatial dimension
The parameter obtaining, the accumulated temperature in such as sowing time to seeding stage, the accumulated temperature in seeding stage to florescence, the accumulated temperature in florescence to period of maturation
Deng, subtracted each other using the daily mean temperature of meteorological observation station acquisition and the biological zero point of winter wheat and obtain daily accumulated temperature, logical
The daily accumulated temperature crossing cumulative difference phenological stage obtains three above accumulated temperature variable, equally adopts and meteorological data interpolation identical
Kriging regression method, obtains the whole accumulated temperature parameter studying spatial variations in area.
U.S.'s thematic mapper of S2 research on utilization area multidate（Thematic Mapper, TM）Image, combined ground is adjusted
The agrotype sample point discovered and seized, using the Decision Tree Algorithm of C5.0, in conjunction with the phenology feature of crops, builds classification
Rule, obtains agrotype scattergram, with 1km grid fit, generates the winter wheat percent purity spatial distribution map of 1km, is
Minimizing non-winter wheat data area assimilates the impact of result, selects the pixel that winter wheat purity is more than 60%, as follow-up county
The basis of domain yardstick Granule weight.
MODIS LAI data in S3 whole period of duration to test block synthesizes in temporal sequence, and each grid cell is given birth to
Become time-serial position.The MRT projection transform instrument providing especially by NASA, the data of Hebei Province is carried out inlaying, throws
Shadow conversion and form conversion, choose 1 to 177 days（Julian date）, totally 23 phase images be overlapped, generate be based on grid cell
Time-serial position.Time serieses MODIS LAI curve is filtered, to eliminate the shortage of data that cloud pollution causes.Specifically
Using coenvelope line（Savitzky-Golay, S-G）Filtering solves the discontinuous problem of leaf area index change procedure, its principle
It is expressed as follows：
In formula：Yj+iRepresent the value in one piece of window on original LAI curve, m is the radius of window, N is convolution number,
The width of window is 2m+1；Yj *Represent filtered LAI value；C represents the filter factor of i-th LAI value.
On the basis of improved to S-G filtering, form coenvelope line filtering method, by LAI timing variations trend by LAI
Value is divided into true value and falsity, with the mode of local circulation iteration, S-G filter value is replaced falsity, and the LAI new with true value synthesis puts down
Sliding curve, progressively matching forms the coenvelope line of LAI time series data.Comprise the following steps that：
（1）Linear interpolation is carried out to the LAI value point having cloud pollution.Because cloud covers or the impact of atmospheric effect, can make
Some pixel values carry larger noise, carry out spatial linear interpolation to these pixels using the credible pixel value of its periphery, generate
Initial LAI data.
（2）Carry out the extraction of LAI timing variations trend with S-G filtering.According to it is assumed that LAI time series data should follow vegetation
The progressive feature of phenology dynamic change, mutation in LAI sequential is considered as because the noise of cloud or air generation causes.Cause
This, smooth to each winter wheat pixel application S-G filtering in research area, obtain LAI seasonal effect in time series variation tendency data
LAItr.
（3）Determine weight Wi of each value point in LAI time serieses, provide basis for estimation for follow-up reconstruction result.Wi
Computing formula be：
Wherein,dmaxMaximum for di.
（4）Generate new LAI time serieses.From above theoretical, containing noisy LAI value than time trend
LAItr value on curve is little, these points is simply replaced and generates new LAI1 timing curve.
（5）LAI1 is carried out with S-G filtering, more filtered data is compared with the LAI value after space interpolation.With
Sample, if the LAI value after space interpolation is less than newly-generated LAI j, is replaced with LAI j, generates new sequence LAI j+1, iteration
It is filtered reconstructing, filter effect is by FkSize evaluating：
LAIik is that i-th pixel filters the LAI value after reconstructing in kth time, and Wi is step（3）Calculated weight, Fk
Judge parameter for the filtered iteration ends of kth time.The condition exiting iteration is Fk-1≥Fk≤Fk+1.
Obtain the cloudless TM data of test block 3 scape, phase is respectively：2009-3-14,2009-5-17,2009-6-10,
Build regression calculation model with corresponding actual measurement LAI respectively.
On March 14th, 2009：Because the middle ten days in March winter wheat is in the stage of turning green, the impact to vegetation index for the exposed soil is relatively
Greatly, the SAVI vegetation index is therefore selected setting up the statistical model between LAI and TM VI to be：
Coefficient of determination R2=0.849.
On May 17th, 2009, the statistical model set up between actual measurement LAI and TM NDVI is：
Coefficient of determination R2=0.742.
On June 10th, 2009, the statistical model set up between actual measurement LAI and TM NDVI is:
Coefficient of determination R2=0.874.
Based on 3 above-mentioned regression calculation models, obtain the spatial distribution map of the 30 meters of winter wheat in research area.
S4 merges ground actual measurement sample point LAI, multidate TM LAI and filtered MODISLAI information, builds two grades of chis
Degree transformation model, generates time serieses rescaling LAI.
Two grades of spatial scaling models include 2 steps, first, by S3 step generate 3 phenological periods actual measurement LAI and
TM VI regression calculation model, 30 meters of formation zone winter wheat TMLAI；Second step, using existing key developmental stages TM LAI with
Ratio coefficient between the S-G MODIS LAI in corresponding period, adjusts other key developmental stages S-G MODIS LAI, little in the winter
The ascent stage of wheat period of duration（Turn green to the jointing stage）, obtain the adjustment LAI in 4 phenological periods（In boot stage, the jointing stage, stand up the phase,
Period of seedling establishment）；Decline stage in winter wheat（Jointing is to the period of maturation）, obtain the adjustment LAI in 4 phenological periods（Boot stage, heading
Phase, florescence, the period of maturation）, it is utilized respectively Logistic curve and be fitted.
Logistic curvilinear equation is as follows：
Wherein, t is MODIS LAI seasonal effect in time series index, and y (t) is t time corresponding LAI value, a and b is matching ginseng
Number, c+d is maximum LAI value, and d is LAI initial value, i.e. first value in LAI time serieses.
The LAI based on the ground sample point key phenological period for the S5, calculates remote sensing observations and the standard deviation of crop modeling simulation, obtains
Obtain remote sensing observations error and the crop modeling error in winter wheat key phenological period.
Model error is considered from 2 parameters optimizing（TDWI and SPAN）, it is calculated analytically acquisition yardstick and repair
The observation error 7 phenological periods for the positive LAI is 0.34；Model error is set to 0.52.
Rescaling LAI and crop modeling simulation LAI that S6 is obtained with S4, and introduce remote sensing observations error and crop modeling
Error, builds four-dimensional variation cost function, minimizes cost function using optimization algorithm, by successive ignition, obtains and optimizes
Crop modeling parameter.
On the basis of S1 carries out crop modeling demarcation, in Experimental Area, run WOFOST crop modeling simulation LAI,
In conjunction with the yardstick correction LAI generating in S4, it is introduced into the remote sensing observations in S5 and modeling error, build four-dimensional variation cost letter
Number, as follows：
Wherein, k represents the number of Optimized model parameter in cost function, xkThe WOFOST parameter of representing optimized is in interval model
The numerical value enclosing, xk0The empirical value of the WOFOST number of parameters of representing optimized, B represents model error, with 2 Optimized model parameters
Error represents, the transposition of T representing matrix, and N represents the number of times of time serieses remote sensing observations data, yiRepresent remote sensing observations LAI, Hi
X () represents WOFOST modeling LAI, Q0Represent the error of remote sensing LAI.
Select the pixel that winter wheat purity is more than 60%, with optimized algorithm, continuous iteration reinitializes 2 WOFOST
Model parameter so that model output LAI and yield change, then in cost function contrast MODIS LAI and
WOFOST model exports under the minimum difference of LAI, obtains the most optimized parameter.SCE_UA algorithm is used for searching in initial parameter space
Rope globally optimal solution is so that cost function Fast Convergent.Can terminate to assimilate when the three below condition of convergence meets one, obtain
The numerical value of parameters optimization.
（1）After continuous 5 circulations, parameter value to be optimized has been retracted to the codomain scope specified；
（2）Target function value cannot improve 0.0001% after 5 circulations；
（3）The number of times of calculation cost function is more than 10000 times.
The optimization crop modeling parameter that S6 obtains is substituted into crop modeling by S7, selects the pixel that winter wheat purity is more than 60%,
Run crop modeling simulation yield by pixel unit, be then aggregated into county domain administration cell, output county domain winter wheat yield.Will be excellent
TDWI and SPAN parameter after change substitutes into WOFOST model running, obtains yield.In conjunction with the winter wheat percentage ratio of pixel, by county domain
Administrative boundaries collect, the winter wheat yield result on output area.
The present invention is filtered to During Growing Period of Winter Wheat time serieses MODISLAI curve by using S-G filtering algorithm,
Eliminate shortage of data and the impact of cloud pollution, generate smooth LAI time-serial position, be that two grades of spatial scaling models carry
Supply high-precision data source.Then time serieses rescaling LAI, effective integration are generated by two grades of spatial scaling models
MODIS LAI seasonal effect in time series tendency information and TM LAI information, improve the precision of LAI, are that assimilation process provides height
The data source of precision.Finally, the LAI after being adjusted by the assimilation and LAI of crop modeling simulation, is improve winter wheat on region and produces
The precision of amount estimation.
The Producer of grain and consumer are required for timely and accurately understanding grain yield information, according to this method, permissible
Obtain yield data in large area in the winter wheat period of maturation, be that national departments concerned carries out the science such as grain feelings judgement, grain regulation and control
Decision-making etc. provides important scientific basis, and can be used as the important evidence of grain trade.
The method of the present invention can be used for the estimation of other crop yields.
Although, above the present invention is described in detail with a general description of the specific embodiments,
On the basis of the present invention, it can be made some modifications or improvements, this will be apparent to those skilled in the art.Cause
This, these modifications or improvements without departing from theon the basis of the spirit of the present invention, belong to the scope of protection of present invention.
Claims (8)
1. a kind of crop yield estimating and measuring method based on the assimilation of spatial scaling data is it is characterised in that comprise the following steps that：
Meteorologic parameter in S1 collection research area, crop parameter, soil parameters and management parameters simultaneously carry out WOFOST crop with it
The single-point yardstick of model is demarcated；Afterwards compartmentalization is carried out based on meteorologic parameter, complete the spatialization of WOFOST crop modeling；
The remote sensing image based on multidate for the S2, the agrotype sample point that combined ground investigation obtains, in conjunction with crops to be measured
Phenology feature, builds classifying ruless, obtains agrotype scattergram and crops percent purity figure to be measured；
S3 is filtered to crops period of duration time serieses MODIS LAI curve to be measured, to eliminate shortage of data and cloud pollution
Impact；Regression calculation model is built based on the actual measurement LAI and multi-temporal vegetation indices TM VI that sample point is surveyed on ground, obtains
The TM LAI in 30 meters of region；
S4 merges ground actual measurement sample point LAI, TM LAI and filtered MODIS LAI information, builds two grades of spatial scaling moulds
Type, generates time serieses rescaling LAI；Two grades of spatial scaling models of described structure are made up of 2 committed steps：
1) select the TM LAI of 3 key developmental stages from 30 meters of crops TM LAI to be measured of S3 formation zone；
2) using the ratio coefficient between the key developmental stages TM LAI and the S-G MODIS LAI in corresponding period of choosing, adjustment
Other key developmental stages S-G MODIS LAI, in ascent stage and the decline stage of crops period of duration to be measured, is utilized respectively
Logistic curve is fitted, and Logistic curvilinear equation is shown in formula (2)：
Wherein, t is MODIS LAI seasonal effect in time series index, and y (t) is t time corresponding LAI value, a and b is fitting parameter, c+d
Maximum LAI value, d is LAI initial value, i.e. first value in LAI time serieses；
The LAI based on the ground sample point key phenological period for the S5, calculates remote sensing observations and the standard deviation of crop modeling simulation, and acquisition is treated
Survey remote sensing observations error and the crop modeling error in crops key phenological period；
S6 simulates, with crop modeling, the time serieses rescaling LAI that LAI and S4 obtains, and introduces remote sensing observations error and crop
Model error, builds four-dimensional variation cost function, minimizes cost function using optimization algorithm, by successive ignition, obtains
Optimize crop modeling parameter；
The optimization crop modeling parameter that S6 obtains is substituted into crop modeling by S7, selects the picture that crops purity to be measured is more than 60%
Unit, runs crop modeling simulation yield by pixel unit, is then aggregated into county domain administration cell, and output county domain crops to be measured are single
Produce, instruct grain-production.
2. the crop yield estimating and measuring method based on the assimilation of spatial scaling data described in claim 1 is it is characterised in that S1 institute
State and compartmentalization is carried out based on meteorologic parameter, the spatialization completing crop modeling refers to reference to remote sensing image, the parameter collected be carried out
The coupling of locus, to the accumulated temperature parameter needed for meteorologic parameter and crop modeling, using Kriging interpolation algorithm, completes to join
Number compartmentalization is demarcated.
3. the crop yield estimating and measuring method based on the assimilation of spatial scaling data described in claim 1 is it is characterised in that S2 institute
The remote sensing image stating multidate is the Landsat TM of multidate.
4. the crop yield estimating and measuring method based on the assimilation of spatial scaling data described in claim 1 is it is characterised in that S2 institute
State and obtain the Decision Tree Algorithm that agrotype scattergram and crops percent purity figure to be measured are using C5.0, in conjunction with agriculture
The phenology feature of crop and spectral signature, build classifying ruless, obtain agrotype scattergram, using 1km grid fit, generate
The crops percent purity spatial distribution map to be measured of 1km, in order to reduce non-crop area to be measured to data assimilation result
Impact, selects the pixel that crops purity to be measured is more than 60%, as the pixel of follow-up County Scale Granule weight.
5. the crop yield estimating and measuring method based on the assimilation of spatial scaling data described in claim 1 is it is characterised in that S3 institute
State and be filtered, using coenvelope line (Savitzky-Golay, S-G) filtering algorithm, its expression formula is shown in formula (1)：
In formula：Yj+iRepresent the value in one piece of window on original LAI curve, m is the radius of window, N is convolution number, window
Width be 2m+1；Yj *Represent filtered LAI value；C represents the filter factor of i-th LAI value.
6. the crop yield estimating and measuring method based on the assimilation of spatial scaling data described in claim 1 is it is characterised in that S5 institute
The remote sensing observations error stating the crops key phenological period to be measured is based on the standard variance surveyed between LAI and yardstick correction LAI
Calculate；Crop modeling source of error described in S5 is in 2 model parameters optimizing, TDWI and SPAN.
7. the crop yield estimating and measuring method based on the assimilation of spatial scaling data described in claim 1 is it is characterised in that S6 institute
State and minimize cost function employing using optimization algorithm
SCE_UA algorithm minimizes to four-dimensional variation cost function, and four-dimensional variation cost function is pressed formula (3) and calculated：
Wherein, k represents the number of Optimized model parameter in cost function, xkThe WOFOST parameter of representing optimized is in interval range
Numerical value, xk0The empirical value of the WOFOST number of parameters of representing optimized, B represents model error, with the error of 2 Optimized model parameters
Represent, the transposition of T representing matrix, N represents the number of times of time serieses remote sensing observations data, yiRepresent remote sensing observations LAI, Hi(x) table
Show WOFOST modeling LAI, Q0Represent the error of remote sensing LAI.
8. the crop yield estimating and measuring method based on the assimilation of spatial scaling data described in claim 1 is it is characterised in that S7 has
Body step is：Select the pixel that crops purity to be measured is more than 60% by continuous iteration, reinitialize 2 crop modelings ginsengs
Number, so that the LAI of crop modeling output and yield constantly change, then contrasts MODIS LAI and work in cost function
Thing model exports under the minimum difference of LAI, can terminate to assimilate when the three below condition of convergence meets one, obtain parameters optimization
Numerical value：
1. after continuous 5 circulations, parameter value to be optimized has been retracted to the codomain scope specified；
2. target function value cannot improve 0.0001% after 5 circulations；
3. the number of times of calculation cost function is more than 10000 times；
Model parameter after optimizing is substituted into crop modeling run, obtain yield, in conjunction with the crops distribution percentage to be measured of pixel
Ratio is collected by county domain Administrative boundaries, the per unit area yield result on output area.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201410155913.0A CN104134095B (en) | 2014-04-17 | 2014-04-17 | Crop yield estimation method based on scale transformation and data assimilation |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201410155913.0A CN104134095B (en) | 2014-04-17 | 2014-04-17 | Crop yield estimation method based on scale transformation and data assimilation |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN104134095A CN104134095A (en) | 2014-11-05 |
| CN104134095B true CN104134095B (en) | 2017-02-15 |
Family
ID=51806768
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201410155913.0A Active CN104134095B (en) | 2014-04-17 | 2014-04-17 | Crop yield estimation method based on scale transformation and data assimilation |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN104134095B (en) |
Families Citing this family (27)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106780323B (en) * | 2016-11-23 | 2020-03-17 | 湖北大学 | Agricultural condition acquisition and real-time updating method and system based on smart phone |
| CN107122739B (en) * | 2017-01-23 | 2020-07-17 | 东北农业大学 | Crop yield estimation model for reconstructing VI time series curve based on Extreme mathematical model |
| CN106951979A (en) * | 2017-02-20 | 2017-07-14 | 中国农业大学 | The crop maturity phase Forecasting Methodology that remote sensing, crop modeling are merged with weather forecast |
| CN107016466B (en) * | 2017-04-11 | 2020-12-22 | 北京林业大学 | Landscape pattern optimization construction method and system |
| CN107274036B (en) * | 2017-07-25 | 2020-04-03 | 中国农业科学院农业信息研究所 | Crop yield prediction method and system |
| CN108509836B (en) * | 2018-01-29 | 2021-10-08 | 中国农业大学 | Crop yield estimation method based on double-polarized synthetic aperture radar and crop model data assimilation |
| CN108304973A (en) * | 2018-02-11 | 2018-07-20 | 中国农业大学 | Area crops maturity period prediction technique based on accumulated temperature, radiation and soil moisture content |
| CN108133298B (en) * | 2018-03-08 | 2022-04-19 | 河南工业大学 | National grain consumption prediction method based on multiple regression model |
| CN108982369B (en) * | 2018-04-28 | 2020-09-01 | 中国农业大学 | Plot scale crop growth monitoring method integrating GF-1WFV and MODIS data |
| CN109033539A (en) * | 2018-07-02 | 2018-12-18 | 河北省科学院地理科学研究所 | The calculation method influenced based on plant growth mechanism model separation key factor pair crop phenology |
| CN108921351A (en) * | 2018-07-06 | 2018-11-30 | 北京兴农丰华科技有限公司 | Crop production forecast method based on trend yield and Meteorological Output |
| CN108984803B (en) * | 2018-10-22 | 2020-09-22 | 北京师范大学 | Method and system for spatializing crop yield |
| CN109800921B (en) * | 2019-01-30 | 2019-12-20 | 北京师范大学 | Regional winter wheat yield estimation method based on remote sensing phenological assimilation and particle swarm optimization |
| CN110008621B (en) * | 2019-04-15 | 2023-01-10 | 中国农业科学院农业资源与农业区划研究所 | Crop model remote sensing assimilation estimation method based on dual stream dependence set square root filtering assimilation algorithm |
| CN110378521B (en) * | 2019-06-28 | 2022-04-22 | 河南农业大学 | Construction and application of yield prediction model of winter wheat in northeast China of Henan province |
| CN110633841B (en) * | 2019-08-13 | 2022-04-01 | 中国农业大学 | Provincial range plot scale data assimilation yield prediction method based on set sampling |
| CN110766308B (en) * | 2019-10-17 | 2020-07-10 | 中国科学院地理科学与资源研究所 | Regional crop yield estimation method based on set assimilation strategy |
| CN110751094B (en) * | 2019-10-21 | 2020-08-28 | 北京师范大学 | Crop yield estimation method based on GEE comprehensive remote sensing image and deep learning method |
| CN113033262B (en) * | 2019-12-25 | 2022-08-05 | 中移(成都)信息通信科技有限公司 | Model training method and crop yield estimation method |
| US11935242B2 (en) | 2020-03-09 | 2024-03-19 | International Business Machines Corporation | Crop yield estimation |
| CN111915096B (en) * | 2020-08-14 | 2021-03-09 | 中国科学院地理科学与资源研究所 | Crop yield early-stage forecasting technology based on crop model, remote sensing data and climate forecasting information |
| CN112052988B (en) * | 2020-08-18 | 2024-03-22 | 中国农业大学 | Crop yield estimation method coupling multi-objective optimization and collection assimilation and application |
| CN114528282A (en) * | 2022-01-20 | 2022-05-24 | 重庆邮电大学 | Improved multi-cloud region MODIS NDVI time sequence reconstruction method based on space-time Kriging |
| CN114529097B (en) * | 2022-02-26 | 2023-01-24 | 黑龙江八一农垦大学 | Multi-scale crop phenological period remote sensing dimensionality reduction prediction method |
| CN115753625B (en) * | 2022-11-02 | 2024-04-19 | 中国农业大学 | Yield estimation method and device for regional crops, electronic equipment and storage medium |
| CN116432859B (en) * | 2023-05-08 | 2023-10-27 | 中山大学 | Crop yield statistical data downscaling method |
| CN117830860A (en) * | 2024-03-06 | 2024-04-05 | 江苏省基础地理信息中心 | Remote sensing automatic extraction method of winter wheat planting structure |
Citations (5)
| 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 |
| CN102323987A (en) * | 2011-09-09 | 2012-01-18 | 北京农业信息技术研究中心 | Crop leaf area index assimilation method |
| CN102651096A (en) * | 2012-04-28 | 2012-08-29 | 中国农业大学 | Method for estimating yield of winter wheat by assimilating characteristics of leaf area index time-sequence curve |
| JP2012203875A (en) * | 2011-03-28 | 2012-10-22 | Tokyo Electric Power Co Inc:The | Yield estimation device and computer program |
| CN103345707A (en) * | 2013-06-04 | 2013-10-09 | 中国科学院遥感与数字地球研究所 | Crop maturation stage remote sensing prediction method based on multi-source remote sensing data |
-
2014
- 2014-04-17 CN CN201410155913.0A patent/CN104134095B/en active Active
Patent Citations (5)
| 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 |
| JP2012203875A (en) * | 2011-03-28 | 2012-10-22 | Tokyo Electric Power Co Inc:The | Yield estimation device and computer program |
| CN102323987A (en) * | 2011-09-09 | 2012-01-18 | 北京农业信息技术研究中心 | Crop leaf area index assimilation method |
| CN102651096A (en) * | 2012-04-28 | 2012-08-29 | 中国农业大学 | Method for estimating yield of winter wheat by assimilating characteristics of leaf area index time-sequence curve |
| CN103345707A (en) * | 2013-06-04 | 2013-10-09 | 中国科学院遥感与数字地球研究所 | Crop maturation stage remote sensing prediction method based on multi-source remote sensing data |
Non-Patent Citations (3)
| Title |
|---|
| Estimating regional winter wheat yield by assimilation of time series of HJ-1 CCD NDVI into WOFOST-ACRM model with Ensemble Kalman Filter;Hongyuan Ma etal.;《Mathematical and Computer Modeling》;20130831;第58卷;第759-770页 * |
| 基于遥感信息与作物模型集合卡尔曼滤波同化地区域冬小麦产量预测;黄健熙 等;《农业工程学报》;20120229;第28卷(第4期);第142-148页 * |
| 条件植被温度指数的四维变分与集合卡尔曼同化方法;王维 等;《农业工程学报》;20111231;第27卷(第12期);第184-190页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN104134095A (en) | 2014-11-05 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN104134095B (en) | Crop yield estimation method based on scale transformation and data assimilation | |
| CN110751094B (en) | Crop yield estimation method based on GEE comprehensive remote sensing image and deep learning method | |
| Huang et al. | Assimilating a synthetic Kalman filter leaf area index series into the WOFOST model to improve regional winter wheat yield estimation | |
| CN102651096B (en) | Method for estimating yield of winter wheat by assimilating characteristics of leaf area index time-sequence curve | |
| CN111104858B (en) | Large-scale crop phenology extraction method based on morphological model method | |
| Ma et al. | Enhancing SWAT with remotely sensed LAI for improved modelling of ecohydrological process in subtropics | |
| Han et al. | Evaluating the impact of groundwater on cotton growth and root zone water balance using Hydrus-1D coupled with a crop growth model | |
| CN107423850B (en) | Regional corn maturity prediction method based on time series LAI curve integral area | |
| CN103345707A (en) | Crop maturation stage remote sensing prediction method based on multi-source remote sensing data | |
| CN105445229A (en) | Nitrogen balance spectroscopy-based wheat spring nitrogenous fertilizer application method, and construction method of nitrogen topdressing amount model thereof | |
| CN109933873A (en) | A kind of method of desertificated area Soil Organic Carbon Density high spatial remote sensing estimation | |
| CN106485002B (en) | In the method for complicated landform climatic province estimation sugarcane potential production | |
| CN103955860A (en) | Regional crop yield estimation method based on ensemble Kalman filter assimilation | |
| CN102346808B (en) | Method for inverting LAI (leaf area index) from HJ-1 satellite data | |
| CN113554232A (en) | Crop yield prediction method and system | |
| CN106446555A (en) | Vegetation change occurrence time detection method based on time series similarity | |
| Dan et al. | Delineating the rice crop activities in Northeast China through regional parametric synthesis using satellite remote sensing time-series data from 2000 to 2015 | |
| Wang et al. | Rice yield estimation based on an NPP model with a changing harvest index | |
| Liu et al. | Winter wheat yield estimation based on assimilated Sentinel-2 images with the CERES-Wheat model | |
| CN108205718A (en) | Production method and system are surveyed in a kind of cereal crops sampling | |
| CN114460015A (en) | Method and model for predicting chlorophyll content of canopy of rice in cold region | |
| WO2024055784A1 (en) | Recommendation model generation method and apparatus, and nitrogen fertilizer application amount recommendation method and apparatus | |
| CN117455062A (en) | Crop yield prediction algorithm based on multi-source heterogeneous agricultural data | |
| Chen et al. | Variation of gross primary productivity dominated by leaf area index in significantly greening area | |
| Hu et al. | Retrieval of photosynthetic capability for yield gap attribution in maize via model-data fusion |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| C14 | Grant of patent or utility model | ||
| GR01 | Patent grant |