CN108959661B - 3S technology-based soil nutrient grade classification map generation method and precision evaluation method thereof - Google Patents

3S technology-based soil nutrient grade classification map generation method and precision evaluation method thereof Download PDF

Info

Publication number
CN108959661B
CN108959661B CN201810946788.3A CN201810946788A CN108959661B CN 108959661 B CN108959661 B CN 108959661B CN 201810946788 A CN201810946788 A CN 201810946788A CN 108959661 B CN108959661 B CN 108959661B
Authority
CN
China
Prior art keywords
classification
soil
land
soil nutrient
map
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201810946788.3A
Other languages
Chinese (zh)
Other versions
CN108959661A (en
Inventor
赵晋陵
董莹莹
黄林生
梁栋
徐超
张东彦
阮莉敏
翁士状
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Anhui University
Institute of Remote Sensing and Digital Earth of CAS
Original Assignee
Anhui University
Institute of Remote Sensing and Digital Earth of CAS
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Anhui University, Institute of Remote Sensing and Digital Earth of CAS filed Critical Anhui University
Priority to CN201810946788.3A priority Critical patent/CN108959661B/en
Publication of CN108959661A publication Critical patent/CN108959661A/en
Application granted granted Critical
Publication of CN108959661B publication Critical patent/CN108959661B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06NCOMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00Computing arrangements based on biological models
    • G06N3/02Neural networks
    • G06N3/08Learning methods
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06QINFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q10/00Administration; Management
    • G06Q10/06Resources, workflows, human or project management; Enterprise or organisation planning; Enterprise or organisation modelling
    • G06Q10/063Operations research, analysis or management
    • G06Q10/0639Performance analysis of employees; Performance analysis of enterprise or organisation operations

Abstract

The invention relates to a method for generating a grade classification map of soil nutrients based on a 3S technology and a precision evaluation method thereof, which overcome the defect that the grade classification of soil nutrients in a large spatial range is difficult to realize compared with the prior art. The invention comprises the following steps: generating a soil nutrient data space grid graph; setting a classification map of soil nutrient grades; extracting cultivated land plots; and generating a soil nutrient grade classification map of a plot scale. Based on the 3S technology, on the basis of extracting cultivated land plots by using field sampling points and RS acquired by a GPS, the spatial distribution characteristics of soil nutrients in the plot scale are analyzed by comprehensively using GIS spatial analysis and land statistics functions, and soil nutrient grade division, spatial mapping and precision evaluation in a large spatial range and on the plot scale are realized.

Description

3S technology-based soil nutrient grade classification map generation method and precision evaluation method thereof
Technical Field
The invention relates to the technical field of remote sensing data analysis, in particular to a soil nutrient grade classification map generation method and a precision evaluation method based on a 3S technology.
Background
The soil quality is an important criterion for determining the quality of cultivated land. The soil quality of cultivated land has important economic value for areas where the economy develops rapidly, but is greatly influenced by the environment and shows the characteristic of metamorphism. The current arable land research has been developed from the direction of only paying attention to quantity to the direction of the weight of quantity, quality and ecological benefit. On the premise of guaranteeing the farmland quantity, how to improve the farmland quality, especially to evaluate the environment quality, has become a hot direction for farmland research. Grading and grading of soil nutrients are of great significance for mastering regional resource conditions and guiding production practices, particularly for soil testing and fertilization.
In the prior art, many scholars propose many technical ideas aiming at the soil nutrient condition. Such as: the hole bin and the like adopt soil samples through GPS positioning, and research the soil fertility level change and the spatial distribution thereof by means of a GIS spatial analysis technology and a Kriging interpolation method; analyzing soil nutrients of cultivated land by means of Sulimacin and the like; road apple et al proposed 2005 characteristics of distribution and variation of farmland topsoil nutrients.
However, most of these studies are based on sampling point data to perform statistical analysis and area difference mapping, but do not scale down to the plot scale, making efficient generalization over a wide range difficult. Therefore, how to form a grade classification map of soil nutrients has become an urgent technical problem to be solved.
Disclosure of Invention
The invention aims to solve the defect that the classification of soil nutrients in a large spatial range is difficult to realize in the prior art, and provides a soil nutrient grade classification map generation method based on a 3S technology and a precision evaluation method thereof to solve the problems.
In order to achieve the purpose, the technical scheme of the invention is as follows:
a method for generating a grade classification map of soil nutrients based on a 3S technology comprises the following steps:
generating a soil nutrient data space grid map, selecting a radial basis function RBF space difference method for carrying out space difference, and generating a grid map of each nutrient;
setting a grading classification map of soil nutrients, calculating to obtain a soil nutrient comprehensive index (SNI) map based on pixels according to a space difference map and a weight coefficient of each nutrient index, and formulating a classification rule according to a value range of the SNI to obtain five classes of grading classification maps of soil nutrients, namely a very high class, a medium class, a low class and a very low class;
extracting cultivated land parcels, namely extracting the cultivated land parcels and spatial distribution thereof by adopting a step-by-step elimination method of layered classification extraction based on multi-temporal remote sensing images of key growth periods of main crops in a research area to generate a Shapefile vector file of the cultivated land parcels;
and generating a soil nutrient grade classification map of a plot scale, taking the farmland plots in a Shapefile format as mask files, and masking the soil nutrient comprehensive evaluation map generated based on the SNI index to form the soil nutrient grade classification map of the plot scale.
The generation of the soil nutrient data space grid map comprises the following steps:
when a soil sample is collected in a farmland, a sub-meter high-precision handheld GPS is used for acquiring longitude and latitude coordinates of a sampling point, the longitude and latitude coordinates are led into an ArcMap space analysis platform, and a point vector file in a Shapefile format is generated;
based on an interpolation method classification tree provided in geostationary analysis, acquiring a soil nutrient value of an unsampled point by using a Radial Basis Function (RBF) difference method;
by utilizing sampling point data, through training a sample by an RBF network, hiding a nonlinear function relation between the soil nutrient value and plane coordinates x and y in a converged network, taking the geographic coordinates of unknown points as network input, and performing simulation prediction on the network to obtain the soil nutrient value of the unknown points, wherein the function expression of the function expression is as follows:
z=f(x,y),
wherein, (x, y) is the geographic coordinate of the sampling point or the prediction point, the RBF network takes x, y as network input and takes the corresponding soil nutrient value as network output, namely the number of nodes of the input layer of the network is 2, and the number of nodes of the output layer is 1;
and (3) performing space difference on each soil nutrient by adopting an RBF difference method on an ArcMap space analysis platform to obtain a corresponding space grid map.
The cultivated land plot extraction method comprises the following steps:
selecting three-scene multi-temporal remote sensing images before planting, in a growth vigorous growing period and after harvesting according to the planting type and planting lunar calendar characteristics of main crops in a research area;
respectively carrying out data preprocessing on each scene image, including radiometric calibration, atmospheric correction and geometric fine correction;
in an eCoginization object-oriented classification software platform, performing hierarchical classification extraction combining a SEATH algorithm and nearest supervision classification on the preprocessed image, and excluding other land types to obtain the cultivated land type;
directly applying an SEATH algorithm to classify the two ground features with the distance close to 2 to generate a Shapefile of the cultivated land parcel;
for the distance between two ground objects is larger than 2, dividing the ground objects into large categories by using a user-defined normalized vegetation index; and subdividing the land types under the large classification, re-analyzing the available classification feature space for the land types by using a nearest supervision classification method, combining the classified large classification and the classified land types together to finish the classification of the cultivated land, and generating a Shapefile file of the cultivated land plots.
The method for generating the grading classification chart of the soil nutrients in the plot scale comprises the following steps:
selecting 4 indexes of soil organic matter SOM, total nitrogen TN or alkaline hydrolysis nitrogen AN, quick-acting phosphorus AP and quick-acting potassium AK;
obtaining the weight coefficient of each evaluation index of soil nutrients;
and (3) solving a comprehensive evaluation index, and solving the nutrient grade of each plot by adopting a comprehensive index constructed based on an addition model, wherein the formula is as follows:
SNI=∑Fi×Wi(i=1,2,3,.......,n),
in the formula: SNI represents the soil nutrient comprehensive index of the plot, FiThe score value, W, of the i-th indexiA weight coefficient indicating an i-th index;
in an ArcMap space analysis platform, pixel values of all nutrient indexes and weight coefficients of all evaluation indexes of soil nutrients are introduced into an SNI formula to obtain a soil nutrient comprehensive evaluation classification diagram, and scales are divided according to value domain grades of all evaluation indexes of the soil nutrients to obtain five classes of soil nutrient grade classification diagrams of high, medium, low and extremely low;
and (4) obtaining a soil nutrient grade classification map of a plot scale by using the extracted farmland plot Shapefile file as a mask.
The SEATH algorithm for classification comprises the following steps:
the characteristics are preferably such that,
the degree of relevance between two land categories on a certain characteristic is judged by adopting the degree of separation, the degree of separation is calculated by the JM distance, so as to evaluate the separability between the land categories, and the formula of the JM distance is shown as follows:
J=2(1-e-B),
in the formula, J represents JM distance, B represents Papanicolaou distance, and the calculation formula is shown as follows:
Figure BDA0001770433290000041
wherein m is1And m2Representing the mean, sigma, of a feature of a representative sample selected from two classes of output1And σ2It represents the variance of a feature in the two class samples. The value range of JM distance is [0, 2 ]]The closer the distance is to 2, the better the separability;
the determination of the characteristic threshold value is carried out,
calculating the optimal threshold values of two ground classes in a certain sample characteristic according to a Gaussian probability distribution formula:
p(x)=p(x|C1)p(C1)p(x|C2)p(C2),
wherein p (x | C)1) Represents land class C1The selected typical sample characteristic value obeys mean value m1Variance is σ1 2Normal distribution of (1), p (x | C)2) Represents land class C2The selected typical sample characteristic value obeys mean value m2Variance is σ2 2Normal ofDistributing;
the formula for the optimal threshold T is as follows:
Figure BDA0001770433290000042
in the formula (I), the compound is shown in the specification,
Figure BDA0001770433290000043
n1and n2Is class C of land1And C2The number of samples selected;
determining which type the object to be determined belongs to according to the selected typical feature sample and the nearest feature space distance constructed by calculation;
selecting a typical ground object sample;
constructing a feature space required by classification;
calculating the distance between the feature used for classification of the unclassified land types and the selected statistical feature of the typical sample by calculating the feature center of each land type, and classifying the unclassified land types into which land type according to the closer distance to the sample of which land type, wherein the calculation formula is as follows:
Figure BDA0001770433290000051
wherein d represents the distance between the sample land class and the land class o to be classified,
Figure BDA0001770433290000052
a feature value representing a class feature f of a typical sample,
Figure BDA0001770433290000053
characteristic value, sigma, representing a characteristic f of the ground class to be classifiedfThe standard deviation of the characteristic f is indicated.
The method also comprises a precision evaluation method for soil nutrient grading, which comprises the following steps:
obtaining corresponding geographic coordinates of the ground by using the partitioned soil nutrient grade data provided by the research area, and constructing a confusion matrix to obtain the total classification accuracy OA, wherein the formula is as follows:
Figure BDA0001770433290000054
in the formula, NaccurateNumber of pixels representing correct classification, which is the sum of diagonals of the confusion matrix, NtotalThe total pixel number is;
constructing a Kappa coefficient, wherein the formula is as follows:
Figure BDA0001770433290000055
in the formula, N represents the total number of pixels in the ground surface real classification; x is the number ofiiRepresents the value on the ith type diagonal in the confusion matrix; k represents the total classification number; x is the number ofi∑Representing the sum corresponding to the column of the ith class; x is the number of∑iAnd the corresponding sum of the rows where the ith class is located is shown.
Advantageous effects
Compared with the prior art, the method for generating the grading classification map of the soil nutrients based on the 3S technology and the precision evaluation method thereof are based on the 3S (GIS, RS remote sensing and GPS) technology, and on the basis of extracting cultivated land plots by using field sampling points and RSs acquired by GPS, the GIS spatial analysis and land statistics functions are comprehensively applied to analyze the spatial distribution characteristics of the soil nutrients in the plot scale, so that the grading classification of the soil nutrients in a large spatial range and on the plot scale, the spatial mapping and the precision evaluation are realized.
The discrete soil nutrient indexes are continuously and spatially expressed, so that the space expansion from the traditional ground random 'point' monitoring to 'surface' analysis is realized, and the grade classification and the space distribution difference of the soil nutrients can be quickly obtained; the method overcomes the one-sidedness of observing and dividing the soil nutrient grade based on the ground fixed point, greatly expands the space monitoring range of the soil nutrient, and reduces the workload and the cost of the traditional fixed point sampling and laboratory analysis and test.
Drawings
FIG. 1 is a sequence diagram of the method of the present invention;
FIG. 2 is a schematic diagram of a technical route for performing layered classification and extraction on cultivated land blocks in the invention;
FIG. 3a is a pseudo color synthesized image using Landsat TM 7, 4, 2 bands;
FIG. 3b is a block diagram space distribution diagram of cultivated land in Beijing City generated by the method of the present invention;
FIG. 4 is a classification chart of soil nutrients in a scale of Beijing city field produced by the method of the present invention.
Detailed Description
So that the manner in which the above recited features of the present invention can be understood and readily understood, a more particular description of the invention, briefly summarized above, may be had by reference to embodiments, some of which are illustrated in the appended drawings, wherein:
as shown in fig. 1, the method for generating a grade classification map of soil nutrients based on 3S technology includes the following steps:
firstly, generating a soil nutrient data space grid graph. And selecting a Radial Basis Function (RBF) space difference method to perform space difference, and generating a grid map of each nutrient. The method comprises the following specific steps:
(1) when soil samples are collected in the farmland, the sub-meter high-precision handheld GPS is used for obtaining longitude and latitude coordinates of sampling points, the longitude and latitude coordinates are led into the ArcMap space analysis platform, and a point vector file in a Shapefile format is generated.
(2) Based on an interpolation method classification tree provided in geostationary analysis, a Radial Basis Function (RBF) difference method is used for obtaining the soil nutrient value of an unsampled point.
Here, the RBF difference method using the radial basis functions is based on a comprehensive consideration of a measure or model of spatial autocorrelation, a level of condition/complexity of the model, a type of interpolation, smoothness of output, processing speed. Radial Basis Function (RBF) interpolation is a combination of a series of exact interpolation; i.e. the interpolation surface has to pass each measured sample value. The RBF interpolation method taking the geographic coordinates as network input only establishes the relation between the spatial geographic coordinates and the soil nutrient values.
By utilizing sampling point data, through training a sample by an RBF network, hiding a nonlinear function relation between the soil nutrient value and plane coordinates x and y in a converged network, taking the geographic coordinates of unknown points as network input, and performing simulation prediction on the network to obtain the soil nutrient value of the unknown points, wherein the function expression of the function expression is as follows:
z=f(x,y),
wherein, (x, y) is the geographic coordinate of the sampling point or the prediction point, the RBF network takes x, y as the network input, and takes the corresponding soil nutrient value as the network output, namely the number of nodes of the input layer of the network is 2, and the number of nodes of the output layer is 1.
(3) And (3) performing space difference on each soil nutrient by adopting an RBF difference method on an ArcMap space analysis platform to obtain a corresponding space grid map.
And secondly, setting a classification map of soil nutrient grades.
According to the prior art method, a soil Nutrient comprehensive index (SNI) (soil Nutrient index) graph based on pixels is obtained through calculation according to a spatial difference graph and a weight coefficient of each Nutrient index, and classification rules are formulated according to the value domain range of the SNI to obtain five classes of soil Nutrient grade classification graphs including extremely high, medium, low and extremely low.
And thirdly, extracting cultivated land plots.
Based on multi-temporal remote sensing images of key growth periods of main crops in a research area, a step-by-step elimination method of layered classification extraction is adopted to extract farmland plots and spatial distribution thereof, and a Shapefile vector file of the farmland plots is generated. The method comprises the following specific steps:
(1) and selecting three-scene multi-temporal remote sensing images before planting, in the growth vigorous growing period and after harvesting according to the planting type and the planting lunar calendar characteristics of main crops in the research area.
(2) And respectively carrying out data preprocessing on each scene image, including radiometric calibration, atmospheric correction and geometric fine correction.
(3) In an eCoginization object-oriented classification software platform, the preprocessed images are subjected to hierarchical classification extraction combining a SEATH algorithm and nearest supervision classification, and the types of cultivated land are obtained by excluding other land types.
As shown in fig. 2, for the distance between two land features close to 2, directly applying the sea algorithm to classify, and generating a Shapefile of the cultivated land parcel; for the distance between two ground objects is larger than 2, dividing the ground objects into large categories by using a user-defined normalized vegetation index; and subdividing the land types under the large classification, re-analyzing the available classification feature space for the land types by using a nearest supervision classification method, combining the classified large classification and the classified land types together to finish the classification of the cultivated land, and generating a Shapefile file of the cultivated land plots.
Here, two major categories, wetland and non-wetland, are taken as examples. For two categories where distances are not particularly large, such as vegetation and non-vegetation, a customized vegetation normalization index can be used to classify non-vegetation (artificial surfaces, others, dry land) and vegetation (grassland, woodland, paddy field); the non-vegetation and vegetation sub-divided land categories can be re-analyzed by using a nearest neighbor supervised classification method for classifying the land categories, and then the classified dry land and paddy field are combined together to complete the classification of the cultivated land.
The SEATH algorithm comprises the following steps of:
A1) the characteristics are preferred.
The degree of separation is used to judge the degree of association between two land classes on a certain characteristic, and the degree of separation is calculated by JM (Jeffries-Matusita) distance, so as to evaluate the separability between the land classes. The formula for the JM distance is as follows:
J=2(1-e-B),
in the formula, J represents JM distance, B represents Papanicolaou distance, and the calculation formula is shown as follows:
Figure BDA0001770433290000081
wherein m is1And m2Representing the mean, sigma, of a feature of a representative sample selected from two classes of output1And σ2It indicates that the samples in both categories are in a particularThe variance of the sign. The value range of JM distance is [0, 2 ]]The closer the distance is to 2, the better the separability.
A2) And determining a characteristic threshold value.
Calculating the optimal threshold values of two ground classes in a certain sample characteristic according to a Gaussian probability distribution formula:
p(x)=p(x|C1)p(C1)p(x|C2)p(C2),
wherein p (x | C)1) Represents land class C1The selected typical sample characteristic value obeys mean value m1Variance is σ1 2Normal distribution of (1), p (x | C)2) Represents land class C2The selected typical sample characteristic value obeys mean value m2Variance is σ2 2Normal distribution of (2);
the formula for the optimal threshold T is as follows:
Figure BDA0001770433290000091
in the formula (I), the compound is shown in the specification,
Figure BDA0001770433290000092
n1and n2Is class C of land1And C2The number of samples selected.
A3) And determining which class the object to be determined belongs to according to the selected typical feature sample and the calculated and constructed nearest characteristic space distance.
A31) Selecting a typical ground object sample;
A32) constructing a feature space required by classification;
A33) calculating the distance between the feature used for classification of the unclassified land types and the selected statistical feature of the typical sample by calculating the feature center of each land type, and classifying the unclassified land types into which land type according to the closer distance to the sample of which land type, wherein the calculation formula is as follows:
Figure BDA0001770433290000093
wherein d represents the distance between the sample land class and the land class o to be classified,
Figure BDA0001770433290000094
a feature value representing a class feature f of a typical sample,
Figure BDA0001770433290000095
characteristic value, sigma, representing a characteristic f of the ground class to be classifiedfThe standard deviation of the characteristic f is indicated.
And fourthly, generating a grading classification chart of soil nutrients in a plot scale.
And taking the farmland plots in Shapefile format as mask files, and masking the soil nutrient comprehensive evaluation graph generated based on the SNI index to form a plot scale soil nutrient grading classification graph. The method comprises the following specific steps:
(1) soil organic matter SOM, total nitrogen TN or alkaline hydrolysis nitrogen AN, quick-acting phosphorus AP and quick-acting potassium AK are selected to have 4 indexes.
(2) And obtaining the weight coefficient of each evaluation index of the soil nutrients. Here, taking the beijing city as an example, as shown in table 1, table 1 is a table of the classification of soil nutrients in the beijing city.
TABLE 1 Beijing City soil nutrient grading basis table
Figure BDA0001770433290000101
Data source: html of 20soil resource management information network of Beijing city, http://202.112.163.254: 8008/new%.
(3) In order to comprehensively evaluate the grading condition of soil nutrients in a certain area, a comprehensive evaluation index is obtained.
And solving the nutrient grade of each land by adopting a comprehensive index constructed based on an addition model, wherein the formula is as follows:
SNI=ΣFi×Wi(i=1,2,3,.......,n),
in the formula: SNI represents the soil nutrient comprehensive index of the plot, FiThe score value, W, of the i-th indexiAnd a weight coefficient representing the i-th index.
(4) And (3) in an ArcMap space analysis platform, substituting the pixel values of all nutrient indexes and the weight coefficients of all evaluation indexes of soil nutrients into an SNI formula to obtain a soil nutrient comprehensive evaluation classification diagram, and dividing scales according to the value domain grades of all evaluation indexes of the soil nutrients to obtain five classes of soil nutrient grade classification diagrams of high, medium, low and extremely low.
(5) And (4) obtaining a soil nutrient grade classification map of a plot scale by using the extracted farmland plot Shapefile file as a mask.
Here, a method for evaluating the accuracy of soil nutrient ranking, which is an accuracy evaluation method for soil nutrient ranking, is also provided with respect to a soil nutrient ranking map, and includes the following steps:
(1) obtaining corresponding geographic coordinates of the ground by using the partitioned soil nutrient grade data provided by the research area, and constructing a confusion matrix to obtain the total classification accuracy OA, wherein the formula is as follows:
Figure BDA0001770433290000111
in the formula, NaccurateNumber of pixels representing correct classification, which is the sum of diagonals of the confusion matrix, NtotalIs the total pixel number.
(2) Constructing a Kappa coefficient, wherein the formula is as follows:
Figure BDA0001770433290000112
in the formula, N represents the total number of pixels in the ground surface real classification; x is the number ofiiRepresents the value on the ith type diagonal in the confusion matrix; k represents the total classification number; x is the number ofi∑Representing the sum corresponding to the column of the ith class; x is the number of∑iAnd the corresponding sum of the rows where the ith class is located is shown.
Taking Beijing as an example, the north of the great plain of North China in Beijing is positioned at 39 degrees 38 'to 40 degrees 51' N,
between 115 deg. 25 'and 117 deg. 20' E. The landform consists of two units, namely a northwest mountain region and a southeast plain, and the general topography is represented as high in the northwest and low in the southeast, and inclines from the northwest to the southeast. The altitude of the plain is generally not more than 100m, most of the plain is 20-60 m, and the altitude of the mountain is generally 1000-1500 m. Beijing belongs to a warm-temperate zone semi-humid climate region, the annual average precipitation is 430.9mm, and the annual average temperature is 13.1 ℃. Wherein, the temperature is coldest in 1 month, and the average temperature is minus 3.9 ℃; in 7 months, the air temperature is hottest and the average air temperature is 26.5 ℃ (Beijing city soil resource management information network). Total city land area 16,410km2. Wherein the plain area is 6338km238.6 percent of the total area of the mountain area, 10,072km2Accounting for 61.4% (the national resource agency of Beijing City). According to the 2010 statistical yearbook, the Beijing municipal administrative district is divided into 14 districts and 2 counties (the end of 2010). Wherein, the gate head ditch area, the Huairou area, the valley area, the dense cloud county and the Yanqing county are defined as ecological conservation development areas and are the counties with the most distributed cultivated land areas. The agricultural land area of the whole city 2008 is 10,959.81km2Cultivated land is 2,316.88
km2Wherein the Yanqing county, the Huairou county and the dense cloud county are the first three counties with the most cultivated land and respectively account for 15.4%, 14.3% and 14.1% of the total area. The main crops in Beijing are wheat, corn and vegetables, and the sowing areas of the wheat, the corn and the vegetables in 2009 are 226,000ha, 151,000ha and 68,000ha respectively. The remote sensing data used for extracting the cultivated land plots in Beijing city is Landsat TM5, and two images (with the track numbers of p123/r32 and p123/r33) are needed for covering the whole research area. Obtaining three scenes with the track number p123/r32 in 2009: day 4, month 15, day 7, month 20 and day 9, month 22. Since the area of the p123/r33 scene in Beijing is very small, only the scene of 9-22 days in 2009 is collected. The Landsat TM image has a spatial resolution of 30m, a coverage of 185km, a revisit period of 16d, and includes 6 multispectral bands (0.45-2.35 μm) and one thermal infrared band (spatial resolution of 120m, 10.40-12.50 μm).
Before cultivated land plot extraction, data preprocessing needs to be respectively carried out on each scene image: radiometric calibration, atmospheric correction, and geometric refinement correction. The geometric correction adopts an automatic correction module Autosync of ERDAS 9.1, and the corrected Landsat ETM + in the same month is used as the reference image, and the RMSE is required to be controlled within 0.5 pixel. Then, digital photogrammetry and remote sensing processing software in ERDAS
Under an LPS (Leica photographic mapping suite) module, geometric corrected Landsat TM and ASTER GDEM elevation data with 30m resolution are simultaneously introduced for orthorectification, so that errors caused by mountain fluctuation in northwest of Beijing are eliminated. After geometric fine correction, the image is led into ENVI software
The FLAASH (Fast Light-of-sight Atmospheric Analysis of Spectral Hypercubes) Atmospheric correction module performs Atmospheric correction. Based on the collected soil sampling points and the grading standard of soil nutrient classification made by a soil fertilizer workstation in Beijing, the four indexes are subjected to spatial difference, and the output resolution is set to be 30m so as to be matched with the resolution of the Landsat TM remote sensing image. And (3) calculating the Spatial distribution map of the soil nutrients in the farmland in the research area in a grid calculator by utilizing a Spatial analysis module (Spatial analysis) of ArcMap. Finally, masking the extracted farmland plots to obtain a spatial distribution map of soil nutrient levels in Beijing City based on the dimensions of the plots, as shown in FIG. 4.
As shown in fig. 3a and 3b, fig. 3a is a pseudo color composite image of the wavelength bands of Landsat TM 7, 4 and 2 in beijing, and fig. 3b is a plot of cultivated land. As can be seen from fig. 3a and 3b, cultivated land in beijing city is mainly distributed in the southeast plains. Wherein, the arable land of the cis-districts, the Tongzhou districts and the great happy districts is distributed most densely and almost extends to the whole administrative district; due to the existence of mountainous terrain in the mountainous area, valley area, Chang Ping area, Miyun county and Yan Qing county, cultivated land is mainly concentrated in the region with flat terrain; compared with other counties with major mountain landforms, the counties in Yanqing counties have more cultivated lands and are mainly distributed in the southwest corner; the open-air area, the lake area, the stone landscape mountain area, the Fengtai area and the urban area have almost no cultivated land distribution.
For soil organic matter index: the content of the mountain area is the highest, and the average value is 23.34 g/kg; the content of the flat valley area is the lowest, and the average value is only 1.78 g/kg; meanwhile, the standard deviation of the soil organic matter in the mountainous area is the largest (23.69), which shows that the spatial dispersion degree of the index is larger than that of other counties; the minimum dispersion is the soil organic matter in the valley region, and the standard deviation is only 0.77. For the total nitrogen index: the highest content of the HuaiRou region is 1.04 g/kg; the content of the flat valley area is the lowest, and the average value is only 0.11 g/kg; the dispersion degree of the Huairou region, the dense cloud county, the Chang Ping region, the cis-meaning region, the Tong Zhou region and the Xing region is similar, the standard deviation is near 0.3, the dispersion degree of the valley region is the minimum, and the standard deviation is only 0.03. Among the fast-acting phosphorus indexes: the maximum content of the HuaiRou area is 50.08 g/kg; the content of the gully area of the door head is minimum, and the average value is only 17.23 g/kg; the degree of dispersion in the Tongzhou region is the most severe, and the standard deviation reaches 50.08; the quick-acting phosphorus index dispersion degree of the gully region of the gate head is minimum, and the standard deviation is 27.13. Among the rapid potassium indicators: the content of the gully area of the gate head is the highest, and the average value reaches 204.00 g/kg; the content of the great happy area is the lowest, and the average value is 110.03 g/kg; the dispersion of the index in the mountainous area is the largest, and the standard deviation reaches 167.60; the dispersion of the HuaiRou region is the smallest with a standard deviation of 63.60.
The foregoing shows and describes the general principles, essential features, and advantages of the invention. It will be understood by those skilled in the art that the present invention is not limited to the embodiments described above, which are merely illustrative of the principles of the invention, but that various changes and modifications may be made without departing from the spirit and scope of the invention, which fall within the scope of the invention as claimed. The scope of the invention is defined by the appended claims and equivalents thereof.

Claims (5)

1. A method for generating a grade classification map of soil nutrients based on a 3S technology is characterized by comprising the following steps:
11) generating a soil nutrient data space grid map, selecting a radial basis function RBF space difference method for carrying out space difference, and generating a grid map of each nutrient;
12) setting a grading classification map of soil nutrients, calculating to obtain a soil nutrient comprehensive index (SNI) map based on pixels according to a space difference map and a weight coefficient of each nutrient index, and formulating a classification rule according to a value range of the SNI to obtain five classes of grading classification maps of soil nutrients, namely a very high class, a medium class, a low class and a very low class;
13) extracting cultivated land parcels, namely extracting the cultivated land parcels and spatial distribution thereof by adopting a step-by-step elimination method of layered classification extraction based on multi-temporal remote sensing images of key growth periods of main crops in a research area to generate a Shapefile vector file of the cultivated land parcels; the cultivated land plot extraction method comprises the following steps:
131) selecting three-scene multi-temporal remote sensing images before planting, in a growth vigorous growing period and after harvesting according to the planting type and planting lunar calendar characteristics of main crops in a research area;
132) respectively carrying out data preprocessing on each scene image, including radiometric calibration, atmospheric correction and geometric fine correction;
133) in an eCoginization object-oriented classification software platform, performing hierarchical classification extraction combining a SEATH algorithm and nearest supervision classification on the preprocessed image, and excluding other land types to obtain the cultivated land type;
directly applying an SEATH algorithm to classify the two ground features with the distance close to 2 to generate a Shapefile of the cultivated land parcel;
for the distance between two ground objects is larger than 2, dividing the ground objects into large categories by using a user-defined normalized vegetation index; subdividing the land types under the large classification, re-analyzing the available classification feature space for the land types by using a nearest supervision classification method, then combining the classified large classification and the classified land types together to finish the classification of the cultivated land, and generating a Shapefile file of the cultivated land parcel;
14) and generating a soil nutrient grade classification map of a plot scale, taking the farmland plots in a Shapefile format as mask files, and masking the soil nutrient comprehensive evaluation map generated based on the SNI index to form the soil nutrient grade classification map of the plot scale.
2. The method for generating the grading classification map of soil nutrients based on the 3S technology as claimed in claim 1, wherein the generation of the spatial grid map of soil nutrients data comprises the following steps:
21) when a soil sample is collected in a farmland, a sub-meter high-precision handheld GPS is used for acquiring longitude and latitude coordinates of a sampling point, the longitude and latitude coordinates are led into an ArcMap space analysis platform, and a point vector file in a Shapefile format is generated;
22) based on an interpolation method classification tree provided in geostationary analysis, acquiring a soil nutrient value of an unsampled point by using a Radial Basis Function (RBF) difference method;
by utilizing sampling point data, through training a sample by an RBF network, hiding a nonlinear function relation between the soil nutrient value and plane coordinates x and y in a converged network, taking the geographic coordinates of unknown points as network input, and performing simulation prediction on the network to obtain the soil nutrient value of the unknown points, wherein the function expression of the function expression is as follows:
z=f(x,y),
wherein, (x, y) is the geographic coordinate of the sampling point or the prediction point, the RBF network takes x, y as network input and takes the corresponding soil nutrient value as network output, namely the number of nodes of the input layer of the network is 2, and the number of nodes of the output layer is 1;
23) and (3) performing space difference on each soil nutrient by adopting an RBF difference method on an ArcMap space analysis platform to obtain a corresponding space grid map.
3. The method for generating the soil nutrient grading map based on the 3S technology as claimed in claim 1, wherein the step of generating the soil nutrient grading map of the land parcel scale comprises the following steps:
31) selecting 4 indexes of soil organic matter SOM, total nitrogen TN or alkaline hydrolysis nitrogen AN, quick-acting phosphorus AP and quick-acting potassium AK;
32) obtaining the weight coefficient of each evaluation index of soil nutrients;
33) and (3) solving a comprehensive evaluation index, and solving the nutrient grade of each plot by adopting a comprehensive index constructed based on an addition model, wherein the formula is as follows:
SNI=∑Fi×Wi(i=1,2,3,.......,n),
in the formula: SNI represents the soil nutrient comprehensive index of the plot, FiThe score value, W, of the i-th indexiA weight coefficient indicating an i-th index;
34) in an ArcMap space analysis platform, pixel values of all nutrient indexes and weight coefficients of all evaluation indexes of soil nutrients are introduced into an SNI formula to obtain a soil nutrient comprehensive evaluation classification diagram, and scales are divided according to value domain grades of all evaluation indexes of the soil nutrients to obtain five classes of soil nutrient grade classification diagrams of high, medium, low and extremely low;
35) and (4) obtaining a soil nutrient grade classification map of a plot scale by using the extracted farmland plot Shapefile file as a mask.
4. The method for generating a classification map of soil nutrients based on 3S technology as claimed in claim 1, wherein the SEATH algorithm for classification includes the following steps:
41) the characteristics are preferably such that,
the degree of relevance between two land categories on a certain characteristic is judged by adopting the degree of separation, the degree of separation is calculated by the JM distance, so as to evaluate the separability between the land categories, and the formula of the JM distance is shown as follows:
J=2(1-e-B),
in the formula, J represents JM distance, B represents Papanicolaou distance, and the calculation formula is shown as follows:
Figure FDA0003117140200000031
wherein m is1And m2Representing the mean, sigma, of a feature of a representative sample selected from two classes of output1And σ2Then the variance of a certain feature in two classes of samples is represented; the value range of JM distance is [0, 2 ]]The closer the distance is to 2, the better the separability;
42) the determination of the characteristic threshold value is carried out,
calculating the optimal threshold values of two ground classes in a certain sample characteristic according to a Gaussian probability distribution formula:
p(x)=p(x|C1)p(C1)p(x|C2)p(C2),
wherein p (x | C)1) Represents land class C1Selected characteristic sample featuresValue obeys a mean value of m1Variance is σ1 2Normal distribution of (1), p (x | C)2) Represents land class C2The selected typical sample characteristic value obeys mean value m2Variance is σ2 2Normal distribution of (2);
the formula for the optimal threshold T is as follows:
Figure FDA0003117140200000041
in the formula (I), the compound is shown in the specification,
Figure FDA0003117140200000042
n1and n2Is class C of land1And C2The number of samples selected;
43) determining which type the object to be determined belongs to according to the selected typical feature sample and the nearest feature space distance constructed by calculation;
431) selecting a typical ground object sample;
432) constructing a feature space required by classification;
433) calculating the distance between the feature used for classification of the unclassified land types and the selected statistical feature of the typical sample by calculating the feature center of each land type, and classifying the unclassified land types into which land type according to the closer distance to the sample of which land type, wherein the calculation formula is as follows:
Figure FDA0003117140200000043
wherein d represents the distance between the sample land class and the land class o to be classified,
Figure FDA0003117140200000044
a feature value representing a class feature f of a typical sample,
Figure FDA0003117140200000045
indicating the place to be classifiedEigenvalues, σ, of class characteristics ffThe standard deviation of the characteristic f is indicated.
5. The method for generating the grading classification map of the soil nutrients based on the 3S technology as claimed in claim 1, further comprising a precision evaluation method for grading the soil nutrients, which comprises the following steps:
51) obtaining corresponding geographic coordinates of the ground by using the partitioned soil nutrient grade data provided by the research area, and constructing a confusion matrix to obtain the total classification accuracy OA, wherein the formula is as follows:
Figure FDA0003117140200000051
in the formula, NaccurateNumber of pixels representing correct classification, which is the sum of diagonals of the confusion matrix, NtotalThe total pixel number is;
52) constructing a Kappa coefficient, wherein the formula is as follows:
Figure FDA0003117140200000052
in the formula, N represents the total number of pixels in the ground surface real classification; x is the number ofiiRepresents the value on the ith type diagonal in the confusion matrix; k represents the total classification number; x is the number ofi∑Representing the sum corresponding to the column of the ith class; x is the number of∑iAnd the corresponding sum of the rows where the ith class is located is shown.
CN201810946788.3A 2018-08-20 2018-08-20 3S technology-based soil nutrient grade classification map generation method and precision evaluation method thereof Active CN108959661B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810946788.3A CN108959661B (en) 2018-08-20 2018-08-20 3S technology-based soil nutrient grade classification map generation method and precision evaluation method thereof

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810946788.3A CN108959661B (en) 2018-08-20 2018-08-20 3S technology-based soil nutrient grade classification map generation method and precision evaluation method thereof

Publications (2)

Publication Number Publication Date
CN108959661A CN108959661A (en) 2018-12-07
CN108959661B true CN108959661B (en) 2021-08-06

Family

ID=64470720

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810946788.3A Active CN108959661B (en) 2018-08-20 2018-08-20 3S technology-based soil nutrient grade classification map generation method and precision evaluation method thereof

Country Status (1)

Country Link
CN (1) CN108959661B (en)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110244018A (en) * 2019-03-28 2019-09-17 国智恒北斗好年景农业科技有限公司 A kind of soil moisture content remote real-time monitoring method based on 3S technology
CN110596008B (en) * 2019-09-06 2020-10-23 中国科学院空天信息创新研究院 Plot-based soil nutrient digital mapping method for agricultural region of Chinese Hongsheng plain
CN110658011A (en) * 2019-11-05 2020-01-07 新疆农业科学院土壤肥料与农业节水研究所(新疆维吾尔自治区新型肥料研究中心) County scale orchard soil quality sampling method
CN111289512B (en) * 2020-02-28 2021-04-13 中国水稻研究所 Rice grain alkali elimination value high-throughput determination method based on deep convolutional neural network
CN115344813B (en) * 2022-08-25 2023-07-11 珠江水利委员会珠江水利科学研究院 Mountain height inversion method based on shadows

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7047135B2 (en) * 2003-01-31 2006-05-16 Deere & Company Method and system of evaluating performance of a crop
CN103930919A (en) * 2011-10-24 2014-07-16 天宝导航有限公司 Agricultural and soil management
CN104459087A (en) * 2014-12-05 2015-03-25 广西壮族自治区林业科学研究院 Forest soil nutrient classification method
CN106888781A (en) * 2017-01-19 2017-06-27 交通运输部科学研究院 The table soil and the method for turf coordinating protection being combined with field investigation based on remote sensing subregion
US9904747B2 (en) * 2015-03-19 2018-02-27 Trimble Inc. Agricultural terrain forming based on soil modeling

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7047135B2 (en) * 2003-01-31 2006-05-16 Deere & Company Method and system of evaluating performance of a crop
CN103930919A (en) * 2011-10-24 2014-07-16 天宝导航有限公司 Agricultural and soil management
CN104459087A (en) * 2014-12-05 2015-03-25 广西壮族自治区林业科学研究院 Forest soil nutrient classification method
US9904747B2 (en) * 2015-03-19 2018-02-27 Trimble Inc. Agricultural terrain forming based on soil modeling
CN106888781A (en) * 2017-01-19 2017-06-27 交通运输部科学研究院 The table soil and the method for turf coordinating protection being combined with field investigation based on remote sensing subregion

Also Published As

Publication number Publication date
CN108959661A (en) 2018-12-07

Similar Documents

Publication Publication Date Title
CN108959661B (en) 3S technology-based soil nutrient grade classification map generation method and precision evaluation method thereof
Koc et al. Evaluating the cooling effects of green infrastructure: A systematic review of methods, indicators and data sources
Bagan et al. Landsat analysis of urban growth: How Tokyo became the world's largest megacity during the last 40 years
Tariq et al. Monitoring land use and land cover changes using geospatial techniques, a case study of Fateh Jang, Attock, Pakistan
Chen et al. A high-resolution monitoring approach of canopy urban heat island using a random forest model and multi-platform observations
Naboureh et al. An integrated object-based image analysis and CA-Markov model approach for modeling land use/land cover trends in the Sarab plain
Zhang et al. An automated early-season method to map winter wheat using time-series Sentinel-2 data: A case study of Shandong, China
Choudhary et al. Random Forest for rice yield mapping and prediction using Sentinel-2 data with Google Earth Engine
Khoshnoodmotlagh et al. Urban morphology detection and it's linking with land surface temperature: A case study for Tehran Metropolis, Iran
Halder et al. An assessment of urban expansion impacts on land transformation of Rajpur-Sonarpur Municipality
Ma et al. Retrieval of high spatial resolution precipitable water vapor maps using heterogeneous earth observation data
Meer et al. Remote sensing application for exploring changes in land-use and land-cover over a district in Northern India
CN112800827A (en) Hyperspectral image classification experimental method
Tang et al. A novel sample selection method for impervious surface area mapping using JL1-3B nighttime light and Sentinel-2 imagery
Garajeh et al. A comparative approach of data-driven split-window algorithms and MODIS products for land surface temperature retrieval
Ahmed et al. GIS-based land suitability mapping for rubber cultivation in Seremban, Malaysia
Rafiq et al. Estimation and validation of remotely sensed land surface temperature in Kashmir valley
Priya et al. Land use land cover representation through supervised machine learning methods: sensitivity on simulation of urban thunderstorms in the east coast of India
Sohail et al. Spatio-temporal analysis of land use dynamics and its potential implications on land surface temperature in lahore district, punjab, pakistan
Sang et al. Farmland distribution dataset of the Yarlung Zangbo–Lhasa–Nyangqu River region of the Tibetan Plateau
Hamzić et al. Structural Characteristics of Patches in the Central Lika Landscape–Application of Spatial and Regression Analysis
Dadgar The effect of land types and consequently land use on soil organic carbon content-case study: Damavand Region of Iran.
Wang et al. Mapping urban land dynamics by automatic generation of ground samples from Globeland30 and random forest classification on the Google Earth Engine
Hamzić et al. Strukturna obilježja uzoraka krajolika Srednje Like–primjena prostorne i regresijske analize
Gessesse et al. Decadal Dynamics of Land Use/Land Cover in the Muga Watershed, Choke Mountain Range, Upper Blue Nile Basin, Ethiopia

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