A multi-time series water eutrophication assessment method based on satellite images, electronic equipment and storage medium

By designing the DSAT_1D_ShuffleNetV2 module and the Shuffle_GRUnet model, combined with the dual-weight logarithmic root mean square error loss function Wwrmse, the problems of model generalization ability and stability in the evaluation of single-period satellite image data were solved, and efficient and accurate assessment of multi-time series water eutrophication was achieved.

CN119625554BActive Publication Date: 2025-09-09HARBIN AEROSPACE STAR DATA SYST TECH CO LTD
View PDF 1 Cites 0 Cited by

Patent Information

Application Number
CN202411663237.8
Authority / Receiving Office
CN · China
Patent Type
Patents(China)
Current Assignee / Owner
Filing Date
2024-11-20
Publication Date
2025-09-09
Estimated Expiration
2044-11-20

AI Technical Summary

Technical Problem

The existing technology for water eutrophication assessment based on single-period satellite image data has poor generalization ability and stability, the accuracy cannot meet application requirements, and it is unable to effectively associate the mutual constraints between data from multiple periods.

Method used

A multi-temporal water eutrophication assessment method based on satellite images was designed. The DSAT_1D_ShuffleNetV2 module was used for feature extraction, and the Shuffle_GRUnet model and the double-weighted logarithmic root mean square error loss function Wwrmse were combined to realize multi-temporal water eutrophication assessment.

Benefits of technology

It has achieved efficient and accurate multi-period water eutrophication assessment, which can monitor the eutrophication status of water bodies faster and more comprehensively, and provide a scientific basis for formulating lake eutrophication control measures.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN119625554B_ABST
    Figure CN119625554B_ABST
Patent Text Reader

Abstract

A multi-time series water eutrophication assessment method, electronic device, and storage medium based on satellite imagery belong to the field of satellite remote sensing quantitative inversion technology. To address the problem of accuracy failing to meet application requirements, the present invention designs a shuffling module and constructs a DSAT_1D_ShuffleNetV2 module, implementing channel shuffling, multi-scale feature extraction, deep and shallow feature fusion, and multi-attention information sharing, effectively mining feature information and more accurately guiding model feature information collection. A dual-weighted logarithmic root mean square error loss function (Wwrmse) is designed to guide the optimization of a multi-time series water eutrophication assessment model. The designed Shuffle_GRUnet model can complete model training related to multi-period data, realizing single-period data prediction and data prediction based on multi-period linkage. The method of the present invention stably and accurately implements multi-time series water eutrophication assessment using satellite remote sensing images.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] The present invention belongs to the technical field of satellite remote sensing quantitative inversion, and in particular relates to a multi-time series water eutrophication assessment method based on satellite images, an electronic device and a storage medium. Background Art

[0002] Eutrophication of lakes has become a global environmental issue. Eutrophication affects aquatic ecosystems worldwide. The accumulation of nitrogen and phosphorus in lakes due to runoff from farmland and the discharge of industrial and domestic wastewater has garnered widespread attention worldwide. However, human impact on aquatic ecosystems continues to increase at an unprecedented rate. Timely, accurate, and comprehensive assessment and understanding of the status of lake eutrophication and monitoring its changing trends are crucial for relevant departments to formulate appropriate countermeasures, mitigate and control lake eutrophication, and fully utilize the various functions of lake water bodies.

[0003] Currently, the conventional method for monitoring water quality involves collecting water samples at sampling points within a lake. These samples are then analyzed in a laboratory and assessed using either a single-parameter evaluation index or a multi-parameter comprehensive evaluation method. While this method can accurately analyze and evaluate numerous water quality indicators, it is time-consuming, labor-intensive, and uneconomical. Furthermore, the number of water samples collected and analyzed is very limited, and the data from these sampling points is only partially representative of the entire lake.

[0004] The emergence and development of remote sensing technology has opened up new avenues for monitoring and studying lake eutrophication. Using satellite or aerial remote sensing information and based on water color remote sensing theory, water quality remote sensing monitoring applies empirical, semi-empirical, or theoretical analysis methods, appropriate remote sensing band data is selected to establish a quantitative remote sensing inversion model for water quality parameters to invert the concentrations of water quality parameters in the water body. Remote sensing methods can be used to conduct large-scale eutrophication surveys and assessments. Traditional statistical relationship models based directly on the relationship between single-scene remote sensing data and measured water quality results can only fuzzily calculate water quality results, resulting in poor monitoring results. Furthermore, the national lake water quality model used in existing lake eutrophication monitoring methods cannot accurately adapt to the water quality monitoring of each lake. Furthermore, this poor adaptability leads to poor monitoring results and unsatisfactory accuracy.

[0005] With the comprehensive development of technologies combining artificial intelligence with remote sensing satellites, higher requirements have been placed on the automation, intelligence, and precision of water eutrophication assessment. Moreover, existing technologies only consider the correlation between single periods and cannot effectively associate the mutual constraints between data from multiple periods, resulting in poor generalization and stability of the model and accuracy that cannot meet application requirements. Summary of the Invention

[0006] The present invention aims to solve the problems of poor model generalization ability and stability and inability to meet application requirements caused by water eutrophication assessment based on single-period satellite image data, and proposes a multi-time series water eutrophication assessment method based on satellite imagery, an electronic device and a storage medium.

[0007] To achieve the above object, the present invention is implemented through the following technical solutions:

[0008] A multi-temporal water eutrophication assessment method based on satellite images includes the following steps:

[0009] S1. Collect paddy field sample data from the first, second, third, and fourth quarters, and use the collected sample data to calculate the nutritional status index of each sample point based on the comprehensive nutritional status index calculation formula, which serves as the target parameter of the model;

[0010] S2. Obtain and preprocess multispectral satellite imagery data that is close in time to the sample data. Calculate the spectral index of each sample point using the preprocessed satellite imagery data based on the spectral index calculation formula, which serves as the input parameter for the model.

[0011] S3. Design shuffle modules shuffle_block1 and shuffle_block2 respectively, and build DSAT_1D_ShuffleNetV2 module based on the designed shuffle modules for feature data extraction;

[0012] S4. Design a dual-weighted logarithmic root mean square error loss function (Wwrmse) to guide the optimization of multi-time series water eutrophication assessment models.

[0013] S5. Based on the DSAT_1D_ShuffleNetV2 module designed in S3 and the Wwrmse loss function designed in step S4, combined with a gated recurrent neural network, design a multi-time series water eutrophication assessment model Shuffle_GRUnet;

[0014] S6. Input the target parameters obtained in step S1 and the input parameters obtained in step S2 into the multi-time series water eutrophication assessment model Shuffle_GRUnet designed in step S5 for model training to obtain the optimal multi-time series water eutrophication assessment model;

[0015] S7. Using the optimal multi-time series water eutrophication assessment model obtained in step S6, perform a multi-time series water eutrophication assessment in the target area.

[0016] Furthermore, the sample data of step S1 includes chlorophyll a data chl_a, total phosphorus data TP, total nitrogen data TN, transparency data SD, and chemical oxygen demand data CODmn.

[0017] Furthermore, the calculation formula of the spectral index in different bands of step S2 is as follows:

[0018] SPI b4,b1 =b4 / b1 (1)

[0019] SPI b4,b2 =b4 / b2 (2)

[0020] SPI b4,b3 =b4 / b3 (3)

[0021] SPI b3,b4,b1 =(b3+b4) / b1 (4)

[0022] SPI b3,b4,b2 =(b3+b4) / b2 (5)

[0023] SPI b1,b2,b4,b3 =(b1+b2+b4) / b3 (6)

[0024] SPI b2,b3,b4,b1 =(b2+b3+b4) / b1 (7)

[0025] SPI b1,b3,b4,b2 =(b1+b3+b4) / b2 (8)

[0026] SPI b3,bb4,b1,b2 =(b3+b4) / (b1+b2) (9)

[0027] Where SPI represents the spectral index, b1 represents the blue band, b2 represents the green band, b3 represents the infrared band, and b4 represents the near-infrared band.

[0028] Furthermore, the specific implementation method of step S3 includes the following steps:

[0029] S3.1. Design the shuffle_block1 module, which consists of 1 channel segmentation module Channel_Split, 5 batch normalization modules BN, 2 improved 1D multi-scale convolution attention MSCA_1D, and 3 1D convolution Conv1D 1,1 , 1 depth 1D convolution DWConv1D 3,1 , 3 Relu activation functions and 1 channel shuffle module Channel_Shuffle. The calculation process expression is as follows:

[0030] [in x1 ,in x2 =Channel_Split(in x ,0.5) (10)

[0031]

[0032] in x2_1 =Conv1D 1,1 (BN(DWConv1D 3,1 (Relu(BN(Conv1D 1,1 (in x2 )))))) (12)

[0033]

[0034] out y1 =Channel_Shuffle(Concat(in x1_1 ,in x2_2 )) (14)

[0035] Among them, x Indicates the input feature information of the shuffle_block1 module, in x1 ,in x2 Respectively represent the channel information after segmentation according to the first channel and the second channel, in x1_1 Indicates in x1 After MSCA_1D operation and in x1 The result of multiplication, in x2_1 Indicates in x2 After Conv1D 1,1 , BN, Relu, DWConv1D 3,1 Feature information after operation, in x2_2 Indicates in x2_1 After BN and MSCA_1D operations and in x2_1 The result of the multiplication is then passed through Conv1D 1,1 , feature information after BN and Relu operations, out y1 Indicates in x1_1 ,in x2_2 After matrix connection and Channel_Shuffle operation, the feature information is obtained. Channel_Split(tz,bl) indicates that the feature information tz is split into channels according to the bl ratio. MSCA_1D(.) indicates the designed 1D multi-scale convolution attention module. BN(.) indicates the batch normalization operation. Relu(.) indicates the activation function operation. Conv1D 1,1 (.) represents a 1D convolution operation with a convolution kernel of 1 and a stride of 1, DWConv1D 3,1(.) represents a depth-1 convolution operation with a convolution kernel of 3 and a stride of 1, Concat(.) represents a matrix concatenation operation, and Channel_Shuffle(.) represents a channel shuffle operation. Represents matrix multiplication;

[0036] S3.2. The designed shuffle_block2 module consists of 6 batch normalization modules BN, 2 improved 1D multi-scale convolution attention MSCA_1D, and 4 1D convolution Conv1D 1,1 , 2 depth 1D convolution DWConv1D 3,2 , 3 Relu activation functions and 1 channel shuffle module Channel_Shuffle. The calculation process expression is as follows:

[0037] in 4_1 =DWConv1D 3,2 (in x4 )(15)

[0038]

[0039] in x4_3 =Conv1D 1,1 (BN(DWConv1D 3,2 (Relu(BN(Conv1D 1,1 (in x4 ))))))(17)

[0040]

[0041] out y2 =Channel_Shuffle(Concat(in x4_2 ,in x4_4 ))(19)

[0042] Among them, x4 Indicates the input feature information of the shuffle_block2 module, in x4_1 Indicates that after DWConv1D 3,2 Feature information after operation, in x4_2 Indicates in x4_1 After BN, MSCA_1D, Conv1D 1,1 The result of the operation is then combined with in x4_1 After multiplication, the feature information after BN and Relu operation is obtained. x4_3 Indicates in x4 After Conv1D 1,1 , BN, Relu, DWConv1D 3,2 、Conv1D1,1 Feature information after operation, inx 4_4 Indicates in x4_3 After BN, MSCA_1D, Conv1D 1,1 The result of the operation is then in x4_3 After multiplication, the feature information after BN and Relu operation is out y2 Indicates in x4_2 ,in x4_4 The feature information after matrix connection and Channel_Shuffle operation, MSCA_1D (.) represents the 1-dimensional multi-scale convolution attention module designed by the present invention, BN (.) represents the batch normalization operation, Relu (.) represents the activation function operation, Conv1D 1,1 (.) represents a 1D convolution operation with a convolution kernel of 1 and a stride of 1, DWConv1D 3,2 (.) represents a depth-1 convolution operation with a kernel size of 3 and a stride of 2. Concat(.) represents a matrix concatenation operation. Channel_Shuffle(.) represents a channel shuffle operation. Represents matrix multiplication;

[0043] S3.3. The designed 1D multi-scale convolution attention module MSCA_1D consists of 4 1D convolutions of different scales and depths DWConv1D, 1 1D convolution Conv1D 1,1 The calculation process is as follows:

[0044] in x5_1 =DWConv1D5(in x5 )(20)

[0045]

[0046] Among them, x5 Represents the input feature information of the 1D multi-scale convolutional attention module MSCA_1D, in x5_1 Indicates in x5 Feature information after DWConv1D5 operation, in x5_2 ,in x5_3 ,in x5_4 Respectively represent in x5_1 After DWConv1D7, DWConv1D 11 、DWConv1D 21 Feature information after operation, out y3 Indicates in x5_1 ,in x5_2 ,in x5_3 ,in x5_4The result after matrix addition and Conv1D1 operation is the same as in x5 Multiplied feature information, DWConv1D q Represents the depth 1D convolution operation of q convolution kernel, q∈[1,5,7,11,21], Conv1D 1,1 The 1 that indicates the convolution kernel is 1 and the step size is 1 is the convolution operation, and Add indicates matrix addition. Represents matrix multiplication;

[0047] S3.4. The calculation process of Channel_Shuffle is as follows:

[0048] S3.4.1. Determine the number of feature map x groups and the number of channels C per group g The calculation formula is as follows:

[0049]

[0050] Among them, C represents the number of channels of the feature map x, G represents the number of groups to be divided, and C g Indicates the number of channels in each group;

[0051] S3.4.2.reshape operation:

[0052] Expand the channel dimension of the feature map x from one dimension to two dimensions, one is the number of groups G, and the other is the number of channels in each group C g ; If the batch size of the feature map is N and the spatial dimension is H×W, the reshape operation can be expressed as:

[0053] x reshaped =x.view(N,G,C g ,H,W) (24)

[0054] S3.4.3.transpose operation:

[0055] The reshaped feature map is transposed in terms of the number of groups and the number of channels in each group, so that the channels in each group can be disrupted. The transposition operation can be expressed as:

[0056] x transposed =x reshaped .transpose(1,2)(25)

[0057] Among them, transpose(1,2) means exchanging the second dimension and the third dimension;

[0058] S3.4.4.flatten operation:

[0059] The transposed feature map xtransposed , flatten back to the original number of channels C to obtain the final channel shuffled feature map. The flatten operation is expressed as:

[0060] x shuffled =x transposed .view(N,C,H,W)(26)

[0061] Among them, view(N,C,H,W) restores the feature map to its original four-dimensional shape;

[0062] S3.5. The designed DSAT_1D_ShuffleNetV2 consists of 1 1D convolution Conv1D 3,2 , 1 1D convolution Conv1D 1,1 , 1 1-dimensional maximum pooling MaxPool 3,2 , 1 global average pooling GlobaPool, 3 designed shuffle_block2, 13 designed shuffle_block1, the calculation process expression is as follows:

[0063] Stage 0 = MaxPool 3,2 (Conv1D 3,2 (input x )) (27)

[0064]

[0065] Stage 2_1 =Concat(interpolate(Stage1,Stage2.shape()),Stage2)(29)

[0066] Stage 3_1 =Concat(interpolate(Stage 2_1 ,Stage3.shape()),Stage3)(30)

[0067] output y =GlobaPool(Conv1D 1,1 (Stage 3_1 )) (31)

[0068] Among them, input x Indicates the input feature information of the designed DSAT_1D_ShuffleNetV2, Stage0 represents input x After Conv1D 3,2 MaxPool 3,2The feature information after the operation, Stage1, Stage2, Stage3 respectively represent the feature information after Stage0 is operated by shuffle_block2 and shuffle_block1, Stage 2_1 It indicates the characteristic information obtained by interpolating the result of Stage1 based on the size of Stage2 and connecting it with the matrix of Stage2. 3_1 Represents Stage 2_1 After interpolation operation based on the Stage3 size, the interpolation result is connected with the Stage3 matrix to obtain the feature information, output y Represents Stage 3_1 After Conv1D 1,1 , feature information of GlobaPool operation, Conv1D 1,1 Represents a 1D convolution operation with a convolution kernel of 1 and a stride of 1, Conv1D 3,2 (.) represents a 1D convolution operation with a kernel size of 3 and a stride of 2. MaxPool 3,2 represents a 1D maximum pooling operation with a convolution kernel of 3 and a stride of 2. shuffle_block2(,v) represents a shuffle_block2 module operation designed v times, v∈[1,1,1], shuffle_block1(,w) represents a shuffle_block1 module operation designed w times, w∈[3,7,3], interpolate(.) represents an interpolation operation, Concat(.) represents a matrix concatenation operation, and GlobaPool(.) represents a global average pooling operation.

[0069] Furthermore, the calculation formula of the dual-weight logarithmic root mean square error loss function Wwrmse designed in step S4 is as follows:

[0070]

[0071] Among them, i represents the i-th period, i.e., the first quarter, the second quarter, the third quarter, and the fourth quarter; n represents the total number of periods; j represents the j-th sample data; m represents the total number of samples; and W i represents the weight of the i-th period, wrmse represents the weight logarithmic root mean square error designed by the present invention, represents the forecast value of period i, represents the true value of period i, represents the average value of the true value in period i, It represents the maximum value of the true value in the i-th period, cov represents the covariance operation, and σ represents the standard deviation operation.

[0072] Furthermore, the specific implementation method of step S5 includes the following steps:

[0073] The Shuffle_GRUnet model designed in S5.1. consists of 1 input layer, 1 DSAT_1D_ShuffleNetV2 layer designed in step S3, 1 gated recurrent neural network layer, and 1 fully connected layer. Each layer contains n components; the input layer X = {1st period, 2nd period, 3rd period, 4th period...nth period}, and each component input parameter contains 9 steps 2 to calculate the spectral index (SPI b4,b1 、SPI b4,b2 、SPI b4,b3 、SPI b3,b4,b1 、SPI b3,b4,b1 、SPI b1,b2,b4,b3 、SPI b2,b3,b4,b1 、SPI b1,b3,b4,b2 、SPI b3,b4,b1,b2 ), input in the format of a three-dimensional tensor, so that DSAT_1D_ShuffleNetV2 and GRU can be used for feature extraction and model training. The three-dimensional tensor used is (number of samples × time steps × number of features), the number of features is 9, and the spectral index stacking method is used to expand it to 224. The insufficient part is supplemented with 0. The calculation process expression is as follows:

[0074]

[0075] Among them, X i represents the input parameters of period i, Indicates the output parameter of period i, namely the comprehensive nutritional status index, DSAT_1D_ShuffleNetV2 i (.) indicates the DSAT_1D_ShuffleNetV2 module operation designed in step 3 in the i-th period, GRU i (.) indicates that the GRU module operation is performed in the i-th period, fc i (.) indicates that the i-th period has undergone a full connection operation;

[0076] S5.2. The referenced GRU consists of an update gate and a reset gate. First, the update gate Z is calculated. t and reset gate r t The value of the input variable x t and the hidden layer result h at the previous moment t-1 The concatenated matrix is ​​input into the update gate after sigmoid nonlinear transformation; the reset gate value will act on h t-1 and the input variable x t The splicing is performed nonlinearly, and then activated by the tanh function to obtain the state of the candidate hidden layer at the current moment By 1-Z ttimes h t-1 Store the information of the previous moment through Z t times Record the information at the current moment and add the two results as the state output h of the hidden layer at the current moment t , and finally the weight matrix of the output layer acts on h t On the top, get the current model output y t , the calculation formula is as follows:

[0077] Z t =sigmoid(W Z [h t-1 ,x t ]) (36)

[0078] r t =sigmoid(W r [h t-1 ,x t ]) (37)

[0079]

[0080] y t =sigmoid(W o h t ) (40)

[0081] Among them, sigmoid represents the sigmoid activation function operation, tanh represents the tanh activation function operation, Z t Indicates that the door state is updated at time t, r t Indicates that the door state is reset at time t, represents the candidate hidden door state at time t, h t represents the hidden layer state at time t, y t represents the model output at time t, x t represents the model input at time t, h t-1 Indicates the hidden layer state at the previous moment, x t , W Z Represents the update gate weight matrix, W r Represents resetting the gate weight matrix, represents the candidate hidden layer weight matrix, W o Represents the output layer weight matrix.

[0082] Furthermore, the implementation method of step S6 is to input the target parameters prepared in step S1 and the input parameters prepared in step S2 into the multi-time series water eutrophication assessment model Shuffle_GRUnet designed in step S5 for model training, and use the Pearson correlation coefficient (r pcc ), mean root mean square error The model accuracy is evaluated to obtain the optimal multi-time series water eutrophication assessment model.

[0083] An electronic device includes a memory and a processor. The memory stores a computer program. When the processor executes the computer program, the steps of the multi-time series water eutrophication assessment method based on satellite images are implemented.

[0084] A computer-readable storage medium stores a computer program, which, when executed by a processor, implements a multi-time series water eutrophication assessment method based on satellite images.

[0085] Beneficial effects of the present invention:

[0086] The present invention discloses a multi-time series water eutrophication assessment method based on satellite imagery. The method designs a DSAT_1D_ShuffleNetV2 feature extraction module, which implements channel shuffling, multi-scale feature extraction, deep and shallow feature fusion, and multi-scale attention information mining. This module achieves multi-dimensional and multi-scale context information acquisition and sharing, effectively retains high-quality information, and more accurately guides model feature information collection. The Shuffle_GRUnet model is designed to implement multi-period-related model training, enabling single-period data prediction and data prediction based on multi-period linkage. A dual-weighted logarithmic root mean square error loss function (Wwrmse) is designed to guide the optimization of the multi-time series water eutrophication assessment model. The method can efficiently and accurately implement multi-period water eutrophication assessment, evaluate and monitor water eutrophication faster, more comprehensively, and more economically over a larger range, comprehensively and timely grasp the eutrophication status of water bodies, and provide a scientific basis for formulating timely control measures for lake eutrophication. BRIEF DESCRIPTION OF THE DRAWINGS

[0087] Figure 1 This is a flow chart of a multi-time series water eutrophication assessment method based on satellite images according to the present invention;

[0088] Figure 2 The sampling point distribution diagram of the present invention, wherein (a) is the sampling point distribution in the first quarter, (b) is the sampling point distribution in the second quarter, (c) is the sampling point distribution in the third quarter, and (d) is the sampling point distribution in the fourth quarter;

[0089] Figure 3 This is the Shuffle_GRUnet model architecture diagram described in the present invention;

[0090] Figure 4 This is the DSAT_1D_ShuffleNetV2 module architecture diagram of the present invention;

[0091] Figure 5 This is the shuffle_block1 module architecture diagram of the present invention;

[0092] Figure 6 This is the shuffle_block2 module architecture diagram of the present invention;

[0093] Figure 7 This is the MSCA_1D module architecture diagram of the present invention;

[0094] Figure 8 This is the Channel_Shuffle module architecture diagram of the present invention;

[0095] Figure 9 These are the water eutrophication results maps described in the present invention, wherein (a) is the water eutrophication results map for the first quarter, (b) is the water eutrophication results map for the second quarter, (c) is the water eutrophication results map for the third quarter, and (d) is the water eutrophication results map for the fourth quarter. DETAILED DESCRIPTION

[0096] In order to make the objectives, technical solutions, and advantages of the present invention more clearly understood, the present invention is further described in detail below with reference to the accompanying drawings and specific embodiments. It should be understood that the specific embodiments described herein are only intended to explain the present invention and are not intended to limit the present invention. That is, the specific embodiments described herein are only some embodiments of the present invention, not all embodiments. Generally, the components of the specific embodiments of the present invention described and illustrated in the drawings herein can be arranged and designed in various different configurations, and the present invention can also have other embodiments.

[0097] Therefore, the following detailed description of the specific embodiments of the present invention provided in the accompanying drawings is not intended to limit the scope of the invention as claimed, but is merely representative of selected specific embodiments of the present invention. All other specific embodiments obtained by those skilled in the art based on the specific embodiments of the present invention without making any creative efforts shall fall within the scope of protection of the present invention.

[0098] In order to further understand the content, features and effects of the present invention, the following specific embodiments are given as examples, and the attached Figure 1 -Attached Figure 9 The detailed instructions are as follows:

[0099] Example 1:

[0100] A multi-temporal water eutrophication assessment method based on satellite images includes the following steps:

[0101] S1. Collect paddy field sample data from the first, second, third, and fourth quarters, and use the collected sample data to calculate the nutritional status index of each sample point based on the comprehensive nutritional status index calculation formula, which serves as the target parameter of the model;

[0102] The sample data of step S1 includes chlorophyll a data chl_a, total phosphorus data TP, total nitrogen data TN, transparency data SD, and chemical oxygen demand data CODmn;

[0103] Furthermore, in step S1, samples were collected from a total of 50 lakes, ponds, and reservoirs, with 150 sampling points collected each quarter, for a total of 600 samples. Water quality sample collection, storage, and transportation strictly adhered to the relevant standards for lake water quality sampling techniques, sample storage, and handling specified in the "Technical Guidelines for Sampling Water Quality in Lakes and Reservoirs" (GB / T 14581-93). A water quality sample collection device was used to collect water quality samples from 0.5 m below the water surface. Water quality samples of different indicators were added with the corresponding fixative, stored at 2-5°C, and delivered to the laboratory within 48 hours. With the exception of transparency (SD), which required field measurement using the Secco plate visual method, the other four indicators were tested in the laboratory using the corresponding methods. Total nitrogen (TN), total phosphorus (TP), and chlorophyll (chl_a) were measured using spectrophotometry, while chemical oxygen demand (CODmn) was measured using the potassium dichromate method.

[0104] Furthermore, the calculation formula of the comprehensive nutritional status index in step S1 is as follows:

[0105]

[0106] Among them, TLI (∑) represents the comprehensive nutritional status index, W k represents the relevant weight of the nutritional status index of the kth parameter, TLI(k) represents the nutritional status index of the kth parameter, and K represents the number of evaluation parameters.

[0107] Taking chlorophyll a (chl_a) as the benchmark parameter, the normalized correlation weight calculation formula of the kth parameter is as follows:

[0108]

[0109] Among them, r fk Represents the correlation coefficient between the kth parameter and the reference parameter chl_a.

[0110] The relevant system statistics are as follows:

[0111] Table 1: Correlation coefficient table of benchmark parameter chl_a

[0112] parameter chl_a TP TN SD CODmn <![CDATA[r fk ]]> 1 0.84 0.82 0.83 0.83 <![CDATA[r fk 2 ]]> 1 0.7056 0.6724 0.6889 0.6889

[0113] The calculation formula of the kth parameter nutritional status index TLI(k) is as follows, where the parameters are chlorophyll a (chl_a), total phosphorus (TP), total nitrogen (TN), transparency (SD), and chemical oxygen demand (CODmn).

[0114] TLI(chl_a)=10×(2.5+1.086lnchl_a); (3)

[0115] TLI(TP)=10×(9.436+1.624lnTP); (4)

[0116] TLI(TN)=10×(5.453+1.694lnTN); (5)

[0117] TLI(SD)=10×(5.118-1.94lnSD); (6)

[0118] TLI(CODmn)=10×(0.109+2.661lnCODmn); (7)

[0119] The unit of chlorophyll a (chl_a) is mg / m 3 The unit of transparency (SD) is m; the units of total phosphorus (TP), total nitrogen (TN) and chemical oxygen demand (CODmn) are all mg / L.

[0120] S2. Obtain and preprocess multispectral satellite imagery data that is close in time to the sample data. Calculate the spectral index of each sample point using the preprocessed satellite imagery data based on the spectral index calculation formula, which serves as the input parameter for the model.

[0121] Furthermore, in step S2, 20 scenes of satellite image data for March, June, September and December are obtained respectively, and the image data are pre-processed by radiometric calibration, atmospheric correction, RPC orthorectification, image fusion and geometric precision correction.

[0122] Furthermore, the spectral index calculation formula based on the pre-processed satellite image data in step S2 is as follows:

[0123] SPI b4,b1 =b4 / b1; (8)

[0124] SPI b4,b2 =b4 / b2; (9)

[0125] SPI b4,b3 =b4 / b3; (10)

[0126] SPI b3,b4,b1 =(b3+b4) / b1; (11)

[0127] SPI b3,b4,b2 =(b3+b4) / b2; (12)

[0128] SPI b1,b2,b4,b3 =(b1+b2+b4) / b3; (13)

[0129] SPI b2,b3,b4,b1 =(b2+b3+b4) / b1; (14)

[0130] SPI b1,b3,b4,b2 =(b1+b3+b4) / b2; (15)

[0131] SPI b3,b4,b1,b2 =(b3+b4) / (b1+b2); (16)

[0132] Among them, SPI represents different spectral indices, b1 represents the blue band, b2 represents the green band, b3 represents the infrared band, and b4 represents the near-infrared band.

[0133] S3. Design shuffle modules shuffle_block1 and shuffle_block2 respectively. Based on the designed shuffle modules, build the DSAT_1D_ShuffleNetV2 module for feature data extraction.

[0134] Furthermore, the specific implementation method of step S3 includes the following steps:

[0135] S3.1. Design the shuffle_block1 module, which consists of 1 channel segmentation module Channel_Split, 5 batch normalization modules BN, 2 improved 1D multi-scale convolution attention MSCA_1D, and 3 1D convolution Conv1D 1,1 , 1 depth 1D convolution DWConv1D 3,1 , 3 Relu activation functions and 1 channel shuffle module Channel_Shuffle. The calculation process expression is as follows:

[0136] [in x1 ,in x2 =Channel_Split(in x ,0.5); (17)

[0137]

[0138] in x2_1 =Conv1D 1,1 (BN(DWConv1D 3,1 (Relu(BN(Conv1D 1,1 (inx2 )))))); (19)

[0139]

[0140] out y1 =Channel_Shuffle(Concat(in x1_1 ,in x2_2 )); (twenty one)

[0141] Among them, x Indicates the input feature information of the shuffle_block1 module, in x1 ,in x2 Respectively represent the channel information after segmentation according to the first channel and the second channel, in x1_1 Indicates in x1 After MSCA_1D operation and in x1 The result of multiplication, in x2_1 Indicates in x2 After Conv1D 1,1 , BN, Relu, DWConv1D 3,1 Feature information after operation, in x2_2 Indicates in x2_1 After BN and MSCA_1D operations and in x2_1 The result of the multiplication is then passed through Conv1D 1,1 , feature information after BN and Relu operations, out y1 Indicates in x1_1 ,in x2_2 After matrix connection and Channel_Shuffle operation, the feature information is obtained. Channel_Split(tz,bl) indicates that the feature information tz is split into channels according to the bl ratio. MSCA_1D(.) indicates the designed 1D multi-scale convolution attention module. BN(.) indicates the batch normalization operation. Relu(.) indicates the activation function operation. Conv1D 1,1 (.) represents a 1D convolution operation with a convolution kernel of 1 and a stride of 1, DWConv1D 3,1 (.) represents a depth-1 convolution operation with a convolution kernel of 3 and a stride of 1, Concat(.) represents a matrix concatenation operation, and Channel_Shuffle(.) represents a channel shuffle operation. Represents matrix multiplication;

[0142] S3.2. The designed shuffle_block2 module consists of 6 batch normalization modules BN, 2 improved 1D multi-scale convolution attention MSCA_1D, and 4 1D convolution Conv1D 1,1, 2 depth 1D convolution DWConv1D 3,2 , 3 Relu activation functions and 1 channel shuffle module Channel_Shuffle. The calculation process expression is as follows:

[0143] in 4_1 =DWConv1D 3,2 (in x4 );(twenty two)

[0144]

[0145] in x4_3 =Conv1D 1,1 (BN(DWConv1D 3,2 (Relu(BN(Conv1D 1,1 (in x4 )))))); (twenty four)

[0146]

[0147] out y2 =Channel_Shuffle(Concat(in x4_2 ,in x4_4 ));(26)

[0148] Among them, x4 Indicates the input feature information of the shuffle_block2 module, in x4_1 Indicates that after DWConv1D 3,2 Feature information after operation, in x4_2 Indicates in x4_1 After BN, MSCA_1D, Conv1D 1,1 The result of the operation is then combined with in x4_1 After multiplication, the feature information after BN and Relu operation is obtained. x4_3 Indicates in x4 After Conv1D 1,1 , BN, Relu, DWConv1D 3,2 、Conv1D 1,1 Feature information after operation, inx 4_4 Indicates in x4_3 After BN, MSCA_1D, Conv1D 1,1 The result of the operation is then in x4_3 After multiplication, the feature information after BN and Relu operation is out y2 Indicates in x4_2 ,in x4_4The feature information after matrix connection and Channel_Shuffle operation, MSCA_1D (.) represents the 1-dimensional multi-scale convolution attention module designed by the present invention, BN (.) represents the batch normalization operation, Relu (.) represents the activation function operation, Conv1D 1,1 (.) represents a 1D convolution operation with a convolution kernel of 1 and a stride of 1, DWConv1D 3,2 (.) represents a depth-1 convolution operation with a kernel size of 3 and a stride of 2. Concat(.) represents a matrix concatenation operation. Channel_Shuffle(.) represents a channel shuffle operation. Represents matrix multiplication;

[0149] S3.3. The designed 1D multi-scale convolution attention module MSCA_1D consists of 4 1D convolutions of different scales and depths DWConv1D, 1 1D convolution Conv1D 1,1 The calculation process is as follows:

[0150] in x5_1 =DWConv1D5(in x5 );(27)

[0151]

[0152]

[0153] Among them, x5 Represents the input feature information of the 1D multi-scale convolutional attention module MSCA_1D, in x5_1 Indicates in x5 Feature information after DWConv1D5 operation, in x5_2 ,in x5_3 ,in x5_4 Respectively represent in x5_1 After DWConv1D7, DWConv1D 11 、DWConv1D 21 Feature information after operation, out y3 Indicates in x5_1 ,in x5_2 ,in x5_3 ,in x5_4 The result after matrix addition and Conv1D1 operation is the same as in x5 Multiplied feature information, DWConv1D q Represents the depth 1D convolution operation of q convolution kernel, q∈[1,5,7,11,21], Conv1D 1,1 The 1 that indicates the convolution kernel is 1 and the step size is 1 is the convolution operation, and Add indicates matrix addition. Represents matrix multiplication;

[0154] S3.4. The calculation process of Channel_Shuffle is as follows:

[0155] S3.4.1. Determine the number of feature map x groups and the number of channels C per group g The calculation formula is as follows:

[0156]

[0157] Among them, C represents the number of channels of the feature map x, G represents the number of groups to be divided, and C g Indicates the number of channels in each group;

[0158] S3.4.2.reshape operation:

[0159] Expand the channel dimension of the feature map x from one dimension to two dimensions, one is the number of groups G, and the other is the number of channels in each group C g ; If the batch size of the feature map is N and the spatial dimension is H×W, the reshape operation can be expressed as:

[0160] x reshaped =x.view(N,G,C g ,H,W); (31)

[0161] S3.4.3.transpose operation:

[0162] The reshaped feature map is transposed in terms of the number of groups and the number of channels in each group, so that the channels in each group can be disrupted. The transposition operation can be expressed as:

[0163] x transposed =x reshaped .transpose(1,2);(32)

[0164] Among them, transpose(1,2) means that the second dimension (number of groups G) and the third dimension (number of channels per group C) are transformed. g ) for exchange;

[0165] S3.4.4.flatten operation:

[0166] The transposed feature map x transposed , flatten back to the original number of channels C, and get the final channel shuffled feature map. The flatten operation can be expressed as:

[0167] x shuffled =x transposed.view(N,C,H,W);(33)

[0168] Among them, view(N,C,H,W) restores the feature map to its original four-dimensional shape;

[0169] S3.5. Based on the designed shuffle_block1 and shuffle_block2 modules, the DSAT_1D_ShuffleNetV2 module is designed for feature data extraction. The specific implementation method includes the following steps:

[0170] The designed DSAT_1D_ShuffleNetV2 consists of 1 1D convolution Conv1D 3,2 , 1 1D convolution Conv1D 1,1 , 1 1-dimensional maximum pooling MaxPool 3,2 , 1 global average pooling GlobaPool, 3 designed shuffle_block2, 13 designed shuffle_block1, the calculation process expression is as follows:

[0171] Stage 0 = MaxPool 3,2 (Conv1D 3,2 (input x )); (34)

[0172]

[0173] Stage 2_1 =Concat(interpolate(Stage1,Stage2.shape()),Stage2); (36)

[0174] Stage 3_1 =Concat(interpolate(Stage 2_1 ,Stage3.shape()),Stage3);(37)

[0175] output y =GlobaPool(Conv1D 1,1 (Stage 3_1 )); (38)

[0176] Among them, inputx represents the input feature information of the designed DSAT_1D_ShuffleNetV2, and Stage0 represents input x After Conv1D 3,2 MaxPool 3,2The feature information after the operation, Stage1, Stage2, Stage3 respectively represent the feature information after Stage0 is operated by shuffle_block2 and shuffle_block1, Stage 2_1 It indicates the characteristic information obtained by interpolating the result of Stage1 based on the size of Stage2 and connecting it with the matrix of Stage2. 3_1 Represents Stage 2_1 After interpolation operation based on the Stage3 size, the interpolation result is connected with the Stage3 matrix to obtain the feature information, output y Represents Stage 3_1 After Conv1D 1,1 , feature information of GlobaPool operation, Conv1D 1,1 Represents a 1D convolution operation with a convolution kernel of 1 and a stride of 1, Conv1D 3,2 (.) represents a 1D convolution operation with a kernel size of 3 and a stride of 2. MaxPool 3,2 represents a 1D maximum pooling operation with a convolution kernel of 3 and a stride of 2. shuffle_block2(,v) represents a shuffle_block2 module operation designed v times, v∈[1,1,1], shuffle_block1(,w) represents a shuffle_block1 module operation designed w times, w∈[3,7,3], interpolate(.) represents an interpolation operation, Concat(.) represents a matrix concatenation operation, and GlobaPool(.) represents a global average pooling operation.

[0177] S4. Design a dual-weighted logarithmic root mean square error loss function (Wwrmse) to guide the optimization of multi-time series water eutrophication assessment models.

[0178] Furthermore, the calculation formula of the dual-weight logarithmic root mean square error loss function Wwrmse designed in step S4 is as follows:

[0179]

[0180] Among them, i represents the i-th period, i.e., the first quarter, the second quarter, the third quarter, and the fourth quarter; n represents the total number of periods; j represents the j-th sample data; m represents the total number of samples; and W i represents the weight of the i-th period, wrmse represents the weight logarithmic root mean square error designed by the present invention, represents the forecast value of period i, represents the true value of period i, represents the average value of the true value in period i, It represents the maximum value of the true value in the i-th period, cov represents the covariance operation, and σ represents the standard deviation operation.

[0181] S5. Based on the DSAT_1D_ShuffleNetV2 module designed in S3 and the Wwrmse loss function designed in step S4, combined with a gated recurrent neural network, design a multi-time series water eutrophication assessment model Shuffle_GRUnet;

[0182] Furthermore, the specific implementation method of step S5 includes the following steps:

[0183] The Shuffle_GRUnet model designed in S5.1. consists of 1 input layer, 1 DSAT_1D_ShuffleNetV2 layer designed in step S3, 1 gated recurrent neural network layer, and 1 fully connected layer. Each layer contains n components; the input layer X = {1st period, 2nd period, 3rd period, 4th period...nth period}, and each component input parameter contains 9 steps 2 to calculate the spectral index (SPI b4,b1 、SPI b4,b2 、SPI b4,b3 、SPI b3,b4,b1 、SPI b3,b4,b1 、SPI b1,b2,b4,b3 、SPI b2,b3,b4,b1 、SPI b1,b3,b4,b2 、SPI b3,b4,b1,b2 ), input in the format of a three-dimensional tensor, so that DSAT_1D_ShuffleNetV2 and GRU can be used for feature extraction and model training. The three-dimensional tensor used is (number of samples × time steps × number of features), the number of features is 9, and the spectral index stacking method is used to expand it to 224. The insufficient part is supplemented with 0. The calculation process expression is as follows:

[0184]

[0185] Among them, X i represents the input parameters of period i, Indicates the output parameter of period i, namely the comprehensive nutritional status index, DSAT_1D_ShuffleNetV2 i (.) indicates the DSAT_1D_ShuffleNetV2 module operation designed in step 3 in the i-th period, GRU i (.) indicates that the GRU module operation is performed in the i-th period, fc i (.) indicates that the i-th period has undergone a full connection operation;

[0186] S5.2. The referenced GRU consists of an update gate and a reset gate. First, the update gate Z is calculated. t and reset gate rt The value of the input variable x t and the hidden layer result h at the previous moment t-1 The concatenated matrix is ​​input into the update gate after sigmoid nonlinear transformation; the reset gate value will act on h t-1 and the input variable x t The splicing is performed nonlinearly, and then activated by the tanh function to obtain the state of the candidate hidden layer at the current moment By 1-Z t times h t-1 Store the information of the previous moment through Z t times Record the information at the current moment and add the two results as the state output h of the hidden layer at the current moment t , and finally the weight matrix of the output layer acts on h t On the top, get the current model output y t , the calculation formula is as follows:

[0187] Z t =sigmoid(W Z [h t-1 ,x t ]); (43)

[0188] r t =sigmoid(W r [h t-1 ,x t ]); (44)

[0189]

[0190]

[0191] y t =sigmoid(W o h t ); (47)

[0192] Among them, sigmoid represents the sigmoid activation function operation, tanh represents the tanh activation function operation, Z t Indicates that the door state is updated at time t, r t Indicates that the door state is reset at time t, represents the candidate hidden door state at time t, h t represents the hidden layer state at time t, y t represents the model output at time t, x t represents the model input at time t, h t-1 Indicates the hidden layer state at the previous moment, x t , W ZRepresents the update gate weight matrix, W r Represents resetting the gate weight matrix, represents the candidate hidden layer weight matrix, W o Represents the output layer weight matrix.

[0193] S6. Input the target parameters obtained in step S1 and the input parameters obtained in step S2 into the multi-time series water eutrophication assessment model Shuffle_GRUnet designed in step S5 for model training to obtain the optimal multi-time series water eutrophication assessment model;

[0194] Furthermore, in step S6, the target parameters prepared in step S1 and the input parameters prepared in step S2 are input into the multi-time series water eutrophication assessment model Shuffle_GRUnet designed in step S5 for model training. pcc ), mean root mean square error The model accuracy was evaluated to obtain the optimal multi-time series water eutrophication assessment model. The calculation formula of the evaluation index is as follows:

[0195]

[0196] Among them, r pcc represents the Pearson correlation coefficient, represents the mean root mean square error, i represents the i-th period, i.e. the first quarter, the second quarter, the third quarter, and the fourth quarter, j represents the j-th sample data, and m represents the total number of samples. represents the forecast value of period i, represents the true value of period i, represents the average value of the true value in period i, cov represents the covariance operation, and σ represents the standard deviation operation.

[0197] S7. Using the optimal multi-time series water eutrophication assessment model obtained in step S6, perform a multi-time series water eutrophication assessment in the target area.

[0198] Furthermore, in step S7, a known nutrient status classification method (such as the evaluation method in the "Surface Water Environmental Quality Evaluation Method (Trial)") is referred to to distinguish the degree of eutrophication of the water body. In this way, based on the spectral characteristic information of the eutrophic water body, a eutrophication segmentation classification index is constructed. The specific classification standards are shown in the following table:

[0199] Table 2: Eutrophication classification standards

[0200]

[0201] In summary: This implementation case fully utilizes the characteristic information of satellite image spectral index, connects the characteristic information of different periods, realizes multi-time series water eutrophication assessment, and solves the problem that water eutrophication assessment based on single-period satellite image data leads to poor generalization ability and stability of the model and the accuracy cannot meet application requirements.

[0202] Example 2:

[0203] An electronic device includes a memory and a processor, wherein the memory stores a computer program, and when the processor executes the computer program, the steps of a multi-time series water eutrophication assessment method based on satellite images described in Example 1 are implemented.

[0204] The computer device of the present invention may include a processor and memory, such as a single-chip microcomputer including a central processing unit. Furthermore, the processor is configured to execute a computer program stored in the memory to implement the steps of the aforementioned multi-time series water eutrophication assessment method based on satellite imagery.

[0205] The processor may be a central processing unit (CPU), other general-purpose processors, digital signal processors (DSP), application-specific integrated circuits (ASIC), field-programmable gate arrays (FPGA), other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, etc. A general-purpose processor may be a microprocessor or any conventional processor.

[0206] The memory may primarily include a program storage area and a data storage area. The program storage area may store an operating system and at least one application required for a function (such as a sound playback function or an image playback function); and the data storage area may store data generated based on the use of the mobile phone (such as audio data, a phone book, etc.). Furthermore, the memory may include high-speed random access memory and non-volatile memory, such as a hard disk, internal memory, a plug-in hard disk, a smart media card (SMC), a secure digital (SD) card, a flash card, at least one disk storage device, a flash memory device, or other volatile solid-state storage device.

[0207] Example 3:

[0208] A computer-readable storage medium stores a computer program, which, when executed by a processor, implements the multi-time series water eutrophication assessment method based on satellite images described in Example 1.

[0209] The computer-readable storage medium of the present invention can be any form of storage medium that can be read by a processor of a computer device, including but not limited to non-volatile memory, volatile memory, ferroelectric memory, etc. The computer-readable storage medium stores a computer program. When the processor of the computer device reads and executes the computer program stored in the memory, the steps of the above-mentioned multi-time series water eutrophication assessment method based on satellite images can be implemented.

[0210] The computer program includes computer program code, which may be in source code form, object code form, executable file, or some intermediate form. The computer-readable medium may include: any entity or device capable of carrying the computer program code, recording medium, USB flash drive, mobile hard disk, magnetic disk, optical disk, computer memory, read-only memory (ROM), random access memory (RAM), electric carrier signal, telecommunication signal, and software distribution medium. It should be noted that the content contained in the computer-readable medium may be appropriately increased or decreased according to the requirements of legislation and patent practice in the jurisdiction. For example, in some jurisdictions, according to legislation and patent practice, computer-readable media do not include electric carrier signals and telecommunication signals.

[0211] It should be noted that relational terms such as "first" and "second" are used only to distinguish one entity or operation from another entity or operation, and do not necessarily require or imply any actual relationship or order between these entities or operations. Moreover, the terms "comprises," "comprising," or any other variations thereof are intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus comprising a series of elements includes not only those elements, but also other elements not explicitly listed, or elements inherent to such process, method, article, or apparatus. In the absence of further limitations, an element defined by the phrase "comprising a ..." does not exclude the presence of additional identical elements in the process, method, article, or apparatus comprising the element.

[0212] Although the present application has been described above with reference to specific embodiments, various modifications may be made thereto and components may be substituted with equivalents without departing from the scope of the present application. In particular, as long as there are no structural conflicts, the various features of the embodiments disclosed herein may be combined with each other in any manner, and the omission of an exhaustive description of these combinations in this specification is solely for the sake of space and resource conservation. Therefore, the present application is not limited to the specific embodiments disclosed herein, but includes all technical solutions within the scope of the claims.

Claims

1. A multi-temporal water eutrophication assessment method based on satellite images, characterized in that: The steps include: S1. Collect paddy field sample data from the first, second, third, and fourth quarters, and calculate the nutrient status index of each sample point based on the formula as the target parameter; S2. Obtain and preprocess multispectral satellite imagery data that is close in time to the sample data. Calculate the spectral index of each sample point using the preprocessed data based on a formula as an input parameter. S3. Design shuffle modules shuffle_block1 and shuffle_block2 respectively, and build DSAT_1D_ShuffleNetV2 module based on the shuffle modules for feature data extraction; Step S3 includes the following steps: S3.1.shuffle_block1 consists of 1 channel segmentation module Channel_Split, 5 batch normalization modules BN, 2 improved 1D multi-scale convolution attention MSCA_1D, 3 1D convolution Conv1D 1,1 , 1 depth 1D convolution DWConv1D 3,1 , 3 Relu activation functions and 1 channel shuffle module Channel_Shuffle; The S3.2.shuffle_block2 module consists of 6 batch normalization modules BN, 2 improved 1D multi-scale convolution attention MSCA_1D, and 4 1D convolution Conv1D 1,1 , 2 depth 1D convolution DWConv1D 3,2 , 3 Relu activation functions and 1 channel shuffle module Channel_Shuffle; S3.

3. The improved 1D multi-scale convolution attention module MSCA_1D consists of four 1D convolutions of different scales and depths DWConv1D, one 1D convolution Conv1D 1,1 composition; S3.5.DSAT_1D_ShuffleNetV2 consists of 1 1D convolution Conv1D 3,2 , 1 1D convolution Conv1D 1,1 , 1 1-dimensional maximum pooling MaxPool 3,2 , 1 global average pooling GlobaPool, 3 shuffle_block2, 13 shuffle_block1; S4. Design a dual-weight logarithmic root mean square error loss function Wwrmse; S5. Design a multi-time series water eutrophication assessment model Shuffle_GRUnet based on DSAT_1D_ShuffleNetV2 and Wwrmse combined with a gated recurrent neural network; The Shuffle_GRUnet model consists of 1 input layer, 1 DSAT_1D_ShuffleNetV2 layer, 1 gated recurrent neural network layer, and 1 fully connected layer; S6. Inputting the target parameters and the input parameters into Shuffle_GRUnet for training to obtain an optimal multi-time series water eutrophication assessment model; S7. Use the optimal multi-time series water eutrophication assessment model to conduct multi-time series water eutrophication assessment in the target area.

2. The multi-time series water eutrophication assessment method based on satellite images according to claim 1 is characterized in that: The sample data in step S1 includes chlorophyll a data chl_a, total phosphorus data TP, total nitrogen data TN, transparency data SD, and chemical oxygen demand data CODmn.

3. The multi-time series water eutrophication assessment method based on satellite images according to claim 2 is characterized in that: The calculation formula of the spectral index in different bands in step S2 is as follows: SPI b4,b1 =b4 / b1 (1) SPI b4,b2 =b4 / b2 (2) SPI b4,b3 =b4 / b3 (3) SPI b3,b4,b1 =(b3+b4) / b1 (4) SPI b3,b4,b2 =(b3+b4) / b2 (5) <h2 style=";text-align:left;direction:ltr">SPI<h2 style=";text-align:left;direction:ltr"> b1,b2,b4,b3 <h2 style=";text-align:left;direction:ltr"> (b1+b2+b4) / b3 (6) <h2 style=";text-align:left;direction:ltr">SPI<h2 style=";text-align:left;direction:ltr"> b2,b3,b4,b1 <h2 style=";text-align:left;direction:ltr"> (b2+b3+b4) / b1 (7) <h2 style=";text-align:left;direction:ltr">SPI<h2 style=";text-align:left;direction:ltr"> b1,b3,b4,b2 <h2 style=";text-align:left;direction:ltr"> (b1+b3+b4) / b2 (8) <h2 style=";text-align:left;direction:ltr">SPI<h2 style=";text-align:left;direction:ltr"> b3,b4,b1,b2 <h2 style=";text-align:left;direction:ltr"> (b3+b4) / (b1+b2) (9) Among them, SPI represents the spectral index, b1 represents the blue band, b2 represents the green band, b3 represents the infrared band, and b4 represents the near-infrared band.

4. The multi-time series water eutrophication assessment method based on satellite images according to claim 3 is characterized in that: The specific implementation method of step S3 includes the following steps: The calculation process expression of the designed shuffle_block1 module is as follows: [in x1 ,in x2 ]=Channel_Split(in x ,0.5) (10) in x2_1 =Conv1D 1,1 (BN(DWConv1D 3,1 (Relu(BN(Conv1D 1,1 (in x2 )))))) (12) out y1 =Channel_Shuffle(Concta(in x1_1 ,in x2_2 )) (14) Among them, x Indicates the input feature information of the shuffle_block1 module, in x1 ,in x2 Respectively represent the channel information after segmentation according to the first channel and the second channel, in x1_1 Indicates in x1 After MSCA_1D operation and in x1 The result of multiplication, in x2_1 Indicates in x2 After Conv1D 1,1 , BN, Relu, DWConv1D 3,1 , BN, Conv1D 1,1 Feature information after operation, in x2_2 Indicates in x2_1 After BN and MSCA_1D operations and in x2_1 The result of the multiplication is then passed through Conv1D 1,1 , feature information after BN and Relu operations, out y1 Indicates in x1_1 ,in x2_2 After matrix connection and Channel_Shuffle operation, the feature information is obtained. Channel_Split(tz,bl) indicates that the feature information tz is split into channels according to the bl ratio. MSCA_1D(.) indicates the designed 1D multi-scale convolution attention module. BN(.) indicates the batch normalization operation. Relu(.) indicates the activation function operation. Conv1D 1,1 (.) represents a 1D convolution operation with a convolution kernel of 1 and a stride of 1, DWConv1D 3,1 (.) represents a depth-1 convolution operation with a convolution kernel of 3 and a stride of 1, Concat(.) represents a matrix concatenation operation, and Channel_Shuffle(.) represents a channel shuffle operation. Represents matrix multiplication; The calculation process expression of the designed shuffle_block2 module is as follows: in 4_1 =DWConv1D 3,2 (in x4 ) (15) in x4_3 =Conv1D 1,1 (BN(DWConv1D 3,2 (Relu(BN(Conv1D 1,1 (in x4 )))))) (17) Among them, x4 Indicates the input feature information of the shuffle_block2 module, in x4_1 Indicates that after DWConv1D 3,2 Feature information after operation, in x4_2 Indicates in x4_1 After BN, MSCA_1D, Conv1D 1,1 The result of the operation is then combined with in x4_1 After multiplication, the feature information after BN and Relu operation is obtained. x4_3 Indicates in x4 After Conv1D 1,1 , BN, Relu, DWConv1D 3,2 , BN, Conv1D 1,1 Feature information after operation, inx 4_4 Indicates in x4_3 After BN, MSCA_1D, Conv1D 1,1 The result of the operation is then in x4_3 After multiplication, the feature information after BN and Relu operation is out y2 Indicates in x4_2 ,in x4_4 After matrix connection and Channel_Shuffle operation, the feature information is obtained. MSCA_1D(.) represents the designed 1D multi-scale convolution attention module, BN(.) represents the batch normalization operation, Relu(.) represents the activation function operation, Conv1D 1,1 (.) represents a 1D convolution operation with a convolution kernel of 1 and a stride of 1, DWConv1D 3,2 (.) represents a depth-1 convolution operation with a kernel size of 3 and a stride of 2. Concat(.) represents a matrix concatenation operation. Channel_Shuffle(.) represents a channel shuffle operation. Represents matrix multiplication; S3.

3. The designed 1D multi-scale convolutional attention module MSCA_1D has the following calculation process expression: in x5_1 =DWConv1D5(in x5 ) (20) Among them, x5 Represents the input feature information of the 1D multi-scale convolutional attention module MSCA_1D, in x5_1 Indicates in x5 Feature information after DWConv1D5 operation, in x5_2 ,in x5_3 ,in x5_4 Respectively represent in x5_1 After DWConv1D7, DWConv1D 11 、DWConv1D 21 Feature information after operation, outy3 represents in x5_1 ,in x5_2 ,in x5_3 ,in x5_4 The result after matrix addition and Conv1D1 operation is the same as in x5 Multiplied feature information, DWConv1D q Represents a depth 1D convolution operation with convolution kernel q, q∈[5,7,11,21], Conv1D 1,1 The 1 that indicates the convolution kernel is 1 and the step size is 1 is the convolution operation, and Add indicates matrix addition. Represents matrix multiplication; S3.

4. The calculation process of Channel_Shuffle is as follows: S3.4.

1. Determine the number of feature map x groups and the number of channels C per group g The calculation formula is as follows: Among them, C represents the number of channels of the feature map x, G represents the number of groups to be divided, and C g Indicates the number of channels in each group; S3.4.2.reshape operation; S3.4.3.transpose operation; S3.4.4.flatten operation; S3.

5. The calculation process of the designed DSAT_1D_ShuffleNetV2 is as follows: Stage0=MaxPool 3,2 (Conv1D 3,2 (input x )) (27) Stage 2_1 =Concat(interpolate(Stage1,Stage2.shape()),Stage2) (29) Stage 3_1 =Concat(interpolate(Stage 2_1 ,Stage3.shape()),Stage3) (30) output y =GlobaPool(Conv1D 1,1 (Stage 3_1 )) (31) Among them, input x Indicates the input feature information of the designed DSAT_1D_ShuffleNetV2, Stage0 represents input x After Conv1D 3,2 MaxPool 3,2 The feature information after the operation, Stage1, Stage2, Stage3 respectively represent the feature information after Stage0 is operated by shuffle_block2 and shuffle_block1, Stage 2_1 It indicates the characteristic information obtained by interpolating the result of Stage1 based on the size of Stage2 and connecting it with the matrix of Stage2. 3_1 Represents Stage 2_1 After interpolation operation based on the Stage3 size, the interpolation result is connected with the Stage3 matrix to obtain the feature information, output y Represents Stage 3_1 After Conv1D 1,1 , feature information of GlobaPool operation, Conv1D 1,1 Represents a 1D convolution operation with a convolution kernel of 1 and a stride of 1, Conv1D 3,2 (.) represents a 1D convolution operation with a kernel size of 3 and a stride of 2. MaxPool 3,2 represents a 1D maximum pooling operation with a convolution kernel of 3 and a stride of 2. shuffle_block2(,v) represents a shuffle_block2 module operation designed v times, v∈[1,1,1], shuffle_block1(,w) represents a shuffle_block1 module operation designed w times, w∈[3,7,3], interpolate(.) represents an interpolation operation, Concat(.) represents a matrix concatenation operation, and GlobaPool(.) represents a global average pooling operation.

5. The multi-time series water eutrophication assessment method based on satellite images according to claim 4 is characterized in that: The calculation formula of the dual-weight logarithmic root mean square error loss function Wwrmse designed in step S4 is as follows: Among them, i represents the i-th period, i.e., the first quarter, the second quarter, the third quarter, and the fourth quarter, n^ represents the total number of periods, j represents the j-th sample data, m represents the total number of samples, and W i represents the weight of the i-th period, wrmse represents the weighted logarithmic root mean square error, represents the forecast value of period i, represents the true value of period i, represents the average value of the true value in period i, It represents the maximum value of the true value in the i-th period, cov represents the covariance operation, and σ represents the standard deviation operation.

6. The multi-time series water eutrophication assessment method based on satellite images according to claim 5 is characterized in that: The specific implementation method of step S5 includes the following steps: The Shuffle_GRUnet model designed in step S5.1 consists of an input layer, a DSAT_1D_ShuffleNetV2 layer designed in step S3, a gated recurrent neural network layer, and a fully connected layer. Each layer contains n components; the input layer X = {1st period, 2nd period, 3rd period, 4th period...nth period}, and the input parameters of each component contain 9 spectral indices (SPI) calculated in step 2. b4,b1 、SPI b4,b2 、SPI b4,b3 、SPI b3,b4,b1 、SPI b3,b4,b1 、SPI b1,b2,b4,b3 、SPI b2,b3,b4,b1 、SPI b1,b3,b4,b2 、SPI b3,b4,b1,b2 ), input in the format of a three-dimensional tensor, so that DSAT_1D_ShuffleNetV2 and GRU can be used for feature extraction and model training. The three-dimensional tensor used is: number of samples × time steps × number of features. The number of features is 9, which is expanded to 224 by spectral index stacking. The insufficient part is supplemented by 0. The calculation process expression is as follows: Among them, X i represents the input parameters of period i, Indicates the output parameter of period i, namely the comprehensive nutritional status index, DSAT_1D_ShuffleNetV2 i (.) indicates the DSAT_1D_SHuffleNetV2 module operation designed in step 3 in the i-th period, GRU i (.) indicates that the GRU module operation is performed in the i-th period, fc i (.) indicates a full connection operation; S5.

2. The referenced GRU consists of an update gate and a reset gate. First, the update gate Z is calculated. t and reset gate r t The value of the input variable x t and the hidden layer result h at the previous moment t-1 The concatenated matrix is ​​input into the update gate after sigmoid nonlinear transformation; the reset gate value will act on h t-1 and the input variable x t The concatenation is performed nonlinearly, and then activated by the tanh function to obtain the state of the candidate hidden layer at the current moment. By 1-Z t times h t-1 Store the information of the previous moment through Z t times Record the information at the current moment and add the two results as the state output h of the hidden layer at the current moment t , and finally the weight matrix of the output layer acts on h t On the top, get the current model output y t , the calculation formula is as follows: Z t =sigmoid(W Z [h t-1 ,x t ]) (36) r t =sigmoid(W r [h t-1 ,x t ]) (37) y t =sigmoid(W o h t ) (40) Among them, sigmoid represents the sigmoid activation function operation, tanh represents the tanh activation function operation, Z t Indicates that the door state is updated at time t, r t Indicates that the door state is reset at time t, represents the candidate hidden door state at time t, h t represents the hidden layer state at time t, y t represents the model output at time t, x t represents the model input at time t, h t-1 Indicates the hidden layer state at the previous moment, x t , W Z Represents the update gate weight matrix, W r Represents resetting the gate weight matrix, Represents the candidate hidden layer weight matrix, W o Represents the output layer weight matrix.

7. The multi-time series water eutrophication assessment method based on satellite images according to claim 6 is characterized in that: The implementation method of step S6 is to input the target parameters prepared in step S1 and the input parameters prepared in step S2 into the multi-time series water eutrophication assessment model Shuffle_GRUnet designed in step S5 for model training, and use the Pearson correlation coefficient r pcc , mean root mean square error The model accuracy is evaluated to obtain the optimal multi-time series water eutrophication assessment model.

8. An electronic device, characterized in that: The method comprises a memory and a processor, wherein the memory stores a computer program, and when the processor executes the computer program, the method realizes the steps of a multi-time series water eutrophication assessment method based on satellite images according to any one of claims 1 to 7.

9. A computer-readable storage medium having a computer program stored thereon, characterized in that: When the computer program is executed by a processor, the multi-time series water eutrophication assessment method based on satellite images according to any one of claims 1 to 7 is implemented.

Citation Information

Patent Citations

  • Water quality prediction method based on DeepTCN-GRU deep learning model

    CN118114820A