Chla monitoring method based on visible-near infrared spectrum and machine learning

Abstract

The invention discloses a Chla monitoring method based on visible-near infrared spectrum and machine learning, which is a modeling method for performing chlorophyll a concentration parameter inversion by using visible-near infrared hyperspectral data, and combines primary screening of spectral characteristic wave bands, training sample expansion based on GAN, spectral characteristic wave band fine screening based on CARS and a regression modeling technology. The problems of high dimension of hyperspectral data and insufficient samples are solved, the overall Chl-a modeling precision is improved, and the method is suitable for water eutrophication monitoring, marine ecological assessment and environment remote sensing application.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of remote sensing, and particularly relates to a Chla monitoring method based on visible-near infrared spectrum and machine learning. BACKGROUND

[0002] Water accounts for about 71% of the earth's surface, and plays an irreplaceable role in various cycles of natural ecological environment, human survival and development. Due to direct or indirect discharge of human life and industrial wastewater, some water bodies are facing eutrophication and pollution of toxic and harmful substances. Among them, water chlorophyll a (Chla) is an important parameter index for indicating the water quality of marine, river, lake and other water bodies, and it is of great significance and value to quickly and quickly obtain Chla for water eutrophication monitoring, pollutant tracing and other work.

[0003] The traditional water Chla monitoring adopts the laboratory analysis method, and its process is field sampling and laboratory analysis. For example, the laboratory determination Chla concentration analysis method includes spectrophotometry and fluorescence spectroscopy. The above laboratory determination Chla process is complicated, time-consuming, and cannot quickly obtain the required Chla value of the point, and the timeliness is poor. The online instrument measurement method is a method for directly measuring Chla in real time through a sensor, and the timeliness is good, but the current sensor is expensive and needs frequent maintenance, and it measures the Chla of a specific point of the water body. Since the Chla value in the water body has the characteristic of changing with space distribution, the Chla measured by the online instrument is not representative, so it is generally installed at the water inlet and outlet, river section and other positions, and cannot be widely laid in natural water bodies.

[0004] Water hyperspectral remote sensing monitoring technology is a kind of remote non-contact method for indirectly monitoring water quality by using the optical characteristics of the measured target, and has the ability to monitor a large range of water bodies, and has good timeliness. Hyperspectral data has continuous wave bands (usually hundreds of wave bands) and high resolution (for example, 1nm), which can capture the specific absorption and reflection characteristics of Chl-a in water, such as the absorption peaks in the blue-green light region (about 440nm) and the red light region (about 675nm), and the reflection platform in the near-infrared region. These characteristics make hyperspectral data have a significant advantage in Chl-a inversion. In the prior art, Chl-a inversion modeling methods mainly include empirical model, semi-empirical model and machine learning model. The empirical model such as water band ratio model uses the reflectance ratio of specific wave bands (such as R(700) / R(675)) to estimate Chl-a concentration; the semi-empirical model such as three-band model considers atmospheric correction and water optical characteristics; the machine learning model includes support vector machine (SVM), neural network (NN) and random forest (RF), etc. These models can handle nonlinear relationships and improve inversion accuracy.

[0005] However, the prior art has the following problems: 1. Insufficient feature band screening: There are many redundant and noisy bands in hyperspectral data, which leads to model overfitting or low computational efficiency. Although traditional screening methods such as principal component analysis (PCA) or stepwise regression are effective, they ignore the correlation between bands and the biological optical properties of Chl-a.

[0006] 2. Insufficient sample data: The number of hyperspectral samples collected in the field is limited, and insufficient training samples lead to poor model generalization ability. Existing data augmentation methods such as SMOTE (Synthetic Minority Over-sampling Technique) mainly target classification problems, and have limited effect on continuous spectral data augmentation for regression tasks.

[0007] 3. Poor model robustness: A single model is difficult to cope with the complexity of water optics (such as suspended solids interference and CDOM influence), and lacks a multi-stage optimization process, leading to a decline in inversion accuracy on the validation set.

[0008] 4. Balance between computational efficiency and accuracy: High-dimensional data processing requires efficient algorithms, but existing methods such as genetic algorithm screening bands are time-consuming, and do not combine generative adversarial networks (GAN) for spectral sample generation.

[0009] Therefore, it is necessary to develop a method with higher precision and stronger universality to improve the automation level of water Chl-a monitoring, thereby improving the automation level of water quality and environmental detection. SUMMARY

[0010] To solve the above technical problems, the present application provides a Chl-a monitoring method based on visible-near infrared spectrum and machine learning, which solves the problems of high dimensionality and insufficient samples of hyperspectral data, improves the overall Chl-a modeling accuracy, and is suitable for water eutrophication monitoring, marine ecological assessment and environmental remote sensing applications.

[0011] To achieve the above purpose, the present application provides a Chl-a monitoring method based on visible-near infrared spectrum and machine learning, comprising the following steps: S1: Data acquisition and preprocessing: S1.1: Obtain visible-near infrared hyperspectral data and Chl-a quantitative experimental data samples, and use the visible-near infrared hyperspectral data to screen data samples with complete spectral range of 350-900 nm, data quality flag ≠ 1, and Chl-a ≠ NaN, and obtain data samples D that can be modeled after data cleaning; S1.2: Data partitioning: Take 2 / 3 of the data sample D as the training set X, with a sample size of N1, and take the remaining 1 / 3 of the sample data D as the validation set, with a sample size of N2; S1.3: Spectral data preprocessing: Spectral smoothing is performed using the Savitzky-Golay method, with the window length set to 5 and the polynomial order set to 2. S1.4: Spectral Data Standardization: To eliminate the effects of scale differences, the values ​​of each band of all samples are standardized using the following formula:

[0012] in For band numbering, and These are the minimum and maximum observed values ​​among all samples in this band, respectively; S2: Initial screening using characteristic bands: S2.1: Calculate the Pearson correlation coefficient between each band and Chl-a concentration:

[0013] in: The reflectance of the k-th sample in the j-th band The Chl-a concentration of the k-th sample The mean of all samples in the j-th band (the mean of the training set) The mean of the labels for all training samples S2.2 on the training set The absolute values ​​of the Pearson correlation coefficients with Chl-a are obtained independently for all bands. The top 100 bands with the highest absolute values ​​of correlation coefficients are selected and denoted as B_coarse={λj1,…,λj100} as the initial feature set X1; S3: Training Sample Augmentation S3.1: Define the structure of the GAN generator G as a fully connected layer (FC) + ReLU activation, with 3 layers and 256 hidden units; define the structure of the discriminator D as FC + LeakyReLU, with 3 layers and 128 hidden units; S3.2 defines the loss function, where the loss function of the discriminator (D) is defined as:

[0014] The loss function of the generator (G) is defined as:

[0015] S3.3 Training Process: ①. Extract the coarse band dataset X1 from the training set, whose bands are the B_coarse selected in step S2; ②. Initialize the parameters of generator G and discriminator D: Use the Adam optimizer with a learning rate of 0.0002; ③. Iterative training, where epochs are 200 and batch_size is 32; ④ Generate synthetic sample X2: Using the trained G, input random noise z2 and c2 to generate 3N1 samples, and add these 3N1 samples to the original X1 to form synthetic sample X2; S4: CARS-based spectral characteristic band fine screening: S4.1: Set up Monte Carlo sampling, randomly sample from the synthetic sample X2 as the training set, where the number of sampling times N is 50, and each time m samples are randomly selected, where the number of m is 0.7-0.8 times the number of synthetic samples X2; S4.2 The number of bands retained in the next iteration Determined by the exponential decay function:

[0016] in To determine the number of bands to retain upon termination (e.g., r can be 10-25); S4.3 Adaptive weighting, calculating the weights of each band. And retain the first one according to probability Characteristic bands:

[0017] S4.4 For each CARS-selected band subset, construct a PLS regression model:

[0018] in To select the band matrix, This is the regression coefficient matrix. The RMSECV values ​​are calculated using the held-out sample.

[0019] S4.5 Select the band combination with the smallest RMSECV as the final selected feature band combination, B_fine={λ1,λ2…,λr}, and the training sample set at this time is... ; S5: Regression Modeling S5.1: Random Forest (RF) is chosen as the modeling method, utilizing mainstream implementation tools such as RandomForestRegressor from sklearn.ensemble. The main parameters and their recommended values ​​are as follows: number of decision trees (n_estimators) is 100–500, initially set to 200; maximum depth of each tree (max_depth) is 10–20, initially set to 15; minimum number of samples per leaf node (min_samples_leaf) is 2–5, initially set to 3; maximum number of features considered per split node (max_features) is 'sqrt' (i.e., ...). ,in This represents the number of spectral bands. S5.2 Cross-validation: Use GridSearchCV to find the optimal combination of parameters such as max_depth, min_samples_leaf, and max_features; S5.3: Model training: Using the training set X3, fit the optimal random forest regressor to obtain the final model. S5.4: Model Performance Evaluation: The root mean squared error (RMSE) and the coefficient of determination (R²) are used to evaluate model performance, as shown in the following formulas:

[0020]

[0021] in: To verify the observations of the i-th sample Chl-a in the set; The predicted value of the i-th sample Chl-a is obtained using the established RF model; The value of Chl-a is the average of all samples in the validation set; where i is the number of data points in the validation set.

[0022] Preferably, the visible-near-infrared hyperspectral data in step S1.1 is the GLORIA spectral reflectance dataset, which includes spectral reflectance data, data quality label data, and water quality data. The GLORIA spectral reflectance dataset is downloaded from the PANGAEA Data Publisher for Earth & Environmental Science website (https: / / doi.pangaea.de / 10.1594 / PANGAEA.948492).

[0023] Preferably, the iterative training in step S3.3 ④ includes the following steps: a: Sample true samples x0 (spectral data in training set X1) and corresponding chlorophyll c0 (Chl-a data in training set X1) from training set X1. b: Generate synthetic sample G(z1|c1), where z1 is random noise data with values ​​from 0 to 1, and c1 is random Chl-a data; c: Update the weights of discriminator D to maximize the discriminative loss; d: Update the weights of generator G to minimize the generation loss.

[0024] The beneficial effects of this invention are as follows: 1. Innovative Multi-Stage Integrated Framework: This invention integrates correlation coefficient-based band coarse screening, GAN training sample expansion, CARS band fine screening, and RF regression modeling into a closed-loop process, solving the problems of high dimensionality and insufficient sample size in hyperspectral data. Traditional methods are mostly single-stage, while this invention's multi-stage optimization improves the overall Chl-a modeling accuracy.

[0025] 2. Innovative robustness for water remote sensing: The method considers the optical complexity of water bodies, mitigates noise interference through data augmentation, and has strong nonlinear modeling capabilities through RF modeling. Compared with traditional models, it is more applicable to different water body types and has universality.

[0026] 3. Innovations in computational efficiency and practicality: Coarse screening reduces computational load, GAN is trained in parallel, CARS iteratively achieves high efficiency, and the entire process can be completed on a standard PC, facilitating real-time application of satellite hyperspectral data.

[0027] 4. Applicable to eutrophication monitoring, marine ecological assessment and environmental remote sensing applications. Attached Figure Description

[0028] Figure 1 This is a flowchart of the overall algorithm of the present invention.

[0029] Figure 2 Data sample in step S1.2 of Example 1 D A map showing the global distribution of [the data / locations].

[0030] Figure 3 Data sample in step S1.2 of Example 1 D The spectral curve of each data sample.

[0031] Figure 4 This is a correlation coefficient spectral curve from step S2 of Example 1.

[0032] Figure 5This is a comparison diagram of the average spectrum of the generated sample obtained by GAN expansion in step S3.3 of Example 1 and the average true spectrum in the preliminary feature set X1.

[0033] Figure 6 This is a graph showing the relationship between the number of band subsets selected and the number of iterations in each iteration of the CARS algorithm in step S4 of Example 1.

[0034] Figure 7 This is a bar chart showing the band selection results and the importance of characteristic bands in random forest modeling in step S4 of Example 1.

[0035] Figure 8 This is a comparison chart of the prediction accuracy of the present invention and traditional random forest modeling that only uses RF modeling. Detailed Implementation

[0036] The technical solution of the present invention will be further explained and described below with reference to the accompanying drawings and specific embodiments. It is worth noting that the following embodiments are only preferred embodiments of the present invention and should not be construed as limiting the present invention. The scope of protection of the present invention should be determined by the contents of the claims. Modifications and substitutions made by those skilled in the art to the technical solution of the present invention without creative effort all fall within the scope of protection of the present invention.

[0037] Example 1 S1: Data Acquisition and Preprocessing S1.1: Download the GLORIA global water remote sensing dataset from the PANGAEA Data Publisher for Earth & Environmental Science website (https: / / doi.pangaea.de / 10.1594 / PANGAEA.948492). This dataset contains visible-near-infrared spectral reflectance data, data quality labeling data, and Chl-a water quality data for lakes, rivers, and oceans worldwide (GLORIA - A globally representative hyperspectral in situ dataset for optical sensing of water quality[J]. Scientific Data, 2023, 10(1).DOI:10.1038 / s41597-023-01973-y.). S1.2: Using all data in the GLORIA dataset, select data samples with a complete spectral range of 350-900nm, data quality set with flag ≠ 1, and chlorophyll a Chl-a ≠ NaN. After data cleaning, define the data samples suitable for modeling as follows:D ( Figures 2-3 ) S1.3: Data Sample D Two-thirds of the data is used as the training set X for the model, with a sample size of 197; the remaining one-third of the data is used as the validation set T, with a sample size of 99. S1.4: Use the Savitzky-Golay method for spectral smoothing, with the window length set to 5 and the polynomial order set to 2; S1.5: To eliminate the impact of scale differences, the values ​​of each band of all samples are standardized, as shown in the following formula:

[0038] in For band numbering, and These are the minimum and maximum observed values ​​among all samples in this band, respectively; S2: Initial screening using characteristic bands: S2.1: Take the standardized data for each band and calculate the Pearson correlation coefficient with the Chl-a concentration data in the training set sample X:

[0039] in: The reflectance of the k-th sample in the j-th band. Chl-a concentration of the kth sample : The mean of all samples in band j (mean of the training set) Mean of labels for all training samples; S2.2: Take all bands (551 bands) from the training set X, and calculate the absolute value of their Pearson correlation coefficients with Chl-a. The top 100 bands with the highest absolute values ​​of correlation coefficients are selected and denoted as B_coarse={λj1,…,λj100}, which serve as the initial feature set X1. Figure 4 ); S3: GAN-based training sample augmentation; S3.1: Define the structure of the GAN generator G as a fully connected layer (FC) + ReLU activation, with 3 layers and 256 hidden units; define the structure of the discriminator D as FC + LeakyReLU, with 3 layers and 128 hidden units; S3.2 defines the loss function, where the loss function of the discriminator (D) is defined as:

[0040] The loss function of the generator (G) is defined as: ; S3.3 Training Process: ① Extract the preliminary feature set X1; ② The parameters of the generator G and discriminator D are initialized using the Adam optimizer, with a learning rate of 0.0002; ③ Set epochs=200 and batch_size=32 to perform iterative training: a: In each iteration, 32 true spectral samples x0 and their corresponding c0 (Chla) values ​​are randomly selected from X1 and input into the initialized generator G; b: Generate a synthetic sample G(z1|c1) using the generator G defined in step S3.1; where The generated 0-1 random noise data is used to simulate the spectrum. The simulated Chl-a value corresponding to each simulated spectral sample; c: Update the weights of discriminator D to maximize the discriminative loss; d: Update the weights of generator G to minimize the generation loss; ④ Using the trained G, input random noise z2 and c2 to generate 3N1 samples. Add these 3N1 generated samples to X1 to obtain the synthetic sample X2 after GAN sample augmentation. Figure 5 ); S4: Selection of spectral characteristic bands based on CARS; S4.1: The Monte Carlo sampling method is used to sample from the synthetic sample X2, where the number of sampling times N is 50, and 158 samples are randomly selected each time to form the training set; the remaining 39 samples are used as the hold-out samples. S4.2: Determining the first exponential decay function The number of bands retained in the next iteration The formula is as follows: , in The number of bands to retain upon termination is set to 0.95. The value of k ranges from 1 to 50. S4.3 Adaptive weighting, calculating the weights of each band. And retain the first one according to probability Characteristic bands ( Figure 6 ): , in Let these be the PLS regression coefficients, and each iteration of training will generate a... ; S4.4 uses CARS after each training iteration to select a subset of bands, i.e., the first... Using characteristic bands, construct a PLS regression model Y:

[0041] in To select the band matrix, Here is the regression coefficient matrix, and E is the residual; the RMSECV value is calculated using the held-out sample.

[0042] S4.5 Selects the band combination with the smallest RMSECV (707nm, 708nm, 711nm, 714nm) as the final selected characteristic band combination, B_fine={707nm, 708nm, 711nm, 714nm}. Figure 7 ), which is the training set sample X3; S5: Regression modeling based on RF; S5.1 uses RandomForestRegressor (RF) as the modeling method, calling mainstream implementation tools such as RandomForestRegressor under sklearn.ensemble. The main parameters and their recommended values ​​are as follows: number of decision trees n_estimators (recommended value: 200), maximum depth of each tree max_depth (recommended value: 15); minimum number of samples per leaf node min_samples_leaf (recommended value: 3); maximum number of features considered per split node max_features, with a recommended value of 2 (of which the number of bands after filtering is 4). S5.2 Model Training: Input the training set X3 into the RF prediction model, use 5-fold cross-validation, and use GridSearchCV to find the optimal combination of parameters: n_estimators=89, max_depth=7, min_samples_leaf=2, max_features=1, to obtain the final modeling model—the Random Forest model. S5.3 Model Performance Evaluation: The root mean squared error (RMSE) and the coefficient of determination (R²) are used to evaluate model performance. The formula for the root mean squared error is as follows: , The formula for the coefficient of determination is: , Use a validation set to evaluate model accuracy, such as Figure 8As shown, the results show that the root mean square error (RMSE) of the model is 11.891 mg / L and the coefficient of determination (R²) is 0.886, indicating that the model has high accuracy. Compared with the model that only uses RF for Chl-a prediction (Leo B. Random Forests[J].Machine Learning, 2001, 45(1):5-32.), the RMSE is reduced by 2.226, indicating that the model has high accuracy.

Claims

1. A Chla monitoring method based on visible-near-infrared spectroscopy and machine learning, characterized in that: Includes the following steps: S1: Data acquisition and preprocessing: Collect visible-near-infrared hyperspectral data and Chl-a quantitative experimental data samples within the monitoring area, and after screening, divide them into training set samples X and validation set samples T. After smoothing, the bands are standardized. S2: Initial screening of characteristic bands: Calculate the Pearson correlation coefficient between the band and Chl-a concentration to screen the preliminary feature set X1; S3: Training sample augmentation: Based on the GAN generator and discriminator, define the loss function, use the initial feature set X1 to train, and generate synthetic samples X2; S4: Feature band fine screening: After sampling the synthetic sample X2 based on CARS, a PLS regression model is constructed through exponential decay function and adaptive weighting, and the band combination with the smallest RMSECV is selected as the feature band combination X3 for fine screening. S5: Regression Modeling: RF modeling is used, and the final model is determined by training with feature band combination X3, and performance is evaluated.

2. The Chla monitoring method based on visible-near-infrared spectroscopy and machine learning according to claim 1, characterized in that: The visible-near-infrared hyperspectral data mentioned in step S1 includes spectral reflectance data, data quality labeling data, and water quality data; the smoothing method is Savitzky-Golay; the band normalization formula is... ; In the formula It is the minimum observed value among all samples in this band. This is the maximum observed value among all samples in this band.

3. The Chla monitoring method based on visible-near-infrared spectroscopy and machine learning according to claim 1, characterized in that: The filtered data samples described in step S1 have a complete spectral range of 350-900nm, data quality concentration flag≠1, and chlorophyll a Chl-a≠NaN.

4. The Chla monitoring method based on visible-near-infrared spectroscopy and machine learning according to claim 1, characterized in that: The formula for calculating the Pearson correlation coefficient in step S2 is as follows: , In the formula The reflectance of the k-th sample in the j-th band. Chl-a concentration of the kth sample : The mean of all samples in band j (mean of the training set) Mean of labels for all training samples.

5. The Chla monitoring method based on visible-near-infrared spectroscopy and machine learning according to claim 1, characterized in that: The screening method described in step S2 is to independently calculate the absolute value of the Pearson correlation coefficient for all bands in the training set X, and select the top 100 bands as the preliminary feature set X1.

6. The Chla monitoring method based on visible-near-infrared spectroscopy and machine learning according to claim 1, characterized in that: The loss function of the GAN generator described in step S3 is: ; The loss function of the GAN discriminator is: .

7. The Chla monitoring method based on visible-near-infrared spectroscopy and machine learning according to claim 1, characterized in that: The training described in step S3 includes the following steps: ① Initialize the parameters of the generator and discriminator; ② Iterative training: a: Sample the training set to obtain the true spectrum x0 and the corresponding chlorophyll c0; b: Generate a synthetic sample G; c: Update the weights of discriminator D to maximize the discriminative loss; d: Update the weights of generator G to minimize the generation loss; e: Using the trained G, input random noise z2 and c2 to generate 3N1 samples, and add these 3N1 samples to the original X1 to form the synthetic sample X2.

8. A Chla monitoring method based on visible-near-infrared spectroscopy and machine learning according to claim 1, characterized in that: The sampling method described in step S4 is Monte Carlo sampling, in which m samples are randomly selected during sampling, and the value of m is 0.7-0.8 times the number of synthetic samples X2.

9. The Chla monitoring method based on visible-near-infrared spectroscopy and machine learning according to claim 1, characterized in that: The exponential decay function mentioned in step S4 is In the formula The number of bands to retain upon termination is 10-25. The formula for adaptive weighting is: In the formula The PLS regression coefficients are denoted as , where the adaptive weighting retains the coefficients before and after the regression. One characteristic band; The PLS regression model is as follows: In the formula To select the band matrix, This is the regression coefficient matrix; The formula for calculating RMSECV is as follows: .

10. The Chla monitoring method based on visible-near-infrared spectroscopy and machine learning according to claim 1, characterized in that: In the RF modeling described in step S5, the number of decision trees is 100-500, the maximum depth of each tree is 10-20, the minimum number of samples per leaf node is 2-5, and the maximum number of features considered for a single split node is [missing information]. In the formula The number of spectral bands; the performance evaluation indicators include root mean square error (RMSE) and coefficient of determination (R²). 2 .

Images