A remote sensing identification and supervision method, medium and system for human activity changes in marine protected areas

Abstract

The present application provides a kind of remote sensing identification supervision method, medium and system for marine protected area human activity change, belong to marine protection technical field, the present application is by constructing physical drive atmospheric correction framework output off-water radiance, two types of features are input into cross-modal cycle distillation model and output human activity target detection, segmentation and category identification result, then utilize marine dynamics numerical model and ensemble Kalman filter assimilation framework, the natural dynamic prediction field is stripped from multi-temporal remote sensing observation and output residual signal, finally, the regularized optimal transport multi-temporal distribution alignment algorithm is used to decompose the transport planning matrix of residual signal into cyclic component and net transfer component, to output human activity change intensity index map with net transfer component, combined with artificial intelligence target identification result generates supervision report, solves the technical problem that natural dynamic process and human activity change are difficult to effectively decoupled, leading to high false detection rate of marine protected area human activity change identification.

Description

Technical Field

[0001] This invention belongs to the field of marine conservation technology, and specifically relates to a remote sensing identification and monitoring method, medium and system for changes in human activities in marine protected areas. Background Technology

[0002] Monitoring human activities in marine protected areas is a crucial means of marine ecological protection. Current technologies typically employ single-modal optical remote sensing imagery combined with traditional change detection algorithms to identify and monitor targets such as ships, aquaculture cages, and oil and gas platforms within the protected areas. When optical imagery is obscured by clouds or fog, some solutions incorporate synthetic aperture radar (SAR) imagery as a supplement; however, the fusion of these two modalities often involves simple feature overlay and lacks cross-modal knowledge transfer mechanisms. In the change detection phase, mainstream methods use pixel differences or ratios from multi-temporal images as change indicators, relying on threshold segmentation to extract changed areas.

[0003] While the aforementioned methods have some applicability in open marine environments, they exhibit significant shortcomings in the complex scenario of marine protected areas. Marine protected areas simultaneously experience strong natural dynamic processes such as tidal fluctuations, ocean current disturbances, and periodic changes in seabed topography. These natural processes also manifest as spatiotemporal variations in radiance values ​​in multi-temporal remote sensing imagery. Under extreme natural disturbance conditions such as typhoons and spring tides, the magnitude of radiance changes caused by natural dynamics far exceeds those caused by human activities. Pixel difference-based change detection methods cannot separate these two types of signals at the fundamental level, leading to a large number of natural dynamic processes being misidentified as changes caused by human activities, severely reducing the reliability of change detection results. In other words, existing technologies suffer from the technical problem of failing to effectively decouple natural dynamic processes from changes caused by human activities, resulting in a high false positive rate for identifying changes in human activities within marine protected areas. Summary of the Invention

[0004] In view of this, the present invention provides a remote sensing identification and monitoring method, medium and system for changes in human activities in marine protected areas, which can solve the technical problem in the prior art that the natural dynamic process and changes in human activities are difficult to decouple effectively, resulting in a high false detection rate in the identification of changes in human activities in marine protected areas.

[0005] The present invention is implemented as follows: The first aspect of the present invention provides a remote sensing identification and monitoring method for changes in human activities in marine protected areas, comprising the following steps:

[0006] Multi-source remote sensing data of marine protected areas were acquired, and Cloude-Pottier polarization decomposition was performed on the fully polarimetric synthetic aperture radar imagery to extract scattering entropy feature maps, scattering angle feature maps, and inverse entropy feature maps.

[0007] Hybrid pixel decomposition was performed on hyperspectral optical images using a Bayesian variational nonnegative matrix decomposition tensor network algorithm. A physical-driven atmospheric correction framework was constructed by combining shortwave infrared band data, and the water radiance inversion results were output.

[0008] The scattering entropy feature map, scattering angle feature map, inverse entropy feature map and water radiance inversion results are input into the cross-modal cyclic distillation model, and the outputs target detection results, semantic segmentation results and human activity category recognition results are output.

[0009] A framework for assimilating ocean dynamics numerical models and remote sensing observation data is constructed. The natural dynamic prediction field and the multi-temporal remote sensing observation field are assimilated by ensemble Kalman filtering to output the residual signal.

[0010] The residual signal is processed by the regularized optimal transmission multi-temporal distribution alignment algorithm, which decomposes the transmission plan matrix into cyclic components and net transfer components, and outputs a human activity change intensity index map.

[0011] Based on the intensity index map of human activity changes, target detection results, and semantic segmentation results, a monitoring report on human activity changes in marine protected areas is generated.

[0012] Specifically, the Cloude-Pottier polarization decomposition involves performing eigenvalue decomposition on the coherence matrix of a fully polarimetric synthetic aperture radar image. The scattering entropy characterizes the degree of scattering randomness, the scattering angle characterizes the type of dominant scattering mechanism, and the anti-entropy characterizes the relative strength of the secondary dominant scattering mechanism and the tertiary scattering mechanism.

[0013] Specifically, the physical-driven atmospheric correction framework couples the radiative transfer equation with a water body bio-optical model, uses shortwave infrared data to estimate aerosol optical thickness and aerosol type under the assumption of near-zero water radiance, and then substitutes the aerosol optical thickness and aerosol type into the radiative transfer equation to invert the water radiance.

[0014] Specifically, the Bayesian variational nonnegative matrix factorization tensor network algorithm organizes hyperspectral optical image data into a spatial-spectral-temporal three-dimensional tensor, applies a Dirichlet distribution prior based on marine optical physics to the endmember matrix, applies a Markov random field prior to the abundance tensor, and performs variational inference through mean-field approximation decomposition of the posterior.

[0015] In the Bayesian variational nonnegative matrix factorization tensor network algorithm, the tensor network part uses a matrix multiplication state structure to compress the three-dimensional abundance tensor. The key dimension value of the matrix multiplication state is determined by multiple sets of experiments on the reconstruction accuracy and memory usage of hyperspectral optical images of different scenes.

[0016] The cross-modal cyclic distillation model consists of three parts: a fully polarimetric synthetic aperture radar polarimetric feature stream, a hyperspectral optical multispectral feature stream, and a cross-modal interactive kernel. The cross-modal interactive kernel implements a soft-label cyclic distillation mechanism, and the two feature streams complete mutual supervision iteration within each mini-batch training.

[0017] The fusion layer of the cross-modal cyclic distillation model employs attention-gated dynamic weighted fusion. The sum of the fusion weights of the fully polarimetric synthetic aperture radar polarimetric feature stream and the hyperspectral optical multispectral feature stream is 1. The gating weights are jointly determined by the Monte Carlo random inactivation uncertainty estimates of the two feature streams.

[0018] The fusion weight is determined by a fusion weight adjustment function. The fusion weight adjustment function takes the number of Monte Carlo random deactivation samples, the current small batch loss of the fully polarimetric synthetic aperture radar polarimetric feature flow, and the current small batch loss of the hyperspectral optical multispectral feature flow as inputs, calculates the fusion adjustment index, and adjusts the number of Monte Carlo random deactivation samples and the fusion weight according to the relationship between the fusion adjustment index and the upper and lower adjustment thresholds.

[0019] Specifically, the ensemble Kalman filter assimilation involves constructing the covariance matrix between the model prediction set and the observation set, updating the state estimate by minimizing the prediction error covariance, fusing the natural dynamic prediction field with the multi-temporal remote sensing observation field, and outputting the residual signal after removing the natural dynamic process.

[0020] Specifically, the regularized optimal transmission multi-temporal distribution alignment algorithm constructs empirical probability distributions for the dual-temporal residual signals and uses the Sinkhorn-Knopp algorithm to solve for the entropy regularized optimal transmission plan matrix. The regularization strength is determined by the adaptive estimation of the naturally varying reference area within the scene.

[0021] The transmission cost matrix of the transmission plan matrix is ​​defined by the characteristic distance of the physical constraints. The Sinkhorn-Knopp algorithm uses logarithmic domain stabilization calculation and parallelization of the graphics processor to accelerate convergence.

[0022] The output head of the cross-modal cyclic distillation model is designed as a multi-task branch, which simultaneously outputs target detection results, semantic segmentation results, and human activity category recognition results. The gradients of the three task branches are backpropagated in a balanced manner through gradient normalization technology.

[0023] The training of the cross-modal cyclic distillation model is divided into two stages: the first stage performs single-stream pre-training on the polarization feature stream of the fully polarized synthetic aperture radar and the multispectral feature stream of hyperspectral optics; the second stage enables the soft-label cyclic distillation mechanism of the cross-modal interactive kernel, and the soft-label distillation loss is calculated using KL divergence.

[0024] In this process, after receiving the water-leaving radiance inversion results, the hyperspectral optical multispectral characteristic flow uses a spectral attention module to dynamically weight each band and adaptively suppresses the band response contaminated by water optical effects based on the sea state estimate. The sea state estimate is jointly calculated from the aerosol optical thickness and the water-leaving radiance inversion results.

[0025] The key dimension of the matrix product state ranges from 8 to 64, with an experimental step size of 8, and the key dimension value is selected when the reconstruction error tends to level off. The distillation temperature coefficient in the soft-label cyclic distillation mechanism ranges from 2 to 8, and is determined by ablation experiments on different temperature coefficient values ​​on the validation set. The upper adjustment threshold is 0.6, and the lower adjustment threshold is 0.3, which is determined by grid search experiments on different threshold combinations on the training and validation sets with the goal of minimizing the uncertainty of the fusion layer output.

[0026] A second aspect of the present invention provides a computer-readable storage medium storing program instructions, which, when executed in a computer, are used to perform the aforementioned remote sensing identification and monitoring method for changes in human activities in marine protected areas.

[0027] A third aspect of the present invention provides a remote sensing identification and monitoring system for changes in human activities in marine protected areas, comprising the aforementioned computer-readable storage medium, wherein the system is a computer, the computer-readable storage medium is disposed within the system, and the system is provided with a microprocessor for executing program instructions stored in the computer-readable storage medium.

[0028] This invention employs a framework for assimilating ocean dynamics numerical models with remote sensing observation data. It integrates the predicted results of tides, currents, and seabed topography changes from regional ocean simulation systems or finite-volume coastal ocean models with multi-temporal remote sensing observation fields using ensemble Kalman filtering. The output is a residual signal after removing natural dynamic processes, thus separating the natural dynamic components from remote sensing observations at a mechanistic level. Building upon this, the invention further employs a regularized optimal transmission multi-temporal distribution alignment algorithm, elevating the residual signal change detection problem from pixel difference comparison to a probabilistic distribution-level quality transmission analysis. The Sinkhorn-Knopp algorithm decomposes the transmission plan matrix into cyclic and net transfer components, further decoupling residual natural fluctuations from anthropogenic changes at an information theory level. The transmission cost of the net transfer component is output as an indicator of the intensity of anthropogenic changes. The physically constrained transmission cost matrix makes the algorithm insensitive to radiated noise and atmospheric residuals, effectively suppressing false detection rates under conditions of similar intensity of natural disturbances after typhoons. In summary, this invention solves the technical problem mentioned in the background art, which is that the difficulty in effectively decoupling natural dynamic processes from changes in human activities leads to a high false detection rate in the identification of changes in human activities in marine protected areas. Attached Figure Description

[0029] Figure 1 This is a flowchart of the method of the present invention.

[0030] Figure 2 This is a schematic diagram of the cross-modal cyclic distillation model.

[0031] Figure 3 This is a schematic diagram of the Kalman filter assimilation process for the set of natural dynamic prediction fields and multi-temporal remote sensing observation fields.

[0032] Figure 4 This is a schematic diagram of the decomposition of the cyclic component and net transfer component of the transmission plan matrix.

[0033] Figure 5 A schematic diagram of the visualization results of the intensity index of changes in human activities.

[0034] Figure 6 This is a schematic diagram comparing the spatial distribution of the scattering entropy feature map, scattering angle feature map, and anti-entropy feature map.

[0035] Figure 7 This is a schematic diagram of the endmember abundance map of mixed pixel decomposition.

[0036] Figure 8 This is a spatial distribution map of changes in human activity in the area reported by the regulatory authorities.

[0037] Figure 9 This is a schematic diagram comparing the target detection results of single-stream inference mode and dual-stream inference mode.

[0038] Figure 10 This is a comparative diagram of the intensity of human activity changes under different meteorological conditions. Detailed Implementation

[0039] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions in the embodiments of the present invention will be clearly and completely described below.

[0040] like Figure 1 The diagram shown is a flowchart of a remote sensing identification and monitoring method for changes in human activities in marine protected areas, provided by the first aspect of this invention. This method includes the following steps:

[0041] S01. Acquire multi-source remote sensing data of marine protected areas, perform Cloude-Pottier polarization decomposition on fully polarimetric synthetic aperture radar images, and extract scattering entropy, scattering angle and anti-entropy feature maps.

[0042] S02. The Bayesian variational nonnegative matrix decomposition tensor network algorithm is used to perform hybrid pixel decomposition on the hyperspectral optical image. Combined with shortwave infrared band data, a physical-driven atmospheric correction framework is constructed to output the water radiance inversion results.

[0043] S03. Input the scattering entropy feature map, scattering angle feature map, anti-entropy feature map and water radiance inversion results into the cross-modal cyclic distillation model, and output the target detection results, semantic segmentation results and human activity category recognition results.

[0044] S04. Construct a framework for the assimilation of ocean dynamics numerical models and remote sensing observation data, and perform ensemble Kalman filtering assimilation of natural dynamic prediction fields and multi-temporal remote sensing observation fields to output residual signals.

[0045] S05. The residual signal is processed by the regularized optimal transmission multi-temporal distribution alignment algorithm, which decomposes the transmission plan matrix into cyclic components and net transfer components, and outputs a human activity change intensity index map.

[0046] S06. Based on the intensity index map of human activity changes, target detection results, and semantic segmentation results, generate a monitoring report on human activity changes in marine protected areas.

[0047] The fully polarimetric synthetic aperture radar imagery originates from a spaceborne or airborne fully polarimetric synthetic aperture radar sensor, the hyperspectral optical imagery originates from a spaceborne hyperspectral sensor, and the shortwave infrared band data originates from a shortwave infrared detector on the same platform or orbit as the hyperspectral optical imagery.

[0048] The Cloude-Pottier polarization decomposition decomposes the coherence matrix of a fully polarimetric synthetic aperture radar image into eigenvalues, and the scattering entropy... Characterizing the degree of randomness in scattering, the value range is: ; Scattering angle Characterizes the type of dominant scattering mechanism, with a value range of . Anti-entropy The relative intensity of the secondary dominant scattering mechanism and the tertiary scattering mechanism is characterized, with values ​​ranging from [value range missing]. These three elements constitute the scattering entropy feature map, the scattering angle feature map, and the anti-entropy feature map, which serve as the input for step S03.

[0049] The physical-driven atmospheric correction framework couples the radiative transfer equation with a water body bio-optical model. It uses shortwave infrared data to estimate aerosol optical thickness and aerosol type under the assumption of near-zero water radiance, and then substitutes the aerosol optical thickness and aerosol type into the radiative transfer equation to invert the water radiance. Nearshore mixed pixels use a pixel linear unmixing model to separate the water body radiation component from the land radiation component. Solar flare contamination pixels are identified and removed using polarization degree images from polarization remote sensing data. The aerosol type discrimination threshold is determined by experimental analysis of 30-50 historical cloudless images of the study area and iterative optimization 3-5 times.

[0050] The Bayesian variational nonnegative matrix factorization tensor network algorithm organizes hyperspectral optical image data into a spatial-spectral-temporal three-dimensional tensor. It applies a Dirichlet distribution prior based on ocean optical physics to the endmember matrix and a Markov random field prior to the abundance tensor. Variational inference uses mean-field approximation to decompose the posterior and utilizes natural gradient updates to accelerate variational parameter convergence. A sparsity-promoted truncated normal prior automatically closes redundant endmembers to automatically determine the number of effective endmembers. The tensor network part uses a matrix multiplication state structure to compress the three-dimensional abundance tensor, reducing storage complexity from... Down to ,in Let be the size of the tensor space. The key dimension of the matrix product state. The range of values ​​is The accuracy and memory usage of reconstruction were determined through multiple experiments on hyperspectral optical images of different scenes. The experimental step size was 8, and the step size was selected when the reconstruction error tended to level off. Values; combined with the GPU tensor shrinking operation library, efficient batch inference is achieved, and the mixed pixel decomposition endmember abundance map is output.

[0051] The Bayesian variational nonnegative matrix factorization tensor network algorithm brings the following technical effects to the solution: the radiation mixing problem caused by the coexistence of aquaculture areas, reefs, and water bodies in nearshore hyperspectral optical images is resolved by probabilistic unmixing under physical constraints; the automatic determination mechanism of the number of endmembers eliminates the subjective error introduced by manually setting hyperparameters; and the uncertainty map output by the Bayesian framework provides confidence weights for subsequent change detection, making the mixed pixel decomposition results more generalizable and interpretable in marine protected areas with variable optical environments.

[0052] The cross-modal cyclic distillation model is based on deep learning and consists of three parts: a fully polarimetric synthetic aperture radar polarimetric feature flow, a hyperspectral optical multispectral feature flow, and a cross-modal interactive kernel.

[0053] The polarization feature flow received scattering entropy feature map, scattering angle feature map, and anti-entropy feature map of the fully polarized synthetic aperture radar are concatenated with the statistical feature map extracted by the convolutional neural network at the channel level to form a fully polarized synthetic aperture radar fused feature map.

[0054] The hyperspectral optical multispectral feature stream receives the water radiance inversion results, uses a spectral attention module to dynamically weight each band, and adaptively suppresses the band response contaminated by water optical effects based on the sea state estimate to form an optical fusion feature map; the sea state estimate is jointly calculated from the aerosol optical thickness output by the physical-driven atmospheric correction framework in step S02 and the water radiance inversion results.

[0055] The cross-modal interactive kernel implements a soft-label cyclic distillation mechanism: the softmax output of the intermediate layer of the fully polarimetric synthetic aperture radar polarimetric feature stream is used as a soft supervision signal and input into the corresponding layer of the hyperspectral optical multispectral feature stream. The softmax output of the intermediate layer of the hyperspectral optical multispectral feature stream is used as a soft supervision signal and input into the corresponding layer of the fully polarimetric synthetic aperture radar polarimetric feature stream. The two streams complete mutual supervision iteration within each mini-batch training. When data for a certain modality is missing, the network switches to single-stream inference mode, and the feature representation of the missing modality is reconstructed from the complete modality through the distillation path.

[0056] The fusion layer employs attention-gated dynamic weighted fusion, with the fusion weights of the polarization feature flows of the fully polarimetric synthetic aperture radar being: The fusion weights of the hyperspectral optical multispectral feature flows are: Both satisfy The gating weights are jointly determined by the Monte Carlo random inactivation uncertainty estimates of the two feature flows, and the fusion weights are automatically reduced for feature flows with high uncertainty.

[0057] The output head is designed as a multi-task branch, simultaneously outputting object detection results, semantic segmentation results, and human activity category recognition results; the gradients of the three task branches are balanced backpropagated using gradient normalization technology; specific task branches are selectively activated during inference to reduce computational overhead.

[0058] The network's CUDA stream allocation strategy assigns the fully polarimetric synthetic aperture radar polarimetric feature stream and the hyperspectral optical multispectral feature stream to two independent CUDA streams for parallel forward propagation. The soft tag transfer operation for cross-modal interaction kernels is assigned to a third CUDA stream responsible for cross-stream data synchronization. Memory allocation for convolutional layer weights is hierarchically managed according to feature map resolution; high-resolution shallow feature maps use half-precision floating-point storage, while low-resolution deep feature maps use full-precision floating-point storage. The number of samples for Monte Carlo random deactivation is also specified. With fusion weight The functional relationship between them is determined by the fusion weight adjustment function.

[0059] The fusion weight adjustment function uses the number of Monte Carlo random deactivation samples. (Unit: times), Current small-batch loss of polarization characteristic flow of fully polarized synthetic aperture radar (Dimensionless) and hyperspectral optical multispectral characteristic flow current small batch loss Calculate the fusion adjustment index using (dimensionless) as input. ;when At that time, the number of Monte Carlo random deactivation samples of the polarization characteristic flow of the fully polarized synthetic aperture radar was adjusted to... and will Lower ;when At that time, keep and Unchanged; when At that time, the number of Monte Carlo random deactivation samples of the hyperspectral optical multispectral characteristic stream was adjusted to and will Lower After each adjustment, renormalization is performed. The thresholds of 0.6 and 0.3 were determined by conducting grid search experiments on different threshold combinations on the training and validation sets, and iterating for 3 to 5 rounds with the goal of minimizing the uncertainty of the fusion layer output.

[0060] The specific steps for establishing the training dataset for the cross-modal cyclic distillation model include: collecting multi-temporal fully polarimetric synthetic aperture radar images and concurrent hyperspectral optical images covering marine protected areas; performing preprocessing such as image registration, radiometric calibration, and atmospheric correction; manually labeling human activity targets such as ships, aquaculture cages, and oil and gas platforms, as well as natural backgrounds such as waves, islands, reefs, and algal blooms, to establish target detection labels, semantic segmentation labels, and human activity category identification labels; retaining samples obscured by clouds and fog to construct single-modal missing training samples; and using random pruning, polarization rotation, and brightness perturbation for data augmentation.

[0061] The training steps of the cross-modal cyclic distillation model specifically include: In the first stage, single-stream pre-training is performed on the fully polarimetric synthetic aperture radar polarimetric feature stream and the hyperspectral optical multispectral feature stream, respectively, with the loss function being the sum of the cross-entropy loss and IoU loss of each task branch; In the second stage, the soft-label cyclic distillation mechanism with a cross-modal interactive kernel is enabled, and the soft-label distillation loss is calculated using KL divergence, with the distillation temperature coefficient... The range of values ​​is By applying different validation sets The values ​​were determined through ablation experiments; the gradient normalization coefficient was calculated by inversely proportional to the gradient norm of each task branch on the current mini-batch; the optimizer used was AdamW, and the initial learning rate ranged from [value missing]. The value is determined by the learning rate grid search experiment on the validation set loss.

[0062] The technical effects of the cross-modal cyclic distillation model are as follows: the cross-modal cyclic distillation mechanism of dual-flow polarimetric synthetic aperture radar and optical imagery enables the two modal branches to transfer intermediate layer knowledge to each other during training; the embedding of physically differentiable polarimetric decomposition into the polarimetric feature flow of fully polarimetric synthetic aperture radar allows the network to retain the physical interpretability of the scattering mechanism; the spectral attention module adaptively suppresses the optical pollution band of water bodies; and the multi-task gradient equalization mechanism avoids a single task dominating the training direction. Overall, the recognition accuracy and robustness of human targets in marine protected areas under marine optical and radar backgrounds are improved.

[0063] The marine dynamics numerical model is a regional marine simulation system or a finite volume coastal marine model, used to numerically simulate tides, currents and seabed topography changes, and output a natural dynamic prediction field. The ensemble Kalman filter assimilation constructs the covariance matrix between the model prediction set and the observation set, updates the state estimate by minimizing the prediction error covariance, and fuses the natural dynamic prediction field with the multi-temporal remote sensing observation field to output the residual signal after removing the natural dynamic process.

[0064] The regularized optimal transmission multi-temporal distribution alignment algorithm reconstructs the problem of change detection between multi-temporal residual signals into an optimal transmission problem between two probability distributions. The transmission cost matrix is ​​defined by the feature distance of physical constraints. Empirical probability distributions are constructed for the dual-temporal residual signals, and the entropy regularized optimal transmission plan matrix is ​​solved using the Sinkhorn-Knopp algorithm. The regularization strength is determined by the adaptive estimation of the natural change reference area within the scene. The transmission plan matrix is ​​decomposed into a cyclic component and a net transfer component. The cyclic component corresponds to natural fluctuations, and the net transfer component corresponds to changes caused by human activities. The transmission cost of the net transfer component is used as the output index of the intensity of changes caused by human activities. The Sinkhorn-Knopp algorithm uses log-domain stabilization calculation and GPU parallelization to accelerate convergence. The regularization strength parameter is determined by selecting a known natural dynamic reference area in each image, statistically analyzing the changes in the feature distribution of the reference area, and using the mean and standard deviation of the transmission cost of the reference area as the upper and lower bounds of the regularization strength, calibrated through 20 to 50 iterations.

[0065] The regularized optimal transmission multi-temporal distribution alignment algorithm brings the following technical effects to the scheme: it upgrades change detection from pixel difference comparison to quality transmission analysis at the probability distribution level; the physically constrained transmission cost matrix makes the algorithm insensitive to radiated noise and atmospheric residuals; the separation of cyclic components and net transfer components decouples natural fluctuations from human activity changes from an information theory perspective, effectively suppressing the false detection rate of human activity change detection under the background of strong natural disturbances after typhoons; and the output human activity change intensity index map has physically interpretable quantitative significance.

[0066] In the soft-label circulating distillation mechanism, the distillation temperature coefficient is... Controlling the smoothness of the softmax output, The larger the value, the smoother the soft label and the more gentle the cross-modal knowledge transfer; the mean-field approximation decomposition posterior index in variational inference approximates the joint posterior distribution as the product of the independent factor distributions of each parameter, thereby reducing the computational cost of variational inference; the Markov random field prior imposes a spatial neighborhood consistency constraint on the abundance tensor, making the endmember abundance distribution of adjacent pixels tend to be smooth; the IoU loss is the negative logarithmic loss of the intersection-union ratio of the predicted box and the ground truth box in target detection, used to optimize the target localization accuracy; the gradient normalization technique scales the gradient vectors of each task branch to the same order of magnitude, avoiding the gradient of a certain task branch from dominating the backpropagation direction due to differences in the magnitude of loss; the key dimension of the matrix multiplication state. The control tensor network balances the expressive power of high-dimensional abundance tensors with storage overhead; the Monte Carlo random inactivation estimates the uncertainty of network predictions by randomly discarding neurons multiple times during the inference phase and taking the variance of multiple outputs; the natural dynamic prediction field refers to the numerical simulation output of regional ocean simulation systems or finite volume coastal ocean models on tides, currents, and seabed topography changes; the sea state estimate refers to the quantitative index characterizing the current ocean optical state obtained by jointly calculating the aerosol optical thickness and water-leaving radiance inversion results, including turbidity estimates and chlorophyll concentration estimates.

[0067] The specific implementation of step S01 is as follows: First, acquire fully polarimetric synthetic aperture radar (MPAR) images covering the marine protected area from a spaceborne or airborne MPAR sensor, and acquire hyperspectral optical images of the same area from a spaceborne hyperspectral sensor. Then, acquire short-wave infrared band data from a short-wave infrared detector on the same platform or orbit as the hyperspectral optical images. Perform registration and radiometric calibration preprocessing on the above multi-source remote sensing data to ensure consistency in spatial location and radiometric magnitude between different data sources. Subsequently, perform Cloude-Pottier polarization decomposition on the MPAR images. Specifically, perform eigenvalue decomposition on the polarization coherence matrix of each pixel, and calculate the scattering entropy from the eigenvalues ​​and eigenvectors respectively. Scattering angle and anti-entropy , where scattering entropy The scattering angle characterizes the randomness of the scattering process. Characterizing the type of dominant scattering mechanism, anti-entropy The relative strengths of the secondary dominant scattering mechanism and the tertiary scattering mechanism are characterized, and the three are used to generate scattering entropy feature maps, scattering angle feature maps, and anti-entropy feature maps, which serve as inputs for subsequent steps.

[0068] The specific implementation of step S02 is as follows: First, a hybrid pixel decomposition algorithm based on Bayesian variational nonnegative matrix factorization tensor network is used to decompose the hyperspectral optical image, organizing the image data into a spatial-spectral-temporal three-dimensional tensor. A Dirichlet distribution prior based on marine optical physics is applied to the endmember matrix, and a Markov random field prior is applied to the abundance tensor to constrain spatial neighborhood consistency. The joint posterior is decomposed into the product of independent factor distributions of each parameter through mean field approximation, and the convergence of variational parameters is accelerated by updating with natural gradients. The sparsity-promoted truncated normal prior automatically closes redundant endmembers to determine the number of effective endmembers. The tensor network part uses a matrix multiplication state structure to compress the three-dimensional abundance tensor, and the key dimension... The reference value range is 8 to 64, the experimental step size is 8, and the value is selected when the reconstruction error decreases gradually. The value will reduce the storage complexity from Down to The system outputs the endmember abundance map of the mixed pixel decomposition. Based on this, a physical-driven atmospheric correction framework is constructed, coupling the radiative transfer equation with the aquatic bio-optical model. Using shortwave infrared data, aerosol optical thickness and aerosol type are estimated under the assumption of near-zero water radiance. The aerosol type discrimination threshold is determined by experimental analysis and iterative optimization 3 to 5 times on 30 to 50 historical cloudless images of the study area. Then, the aerosol optical thickness and aerosol type are substituted into the radiative transfer equation to invert the water radiance. Nearshore mixed pixels are separated from the water radiative component and the land radiative component using a pixel linear unmixing model. Solar flare contamination pixels are identified and removed using polarization degree images from polarization remote sensing data. Finally, the water radiance inversion results are output.

[0069] The specific implementation of step S03 is as follows: The cross-modal cyclic distillation model consists of three parts: a fully polarimetric synthetic aperture radar (SAR) polarimetric feature stream, a hyperspectral optical multispectral feature stream, and a cross-modal interactive kernel. The SAR polarimetric feature stream receives scattering entropy feature maps, scattering angle feature maps, and inverse entropy feature maps, and concatenates these with the statistical feature maps extracted by the convolutional neural network at the channel level to form a SAR fusion feature map. The hyperspectral optical multispectral feature stream receives the water radiance inversion results and uses a spectral attention module to dynamically weight each band. Based on the sea state estimate, it adaptively suppresses the band response polluted by water optical effects. The sea state estimate is jointly calculated from the aerosol optical thickness and the water radiance inversion results, including turbidity and chlorophyll concentration estimates. The cross-modal interactive kernel implements a soft-label cyclic distillation mechanism. The softmax outputs of the intermediate layers of the two feature streams serve as soft supervision signals, completing mutual supervision iterations within each mini-batch training. When data for a certain modality is missing, the network switches to a single-stream inference mode. The fusion layer employs attention-gated dynamic weighted fusion, with the gating weights determined by a fusion weight adjustment function, which uses the number of Monte Carlo random deactivation samples as a parameter. The current mini-batch loss of the two feature streams is used as input to calculate the fusion adjustment index. ,when The polarization characteristic flow of the fully polarized synthetic aperture radar will be... Adjusted to And Lower, when When performing symmetric operations on hyperspectral optical multispectral feature streams, when The time remains constant, and renormalization is performed after each adjustment. The output head is designed with multiple task branches, simultaneously outputting object detection results, semantic segmentation results, and human activity category recognition results. The gradients of the three task branches are backpropagated in a balanced manner using gradient normalization technology.

[0070] The specific implementation of step S04 is as follows: A marine dynamics numerical model is constructed using a regional marine simulation system or a finite-volume coastal marine model as the core. Seabed topography data, meteorological forcing fields, and boundary conditions of the study area are input. Numerical simulations of tides, currents, and seabed topography changes are performed, and a natural dynamic prediction field is output. Subsequently, a remote sensing data assimilation framework is constructed. The covariance matrix of the model prediction set and the multi-temporal remote sensing observation set is constructed using the ensemble Kalman filter method. The state estimate is updated based on minimizing the prediction error covariance. The natural dynamic prediction field is fused with the multi-temporal remote sensing observation field, and a residual signal after removing the natural dynamic processes is output. Theoretically, this residual signal has removed the main contributions of natural dynamics such as tides and currents, thus improving the signal-to-noise ratio of anthropogenic changes in the input signal for subsequent change detection.

[0071] The specific implementation of step S05 is as follows: Empirical probability distributions are constructed for the dual-phase residual signals. The transmission cost matrix is ​​defined using the characteristic distance of physical constraints. The Sinkhorn-Knopp algorithm is used for stabilization calculations in the logarithmic domain, and convergence is accelerated through parallelization using a graphics processor. The entropy-regularized optimal transmission plan matrix is ​​solved. The regularization strength is adaptively determined by the statistical characteristics of the natural variation reference area within the scene. The mean and standard deviation of the transmission cost in the reference area are used as the upper and lower bounds of the regularization strength, calibrated through 20 to 50 iterations. The transmission plan matrix is ​​decomposed into cyclic components and net transfer components. The cyclic components correspond to residual natural fluctuations, and the net transfer components correspond to changes caused by human activities. The transmission cost of the net transfer components is output as an index of the intensity of changes caused by human activities.

[0072] The specific implementation of step S06 is as follows: the human activity change intensity index map output in step S05 is spatially superimposed and semantically correlated with the target detection results, semantic segmentation results and human activity category identification results output in step S03. For areas where the change intensity exceeds the judgment threshold, the human activity category is confirmed in combination with the target identification results. Finally, a monitoring report on human activity changes in marine protected areas is generated. The report includes a spatial distribution map of the changed area, statistics on human activity categories and quantitative indicators of change intensity.

[0073] It should be noted that the key technologies of this invention include: the Bayesian variational nonnegative matrix factorization tensor network algorithm solves the radiation mixing problem caused by the coexistence of multiple types of land cover in nearshore hyperspectral imagery through probabilistic demixing and automatic endmember number determination; its output uncertainty map provides confidence weights for subsequent change detection; the cross-modal cyclic distillation mechanism enables the two modal branches of fully polarimetric synthetic aperture radar and optics to transfer intermediate layer knowledge during training; the embedding of physically differentiable polarimetric decomposition into the polarimetric feature flow of fully polarimetric synthetic aperture radar preserves the interpretability of the scattering mechanism; and the multi-task gradient equalization mechanism avoids single-task dominating training, thus improving the overall efficiency. Robustness of human target identification in complex marine environments; a two-level decoupled system is formed by regularized optimal transmission multi-temporal distribution alignment algorithm and ensemble Kalman filter assimilation framework. Numerical model assimilation removes deterministic natural dynamic components that can be physically modeled, and optimal transmission decomposition further removes residual random natural fluctuations. The two work together to enable the effective extraction of human activity change signals under strong natural disturbances. From the signal generation mechanism level, it completely solves the inherent defect of pixel difference methods that cannot distinguish between two types of changes. The three key technologies support each other in their respective links and jointly realize highly reliable remote sensing identification and monitoring of human activity changes in marine protected areas.

[0074] It should be noted that in marine protected areas where there are numerous hyperspectral image bands, complex nearshore land cover types, and overlapping radiation, the mixed pixel problem can severely affect the accuracy of water optical parameter inversion, thereby reducing the reliability of subsequent target identification and change detection. This technical problem arises because the spatial resolution of hyperspectral images often cannot perfectly match the boundaries of land cover. A single pixel often contains radiation contributions from multiple land cover types, such as aquaculture areas, reefs, and water bodies. Traditional non-negative matrix factorization methods require manually pre-setting the number of endmembers; if the pre-set value does not match the actual number of endmembers, systematic errors are introduced. Furthermore, the marine optical environment changes with time and sea conditions, leading to instability in the spectral characteristics of endmembers, further exacerbating the difficulty of unmixing. The usual solution to this technical problem is to use classical non-negative matrix factorization or linear unmixing models, relying on manual experience to set the number of endmembers and select reference endmember spectra. However, such methods lack explicit modeling of endmember spectral uncertainties, resulting in insufficient generalization ability in marine protected areas with variable optical environments. Furthermore, they fail to provide confidence information for unmixing results in subsequent processing steps, leading to the inability of the change detection stage to distinguish between pseudo-changes caused by unmixing errors and genuine human activity changes. This invention effectively solves this technical problem. The Bayesian variational nonnegative matrix factorization tensor network algorithm applies a Dirichlet distribution prior based on marine optical physics to the endmember matrix and a Markov random field prior to the abundance tensor. Under physical constraints, it probabilistically models the unmixing process. The sparsity-promoted truncated normal prior automatically closes redundant endmembers during variational inference, fundamentally eliminating subjective errors introduced by manually setting the number of endmembers. The uncertainty map output by the Bayesian framework provides pixel-by-pixel confidence weights for subsequent change detection, enabling change detection to automatically reduce decision confidence in high uncertainty regions, thereby suppressing false detections caused by unmixing error propagation. Overall, this effectively overcomes the mixed pixel problem under variable marine optical environments.

[0075] A second aspect of the present invention provides a computer-readable storage medium storing program instructions, which, when executed in a computer, are used to perform the aforementioned remote sensing identification and monitoring method for changes in human activities in marine protected areas.

[0076] A third aspect of the present invention provides a remote sensing identification and monitoring system for changes in human activities in marine protected areas, comprising the aforementioned computer-readable storage medium. The system can be any one of a computer, a server, or a microcontroller. The computer-readable storage medium is disposed within the system, and the system is provided with a microprocessor that executes the program instructions stored in the computer-readable storage medium.

[0077] Specifically, the principle of this invention is:

[0078] The fundamental reason why this invention can solve the above-mentioned technical problems is that its technical approach starts from the signal generation mechanism and constructs a two-level decoupling system, rather than relying on statistical thresholds to perform post-processing segmentation of mixed signals.

[0079] The first stage of decoupling relies on numerical models of ocean dynamics. Regional ocean simulation systems or finite-volume coastal ocean models can perform numerical simulations of the physical mechanisms of tides, currents, and seabed topography changes, outputting a natural dynamic prediction field. The ensemble Kalman filter assimilation method constructs the covariance matrix of the model prediction set and the remote sensing observation set, updating the state estimate by minimizing the prediction error covariance, thus fusing the natural dynamic prediction field with the multi-temporal remote sensing observation field and outputting a residual signal. This residual signal theoretically eliminates the main contributions of periodic natural dynamics such as tides and currents, significantly improving the signal-to-noise ratio of anthropogenic changes in the input signal for subsequent change detection.

[0080] The second-level decoupling relies on a regularized optimal transmission multi-temporal distribution alignment algorithm. This algorithm models the dual-temporal residual signals as empirical probability distributions and uses the Sinkhorn-Knopp algorithm to solve for the entropy-regularized optimal transmission plan matrix. The core advantage of the optimal transmission framework lies in its transmission cost matrix, which is defined by the characteristic distance of physical constraints, exhibiting inherent robustness to radiated noise and atmospheric correction residuals. After decomposing the transmission plan matrix into cyclic and net transfer components, the cyclic component corresponds to residual natural fluctuations, while the net transfer component corresponds to changes caused by human activities. The regularization strength is adaptively determined by the statistical characteristics of the natural variation reference area within the scene, enabling the algorithm to automatically adjust the separation boundary based on the natural dynamic intensity of the current scene, rather than relying on a fixed threshold.

[0081] The two-stage decoupling mechanism is logically complementary: numerical model assimilation removes deterministic natural dynamic components that can be modeled by the physical model, while optimal transport decomposition further removes the random natural fluctuation components remaining after numerical model assimilation. Their synergy enables the effective extraction of signals from human-induced changes even under strong natural disturbances, fundamentally solving the inherent limitation of pixel-difference methods in distinguishing between the two types of changes at the signal level.

[0082] The following provides a specific embodiment 1 of the present invention, and the specific implementation of each step in this embodiment 1 is described in detail below.

[0083] The specific implementation method of step S01 is as follows.

[0084] Acquire fully polarimetric synthetic aperture radar (MPAR) images of marine protected areas collected by spaceborne or airborne MPAR sensors, and analyze each pixel. coherence matrix Eigenvalue decomposition is performed, specifically as follows:

[0085] ;

[0086] In the formula, For the first There are eigenvalues ​​that satisfy... , The corresponding unit eigenvector is represented as a column vector. , for The first component has a value range of . , for The conjugate transpose of . Here, the eigenvalue index takes the value 1, 2, or 3. Scattering entropy The calculation formula is as follows:

[0087] ;

[0088] In the formula, For the first The normalized probability of each eigenvalue is dimensionless. , A value approaching 0 indicates a single scattering mechanism dominating, while a value approaching 1 indicates completely random scattering. Scattering angle The calculation formula is as follows:

[0089] ;

[0090] In the formula, For the first The angle corresponding to the magnitude of the first component of each eigenvector. This represents the modulo operation. , near Corresponding to surface scattering, close to Corresponding to secondary scattering. Anti-entropy. The calculation formula is as follows:

[0091] ;

[0092] In the formula, , A value close to 1 indicates that the secondary dominant scattering mechanism is much stronger than the tertiary scattering mechanism. These three values ​​are used to generate scattering entropy feature maps, scattering angle feature maps, and anti-entropy feature maps, which are then used as inputs for step S03.

[0093] The specific implementation method of step S02 is as follows.

[0094] For spaceborne hyperspectral optical images, a Bayesian variational nonnegative matrix decomposition tensor network algorithm is used for hybrid pixel decomposition. The hyperspectral optical images are organized into a spatial-spectral-temporal three-dimensional tensor. ,in For the number of rows, For column numbers, Let be the number of bands. The endmember matrix is ​​denoted as . , The abundance tensor is denoted as the number of endmembers. The reconstructed relationships are as follows:

[0095] ;

[0096] In the formula, This represents the modular product operation along the third dimension. , All units are radiance units. The counterpart matrix Imposing a Dirichlet distribution prior based on ocean optical physics ,in The Dirichlet concentration parameter has an empirical value of 0.1; for the abundance tensor... We apply a Markov random field prior to smooth out the abundance of adjacent pixels, with the following constraints:

[0097] ;

[0098] In the formula, Represents spatial neighborhood cell pairs, , These are the cell row and column indices, , Use the row and column indexes of its spatial neighborhood cells. For pixels Abundance vector at location (dimensionless). The smoothness intensity coefficient (dimensionless) has an empirical value of 0.01 to 0.1. The Euclidean norm is used. Variational inference employs a mean-field approximation, approximating the joint posterior distribution as the product of independent factors of each parameter. Accelerate convergence by using natural gradient updates, and provide evidence for the lower bound. The formula for (dimensionless) is as follows:

[0099] ;

[0100] In the formula, Represents the variational distribution The expected value of each term is a dimensionless logarithmic probability expectation. For the log-likelihood of the observed data, assuming the observed noise follows a Gaussian distribution, then:

[0101] ;

[0102] In the formula, To observe the noise variance, the dimension is... , The dimensions are The result of dividing the two is dimensionless. The total number of elements in the tensor. Redundant endmembers are automatically closed through a sparsity-promoted truncated normal prior; an endmember is automatically removed when its abundance variance converges to near zero. Three-dimensional abundance tensor. Using matrix product state structure compression, the storage complexity is reduced from... Down to ,in Let be the size of the tensor space. As the key dimension, The accuracy and memory usage of reconstruction were determined through multiple experiments on hyperspectral optical images of different scenes. The experimental step size was 8, and the step size was selected when the reconstruction error tended to level off. value.

[0103] The physics-driven atmospheric correction framework couples the radiative transfer equation with a water bio-optical model to estimate aerosol optical thickness under the assumption of near-zero water radiance in the shortwave infrared band. (Dimensionless) and aerosol type, then substitute into the radiative transfer equation to invert the water radiance. (unit: Nearshore mixed pixels use a linear unmixing model to separate the water radiation component from the land radiation component. Apparent radiance of band pixels (unit: Decomposed into:

[0104] ;

[0105] In the formula, The proportion of the radiation component of the water body (dimensionless). The proportion of the land radiation component (dimensionless) satisfies and , For the first Water radiance of the band (unit: ), For the first Land radiance in band (unit: ), For the first Band radiation error term (unit: The dimensions on both sides of the equal sign are... Solar flare contamination pixels were identified and removed using polarization images; the aerosol type discrimination threshold was determined through experimental analysis of 30–50 historical cloudless images of the study area, and iterative optimization was performed 3–5 times.

[0106] Sea state estimates are derived from aerosol optical thickness. With water radiance Joint calculation, including turbidity estimates Compared with chlorophyll concentration estimates The formula for estimating turbidity is as follows:

[0107] ;

[0108] In the formula, Red band water radiance (unit: ), Green band water radiance (unit: ), It is a dimensionless ratio. The optical thickness of the aerosol is dimensionless. and Here is the dimensionless adjustment coefficient, all of which are NTU, so that... The dimension is NTU. and The data were obtained by fitting measured water quality data and remote sensing data of the study area.

[0109] The specific implementation method of step S03 is as follows.

[0110] The cross-modal cyclic distillation model consists of a fully polarimetric synthetic aperture radar (SAR) polarimetric feature flow, a hyperspectral optical multispectral feature flow, and a cross-modal interactive kernel. The SAR polarimetric feature flow concatenates the scattering entropy feature map, scattering angle feature map, and inverse entropy feature map with the statistical feature map extracted by the convolutional neural network at the channel level to form a fused SAR feature map. (Dimensionless feature tensor). The hyperspectral optical multispectral feature current receiver retrieves the water radiance inversion results, and uses a spectral attention module to dynamically weight each band. Band weight From aerosol optical thickness With water radiance The jointly calculated sea state estimate is driven by the following formula:

[0111] ;

[0112] In the formula, For the first Sensitivity coefficient of turbidity to band (dimension: The data was obtained through regression analysis of multiple measured data points in the study area. For band summation index, For dimensionless normalized weights, The product is dimensionless, ensuring the dimensionless nature of the independent variable in the exponential function. After adaptively suppressing bands contaminated by water optical effects, a multispectral feature map is formed using hyperspectral optical multispectral characteristic flow. (Dimensionless characteristic tensor).

[0113] Cross-modal interactive kernels enable soft-label cyclic distillation, and the polarization characteristic flow of fully polarized synthetic aperture radar is achieved. The layer logit output is Hyperspectral optical multispectral characteristic flow The layer logit output is , The network layer number is used. The soft label distillation loss is expressed using KL divergence, as shown in the following formula:

[0114] ;

[0115] In the formula, The temperature coefficient is The softening softmax function, By different on the validation set The value was determined by ablation experiments. Let KL divergence be a dimensionless property. This is the temperature compensation coefficient (dimensionless). This represents dimensionless distillation loss. The fusion layer employs attention-gated dynamic weighting, and the fusion features... The calculation formula is as follows:

[0116] ;

[0117] In the formula, The fusion weights for the polarization characteristic flow of the fully polarized synthetic aperture radar are dimensionless. The fusion weights for hyperspectral optical multispectral feature flows (dimensionless) satisfy the following conditions: , , , All dimensions are dimensionless feature tensors, and the gating weights are determined by the Monte Carlo random deactivation uncertainty estimate. Monte Carlo random deactivation is applied to a certain feature flow. Sub-random reasoning, uncertainty estimation To output the variance multiple times, the formula is as follows:

[0118] ;

[0119] In the formula, For the first The output of the sub-random inactivation inference (dimensionless). for The mean of the outputs (dimensionless). For the index of the number of random samples, This is a dimensionless uncertainty estimate; high-uncertainty characteristic flows automatically have reduced fusion weights. Fusion adjustment index. The calculation formula is as follows:

[0120] ;

[0121] In the formula, The current small-batch loss (dimensionless) of the polarization characteristic flow of a fully polarized synthetic aperture radar. The current mini-batch loss (dimensionless) for hyperspectral optical multispectral feature flow. The number of random inactivation samples in Monte Carlo. It is a dimensionless attenuation factor. This is a dimensionless fusion adjustment index. When... hour, Adjusted to and Lowered by 0.05 to 0.10; when hour, and Remain unchanged; when hour, Adjusted to and Adjust by 0.05–0.10; renormalize after each adjustment. Gradient normalization coefficients for multi-task branches Calculate using the following formula:

[0122] ;

[0123] In the formula, For the first The gradient vector (dimensionless) of each task branch on the current mini-batch. for The Euclidean norm (dimensionless). This is the mean (dimensionless) of the gradient norm of each task branch. This is the task branch number, with values ​​of 1, 2, and 3. These are dimensionless normalized coefficients. Total training loss. The formula is as follows:

[0124] ;

[0125] In the formula, For the first Cross-entropy loss of each task branch (dimensionless). For the first The intersection-union ratio (IUU) loss (dimensionless) for each task branch is specifically: ,in For the first Intersection over Union (IoU) ratio (dimensionless) between predicted bounding boxes and ground truth bounding boxes for task branches. This is the distillation loss weighting coefficient (dimensionless), with an empirical value of 0.5–1.0. The total loss is dimensionless. The optimizer employs an adaptive moment estimation optimization algorithm with weighted decay, and the initial learning rate ranges from [value missing]. The value is determined by the learning rate grid search experiment on the validation set loss.

[0126] The specific implementation method of step S04 is as follows.

[0127] Construct a regional ocean simulation system or a finite-volume coastal ocean model to numerically simulate tides, currents, and changes in seabed topography, and output a natural dynamic prediction field. Ensemble Kalman filtering assimilation constructs the covariance matrix between the prediction set and the observation set, using minimizing the prediction error covariance as the criterion, to integrate the natural dynamic prediction field with the multi-temporal remote sensing observation field. Fusion, outputting the residual signal after removing the natural dynamic process. ,in For the set of Kalman filter analysis fields, For phase index, , , , The dimensions are consistent, and are the same as those of remote sensing observations.

[0128] The specific implementation method of step S05 is as follows.

[0129] For dual-phase residual signals and Construct empirical probability distributions respectively and (All are dimensionless), among which and Given two time phases as time indices, the change detection problem is reconstructed into an optimal transmission problem, with the transmission cost matrix... Defined by the feature distance of physical constraints, where for The number of samples, for Number of samples, elements for No. One sample and No. Physical feature distance (dimensionless normalized distance) between samples. Entropy regularization of the optimal transport plan matrix. This was obtained by solving the following problem:

[0130] ;

[0131] In the formula, For and For a set of transmission plans distributed at the edge, This is the matrix inner product (dimensionless). The Shannon entropy (dimensionless) of the transmission plan. For the transmission plan matrix Line 1 Column elements (dimensionless). The regularization intensity parameter (dimensionless) is determined by adaptive estimation of the reference area of ​​natural variations within the scene. It is determined through 20–50 iterations of experimental calibration, using the mean and standard deviation of the transmission cost of the reference area as upper and lower bounds. The dimensionless probabilistic mass transfer weight matrix is ​​used. The above problem employs the Sinkhorn-Knoop algorithm for stabilization in the logarithmic domain, and leverages GPU parallelization to accelerate convergence. Transfer plan matrix. Decomposed into cyclic components With net transfer component The formula is as follows:

[0132] ;

[0133] In the formula, the cyclic component for Its transpose's symmetrical component (dimensionless) corresponds to a bidirectional reversible natural fluctuation; net transfer component The antisymmetric part (dimensionless) corresponds to unidirectional and irreversible changes caused by human activity. for The transpose matrix of the vector. The transmission cost of the net transferred components is used as the output index of the intensity of human activity change, as shown in the following formula:

[0134] ;

[0135] In the formula, This is a dimensionless normalized physical feature distance matrix. The dimensionless net transfer component. It is a dimensionless index of the intensity of changes in human activities.

[0136] The specific implementation method of step S06 is as follows.

[0137] Based on the human activity change intensity index map, target detection results, and semantic segmentation results, combined with the human activity category identification results output in step S03, the multi-layer results are overlaid and analyzed according to spatial location to generate a monitoring report on human activity changes in marine protected areas. The report includes a distribution map of human activity types, a change intensity classification map, and early warning information for key areas.

[0138] To better understand and implement this invention, the following is a specific application scenario of this invention, Example 2:

[0139] To verify the effectiveness of the invention, technicians constructed a test scenario for a marine protected area, covering an area of ​​approximately 1200 square kilometers. The protected area contains human-caused targets such as aquaculture cages, patrol boats, and subsea pipelines, and is also subject to significant natural tidal and ocean current dynamics due to the monsoon influence. The fully polarimetric synthetic aperture radar (SAR) imagery used in the test was acquired by a spaceborne sensor, with an imaging interval of 16 days, fully polarized mode, and a spatial resolution of 5m. The hyperspectral optical imagery has 242 bands, covering a spectral range of 400 to 2500 nm, with a spatial resolution of 30m. Shortwave infrared data was acquired concurrently with the hyperspectral optical imagery. The dataset covers three typical meteorological conditions: clear sky, thin clouds, and post-typhoon conditions. The post-typhoon condition is a key scenario for verifying the ability to decouple natural dynamics.

[0140] First, step S01 is executed to perform Cloude-Pottier polarization decomposition on the fully polarimetric synthetic aperture radar image, extracting the scattering entropy feature map, scattering angle feature map, and inverse entropy feature map, such as... Figure 2 The diagram shows the structure of a cross-modal cyclic distillation model, illustrating how two feature flows achieve soft-label cyclic distillation through a cross-modal interaction kernel. In the aquaculture cage region, the scattering entropy... The values ​​are concentrated between 0.6 and 0.8, reflecting the composite scattering characteristics of the cage structure for radar waves; in open water areas, the scattering entropy... A value below 0.3 indicates a simple scattering mechanism dominated by Bragg scattering.

[0141] Subsequently, step S02 is executed, and the hyperspectral optical image is decomposed into hybrid pixels using a Bayesian variational nonnegative matrix factorization tensor network algorithm, with the matrix product state bond dimension... Through multiple sets of experiments The time reconstruction error decreases gradually, therefore, select... As the final value, the corresponding storage complexity is... Down to The aerosol type discrimination threshold was determined through experimental analysis and iterative optimization four times using 45 historical cloudless images of the protected area. The water radiance inversion results showed that the estimated chlorophyll concentration in the aquaculture area was significantly higher than that in the surrounding open sea, and the estimated turbidity was significantly higher in the nearshore reef area. This result is consistent with historical water quality monitoring records for this sea area. The number of endmembers for mixed pixel decomposition was automatically determined to be four using a sparse-promoted truncated normal prior, corresponding to four types of endmembers: water, aquaculture substrate, reef, and terrestrial vegetation, without the need for manual pre-setting.

[0142] Execute step S03, feeding the outputs of steps S01 and S02 into the cross-modal cyclic distillation model. The model training employs a two-stage strategy, and the distillation temperature coefficient is determined in the validation set ablation experiments. The initial learning rate is set to The fusion weight adjustment function dynamically adjusts the number of Monte Carlo random deactivation samples and the fusion weights during training based on the mini-batch loss of the two feature streams. The upper adjustment threshold of 0.6 and the lower adjustment threshold of 0.3 were determined through grid search experiments. Table 1 shows the target detection results, semantic segmentation results, and human activity category recognition results, listing the detection count for each human activity category.

[0143] Table 1. Statistics on the number of tests conducted by each activity category

[0144]

[0145] Execute step S04, using the regional ocean simulation system as the core, input the seabed topography data and the meteorological forcing field during the typhoon's passage, perform numerical simulation of tides and ocean currents, and output the natural dynamic prediction field, such as... Figure 3 The diagram illustrates the ensemble Kalman filter assimilation process between the natural dynamic prediction field and the multi-temporal remote sensing observation field. The ensemble Kalman filter assimilation uses a prediction ensemble with a size of 50 to fuse the natural dynamic prediction field and the multi-temporal remote sensing observation field, outputting a residual signal. In the post-typhoon scenario, approximately 62% of the variation area in the original multi-temporal image difference map originates from changes in sea surface scattering characteristics caused by the typhoon. After ensemble Kalman filter assimilation, this natural dynamic contribution to the residual signal is effectively suppressed.

[0146] Execute step S05, and process the residual signal using a regularized optimal transmission multi-temporal distribution alignment algorithm, such as... Figure 4 The diagram illustrates the decomposition of the transmission plan matrix into cyclic and net transfer components. The cyclic component corresponds to residual natural fluctuations, while the net transfer component corresponds to changes caused by human activities. The regularization intensity was determined by statistically analyzing the characteristic distribution changes of three known natural dynamic reference zones within the scene, and calibrated through 35 iterative experiments. The Sinkhorn-Knopp algorithm performs stabilization calculations in the logarithmic domain and is accelerated through parallelization using a graphics processor, achieving a convergence iteration count of approximately 200. The output results of the human activity change intensity index map show that areas with higher change intensity are concentrated at the expansion boundaries of aquaculture cages and the locations of newly added illegal structures, such as... Figure 5 The image shows a visualization of the intensity index of changes in human activities.

[0147] like Figure 6The image shows a spatial distribution comparison of the scattering entropy feature map, scattering angle feature map, and anti-entropy feature map. It reveals the differences in the three polarization feature dimensions among different man-made targets, verifying the effectiveness of Cloude-Pottier polarization decomposition in target differentiation. Figure 7 The image shows the abundance map of endmembers in the mixed pixel decomposition. The spatial distribution of the abundance of the four types of endmembers clearly reflects the spatial pattern of water bodies, aquaculture substrates, reefs and terrestrial vegetation within the protected area, which is consistent with the results of historical land cover surveys.

[0148] The results from steps S03 to S05 are spatially overlaid and combined with the results of human activity category identification to generate a regulatory report. The report marks the spatial location and quantitative values ​​of the intensity of change of six suspected illegal structures. For example... Figure 8 The image shows a spatial distribution map of areas affected by changes in human activity as reported in the regulatory report. Different color depths correspond to different levels of change intensity.

[0149] To further evaluate the performance of the cross-modal cyclic distillation model under single-mode missing conditions, technicians simulated a scenario where hyperspectral optical images were completely obscured by clouds and fog. In this case, the model switched to single-flow inference mode, relying solely on the fully polarimetric synthetic aperture radar polarimetric feature flow to reconstruct the feature representation of the missing mode through the distillation path. The identification results for aquaculture cages and patrol boats still maintained high reliability. Figure 9 The diagram shows a comparison of target detection results between single-stream and dual-stream inference modes. Figure 10 The diagram shows a comparison of the intensity index of human activity changes under different meteorological conditions. The results of the conditions after the typhoon passed are highly consistent with those of the clear sky conditions, which verifies the suppression effect of the two-level decoupling system on the background of strong natural disturbances. This is consistent with the quantitative statistical results in Table 2.

[0150] Table 2. Statistical Table of Human Activity Intensity Indicators under Different Meteorological Conditions

[0151]

[0152] Compared with traditional methods for detecting pixel difference changes, this invention constructs a two-level decoupling system based on the signal generation mechanism. The numerical model assimilation directly introduces the physical mechanisms of natural dynamics such as tides and ocean currents into the signal separation process, so that the natural dynamic components are removed in the residual signal output stage, rather than relying on post-processing thresholds to segment the mixed signal. The regularized optimal transmission framework performs quality transmission analysis on residual fluctuations from the probability distribution level. The physically constrained transmission cost matrix gives the algorithm inherent robustness to radiated noise. The separation of cyclic components and net transfer components realizes the essential decoupling of natural fluctuations and changes caused by human activities from the information theory level, so that the monitoring results still have high reliability under extreme natural disturbance conditions.

[0153] It should be noted that the variables involved in this invention are explained in detail in Tables 3 and 4.

[0154] Table 3. Variable Explanation Table (Part 1)

[0155]

[0156] Table 4. Variable Explanation Table (Part Two)

[0157]

[0158] The above description is merely a specific embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any changes or substitutions that can be easily conceived by those skilled in the art within the scope of the technology disclosed in the present invention should be included within the scope of protection of the present invention.

Claims

1. A remote sensing identification and monitoring method for changes in human activities in marine protected areas, characterized in that, Includes the following steps: Multi-source remote sensing data of marine protected areas were acquired, and Cloude-Pottier polarization decomposition was performed on the fully polarimetric synthetic aperture radar imagery to extract scattering entropy feature maps, scattering angle feature maps, and inverse entropy feature maps. Hybrid pixel decomposition was performed on hyperspectral optical images using a Bayesian variational nonnegative matrix decomposition tensor network algorithm. A physical-driven atmospheric correction framework was constructed by combining shortwave infrared band data, and the water radiance inversion results were output. The scattering entropy feature map, scattering angle feature map, inverse entropy feature map and water radiance inversion results are input into the cross-modal cyclic distillation model, and the outputs target detection results, semantic segmentation results and human activity category recognition results are output. A framework for assimilating ocean dynamics numerical models and remote sensing observation data is constructed. The natural dynamic prediction field and the multi-temporal remote sensing observation field are assimilated by ensemble Kalman filtering to output the residual signal. The residual signal is processed by the regularized optimal transmission multi-temporal distribution alignment algorithm, which decomposes the transmission plan matrix into cyclic components and net transfer components, and outputs a human activity change intensity index map. Based on the intensity index map of human activity changes, target detection results, and semantic segmentation results, a monitoring report on human activity changes in marine protected areas is generated.

2. The remote sensing identification and monitoring method according to claim 1, characterized in that, The Cloude-Pottier polarization decomposition specifically involves decomposing the coherence matrix of a fully polarimetric synthetic aperture radar image into eigenvalues. The scattering entropy characterizes the degree of scattering randomness, the scattering angle characterizes the type of dominant scattering mechanism, and the anti-entropy characterizes the relative strength of the secondary dominant scattering mechanism and the tertiary scattering mechanism.

3. The remote sensing identification and monitoring method according to claim 2, characterized in that, The physical-driven atmospheric correction framework specifically couples the radiative transfer equation with a water body bio-optical model, uses shortwave infrared data to estimate aerosol optical thickness and aerosol type under the assumption of near-zero water radiance, and then substitutes the aerosol optical thickness and aerosol type into the radiative transfer equation to invert water radiance.

4. The remote sensing identification and monitoring method according to claim 3, characterized in that, The Bayesian variational nonnegative matrix factorization tensor network algorithm specifically organizes hyperspectral optical image data into a spatial-spectral-temporal three-dimensional tensor, applies a Dirichlet distribution prior based on marine optical physics to the endmember matrix, applies a Markov random field prior to the abundance tensor, and performs variational inference through mean-field approximation decomposition of the posterior.

5. The remote sensing identification and monitoring method according to claim 4, characterized in that, In the Bayesian variational nonnegative matrix factorization tensor network algorithm, the tensor network part uses a matrix multiplication state structure to compress the three-dimensional abundance tensor. The key dimension value of the matrix multiplication state is determined by multiple sets of experiments on the reconstruction accuracy and memory usage of hyperspectral optical images of different scenes.

6. The remote sensing identification and monitoring method according to claim 5, characterized in that, The cross-modal cyclic distillation model consists of three parts: a fully polarimetric synthetic aperture radar polarimetric feature stream, a hyperspectral optical multispectral feature stream, and a cross-modal interactive kernel. The cross-modal interactive kernel implements a soft-label cyclic distillation mechanism, and the two feature streams complete mutual supervised iteration within each mini-batch training.

7. The remote sensing identification and monitoring method according to claim 6, characterized in that, The fusion layer of the cross-modal cyclic distillation model adopts attention-gated dynamic weighted fusion. The sum of the fusion weights of the fully polarimetric synthetic aperture radar polarimetric feature stream and the hyperspectral optical multispectral feature stream is 1. The gating weights are jointly determined by the Monte Carlo random inactivation uncertainty estimates of the two feature streams.

8. The remote sensing identification and monitoring method according to claim 7, characterized in that, The fusion weight is determined by the fusion weight adjustment function, which takes the number of Monte Carlo random deactivation samples, the current small batch loss of the fully polarimetric synthetic aperture radar polarimetric feature flow, and the current small batch loss of the hyperspectral optical multispectral feature flow as inputs, calculates the fusion adjustment index, and adjusts the number of Monte Carlo random deactivation samples and the fusion weight according to the relationship between the fusion adjustment index and the upper and lower adjustment thresholds.

9. A computer-readable storage medium, characterized in that, The computer-readable storage medium stores program instructions, which, when executed in a computer, are used to perform a remote sensing identification and monitoring method for changes in human activities in marine protected areas as described in any one of claims 1-8.

10. A remote sensing identification and monitoring system for changes in human activities in marine protected areas, characterized in that, The system comprises the computer-readable storage medium of claim 9, wherein the system is a computer, the computer-readable storage medium is disposed within the system, and the system is provided with a microprocessor that executes program instructions stored in the computer-readable storage medium.

Images