CN110909652B - Monthly-scale dynamic extraction method of crop planting structure based on texture feature optimization - Google Patents

Abstract

The invention discloses a method for dynamically extracting a monthly scale of a crop planting structure with optimal textural features, which comprises the following steps of: s1: determining the space range of an analysis area, preparing data, collecting a time series satellite remote sensing data set not larger than a month scale, uniformly processing the time series satellite remote sensing data set into month scale data in time, and simultaneously completing pre-acquisition of sample data in a research area; s2: utilizing the preprocessed monthly scale satellite remote sensing image data to calculate texture features of the image based on the gray level co-occurrence matrix; s3: calculating the mean value and the variance of different texture characteristic quantities based on the actually measured samples, and calculating the differentiable capacity of the texture characteristic quantities among different samples; the method solves the problems that the traditional method for acquiring the crop planting structure information neglects the screening of the classification characteristic quantity and increases the time complexity.

Description

Monthly-scale dynamic extraction method of crop planting structure based on texture feature optimization

technical field

The invention relates to the field of remote sensing planting structure monitoring, in particular to a monthly-scale dynamic extraction method for crop planting structures with optimized texture features.

Background technique

Crop planting structure reflects the utilization of agricultural production resources by human agricultural production in a spatial range. It is an important information for studying crop types, quantity structure and spatial distribution characteristics, and is also the basis for adjusting and optimizing crop structure and fine management of agricultural water use. The traditional method of obtaining crop planting structure information mainly classifies images based on the texture information calculated by the gray level co-occurrence matrix, but ignores the screening of classification features, which increases the time complexity.

Texture features are widely used in the classification process of high-resolution images. High-resolution satellite images have quickly occupied the satellite application market of urban planning, regional environmental monitoring and evaluation with their clear and delicate ground object structure, shape and texture information. For example, Dekker R J. 2003 established texture features based on SAR to analyze urban building areas, and found that there is a certain influence between the texture feature parameters of crops, and the errors caused by the classification accuracy need to be removed in the identification process; Chen Junying et al. 2007 based on The spectral and texture features of IKONOS images act on the decision tree classification and identification process of agricultural vegetation, and it is found that the classification accuracy of only spectral features can reach 83%, and the accuracy of texture features can be improved to 91%; Elmqvist B 2008, ZhouW 2009, etc. use high spatial resolution land use mapping using high-rate images, and verifying the analysis results in the field; Ye Shiping 2008 obtained Quick Bird image texture features based on the gray-scale co-occurrence matrix calculation principle, so that the extraction accuracy of built-up areas was as high as 93%; Hou Xue et al. 2013 used SPOT-NDVI time series analysis The research results are consistent with the existing results; Liu Kebao et al. 2014 used RapidEye 5-meter resolution remote sensing data to extract crops, and found that 5-meter resolution can reduce a large amount of errors in the classification process and improve the classification accuracy. In 2015, Liu Guodong et al. completed a remote sensing sampling survey based on the texture of high-scoring 2/8-meter images, and used crop phenology and multi-temporal images to complete crop extraction, and the final accuracy was better than 80%; Piazza GA et al. Research on resolution shadows, complete forest mapping, classify and compare based on different classifiers, get the best classifier and apply it to large-scale identification. Song Qian et al. 2016, Zhang Chao et al. 2016 extracted texture features of specific crops based on Gaofen-1 PMS, and realized the change map of crop planting area through object-oriented classification; Zheng Lijuan 2017 was based on the application of Gaofen series satellite analysis in agriculture, and The advantages of Gaofen-1 and Gaofen-6 in the extraction process are put forward. However, these do not calculate and filter the texture features, which leads to too many classification features and increases the classification time complexity; at the same time, the correlation between the features is ignored, which affects the final classification accuracy, and it is difficult to realize the application of remote sensing crop classification. .

SUMMARY OF THE INVENTION

In order to solve the problems existing in the prior art, the present invention provides a monthly-scale dynamic extraction method of crop planting structure with optimized texture features, which solves the problem that the traditional method for obtaining crop planting structure information ignores the screening of classification feature quantities, which increases the time complexity. degree issue.

The technical scheme adopted in the present invention is that the monthly-scale dynamic extraction method of crop planting structure with preferred texture features comprises the following steps:

S1: Determine the spatial scope of the analysis area and prepare the data, collect time series satellite remote sensing data sets not larger than the monthly scale, and uniformly process it into monthly scale data in time, and complete the pre-acquisition of sample data in the study area;

S2: According to the preprocessed satellite remote sensing image data, the image texture features are calculated based on the gray level co-occurrence matrix, and eight feature quantities are used to describe the texture characteristics of crops;

S3: Calculate the mean and variance of different texture feature quantities based on the measured samples, and calculate the distinguishability of texture feature quantities between different samples;

S4: Establish an optimal formula based on the distinguishability of each feature quantity, and use the formula to determine the optimal number of texture feature quantities participating in the classification, and construct it as a new image;

S5: Use the random forest classifier to finely identify the types of crops in the study area, realize the fine management of crops, generate a thematic map of the spatiotemporal distribution of crops in a complete time series, and verify the accuracy.

Preferably, S1 includes the following steps:

S11: According to the location and extent of the study area, select GF-1WFV data with high temporal resolution and high spatial resolution. If the data source cannot be completely covered, consider using sention-2, Gaofen-2, landsat8 or HJ-1A/B instead, while investigating the types of crops and their respective growth phenology within the scope of the examples;

S12: Perform remote sensing image processing on the collected data. If there is alternative data, it needs to be resampled to a unified spatial resolution;

S13: The representativeness, typicality and timeliness of sample collection should be considered. The study area is divided into n areas with the same area by establishing a regular grid, and different crop samples are selected in each area.

Preferably, S2 includes the following steps:

S21: Calculate the amount of texture feature information based on the gray-level co-occurrence matrix, and count the gray-level correlation coefficient between two pixels at a certain distance according to the gray-level co-occurrence matrix GLCM, and its expression is:

p(i,j)=[p(i,j,d,θ)]

Among them, P(i,j) is the frequency of the same pixel pair when the distance and direction are determined; d is the distance from the pixel point, and the angle of the line vector between the two pixels is θ, usually θ takes 0°, 45° , 90° and 135°;

S22: When calculating the texture by using the grayscale co-occurrence matrix, select eight feature quantities to characterize the characteristics of the texture:

Average: It reflects the regularity of the average gray level and texture in the window. The calculation formula is:

Variance: It reflects the degree of deviation of the matrix elements from the mean value and the size of the grayscale change. The calculation formula is:

where μ is the mean of p(i,j),

Contrast: It reflects the clarity of the image and the depth of the texture grooves. The calculation formula is:

Inverse gap: It reflects the smoothness of the image distribution and is a measure of the uniformity of the image. Its calculation formula is:

Difference degree: It is used to detect the difference degree of the image, and its calculation formula is:

The amount of information contained in the image: it measures the randomness of the image texture, which is a characteristic parameter to measure the randomness of the gray level distribution, and represents the degree of confusion of the gray level of the image. The calculation formula is:

Image gray distribution uniformity: It reflects the image gray distribution uniformity and texture thickness. The calculation formula is:

Similarity: It reflects the similarity of the elements of the spatial grayscale co-occurrence matrix in the row or column direction. The calculation formula is:

where:

Preferably, S3 includes the following steps:

S31: Count the mean and variance of different samples on each feature quantity, and the calculation formula is:

S32: Calculate the degree of separability of the constructed samples based on the Bavarian distance, and calculate the distinguishability between different samples for different texture feature quantities. The calculation formula is:

In the formula, in the formula, μ is the mean value of the same feature between two different categories on the image, and σ is the standard deviation of the same feature between two different categories.

Preferably, S4 includes the following steps:

S41: Calculate the overall separation ability of each feature quantity in different crop categories, and the specific calculation formula is:

In the formula, D _ij is the degree of separability among different samples;

S42: Sort the feature quantities in descending order of the distinguishing ability value corresponding to each feature quantity;

S43: Calculate the cumulative sum of the separability of each feature quantity, and the calculation formula is:

S44: Construct a corresponding function expression based on the cumulative separation value corresponding to the texture feature quantity and the corresponding precision, and find that it satisfies the basic law of the logarithmic function, so the derivative corresponding to the function is obtained, that is, the slope of each point, and the calculation formula is: :

y=0.0625ln(x)+0.642

R ² =0.9191

Among them, x represents the cumulative D value corresponding to each feature quantity;

S45: Calculate the difference between the slopes. If the difference is less than 0.0001, the calculated points correspond to the optimum number of texture features, and output the results. The specific preferred calculation formula is:

D=|y' _i -y' _i-1 |

S46: Combine the filtered texture feature quantities to form a new image.

Preferably, S5 includes the following steps:

S51: In the unit of month, identify the planting structure in the study area, and generate a thematic map of the temporal and spatial distribution of crops in a complete time series;

S52: Validate the result according to the validation sample, and obtain the overall classification accuracy and Kappa coefficient.

The beneficial effects of the monthly-scale dynamic extraction method of crop planting structure with the preferred texture feature of the present invention are as follows:

1. The present invention optimizes and ultimately serves the identification of crop planting structures based on image textures, which can not only effectively identify crop planting structures, but also save time and reduce computer workload.

2. The technology has the characteristics of fast calculation and strong applicability, which effectively improves the classification limitations of medium and high-resolution data, improves the classification rate, and improves the accuracy and efficiency of Gaofen-1 in classification. Business promotion is of great significance.

Description of drawings

FIG. 1 is a general flow chart of the monthly-scale dynamic extraction method of crop planting structure with preferred texture features of the present invention.

FIG. 2 is a step-by-step flow chart of S1 of the monthly-scale dynamic extraction method of crop planting structure with preferred texture features of the present invention.

FIG. 3 is a step-by-step flow chart of S2 of the monthly-scale dynamic extraction method of crop planting structure with preferred texture features of the present invention.

FIG. 4 is a step-by-step flow chart of S3 of the monthly-scale dynamic extraction method of crop planting structure with preferred texture features of the present invention.

FIG. 5 is a step-by-step flow chart of S4 of the monthly-scale dynamic extraction method of crop planting structure with preferred texture features of the present invention.

FIG. 6 is a step-by-step flow chart of S5 of the monthly-scale dynamic extraction method of crop planting structure with preferred texture features of the present invention.

Detailed ways

The specific embodiments of the present invention are described below to facilitate those skilled in the art to understand the present invention, but it should be clear that the present invention is not limited to the scope of the specific embodiments. For those of ordinary skill in the art, as long as various changes Such changes are obvious within the spirit and scope of the present invention as defined and determined by the appended claims, and all inventions and creations utilizing the inventive concept are within the scope of protection.

As shown in Figure 1, a monthly-scale dynamic extraction method of crop planting structure includes the following steps:

S1: Determine the spatial scope of the analysis area and prepare the data, collect time series satellite remote sensing data sets not larger than the monthly scale, and uniformly process it into monthly scale data in time, and complete the pre-acquisition of sample data in the study area;

S2: According to the preprocessed satellite remote sensing image data, the image texture features are calculated based on the gray level co-occurrence matrix, and eight feature quantities are used to describe the texture characteristics of crops;

S3: Calculate the mean and variance of different texture feature quantities based on the measured samples, and calculate the distinguishing ability of texture feature quantities between different samples, and finally evaluate the distinguishing ability of different feature quantities between samples;

S4: Establish an optimal formula based on the distinguishability of each feature quantity, and use the formula to determine the optimal number of texture feature quantities participating in the classification, and construct it as a new image;

S5: Use the random forest classifier to finely identify the types of crops in the study area, realize the fine management of crops, generate a thematic map of the spatiotemporal distribution of crops in a complete time series, and verify the accuracy.

As shown in Figure 2, S1 of this embodiment includes the following steps:

S11: According to the location and scope of the study area, select the GF-1WFV data independently developed by my country with high temporal resolution and high spatial resolution. If the data source cannot be completely covered, consider using sention-2, with a high score of 2 No., landsat8 or HJ-1A/B instead, while investigating the crop types within the scope of the examples and their respective growth phenology;

S12: Perform remote sensing image processing on the collected data. If there is alternative data, it needs to be resampled to a unified spatial resolution;

S13: The representativeness, typicality, and timeliness of sample collection should be considered. The study area is divided into n areas with the same area by establishing a regular grid, and crop samples are selected in each area.

As shown in Figure 3, S2 of this embodiment includes the following steps:

S21: Calculate the amount of texture feature information based on the gray-level co-occurrence matrix, and count the gray-level correlation coefficient between two pixels at a certain distance according to the gray-level co-occurrence matrix GLCM, which represents the probability distribution of repeated gray levels, and its expression is:

p(i,j)=[p(i,j,d,θ)]

Among them, P(i,j) is the frequency of the same pixel pair when the distance and direction are determined; d is the distance from the pixel point, and the angle of the line vector between the two pixels is θ, usually θ takes 0°, 45° , 90° and 135°;

When using the gray level co-occurrence matrix to calculate the texture, eight commonly used feature quantities are selected to characterize the characteristics of the texture.

Mean: It reflects the average gray level and the regularity of the texture in the window.

Variance: It reflects the degree to which the matrix elements deviate from the mean and the size of the grayscale change. When the grayscale changes greatly, the larger the deviation from the mean value is, the larger the value is.

where μ is the mean of p(i,j). This value is larger when the grayscale changes in the image are larger.

Contrast: Reflects the clarity of the image and the depth of the texture grooves. The deeper the texture groove, the larger its value. The larger the value, the higher the clarity of the texture.

Homogeneity: Also known as the inverse gap, it reflects the smoothness of the image distribution and is a measure of the uniformity of the image. If the local gray level of the image is uniform, the value of the inverse difference moment is larger.

Dissimilarity: Used to detect the dissimilarity of images. If the difference in the local area varies greatly, the value will be larger.

Entropy: Measure the randomness of image texture, that is, the amount of information contained in the image. It is a characteristic parameter that measures the randomness of gray level distribution, and represents the degree of confusion in the gray level of the image. The greater the entropy, the greater the class uncertainty of the sample. On the contrary, if the gray level in the image is uniform, the direct value is small.

Angular Second Moment: It reflects the uniformity of the grayscale distribution of the image and the thickness of the texture, and is a measure of the uniformity of the grayscale distribution of the image. When the distribution of elements in the gray-level co-occurrence matrix is more concentrated near the main diagonal, it means that the gray-level distribution of the image in the local area is relatively uniform, the image has a thicker texture, and the value of the second-order moment of the angle is correspondingly large.

Correlation: It reflects the similarity of the elements of the spatial grayscale co-occurrence matrix in the row or column direction. If the matrix elements are uniform and equal, the correlation degree is high and the value is large; on the contrary, if the difference between the elements is large, the correlation is small.

where:

As shown in Figure 4, S3 of this embodiment includes the following steps:

S31: Count the mean and variance of different samples on each feature quantity, and the calculation formula is:

S32: Calculate the degree of separability of the constructed samples based on the Bhattacharyya distance, and calculate the distinguishability between different samples for different texture features. The calculation formula is:

In the formula, in the formula, μ is the mean value of the same feature between two different categories on the image, and σ is the standard deviation of the same feature between two different categories.

As shown in Figure 5, S4 of this embodiment includes the following steps:

S41: Calculate the overall separation ability of each feature quantity in different crop categories, and the specific calculation formula is as follows:

In the formula, D _ij is the degree of separability among different samples;

S42: Sort the feature quantities in descending order of the distinguishing ability value corresponding to each feature quantity;

S43: Calculate the cumulative sum of the separability of each feature quantity, and the calculation formula is:

S44: Construct a corresponding function expression based on the cumulative separation value corresponding to the texture feature quantity and the corresponding precision, and find that it satisfies the basic law of the logarithmic function, so the derivative corresponding to the function is obtained, that is, the slope of each point, as follows:

y=0.0625ln(x)+0.642

R ² =0.9191

Among them, x represents the cumulative D value corresponding to each feature quantity;

S45: Calculate the difference value of the slope, if the difference value is less than 0.0001, the corresponding point is the optimal number of the obtained texture feature quantity. and output the result. The specific optimal calculation formula is obtained as follows:

D=|y' _i -y' _i-1 |

S46: Combine the filtered texture feature quantities to form a new image.

As shown in Figure 6, S5 of this embodiment includes the following steps:

S51: Identify the ground objects in the study area according to the random forest principle;

S52: Set various parameters required by the classifier, input classification samples to identify and classify the crop planting structure in the research area, and complete the dynamic identification.

S53: Identify crops in the study area, and generate a thematic map of the spatial and temporal distribution of crops in a complete time series;

S54: Validate the result according to the validation sample, and obtain the overall classification accuracy and Kappa coefficient.

When this embodiment is implemented, the method proposed by the present invention is based on the principle that different crops have different texture features on the image. First, the types of crops in the experimental area are analyzed and summarized, and then based on the separability and accuracy of texture feature quantities. According to the calculation principle, the optimal texture feature quantity is screened, and then a new image is formed. Finally, the random forest classification method is used to identify the crop planting structure information. The technical solution has the characteristics of simplicity, effectiveness and strong applicability, and can quickly and accurately obtain a wide range of crop spatial distribution information, and improves the classification rate and classification accuracy caused by too many features involved in the classification in the traditional method. It improves the computational efficiency and accuracy of remote sensing monitoring of crop planting structure, and is helpful for the commercialization of remote sensing technology to monitor crop planting structure.

Claims (1)

1. The method for dynamically extracting the monthly scale of the crop planting structure with the optimized textural features is characterized by comprising the following steps of:

s1: determining the space range of an analysis area, preparing data, collecting a time series satellite remote sensing data set not larger than a month scale, uniformly processing the time series satellite remote sensing data set into month scale data in time, and simultaneously completing pre-acquisition of sample data in a research area;

s2: calculating image texture characteristics based on a gray level co-occurrence matrix according to the preprocessed satellite remote sensing image data, and describing the texture characteristics of crops by using eight characteristic quantities;

s3: calculating the mean value and the variance of different texture characteristic quantities based on the actually measured samples, and calculating the differentiable capacity of the texture characteristic quantities among different samples;

s4: establishing an optimal calculation formula based on the distinguishing capability of each characteristic quantity, determining the optimal number of the texture characteristic quantities participating in classification by using the formula, and constructing a new image;

s5: finely identifying the crop types in the research area by using a random forest classifier, realizing fine management of the crops, generating a space-time distribution thematic map of the crops with a complete time sequence and verifying the precision;

the S1 comprises the following steps:

s11: selecting GF-1WFV data with high time resolution and high spatial resolution according to the position and range of the research area, and considering using the position-2, high mark two, landsat8 or HJ-1A/B to replace if the data source can not completely cover, and simultaneously investigating the type of crops and the respective growth phenological period;

s12: processing the collected data by remote sensing images, and if substitute data appears, resampling the uniform spatial resolution;

s13: the representativeness, the typicality and the timeliness of the samples need to be considered for collecting the samples, a regular grid is established to divide a research area into n areas with the same area, and different crop samples are selected in each area;

the S2 comprises the following steps:

s21: calculating texture characteristic information quantity based on a gray level co-occurrence matrix, and counting a gray level correlation coefficient between two pixel points at a certain distance according to a gray level co-occurrence matrix GLCM, wherein the expression is as follows:

wherein,

is the frequency at which the same pixel pair occurs with distance and direction determinations; i, j represents the spatial position of the pixel;

the distance from a pixel point is the angle of a connecting line vector of two pixels

，

Taking 0 degree, 45 degrees, 90 degrees and 135 degrees;

s22, when the texture is calculated by utilizing the gray level co-occurrence matrix, eight characteristic quantities are selected to characterize the characteristics of the texture:

average value: the calculation formula of the regular degree of the gray level average value and the texture in the window is as follows:

wherein L represents the number of pixels;

variance: the degree of matrix element deviation from the mean value and the gray scale change are reflected, and the calculation formula is as follows:

wherein μ is the mean of p (i, j);

and (3) comparison: the degree of the definition of the image and the depth of the texture groove is reflected, and the calculation formula is as follows:

reverse difference: the smoothness of the image distribution is reflected, and is a measure of the uniformity of the image, and the calculation formula is as follows:

degree of difference: detecting the difference degree of the images;

amount of information contained in the image: the randomness of the image texture is measured, the characteristic parameters of the randomness of the gray level distribution are measured, the chaotic degree of the gray level of the image is represented, and the calculation formula is as follows:

uniformity of gray level distribution of the image: the image gray level distribution uniformity and texture thickness are reflected, and the calculation formula is as follows:

degree of similarity: the similarity degree of the elements of the space gray level co-occurrence matrix in the row or column direction is reflected, and the calculation formula is as follows:

wherein,

represents the mean value in the row direction;

indicates all in the column directionA value;

represents the standard deviation in the row direction;

represents the standard deviation in the column direction;

the S3 comprises the following steps:

s31: and (3) counting the mean value and the variance of different samples on each characteristic quantity, wherein the calculation formula is as follows:

s32: the method comprises the following steps of constructing separable degree calculation of samples based on the Pasteur distance, and calculating the differentiable capacity between different samples aiming at different texture characteristic quantities, wherein the calculation formula is as follows:

in the formula, mu _a Is the mean value, mu, of the same feature in the image between the a classes _b Mean, σ, between b classes for the same feature of the image _a Standard deviation, σ, between a classes for the same feature _b Standard deviation between categories b for the same feature;

the S4 comprises the following steps:

s41: calculating the total separation capacity of each characteristic quantity in different crop categories, wherein the specific calculation formula is as follows:

in the formula, D _ij The degree of separability between different samples;

s42: sorting the characteristic quantities according to the sequence of the distinguishing capacity values corresponding to the characteristic quantities from large to small;

s43: and calculating the cumulative sum of the separability of the characteristic quantities, wherein the calculation formula is as follows:

s44: constructing a corresponding function expression based on the accumulated separation value corresponding to the texture characteristic quantity and the corresponding precision, finding that the function expression meets the basic rule of a logarithmic function, and acquiring a derivative corresponding to the function, namely the slope of each point, wherein the calculation formula is as follows:

y = 0.0625ln(x) + 0.642

R² = 0.9191

wherein x represents the accumulated D value corresponding to each characteristic quantity; y represents the precision corresponding to the accumulated D value;

is the first partial derivative of y;

s45: and calculating the difference of the slopes, if the difference is less than 0.0001, the calculated point corresponds to the optimal number of the texture feature quantity, outputting the result, and acquiring a specific optimal calculation formula as follows:

s46: combining the screened texture characteristic quantities to form a new image;

the S5 comprises the following steps:

s51: recognizing the ground features of the research area according to a random forest principle;

s52: setting various parameters required by a classifier, inputting classification samples to identify and classify crop planting structures in a research area, and completing dynamic identification;

s53: identifying crops in a research area to generate a complete time sequence crop space-time distribution thematic map;

s54: and verifying the result according to the verification sample to obtain the overall classification precision and the Kappa coefficient.

Images