CN113537108B - Crop chlorophyll content inversion method based on unmanned aerial vehicle hyperspectral remote sensing - Google Patents
Crop chlorophyll content inversion method based on unmanned aerial vehicle hyperspectral remote sensing Download PDFInfo
- Publication number
- CN113537108B CN113537108B CN202110837401.2A CN202110837401A CN113537108B CN 113537108 B CN113537108 B CN 113537108B CN 202110837401 A CN202110837401 A CN 202110837401A CN 113537108 B CN113537108 B CN 113537108B
- Authority
- CN
- China
- Prior art keywords
- cluster
- vector
- field area
- field
- reflectivity
- 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
- 229930002875 chlorophyll Natural products 0.000 title claims abstract description 68
- 235000019804 chlorophyll Nutrition 0.000 title claims abstract description 68
- ATNHDLDRLWWWCB-AENOIHSZSA-M chlorophyll a Chemical compound C1([C@@H](C(=O)OC)C(=O)C2=C3C)=C2N2C3=CC(C(CC)=C3C)=[N+]4C3=CC3=C(C=C)C(C)=C5N3[Mg-2]42[N+]2=C1[C@@H](CCC(=O)OC\C=C(/C)CCC[C@H](C)CCC[C@H](C)CCCC(C)C)[C@H](C)C2=C5 ATNHDLDRLWWWCB-AENOIHSZSA-M 0.000 title claims abstract description 68
- 238000000034 method Methods 0.000 title claims abstract description 32
- 230000002776 aggregation Effects 0.000 claims abstract description 19
- 238000004220 aggregation Methods 0.000 claims abstract description 19
- 239000011159 matrix material Substances 0.000 claims abstract description 19
- 238000005259 measurement Methods 0.000 claims abstract description 19
- 230000003595 spectral effect Effects 0.000 claims abstract description 19
- 238000012216 screening Methods 0.000 claims abstract description 16
- 238000012549 training Methods 0.000 claims abstract description 9
- 238000007637 random forest analysis Methods 0.000 claims abstract description 8
- 239000013598 vector Substances 0.000 claims description 136
- 238000002310 reflectometry Methods 0.000 claims description 71
- 238000012360 testing method Methods 0.000 claims description 25
- 238000003066 decision tree Methods 0.000 claims description 23
- 238000004364 calculation method Methods 0.000 claims description 12
- 238000004422 calculation algorithm Methods 0.000 claims description 7
- 239000000284 extract Substances 0.000 claims description 2
- 238000012544 monitoring process Methods 0.000 abstract description 14
- 238000001228 spectrum Methods 0.000 abstract description 5
- 230000009286 beneficial effect Effects 0.000 abstract 1
- 241000209140 Triticum Species 0.000 description 8
- 235000021307 Triticum Nutrition 0.000 description 8
- 238000005070 sampling Methods 0.000 description 5
- IJGRMHOSHXDMSA-UHFFFAOYSA-N Atomic nitrogen Chemical compound N#N IJGRMHOSHXDMSA-UHFFFAOYSA-N 0.000 description 2
- 230000001066 destructive effect Effects 0.000 description 2
- 230000007774 longterm Effects 0.000 description 2
- 238000007781 pre-processing Methods 0.000 description 2
- 239000002689 soil Substances 0.000 description 2
- 238000012795 verification Methods 0.000 description 2
- 240000002768 Alpinia galanga Species 0.000 description 1
- 235000006887 Alpinia galanga Nutrition 0.000 description 1
- 241000196324 Embryophyta Species 0.000 description 1
- 230000009418 agronomic effect Effects 0.000 description 1
- 238000012937 correction Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 230000002349 favourable effect Effects 0.000 description 1
- 230000035558 fertility Effects 0.000 description 1
- 238000010801 machine learning Methods 0.000 description 1
- 230000004060 metabolic process Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 229910052757 nitrogen Inorganic materials 0.000 description 1
- 230000029553 photosynthesis Effects 0.000 description 1
- 238000010672 photosynthesis Methods 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
- G06F18/214—Generating training patterns; Bootstrap methods, e.g. bagging or boosting
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/23—Clustering techniques
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/24—Classification techniques
- G06F18/243—Classification techniques relating to the number of classes
- G06F18/24323—Tree-organised classifiers
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A40/00—Adaptation technologies in agriculture, forestry, livestock or agroalimentary production
- Y02A40/10—Adaptation technologies in agriculture, forestry, livestock or agroalimentary production in agriculture
Landscapes
- Engineering & Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Theoretical Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Artificial Intelligence (AREA)
- Evolutionary Biology (AREA)
- Evolutionary Computation (AREA)
- Physics & Mathematics (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Investigating Or Analysing Materials By Optical Means (AREA)
Abstract
The invention discloses a crop chlorophyll content inversion method based on unmanned aerial vehicle hyperspectral remote sensing, which divides a field into a plurality of field areas to obtain chlorophyll content measurement values of each field area; acquiring spectrum information of each field area by utilizing a hyperspectral remote sensing image of the unmanned aerial vehicle; clustering the field areas based on the spectral information of the field areas to obtain various aggregation sets; screening characteristic wave bands of each cluster set to construct a characteristic matrix of each cluster set; and respectively corresponding the feature matrix of each aggregation set with the chlorophyll content measurement value of the corresponding field area, and training to construct a random forest regression model, namely an inversion model. The inversion method is beneficial to monitoring the chlorophyll content of crops in a large area, does not damage the crops, provides a novel and reliable method for monitoring the chlorophyll content of the leaves of the crops, and has higher practical value.
Description
Technical Field
The invention relates to the technical field of processing and application of remote sensing information, in particular to a crop chlorophyll content inversion method based on hyperspectral remote sensing of an unmanned aerial vehicle.
Background
In traditional agriculture, crop growth parameter monitoring is mainly achieved through manual destructive sampling and manual measurement to obtain related agronomic indexes. The traditional method is difficult to apply to large-area crop growth monitoring, and the traditional method has large workload. However, during the growth of crops and the photosynthesis and physiological metabolism of crops, the growth vigor of crops comprises physiological and biochemical parameters, such as leaf chlorophyll content, leaf area index, plant nitrogen content and the like, which are closely related to the spectral characteristics of canopy. Therefore, the remote sensing means can be utilized to acquire the spectrum information of crops in real time, and the physiological and biochemical parameters and the like of the crops can be indirectly deduced. Satellite remote sensing can be monitored in a large range and has low cost, but is easily affected by cloud layers, has long operation period and lacks the space-time resolution required by accurate agriculture. In recent years, unmanned aerial vehicles are widely applied to various fields of agricultural monitoring at home and abroad, and the unmanned aerial vehicle is flexible in operation, high in data acquisition speed and high in space-time resolution, and can be applied to complex farmland environments.
Disclosure of Invention
In order to overcome the defects in the prior art, the invention provides the crop chlorophyll content inversion method based on the hyperspectral remote sensing of the unmanned aerial vehicle, which is favorable for monitoring the chlorophyll content of crops in a large area, does not damage the crops, provides a novel and reliable method for monitoring the chlorophyll content of the blades of the crops, and has higher practical value.
In order to achieve the above purpose, the present invention adopts the following technical scheme, including:
a crop chlorophyll content inversion method based on unmanned aerial vehicle hyperspectral remote sensing comprises the following steps:
s1, measuring chlorophyll content of leaves in a field, dividing the field into M field areas, and obtaining a chlorophyll content measurement value of each field area;
constructing a tag set Y, Y= { Y according to chlorophyll content measurement values of each field area m }，m＝1,2,…,M，y m Chlorophyll content measurement for the mth field;
s2, acquiring hyperspectral remote sensing images of fields by using an unmanned aerial vehicle, and obtaining the average value of the reflectivity of each field area on each wave band on N wave bands;
according to the average value of the reflectivity of each field area on each wave band, a sample set X is constructed, for the m field area at the n-th wave band lambda n The average value of the reflectivity;
s3, clustering vectors of field areas in the sample set X to obtain J clustering sets;
each cluster set comprises vectors of a plurality of field areas; the vector formed by the average value of the reflectivity of the field area on each wave band is the vector of the field area;
wherein the j-th cluster set isI.e. < ->
j＝1,2,…,J；
K j For the j-th cluster set C j The total number of vectors for the field area;
for the j-th cluster set C j Vectors of the kth field region, i.e. the jth cluster center z j Vector of the mean value of the reflectivity of the kth field over the individual wavelength bands,/->I.e. < -> For the j-th cluster set C j The average value of the reflectivity of the kth field area on the nth wave band;
s4, screening characteristic wave bands from each cluster set respectively, and screening D characteristic wave bands from each cluster set respectively; respectively constructing a feature matrix of each cluster set according to the average value of the reflectivity of each field area in each cluster set on each feature wave band;
wherein the j-th cluster set C j Is given by (a)
For the j-th cluster set C j The kth field region in (a) is in the (d) th characteristic band +.>The average value of the reflectivity;
s5, respectively corresponding the feature matrix of each aggregation set to the chlorophyll content measurement value of the corresponding field area in the tag set Y, and training to obtain an inversion model;
the input of the inversion model is a vector formed by the average value of the reflectivity of a certain field area on each wave band, namely the input is the vector of the certain field area, and the output of the inversion model is the inversion value of the chlorophyll content of the field area.
In step S2, the mth field area is at the nth band lambda n On the average of the reflectivityThe calculation mode of (2) is as follows:
wherein,for the nth band lambda n The reflectivity value of the a pixel point in the m field area; n=1, 2, …, N; m=1, 2, …, M; a=1, 2, …, a m ，A m Is the total number of pixel points in the mth field area.
In step S3, the vectors of the field areas in the sample set X are clustered, specifically as follows:
s31, randomly selecting vectors of J field areas in a sample set X, wherein J is less than M; taking the vectors of the J field areas as initial values of clustering centers in a clustering algorithm to form a clustering center set Z; taking the spectral angular distance as a clustering distance parameter in a clustering algorithm;
s32, respectively calculating the spectral angular distance between the vector of each field area in the sample set X and each cluster center in the cluster center set Z, respectively clustering the vector of each field area in the sample set X onto the cluster center with the smallest spectral angular distance to obtain the vector of the field area clustered by each cluster center, wherein the vector of the field area clustered by each cluster center forms a cluster set corresponding to the cluster center, and J cluster sets are obtained;
vector X of the mth field in sample set X m With the j-th cluster center Z in the cluster center set Z j Spectral angular distance between SAD (x) m ,z j ) The calculation method is as follows:
wherein SAD (x) m ,z j ) Vector x for the mth field in the sample set m With the j-th cluster center Z in the cluster center set Z j Spectral angular distance between;
x m the vector of the m field area in the sample set X is the vector formed by the average value of the reflectivity of the m field area in the sample set X on each wave band,i.e. < -> For the mth field region in sample set X at the nth band lambda n The average value of the reflectivity; m=1, 2, …, M;
z j is the j-th cluster center Z in the cluster center set Z j Namely, gatherVector formed by average value of reflectivity of jth field area in class center set Z on each wave bandI.e. < -> For the jth field area in the cluster center set Z at the nth wave band lambda n The average value of the reflectivity; j=1, 2, …, J;
s33, respectively calculating average vectors of all the aggregation sets;
wherein the j-th cluster set C j Is the average vector delta of (2) j The calculation mode of (2) is as follows:
s34, updating the clustering centers of the clustering sets according to the average vectors of the clustering sets respectively to obtain an optimal clustering center, and obtaining an optimal clustering set corresponding to the optimal clustering center;
wherein, the j-th cluster set C is calculated j Is of the center vector delta of (2) j With the jth cluster center z j If the difference is smaller than or equal to the set threshold T, the j-th cluster center z j Without updating, the jth cluster center z j I.e. the optimal cluster center, the jth cluster set C j The optimal aggregate set is obtained; if the difference is greater than the set threshold value T, the jth clustering center z is subjected to j Updating the j-th cluster center z j Updated to the value of the j-th cluster set C j Is of the center vector delta of (2) j Re-clustering, i.e. re-executing steps S32-S34 until the optimal j-th cluster center z is obtained j The optimal jth cluster center z j Corresponding collection C j And the optimal aggregate set is obtained.
In step S4, feature band screening is performed on the jth cluster set, and the specific manner is as follows:
s41, initializing d=1, randomly selecting a band λ in the j-th cluster set ω The wave band lambda ω Vector of (3)The start vector h as the d-th projection d The method comprises the steps of carrying out a first treatment on the surface of the Band lambda in the j-th cluster set ω The vector formed by the average value of the reflectivity of each field area is the wave band lambda ω Vector of-> I.e. < -> For the j-th cluster set C j The kth field area in the middle is in the wave band lambda ω The average value of the reflectivity;
s42, performing the d-th projection on the j-th cluster set, and respectively calculating the initial vector h of the d-th projection d Projection vectors of each wave band in the j-th cluster set are projected;
wherein the start vector h of the d-th projection d For the nth band lambda in the jth cluster set n Vector of (3)Projection vector of +.>The calculation method is as follows:
wherein n=1, 2, …, N and expressed +.>Is a transpose of (2); h is a d T H is expressed as d Is a transpose of (2);
the start vector h for the d-th projection d For the j-th cluster set C j N-th band lambda of n Vector of->Is a projection vector of (a);
the j-th cluster set C j N-th band lambda of n Vector of (3)The method comprises the following steps: the j-th cluster set C j N-th band lambda of n Vector of the average value of the reflectivity of each field area +.>I.e. < -> For the j-th cluster set C j The average value of the reflectivity of the kth field area on the nth wave band; n=1, 2, …, N;
s43, recording the wave band corresponding to the largest projection vector, wherein the wave band is the characteristic wave band of the (d) projection of the (j) th aggregation set
S44, carrying out next projection, namely, d+1st projection, on the jth cluster set, and carrying out characteristic wave band of the d projectionIs the reflectance mean vector x of (2) d As the start vector h of the (d+1) th projection d+1 According to the steps S42-S43, d+1th projection is carried out to obtain characteristic wave band of d+1th projection +.>
S45, if d+1=D, the next projection is not performed, and the screening of the D characteristic wave bands of the j-th aggregation set is completed; if d+1 is less than D, continuing to project the j-th cluster set for the next time in a mode of step S44;
s46, from the j-th cluster set C j After screening D characteristic wave bands, according to the j-th aggregation set C j The average value of the reflectivity of each field area on each characteristic wave band is used for constructing a j-th cluster set C j Feature matrix of (a)
In step S5, the training method of the inversion model is as follows:
s51, the feature matrix of each cluster is equivalent to a decision tree, the feature wave bands of the cluster are the feature wave bands of the corresponding decision tree, and the feature matrix of each cluster is respectively corresponding to the tag set Y;
s52, training and generating each decision tree based on the feature matrix of each aggregation set and the chlorophyll content measurement value of the corresponding field area in the tag set Y, integrating a plurality of decision trees, and establishing a random forest regression model, namely an inversion model.
When the inversion model is utilized for prediction, a test sample is input, the test sample is a vector formed by the average value of the reflectivities of a certain field area on each wave band, each decision tree in the random forest regression model extracts the average value of the reflectivities of the corresponding characteristic wave bands in the test sample, each decision tree predicts according to the average value of the reflectivities of the corresponding characteristic wave bands, each decision tree outputs a chlorophyll content predicted value of the field area, an average value is calculated for the chlorophyll content predicted value of the field area output by each decision tree, and the average value is used as a final inversion value of the chlorophyll content of the field area.
N＝176。
The invention has the advantages that:
(1) Traditional crop chlorophyll content monitoring is achieved through destructive sampling, and monitoring of a single sampling point can be achieved only, so that the monitoring method is difficult to apply to large-area crop monitoring. According to the invention, the unmanned aerial vehicle hyperspectral remote sensing technology is introduced into the crop chlorophyll content monitoring, so that the breakthrough from point to surface, namely from a single sampling point to a field area is realized, the monitoring of the chlorophyll content of the crops in a large area is facilitated, and the damage to the crops is avoided.
(2) According to the invention, the spectral angular distance is utilized to perform spectral clustering on the field area, so that the characteristics contained in the spectrum and the difference among the spectra can be shown to the greatest extent after clustering, and the method is different from the conventional integral modeling. The invention builds the inversion model by clustering and modeling, and can lead the accuracy of the random forest regression model, namely the inversion model, to be higher based on the idea of differential modeling.
(3) According to the invention, characteristic wave bands are screened, modeling dimension is reduced, redundancy and complexity of an inversion model are reduced, and an overfitting phenomenon is avoided.
Drawings
FIG. 1 is a flow chart of the present invention.
Fig. 2 is an acquired hyperspectral remote sensing image of an unmanned aerial vehicle of a test field.
Fig. 3 is a schematic diagram of inversion results of chlorophyll content of a portion of a field area in a test field.
FIG. 4 is a verification accuracy graph of an inversion model.
Detailed Description
The following description of the embodiments of the present invention will be made clearly and completely with reference to the accompanying drawings, in which it is apparent that the embodiments described are only some embodiments of the present invention, but not all embodiments. All other embodiments, which can be made by those skilled in the art based on the embodiments of the invention without making any inventive effort, are intended to be within the scope of the invention.
As shown in fig. 1, the crop chlorophyll content inversion method based on unmanned aerial vehicle hyperspectral remote sensing provided by the invention comprises the following steps:
s1, measuring chlorophyll content of the blades in a test field, dividing the field into M field areas, and obtaining a chlorophyll content measurement value of each field area. Constructing a tag set Y, Y= { Y according to chlorophyll content measurement values of each field area m }，m＝1,2,…,M，y m Chlorophyll content measurement for the mth field.
In step S1, the test field is divided into M field areas, in this embodiment, m=126; and (3) measuring the chlorophyll content of the leaf blade of each field area by using a daily SPAD-502 type handheld chlorophyll meter, wherein the chlorophyll content of the leaf blade of each sampling point is screened out from each field area to be used as a chlorophyll content measurement value of the field area.
S2, acquiring hyperspectral remote sensing images of the test field by using the unmanned aerial vehicle, and obtaining the average value of the reflectivity of each field area on each wave band on N wave bands. According to the average value of the reflectivity of each field area on each wave band, a sample set X is constructed, for the m field area at the n-th wave band lambda n And the average value of the reflectivity.
In the step S2, a Gaiasky-mini2 hyperspectral imager is carried by using a Dajiang M600PRO six-wing unmanned aerial vehicle, and hyperspectral remote sensing images of the wheat in the jointing period of a Mongolian horse store test station long-term positioning test field (116 DEG 37'E,33 DEG 13' N) are obtained. The unmanned aerial vehicle acquires a plurality of hyperspectral remote sensing images of the test field, and the hyperspectral remote sensing images of the test field are spliced by using HiSpectralStitcher software to obtain the complete hyperspectral remote sensing images of the test field. Pre-processing hyperspectral remote sensing images of a test field using specv iew software, the pre-processing comprising: and (5) carrying out radiometric calibration and atmospheric correction, and obtaining reflectivity values of all pixel points in the hyperspectral remote sensing image after pretreatment.
In step S2, the mth field area is at the nth band lambda n On the average of the reflectivityThe calculation mode of (2) is as follows:
wherein,for the nth band lambda n The reflectivity value of the a pixel point in the m field area; n=1, 2, …, N; m=1, 2, …, M; a=1, 2, …, a m ，A m Is the total number of pixel points in the mth field area.
In step S2, an ROI tool in ENVI software may be further used to extract a region of interest from each field region, and on N bands, reflectance values of each pixel point in the region of interest of each field region on each band are respectively derived, and a reflectance average value of the region of interest of each field region on each band is calculated, where the reflectance average value of the region of interest of each field region on each band is used as a reflectance average value of the field region on each band. That is, the mth field region is at the nth band lambda n On the average of the reflectivityThe method comprises the following steps:
wherein,for the nth band lambda n The reflectivity value of the a' th pixel point in the interested area of the m field area; n=1, 2, …, N; m=1, 2, …, M; a ' =1, 2, …, a ' ' m ，A′ m The total number of pixel points in the interested area of the mth field area.
And S3, clustering vectors of field areas in the sample set X to obtain J clustering sets. Each cluster set comprises vectors of a plurality of field areas; the vector formed by the average value of the reflectivity of the field area on each wave band is the vector of the field area.
In step S3, the vectors of the field areas in the sample set X are clustered, specifically as follows:
s31, randomly selecting vectors of J field areas in a sample set X, wherein J is less than M, and taking the vectors of the J field areas as initial values of clustering centers in a clustering algorithm to form a clustering center set Z; and taking the spectral angular distance as a clustering distance parameter in a clustering algorithm.
S32, respectively calculating the spectral angular distances between the vectors of all the field areas in the sample set X and all the clustering centers in the clustering center set Z, respectively clustering the vectors of all the field areas in the sample set X onto the clustering center with the minimum spectral angular distance to obtain the vectors of the field areas clustered by all the clustering centers, and forming the clustering set corresponding to the clustering center by the vectors of the field areas clustered by all the clustering centers.
Wherein, vector X of the mth field area in sample set X m With the j-th cluster center Z in the cluster center set Z j Spectral angular distance between SAD (x) m ,z j ) The calculation method is as follows:
wherein SAD (x) m ,z j ) Vector x for the mth field in the sample set m With the j-th cluster center Z in the cluster center set Z j The spectral angular distance between them.
x m The vector of the m field area in the sample set X is the vector formed by the average value of the reflectivity of the m field area in the sample set X on each wave band,i.e. < -> For the mth field region in sample set X at the nth band lambda n The average value of the reflectivity; m=1, 2, …, M.
z j Is the j-th cluster center Z in the cluster center set Z j Namely, the vector formed by the average value of the reflectivity of the jth field area in the cluster center set Z on each wave band,i.e. < -> For the jth field area in the cluster center set Z at the nth wave band lambda n The average value of the reflectivity; j=1, 2, …, J.
The j-th cluster center Z in the cluster center set Z j The corresponding cluster set is the j-th cluster set C j ，I.e. < -> For the j-th cluster center z j Vectors of the kth field region clustered, i.e. the jth cluster set C j Vectors of the kth field region, i.e. the jth cluster center z j Vector of the average value of the reflectivity of the k-th field region of the cluster over each band, +.>I.e.K j For the j-th cluster set C j Total number of vectors of field area, +.>For the j-th cluster set C j The average value of the reflectivity of the kth field area on the nth band.
When the vector of a certain field area is equal to the clustering center, the spectrum angle distance between the vector and the clustering center is the smallest, so that the clustering center comprises the vector of the certain field area equal to the clustering center, namely the vector is equivalent to the vector comprising the clustering center.
S33, respectively calculating average vectors of the aggregation sets.
Wherein the j-th cluster set C j Is the average vector delta of (2) j The calculation mode of (2) is as follows:
s34, updating the clustering centers of the clustering sets according to the average vectors of the clustering sets respectively to obtain the optimal clustering centers and the clustering sets corresponding to the optimal clustering centers, namely the optimal clustering sets.
Wherein, the j-th cluster set C is calculated j Is directed in the center of (C)Quantity delta j With the jth cluster center z j If the difference is smaller than or equal to the set threshold T, the j-th cluster center z j Without updating, the jth cluster center z j Namely the optimal cluster center, the j-th cluster set C j Namely, the optimal aggregate set; if the difference is greater than the set threshold value T, the jth clustering center z is subjected to j Updating the j-th cluster center z j Updated to the value of the j-th cluster set C j Is of the center vector delta of (2) j Re-clustering, i.e. re-executing steps S32-S34 until the optimal j-th cluster center z is obtained j The optimal jth cluster center z j Corresponding collection C j And the optimal aggregate set is obtained.
S4, screening characteristic wave bands from each cluster set respectively, and screening D characteristic wave bands from each cluster set respectively; and respectively constructing the feature matrixes of each cluster set according to the average value of the reflectivity of each field area in each cluster set on each feature wave band.
In step S4, the optimal jth cluster set C obtained in step S3 is subjected to j The characteristic wave band screening is carried out, and the specific mode is as follows:
s41 initializes d=1, randomly selects a band λ in the j-th cluster set ω The wave band lambda ω Vector of (3)The start vector h as the d-th projection d The method comprises the steps of carrying out a first treatment on the surface of the Band lambda in the j-th cluster set ω The vector formed by the average value of the reflectivity of each field area is the wave band lambda ω Vector of-> I.e. < -> For the j-th cluster set C j The kth field area in the middle is in the wave band lambda ω And the average value of the reflectivity.
S42, for the j-th cluster set C j Performing the d-th projection, and respectively calculating the initial vector h of the d-th projection d For the j-th cluster set C j Projection vectors of respective bands in the (c).
The start vector h of the d-th projection d For the j-th cluster set C j N-th band lambda of n Vector of (3)Projection vector of (a)The calculation method is as follows:
wherein n=1, 2, …, N and expressed +.>Is a transpose of (2); h is a d T H is expressed as d Is a transpose of (a).
The start vector h for the d-th projection d For the j-th cluster set C j N-th band lambda of n Vector of->Is included in the projection vector of (a).
The j-th cluster set C j N-th band lambda of n Vector of (3)The method comprises the following steps: the j-th cluster set C j N-th band lambda of n Vector of the average value of the reflectivity of each field area +.>I.e. < -> For the j-th cluster set C j The average value of the reflectivity of the kth field area on the nth wave band; n=1, 2, …, N.
S43, recording the wave band corresponding to the largest projection vector, wherein the wave band is the j-th aggregation set C j Characteristic wave band of the d-th projection of (2)
S44, for the j-th cluster set C j The next projection, namely the (d+1) th projection, is carried out, and the characteristic wave band of the (d) th projection is obtainedIs the reflectance mean vector x of (2) d As the start vector h of the (d+1) th projection d+1 According to the steps S42-S43, d+1th projection is carried out to obtain characteristic wave band of d+1th projection +.>
S45, if d+1=D, the next projection is not performed, namely the j-th cluster C is completed j Screening D characteristic wave bands; if d+1 < D, then the j-th cluster C is processed in the manner of step S44 j Continuing to go downAnd (5) one projection.
S46, from the j-th cluster set C j After screening D characteristic wave bands, according to the j-th aggregation set C j The average value of the reflectivity of each field area on each characteristic wave band is used for constructing a j-th cluster set C j Feature matrix of (a) For the j-th cluster set C j The kth field region in (a) is in the (d) th characteristic band +.>And the average value of the reflectivity.
S5, respectively corresponding the feature matrix of each aggregation set to the chlorophyll content measurement value of the corresponding field area in the tag set Y, and training to obtain an inversion model;
the input of the inversion model is a vector formed by the average value of the reflectivity of a certain field area on each wave band, namely the input is the vector of the field area, and the output of the inversion model is the inversion value of the chlorophyll content of the field area.
In step S5, the training method of the inversion model is as follows:
s51, the feature matrix of each cluster is equivalent to a decision tree, the feature wave bands of the cluster are the feature wave bands of the corresponding decision tree, and the feature matrix of each cluster is respectively corresponding to the tag set Y;
s52, training to obtain each decision tree based on the feature matrix of each optimal aggregation set and chlorophyll content measurement values of the corresponding field area in the tag set Y, integrating a plurality of decision trees, and establishing a random forest regression model, namely an inversion model.
When the inversion model is utilized for prediction, a test sample is input, the test sample is a vector formed by the average value of the reflectivities of a certain field area on each wave band, each decision tree in the random forest regression model predicts according to the average value of the reflectivities of the corresponding characteristic wave bands in the test sample, each decision tree outputs a chlorophyll content predicted value of the field area, an average value is calculated on the chlorophyll content predicted value of the field area output by each decision tree, and the average value is used as a final inversion value of the chlorophyll content of the field area.
In this embodiment, the test area is in a factory (116 ° 37'e,33 ° 13' n) of the black soil ecological environment of the alpinia galanga in rural areas of agriculture, fig. 2 is an acquired hyperspectral remote sensing image of the unmanned aerial vehicle in the test field, and the test field is in a white dotted line frame. Fig. 3 shows the result of model prediction of chlorophyll content of a part of field block area in a test field, and shows that the total chlorophyll content value of winter wheat in the jointing period ranges from 36 to 53, and the model can realize accurate prediction of chlorophyll content of winter wheat leaves, namely relative chlorophyll content from point to face. FIG. 4 is a verification accuracy graph of an inversion model, and the model predicted SPAD value and the measured SPAD value have strong positive correlation to determine a coefficient R as shown in FIG. 4 2 I.e. the goodness of fit reaches 0.77.
The method is based on hyperspectral remote sensing influence of an unmanned aerial vehicle, fully utilizes rich hyperspectral information of wheat fields, combines measured values of chlorophyll content of wheat leaves in each field area, adopts a pre-clustered RF stepwise regression advanced machine learning algorithm, establishes an inversion model, inverts chlorophyll content of winter wheat leaves in a jointing period, obtains a good inversion result, can intuitively and accurately display spatial distribution of chlorophyll content of winter wheat in a long-term positioning test field and characteristics thereof, provides a novel and reliable method for monitoring chlorophyll information of the wheat leaves and judging influence of soil fertility on chlorophyll content, and has higher practical value.
The above embodiments are merely preferred embodiments of the present invention and are not intended to limit the present invention, and any modifications, equivalent substitutions and improvements made within the spirit and principles of the present invention should be included in the scope of the present invention.
Claims (7)
1. The crop chlorophyll content inversion method based on unmanned aerial vehicle hyperspectral remote sensing is characterized by comprising the following steps of:
s1, measuring chlorophyll content of leaves in a field, dividing the field into M field areas, and obtaining a chlorophyll content measurement value of each field area;
constructing a tag set Y, Y= { Y according to chlorophyll content measurement values of each field area m }，m＝1,2,…,M，y m Chlorophyll content measurement for the mth field;
s2, acquiring hyperspectral remote sensing images of fields by using an unmanned aerial vehicle, and obtaining the average value of the reflectivity of each field area on each wave band on N wave bands;
according to the average value of the reflectivity of each field area on each wave band, a sample set X is constructed,m＝1,2,…,M，n＝1,2,…,N；/>for the m field area at the n-th wave band lambda n The average value of the reflectivity;
s3, clustering vectors of field areas in the sample set X to obtain J clustering sets;
each cluster set comprises vectors of a plurality of field areas; the vector formed by the average value of the reflectivity of the field area on each wave band is the vector of the field area;
wherein the j-th cluster set isI.e. < ->
j＝1,2,…,J；
K j For the j-th cluster set C j In (a) and (b)Total number of vectors for field areas;
for the j-th cluster set C j Vectors of the kth field region, i.e. the jth cluster center z j Vector of the mean value of the reflectivity of the kth field over the individual wavelength bands,/->I.e. < -> For the j-th cluster set C j The average value of the reflectivity of the kth field area on the nth wave band;
s4, screening characteristic wave bands from each cluster set respectively, and screening D characteristic wave bands from each cluster set respectively; respectively constructing a feature matrix of each cluster set according to the average value of the reflectivity of each field area in each cluster set on each feature wave band;
wherein the j-th cluster set C j Is given by (a)
For the j-th cluster set C j The kth field region in (a) is in the (d) th characteristic band +.>The average value of the reflectivity;
s5, respectively corresponding the feature matrix of each aggregation set to the chlorophyll content measurement value of the corresponding field area in the tag set Y, and training to obtain an inversion model;
the input of the inversion model is a vector formed by the average value of the reflectivity of a certain field area on each wave band, namely the input is the vector of the certain field area, and the output of the inversion model is the inversion value of the chlorophyll content of the field area.
2. The method for inverting chlorophyll content of crops based on hyperspectral remote sensing of unmanned aerial vehicle as claimed in claim 1, wherein in step S2, the mth field area is at the nth wave band lambda n On the average of the reflectivityThe calculation mode of (2) is as follows:
wherein,for the nth band lambda n The reflectivity value of the a pixel point in the m field area; n=1, 2, …, N; m=1, 2, …, M; a=1, 2, …, a m ，A m Is the total number of pixel points in the mth field area.
3. The method for inverting the chlorophyll content of crops based on hyperspectral remote sensing of an unmanned aerial vehicle according to claim 1, wherein in the step S3, the vectors of the field areas in the sample set X are clustered, and the specific mode is as follows:
s31, randomly selecting vectors of J field areas in a sample set X, wherein J is less than M; taking the vectors of the J field areas as initial values of clustering centers in a clustering algorithm to form a clustering center set Z; taking the spectral angular distance as a clustering distance parameter in a clustering algorithm;
s32, respectively calculating the spectral angular distance between the vector of each field area in the sample set X and each cluster center in the cluster center set Z, respectively clustering the vector of each field area in the sample set X onto the cluster center with the smallest spectral angular distance to obtain the vector of the field area clustered by each cluster center, wherein the vector of the field area clustered by each cluster center forms a cluster set corresponding to the cluster center, and J cluster sets are obtained;
vector X of the mth field in sample set X m With the j-th cluster center Z in the cluster center set Z j Spectral angular distance between SAD (x) m ,z j ) The calculation method is as follows:
wherein SAD (x) m ,z j ) Vector x for the mth field in the sample set m With the j-th cluster center Z in the cluster center set Z j Spectral angular distance between;
x m the vector of the m field area in the sample set X is the vector formed by the average value of the reflectivity of the m field area in the sample set X on each wave band,i.e. < -> For the mth field region in sample set X at the nth band lambda n The average value of the reflectivity; m=1, 2, …, M;
z j is the j-th cluster center Z in the cluster center set Z j Namely, the vector formed by the average value of the reflectivity of the jth field area in the cluster center set Z on each wave bandI.e. < -> For the jth field area in the cluster center set Z at the nth wave band lambda n The average value of the reflectivity; j=1, 2, …, J;
s33, respectively calculating average vectors of all the aggregation sets;
wherein the j-th cluster set C j Is the average vector delta of (2) j The calculation mode of (2) is as follows:
s34, updating the clustering centers of the clustering sets according to the average vectors of the clustering sets respectively to obtain an optimal clustering center, and obtaining an optimal clustering set corresponding to the optimal clustering center;
wherein, the j-th cluster set C is calculated j Is of the center vector delta of (2) j With the jth cluster center z j If the difference is smaller than or equal to the set threshold T, the j-th cluster center z j Without updating, the jth cluster center z j I.e. the optimal cluster center, the jth cluster set C j The optimal aggregate set is obtained; if the difference is greater than the set threshold value T, the jth clustering center z is subjected to j Updating the j-th cluster center z j Updated to the value of the j-th cluster set C j Is of the center vector delta of (2) j Re-clustering, i.e. re-executing steps S32-S34 until the optimal j-th cluster center z is obtained j The optimal jth cluster center z j Corresponding collection C j And the optimal aggregate set is obtained.
4. The method for inverting the chlorophyll content of crops based on hyperspectral remote sensing of an unmanned aerial vehicle according to claim 1, wherein in the step S4, characteristic wave band screening is performed on a j-th cluster set, and the specific mode is as follows:
s41, initializing d=1, randomly selecting a band λ in the j-th cluster set ω The wave band lambda ω Vector of (3)The start vector h as the d-th projection d The method comprises the steps of carrying out a first treatment on the surface of the Band lambda in the j-th cluster set ω The vector formed by the average value of the reflectivity of each field area is the wave band lambda ω Vector of-> I.e. < -> For the j-th cluster set C j The kth field area in the middle is in the wave band lambda ω The average value of the reflectivity;
s42, performing the d-th projection on the j-th cluster set, and respectively calculating the initial vector h of the d-th projection d Projection vectors of each wave band in the j-th cluster set are projected;
wherein the start vector h of the d-th projection d For the nth band lambda in the jth cluster set n Vector of (3)Projection vector of +.>The calculation method is as follows:
wherein n=1, 2, …, N and expressed +.>Is a transpose of (2); h is a d T H is expressed as d Is a transpose of (2);
the start vector h for the d-th projection d For the j-th cluster set C j N-th band lambda of n Vector of->Is a projection vector of (a);
the j-th cluster set C j N-th band lambda of n Vector of (3)The method comprises the following steps: the j-th cluster set C j N-th band lambda of n Vector of the average value of the reflectivity of each field area +.>I.e. < -> For the j-th cluster set C j The average value of the reflectivity of the kth field area on the nth wave band; n=1, 2, …, N;
s43, recording the wave band corresponding to the largest projection vector, wherein the wave band is the characteristic wave band of the (d) projection of the (j) th aggregation set
S44, carrying out next projection, namely, d+1st projection, on the jth cluster set, and carrying out characteristic wave band of the d projectionIs the reflectance mean vector x of (2) d As the start vector h of the (d+1) th projection d+1 According to the steps S42-S43, d+1th projection is carried out to obtain characteristic wave band of d+1th projection +.>
S45, if d+1=D, the next projection is not performed, and the screening of the D characteristic wave bands of the j-th aggregation set is completed; if d+1 is less than D, continuing to project the j-th cluster set for the next time in a mode of step S44;
s46, from the j-th cluster set C j After screening D characteristic wave bands, according to the j-th aggregation set C j The average value of the reflectivity of each field area on each characteristic wave band is used for constructing a j-th cluster set C j Feature matrix of (a)
5. The method for inverting the chlorophyll content of crops based on hyperspectral remote sensing of an unmanned aerial vehicle according to claim 1, wherein in step S5, the inversion model is trained as follows:
s51, the feature matrix of each cluster is equivalent to a decision tree, the feature wave bands of the cluster are the feature wave bands of the corresponding decision tree, and the feature matrix of each cluster is respectively corresponding to the tag set Y;
s52, training and generating each decision tree based on the feature matrix of each aggregation set and the chlorophyll content measurement value of the corresponding field area in the tag set Y, integrating a plurality of decision trees, and establishing a random forest regression model, namely an inversion model.
6. The method for inverting the chlorophyll content of crops based on hyperspectral remote sensing of an unmanned aerial vehicle according to claim 5, wherein when the inversion model is utilized for prediction, a test sample is input, the test sample is a vector formed by the average value of the reflectivities of a certain field area on each wave band, each decision tree in the random forest regression model extracts the average value of the reflectivities of the corresponding characteristic wave bands in the test sample respectively, each decision tree predicts according to the average value of the reflectivities of the corresponding characteristic wave bands, each decision tree outputs a chlorophyll content predicted value of the field area, an average value is calculated for the chlorophyll content predicted value of the field area output by each decision tree, and the average value is used as the ultimate inversion value of the chlorophyll content of the field area.
7. The method for inverting the chlorophyll content of crops based on hyperspectral remote sensing of unmanned aerial vehicle according to claim 1, wherein n=176.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110837401.2A CN113537108B (en) | 2021-07-23 | 2021-07-23 | Crop chlorophyll content inversion method based on unmanned aerial vehicle hyperspectral remote sensing |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110837401.2A CN113537108B (en) | 2021-07-23 | 2021-07-23 | Crop chlorophyll content inversion method based on unmanned aerial vehicle hyperspectral remote sensing |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113537108A CN113537108A (en) | 2021-10-22 |
| CN113537108B true CN113537108B (en) | 2024-02-23 |
Family
ID=78088858
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110837401.2A Active CN113537108B (en) | 2021-07-23 | 2021-07-23 | Crop chlorophyll content inversion method based on unmanned aerial vehicle hyperspectral remote sensing |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113537108B (en) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114019082B (en) * | 2021-11-19 | 2024-05-14 | 安徽省农业科学院土壤肥料研究所 | Soil organic matter content monitoring method and system |
| CN116089660A (en) * | 2023-03-16 | 2023-05-09 | 西安中创博远网络科技有限公司 | Data storage method of monitoring system |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110398466A (en) * | 2019-08-05 | 2019-11-01 | 北京绿土科技有限公司 | Crop growth state monitoring method based on remote-sensing inversion |
| CN112504979A (en) * | 2020-11-12 | 2021-03-16 | 内蒙古工业大学 | Water body chlorophyll concentration inversion method and device based on hyperspectral remote sensing |
| CN112749627A (en) * | 2020-12-09 | 2021-05-04 | 北京星衡科技有限公司 | Method and device for dynamically monitoring tobacco based on multi-source remote sensing image |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110455722A (en) * | 2019-08-20 | 2019-11-15 | 中国热带农业科学院橡胶研究所 | Rubber tree blade phosphorus content EO-1 hyperion inversion method and system |
-
2021
- 2021-07-23 CN CN202110837401.2A patent/CN113537108B/en active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110398466A (en) * | 2019-08-05 | 2019-11-01 | 北京绿土科技有限公司 | Crop growth state monitoring method based on remote-sensing inversion |
| CN112504979A (en) * | 2020-11-12 | 2021-03-16 | 内蒙古工业大学 | Water body chlorophyll concentration inversion method and device based on hyperspectral remote sensing |
| CN112749627A (en) * | 2020-12-09 | 2021-05-04 | 北京星衡科技有限公司 | Method and device for dynamically monitoring tobacco based on multi-source remote sensing image |
Non-Patent Citations (1)
| Title |
|---|
| 基于正交变换与SPXY样本划分的冬小麦叶绿素诊断;毛博慧;孙红;刘豪杰;张俊逸;李民赞;杨立伟;;农业机械学报;20171225(第S1期);全文 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113537108A (en) | 2021-10-22 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Li et al. | Multi-LUTs method for canopy nitrogen density estimation in winter wheat by field and UAV hyperspectral | |
| CN113537108B (en) | Crop chlorophyll content inversion method based on unmanned aerial vehicle hyperspectral remote sensing | |
| CN112051222A (en) | River and lake water quality monitoring method based on high-resolution satellite image | |
| CN111241912A (en) | Multi-vegetation index rice yield estimation method based on machine learning algorithm | |
| CN113466143B (en) | Soil nutrient inversion method, device, equipment and medium | |
| CN114529097B (en) | Multi-scale crop phenological period remote sensing dimensionality reduction prediction method | |
| CN108520127A (en) | A kind of EO-1 hyperion inversion method of seeds leaf area index | |
| CN115372282B (en) | Farmland soil water content monitoring method based on hyperspectral image of unmanned aerial vehicle | |
| CN112836725A (en) | Weak supervision LSTM recurrent neural network rice field identification method based on time sequence remote sensing data | |
| CN114692991A (en) | Wolfberry yield prediction method and system based on deep learning | |
| CN116912690A (en) | Forest leaf area index inversion acquisition method and system based on data fusion | |
| CN114676636A (en) | Grassland area soil moisture rapid inversion method integrating vegetation and habitat characteristics | |
| CN115169728A (en) | Soil fertility prediction method based on simplified neural network | |
| CN115759524A (en) | Soil productivity grade identification method based on remote sensing image vegetation index | |
| CN114782840A (en) | Real-time wheat phenological period classification method based on unmanned aerial vehicle RGB images | |
| CN112257531B (en) | Remote sensing monitoring method for forest land change based on diversity feature combination | |
| CN116501925B (en) | Method, system, equipment and medium for predicting leaf area index under multi-source remote sensing data | |
| CN116994126A (en) | Crop leaf area index obtaining method and device based on canopy reflectivity spectrum | |
| CN116151454A (en) | Method and system for predicting yield of short-forest linalool essential oil by multispectral unmanned aerial vehicle | |
| CN112991425B (en) | Water area water level extraction method and system and storage medium | |
| CN114692971A (en) | Crop yield prediction method and device based on yield difference | |
| CN114154694A (en) | Method for predicting plant growth state in cotton topping period based on multispectral remote sensing of unmanned aerial vehicle | |
| Rananavare et al. | Monocot Crop Yield Prediction using Sentinal-2 Satellite Data | |
| CN116256334A (en) | Crop leaf net photosynthetic rate inversion method, system, equipment and medium | |
| CN114036461A (en) | Crop total nitrogen content inversion method and system based on hyperspectral remote sensing |
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 |