A ship GNSS spoofing detection method based on timing constraint prior correlation weight

Abstract

The application relates to the technical field of satellite navigation and artificial intelligence, and discloses a ship GNSS spoofing detection method based on a time sequence constraint prior correlation weight, which comprises the following steps: attitude estimation is carried out by using an inertial measurement unit and a compass, and the accelerometer observation noise is adaptively adjusted based on a gravity amplitude deviation; the carrier coordinate system acceleration measured by the IMU is converted to a navigation coordinate system, and a GNSS second-order differential acceleration and an IMU acceleration residual sequence are constructed; a detection algorithm framework based on a Transformer is proposed, the framework firstly introduces a time sequence constraint prior correlation weight matrix, then combines a sequence correlation matrix and the prior correlation matrix, uses a symmetric KL divergence to depict correlation differences, and finally combines the reconstruction error and the correlation differences to obtain an abnormal score, which is used for detecting GNSS spoofing attacks. The application can effectively detect GNSS spoofing attacks against ships and has good detection precision and stability.

Description

Technical Field

[0001] This invention relates to the fields of satellite navigation and artificial intelligence technology, and more particularly to a method for detecting ship GNSS deception based on temporal-constrained prior association weights. Background Technology

[0002] Global Navigation Satellite System (GNSS) is a crucial infrastructure for ship positioning, navigation, and timing, widely used in tasks such as autonomous navigation, track keeping, and berthing / unberthing assistance. In ship electronic navigation systems, GNSS provides critical spatiotemporal reference for route planning, electronic chart display, radar fusion, and warning decision-making; the reliability of its output information directly affects ship navigation safety and the stability of maritime traffic operations.

[0003] However, GNSS signals are weak in power and highly vulnerable to interference and deception attacks. Deception attacks, in particular, induce receivers to output incorrect navigation results by transmitting forged satellite signals. These attacks are highly covert and misleading, potentially causing ship track deviations, course deviations, and berthing errors, and in severe cases, even collisions and groundings. Therefore, researching GNSS deception detection methods suitable for ship scenarios is of great significance.

[0004] Existing spoofing detection methods mainly include signal authentication, signal feature analysis, raw observation detection, and consistency detection with external auxiliary sensors. Among these, methods based on signal authentication and feature analysis require modifications to receiver hardware or signal systems, resulting in high engineering costs. Methods based on raw observations are susceptible to noise, multipath propagation, and environmental variations. Consistency detection methods based on auxiliary sensors are mostly geared towards land-based scenarios, while the motion characteristics of the target vehicle and environmental disturbances differ in marine scenarios. Therefore, the application of existing methods in marine scenarios still requires targeted adaptation.

[0005] In recent years, time-series anomaly detection based on deep learning has provided a new approach for GNSS deception detection. The Transformer architecture has strong temporal dependency modeling capabilities, and the Anomaly Transformer is suitable for multi-source temporal anomaly identification. However, when directly applied to ship GNSS deception detection, it still suffers from problems such as the difficulty of adapting prior association structures to the gradual accumulation of errors and insufficient support for continuous online detection. Therefore, it is necessary to propose a high-precision GNSS deception detection method for maritime scenarios. Summary of the Invention

[0006] To address the aforementioned technical problems, a ship GNSS spoofing detection method based on temporal-constrained prior association weights is proposed. This method aims to solve the problem that GNSS is vulnerable to spoofing attacks during ship navigation, and that existing technologies have low accuracy in detecting GNSS spoofing in marine environments.

[0007] The technical means employed in this invention are as follows: A method for detecting ship GNSS deception based on temporal-constrained prior association weights includes the following steps: S1: The angular velocity measurement value of the gyroscope based on the IMU, the acceleration measurement value of the accelerometer and the bow observation value of the compass are used to perform attitude estimation using the error state extended Kalman filter to obtain the ship attitude estimation result; S2: Using the ship attitude estimation results, calculate the rotation matrix from the navigation coordinate system to the carrier coordinate system; based on the rotation matrix, project the three-axis acceleration measured by the IMU in the carrier coordinate system to the navigation coordinate system to obtain the IMU acceleration observation in the navigation coordinate system; perform gravity compensation on the IMU acceleration observation to obtain the IMU linear acceleration in the navigation coordinate system. S3: Perform second-order difference calculation on the GNSS position sequence output by the receiver to obtain the GNSS acceleration in the navigation coordinate system; then, subtract the GNSS acceleration at the corresponding time from the IMU linear acceleration in the navigation coordinate system, extract the northward acceleration residual and the eastward acceleration residual to form a two-dimensional feature, and truncate it through a multivariate time series sliding window to form the input sequence; S4: For the input sequence, a deception detection model with temporal constraint prior correlation weights is constructed based on Transformer. The deception detection model includes a data encoding module, a multi-head attention mechanism module, and a decoder module. The data encoding module uses one-dimensional convolution and positional encoding to map the input sequence to the latent space, enabling the deception detection model to perceive and capture the temporal information of each data point in the input sequence. The multi-head attention mechanism module extracts the correlation weights between each time step in the input sequence and outputs the sequence correlation matrix and a high-dimensional latent space representation, while constructing a temporal constraint prior correlation matrix that approximates GNSS deception features. The decoder module reconstructs the high-dimensional latent space representation output by the multi-head attention mechanism module back to the original physical feature space through a linear projection layer. S5: Train the deception detection model using the input sequence, then process the navigation data to be detected according to the procedures described in S1 to S3 to obtain the input sequence to be detected; input the input sequence to be detected into the trained deception detection model; perform deception detection by evaluating anomaly scores, the anomaly scores including reconstruction error and association difference, the reconstruction error being determined by the mean square error between the reconstructed sequence output by the deception detection model and the input sequence, the association difference being determined by the symmetric KL divergence between sequence association and prior association under the multi-head attention mechanism, the anomaly score being obtained by combining the reconstruction error and the association difference; compare the anomaly score with a detection threshold adaptively determined based on the statistical characteristics of the initial normal time period, if the anomaly score exceeds the detection threshold, it is marked as the beginning of a GNSS deception attack.

[0008] Furthermore, in S1, during the state prediction stage, the gyroscope is subjected to zero bias compensation and the attitude is updated based on the result of the zero bias compensation; during the measurement update stage, the accelerometer observation vector is normalized to obtain the gravity direction estimate, and an adaptive noise covariance based on the gravity amplitude deviation is introduced to update the gravity direction.

[0009] Furthermore, the implementation process of "introducing adaptive noise covariance based on gravity amplitude deviation for gravity direction update" is as follows: Set minimum standard deviation, maximum standard deviation, and threshold values; the minimum standard deviation corresponds to the background noise level of the normalized directional observation of the accelerometer under stationary or stable navigation conditions; the maximum standard deviation corresponds to the upper limit of the observation uncertainty under dynamic conditions, which is used to significantly reduce the observation weight under normal maneuver conditions; the threshold values ​​are preset with reference to the upper quantile of the normalized amplitude deviation of the ship under typical sea conditions. To effectively suppress the observation weights at high dynamic moments, the adaptive noise standard deviation is dynamically adjusted between the minimum and maximum standard deviations based on the gravity amplitude deviation with a square growth law, and the adaptive standard deviation is calculated, thereby obtaining the adaptive noise covariance matrix. Determine whether the gravity amplitude deviation is less than or equal to the threshold value; if the gravity amplitude deviation is less than or equal to the threshold value, then perform gravity direction update using the adaptive noise covariance matrix; if the gravity amplitude deviation is greater than the threshold value, then skip this gravity direction update.

[0010] Furthermore, in S2, the mathematical expression for "forming the input sequence by truncating through a multivariate time series sliding window" is:

[0011]

[0012] in, This represents the input sequence, with a length of... ; for The two-dimensional acceleration residual vector at time step, whose elements are the northward acceleration residuals. and eastward acceleration residual .

[0013] Furthermore, in S4, the process of "constructing a temporal-constrained prior correlation matrix for approximate GNSS deception features" is as follows: The multi-head attention mechanism module includes multiple parallel attention heads, each with an effective bandwidth at each time step, used to control the diffusion degree of the weights of the prior correlation matrix generated by that attention head at the corresponding time step. Since the input sequence is constructed based on a multivariate time series sliding window, a GNSS spoofing attack will cause changes in the residual data at the end of the multivariate time series sliding window, significantly reducing the cross-time step correlation between the residual data at the end of the window and the residual data at the beginning of the window. To reflect this cross-time step correlation, the prior correlation matrix is ​​designed with a diagonal trend structure that is wider at the beginning and narrower at the end. The diagonal bandwidth at the beginning of the window is larger to accommodate the time correlation of normal acceleration residuals, while the bandwidth at the end gradually shrinks to show the local correlation characteristics of the end samples. The expression for the baseline bandwidth, which monotonically changes with time and exhibits a trend of wider at the beginning and narrower at the end, is defined as follows:

[0014] in, For time step index, The length of the sliding window for the multivariate time series; p Shape parameters that represent how bandwidth changes over time; and This refers to bandwidth hyperparameters. For the first i The baseline bandwidth for each location; To allow the deception detection model to adaptively fine-tune without disrupting the trend, a bandwidth modulation parameter predicted by the linear mapping layer for encoded features is introduced. This parameter, combined with an activation function, constructs a bounded modulation factor. Within a set upper and lower limit, the baseline bandwidth is adaptively fine-tuned to obtain the effective bandwidth for each attention point at different locations. The bounded modulation factor is defined as follows:

[0015]

[0016] in, For attention head index; These are the modulation parameters predicted by the network; and Modulation range; The slope; For the first Height in position Effective bandwidth; To ensure the prior correlation matrix is ​​symmetric, a geometric mean is applied to the effective bandwidth; the paired bandwidths are obtained using the geometric mean, expressed as:

[0017] in, Indicates position Symmetric bandwidth; After obtaining the paired bandwidths, a Gaussian similarity is constructed based on the distance between time points, expressed as:

[0018] in, These are the elements of the prior similarity matrix; To avoid an overly sharp distribution in the prior correlation matrix, a uniform mixture term is introduced, expressed as:

[0019] in, For elements The prior similarity matrix formed by; To achieve a uniform bottom ratio; It is a matrix of all 1s; Diagonal scaling and row normalization are performed on the prior similarity matrix to generate a probabilistic prior association matrix:

[0020]

[0021] in, For the first h Prior association matrix of the size; To scale vectors A diagonal matrix composed of main diagonal elements; for The Middle i Row vectors; for The Middle i Line 1 j Column elements; To The result of normalization is easier to correlate with the sequence. Perform distributional variation measurement; Middle elements Iterative updates are performed as follows:

[0022] in, t This represents the number of iterations.

[0023] Furthermore, in S5, the "training the deception detection model using the input sequence" employs a two-stage alternating adversarial training strategy, the process of which is as follows: For each batch of data, the training process is divided into two stages, with the loss functions for the two stages being:

[0024]

[0025] in, Indicates reconstruction error; The associated difference loss term is calculated using symmetric KL divergence. This is the prior correlation matrix; It is a sequence incidence matrix; It is a weighting coefficient; This indicates that gradient calculation has stopped; This represents the parameters that influence the prior association matrix; This represents the parameters that affect the sequence correlation matrix; In the first stage of the training process, the relevant parameters of sequence association are fixed, and only the network parameters that generate prior associations are updated to minimize the association difference between prior associations and sequence associations, so that the prior association matrix fits the current attention pattern. In the second stage of the training process, the relevant parameters of the prior association are fixed, and only the network parameters of the generated sequence association are updated to maximize the association difference between the prior association and the sequence association, while minimizing the reconstruction error, thereby encouraging the model to capture data-specific association structures. The definition of the reconstruction error is as follows:

[0026] in, Output the feature matrix of the sequence to the model; Input the feature matrix of the sequence into the model; W This represents the length of the sliding window.

[0027] Furthermore, the anomaly score is composed of reconstruction error and correlation difference, defined as follows: During the testing phase, to obtain the time-by-time reconstruction loss of the samples:

[0028] Among them, coefficient This indicates that the average is taken across two acceleration feature dimensions; This means taking the diagonal elements of a matrix to form a vector; Meanwhile, in order to obtain the time-by-time correlation difference weights in the samples, the average is taken in the encoder layer dimension during the testing phase;

[0029] in, L This represents the number of encoder layers. Indicates the first l The correlation difference vector of each coding layer; During the testing phase, weights need to be assigned based on the correlation differences at each time step to reflect the degree of correlation deviation at each time step. The correlation differences are then negative and normalized using the Softmax function to obtain the correlation weight coefficients:

[0030] Multiplying the correlation weight coefficient by the value at the last moment of the reconstructed error sequence yields the GNSS deception anomaly score at the current moment:

[0031] in, Representing vectors The final time component, Representing vectors The final time component.

[0032] Compared with the prior art, the present invention has the following advantages: This invention utilizes an Inertial Measurement Unit (IMU) and a compass for attitude estimation, and adaptively adjusts the accelerometer observation noise based on gravity amplitude deviation. It transforms the carrier coordinate system acceleration measured by the IMU to the navigation coordinate system, constructing a sequence of GNSS second-order differential acceleration and IMU acceleration residuals. A Transformer-based detection algorithm framework is proposed. This framework first introduces a time-constrained prior correlation weight matrix, then combines the sequence correlation matrix and the prior correlation matrix, using symmetric KL divergence to characterize correlation differences. Finally, it combines the reconstruction error with the correlation differences to derive anomaly scores for detecting GNSS spoofing attacks. This invention can effectively detect GNSS spoofing attacks targeting ships, exhibiting good detection accuracy and stability. Attached Figure Description

[0033] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0034] Figure 1 This is an overall flowchart of a ship GNSS deception detection method based on temporal constraint prior association weights according to the present invention. Figure 2 This is an overall architecture diagram of a ship GNSS deception detection method based on temporal constraint prior association weights according to the present invention; Figure 3 This is a diagram of the Transformer neural network architecture based on temporal constraint prior association weights. Detailed Implementation

[0035] It should be noted that, unless otherwise specified, the embodiments and features described in the present invention can be combined with each other. The present invention will now be described in detail with reference to the accompanying drawings and embodiments.

[0036] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. The following description of at least one exemplary embodiment is merely illustrative and is in no way intended to limit the present invention or its application or use. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0037] It should be noted that the terminology used herein is for the purpose of describing particular embodiments only and is not intended to limit the scope of exemplary embodiments according to the invention. As used herein, the singular form is intended to include the plural form as well, unless the context clearly indicates otherwise. Furthermore, it should be understood that when the terms "comprising" and / or "including" are used in this specification, they indicate the presence of features, steps, operations, devices, components, and / or combinations thereof.

[0038] Unless otherwise specifically stated, the relative arrangement, numerical expressions, and values ​​of the components and steps described in these embodiments do not limit the scope of the invention. It should also be understood that, for ease of description, the dimensions of the various parts shown in the drawings are not drawn to actual scale. Techniques, methods, and devices known to those skilled in the art may not be discussed in detail, but where appropriate, such techniques, methods, and devices should be considered part of the specification. In all examples shown and discussed herein, any specific values ​​should be interpreted as merely exemplary and not as limitations. Therefore, other examples of exemplary embodiments may have different values. It should be noted that similar reference numerals and letters in the following figures denote similar items; therefore, once an item is defined in one figure, it need not be further discussed in subsequent figures.

[0039] In the description of this invention, it should be understood that the orientation or positional relationship indicated by directional terms such as "front, back, up, down, left, right", "horizontal, vertical, horizontal" and "top, bottom" is generally based on the orientation or positional relationship shown in the accompanying drawings, and is only for the convenience of describing this invention and simplifying the description. Unless otherwise stated, these directional terms do not indicate or imply that the device or element referred to must have a specific orientation or be constructed and operated in a specific orientation, and therefore should not be construed as a limitation on the scope of protection of this invention. The directional terms "inner" and "outer" refer to the inner and outer contours relative to the outline of each component itself.

[0040] For ease of description, spatial relative terms such as "above," "over," "on the upper surface of," "above," etc., are used herein to describe the spatial positional relationship of a device or feature as shown in the figures to other devices or features. It should be understood that spatial relative terms are intended to encompass different orientations in use or operation besides the orientation of the device as described in the figures. For example, if the device in the figures is inverted, a device described as "above" or "above" other devices or structures would subsequently be positioned as "below" or "under" other devices or structures. Thus, the exemplary term "above" can include both "above" and "below." The device may also be positioned in other different ways (rotated 90 degrees or in other orientations), and the spatial relative descriptions used herein will be interpreted accordingly.

[0041] Furthermore, it should be noted that the use of terms such as "first" and "second" to define components is merely for the purpose of distinguishing the corresponding components. Unless otherwise stated, the above terms have no special meaning and therefore should not be construed as limiting the scope of protection of this invention.

[0042] Example 1: like Figure 1 As shown in the figure, this embodiment provides a ship GNSS spoofing detection method based on temporal constraint prior association weights. The execution architecture diagram of this method is as follows. Figure 2 As shown, this method mainly includes the following steps: S1: Acquire acceleration, angular velocity, and heading data from shipborne navigation sensors and perform ship attitude estimation. Based on angular velocity measurements from the IMU gyroscope, acceleration measurements from the accelerometer, and heading observations from the compass, an error-state extended Kalman filter is used for attitude estimation to obtain the ship attitude estimation result. In the state prediction stage, zero-bias compensation is performed on the gyroscope, and attitude updates are completed based on the results of zero-bias compensation. In the measurement update stage, to eliminate the influence of acceleration amplitude, the accelerometer observation vector is normalized to obtain the gravity direction estimate. An adaptive noise covariance based on gravity amplitude deviation is introduced for gravity direction updates, thereby improving the stability of ship attitude estimation.

[0043] Specifically, the implementation process of gravity direction update based on adaptive noise covariance of gravity amplitude deviation is as follows: A minimum standard deviation, a maximum standard deviation, and a threshold value are set. The minimum standard deviation corresponds to the background noise level of the normalized direction observation of the accelerometer under stationary or stable navigation conditions. The maximum standard deviation corresponds to the upper limit of the observation uncertainty under dynamic conditions, used to significantly reduce the observation weight under normal maneuvering conditions. The threshold value is preset with reference to the upper quantile of the normalized amplitude deviation of the ship under typical sea states. To effectively suppress the observation weight at high dynamic moments, the adaptive noise standard deviation is dynamically adjusted between the minimum and maximum standard deviations based on the gravity amplitude deviation using a square growth law, calculating the adaptive standard deviation, and then obtaining the adaptive noise covariance matrix. It is then determined whether the gravity amplitude deviation is less than or equal to the threshold value: if yes, the gravity direction measurement update is performed using the adaptive noise covariance matrix; if not, the current gravity direction update is skipped.

[0044] Specifically, the shipborne autonomous navigation sensors are first modeled, including gyroscope observation models, compass observation models, and accelerometer observation models; the expression for the gyroscope observation model is:

[0045]

[0046] in, These are gyroscope observations; This represents the actual angular velocity in the carrier coordinate system. The zero-bias term of the gyroscope is modeled as a random walk; and All are Gaussian white noise; This refers to the IMU sampling time.

[0047] The expression for the compass observation model is:

[0048] in, These are compass readings; This represents the ship's actual bow heading. It is Gaussian white noise; This is the sampling time for the compass.

[0049] The expression for the accelerometer observation model is:

[0050] in, These are accelerometer measurements; This represents the actual acceleration in the carrier's coordinate system. This refers to the gravitational acceleration in the carrier's coordinate system. The total observed noise of the accelerometer; considering the ship's low maneuverability and high inertia in low wave class and steady-state navigation, the visible .

[0051] The attitude estimation process employs an error-state extended Kalman filter. To obtain the attitude transfer matrix, the attitude angle error and the gyroscope bias error are estimated, and the error state is defined as follows:

[0052] in, The angle error between the estimated attitude and the true attitude; Let be the zero bias error of the gyroscope. Then the covariance matrix can be defined as... .

[0053] In the state prediction phase, after debiasing the gyroscope, the attitude is updated via quaternion integration. For covariance prediction, the process noise and its discretized approximation are defined as:

[0054]

[0055] in, The white noise variance of the gyroscope angular velocity measurement; The variance of the random walk drive noise for zero bias of the gyroscope; It is a three-dimensional identity matrix; This is the IMU sampling period.

[0056] During the measurement update phase, to eliminate the influence of accelerometer amplitude, only the direction of the gravity vector is predicted, and the accelerometer observation vector is normalized:

[0057] The predicted direction of gravity is:

[0058] in, It is derived from the current attitude quaternion The rotation matrix from the navigation coordinate system to the vehicle coordinate system is calculated.

[0059] To adjust the confidence level of the accelerometer gravity observations during ship motion, the adaptive noise covariance is specified here as follows:

[0060]

[0061] in, and These are the minimum and maximum standard deviations, respectively. Defined as gravity amplitude deviation; It is a set threshold value; This is the accelerometer noise covariance matrix, used in the calculation of the innovation covariance during Kalman filtering. Update the gravity observation vector in real time; otherwise, skip the current gravity update.

[0062] The bow measurement is Output from the compass. The predicted bow heading is:

[0063]

[0064] in, This is the bow vector in the navigation coordinate system; It is the forward unit vector in the carrier coordinate system.

[0065] S2: Perform coordinate system transformation for IMU acceleration based on the ship attitude estimation results. Using the ship attitude estimation results, calculate the rotation matrix from the navigation coordinate system to the carrier coordinate system; based on the rotation matrix, project the three-axis acceleration measured by the IMU in the carrier coordinate system onto the navigation coordinate system to obtain the IMU acceleration observation in the navigation coordinate system; perform gravity compensation on the IMU acceleration observation to obtain the IMU linear acceleration in the navigation coordinate system.

[0066] S3: Extract the acceleration residual sequence and construct a multivariate time series sliding window. Perform second-order difference calculation on the GNSS position sequence output by the receiver to obtain the GNSS acceleration in the navigation coordinate system; then, subtract the GNSS acceleration at the corresponding time from the IMU linear acceleration in the navigation coordinate system to extract the northward and eastward acceleration residuals to form a two-dimensional feature, and truncate it through a multivariate time series sliding window to form the input sequence.

[0067] Specifically, using the current moment as a baseline, multiple consecutive historical data points, including the current moment, are extracted to form a single time window. By sliding the window step by step, the original long-segment time series is transformed into a series of overlapping local time series, thereby constructing a time-by-time input sequence.

[0068]

[0069] in, This represents the input sequence, with a length of... ; for The two-dimensional acceleration residual vector at time step, whose elements are the northward acceleration residuals. and eastward acceleration residual .

[0070] S4: For the input sequence, a deception detection model with temporal constraint prior association weights is constructed based on Transformer. The deception detection model includes a data encoding module, a multi-head attention mechanism module, and a decoder module. The data encoding module uses one-dimensional convolution and positional encoding to map the input sequence to the latent space, enabling the deception detection model to perceive and capture the temporal information of each data point in the input sequence. The multi-head attention mechanism module extracts the association weights between each time step in the input sequence and outputs the sequence association matrix and a high-dimensional latent space representation, while constructing a temporal constraint prior association matrix that approximates GNSS deception features. The decoder module reconstructs the high-dimensional latent space representation output by the multi-head attention mechanism module back to the original physical feature space through a linear projection layer.

[0071] Specifically, the process of "constructing a temporal-constrained prior correlation matrix for approximate GNSS deception features" is as follows: The multi-head attention mechanism module includes multiple parallel attention heads, each with an effective bandwidth at each time step, used to control the diffusion degree of the weights of the prior correlation matrix generated by that attention head at the corresponding time step. Since the input sequence is constructed based on a multivariate time series sliding window, a GNSS spoofing attack will cause changes in the residual data at the end of the multivariate time series sliding window, significantly reducing the cross-time step correlation between the residual data at the end of the window and the residual data at the beginning of the window. To reflect this cross-time step correlation, the prior correlation matrix is ​​designed with a diagonal trend structure that is wider at the beginning and narrower at the end. The diagonal bandwidth at the beginning of the window is larger to accommodate the time correlation of normal acceleration residuals, while the bandwidth at the end gradually shrinks to show the local correlation characteristics of the end samples. The expression for the baseline bandwidth, which monotonically changes with time and exhibits a trend of wider at the beginning and narrower at the end, is defined as follows:

[0072] in, For time step index, The length of the sliding window for the multivariate time series; p Shape parameters that represent how bandwidth changes over time; and This refers to bandwidth hyperparameters. For the first i The baseline bandwidth for each location; To allow the deception detection model to adaptively fine-tune without disrupting the trend, a bandwidth modulation parameter predicted by the linear mapping layer for encoded features is introduced. This parameter, combined with an activation function, constructs a bounded modulation factor. Within a set upper and lower limit, the baseline bandwidth is adaptively fine-tuned to obtain the effective bandwidth for each attention point at different locations. The bounded modulation factor is defined as follows:

[0073]

[0074] in, For attention head index; These are the modulation parameters predicted by the network; and Modulation range; The slope; For the first Height in position Effective bandwidth; To ensure the prior correlation matrix is ​​symmetric, a geometric mean is applied to the effective bandwidth; the paired bandwidths are obtained using the geometric mean, expressed as:

[0075] in, Indicates position Symmetric bandwidth; After obtaining the paired bandwidths, a Gaussian similarity is constructed based on the distance between time points, expressed as:

[0076] in, These are the elements of the prior similarity matrix; To avoid an overly sharp distribution in the prior correlation matrix, a uniform mixture term is introduced, expressed as:

[0077] in, For elements The prior similarity matrix formed by; To achieve a uniform bottom ratio; It is a matrix of all 1s; Diagonal scaling and row normalization are performed on the prior similarity matrix to generate a probabilistic prior association matrix:

[0078]

[0079] in, For the first h Prior association matrix of the size; To scale vectors A diagonal matrix composed of main diagonal elements; for The Middle i Row vectors; for The Middle i Line 1 j Column elements; To The result of normalization is easier to correlate with the sequence. Perform distributional variation measurement; Middle elements Iterative updates are performed as follows:

[0080] in, t This represents the number of iterations.

[0081] S5: Training and Evaluation Phase. During the training phase, a deception detection model is trained using an input sequence comprised of longitude and latitude calculated by the GNSS receiver during normal ship navigation, ship's heading output from the compass, and acceleration and angular velocity data output from the six-axis IMU. The input sequence to be detected is then input into the trained deception detection model. Deception detection is performed by evaluating anomaly scores, which include reconstruction error and correlation difference. The reconstruction error is determined by the mean square error between the reconstructed sequence output by the deception detection model and the input sequence. The correlation difference is determined by the symmetric KL divergence between sequence correlation and prior correlation under a multi-head attention mechanism. The anomaly score is obtained by combining these two factors.

[0082] Specifically, the "training of the deception detection model using the input sequence" employs a two-stage alternating adversarial training strategy, the process of which is as follows: For each batch of data, the training process is divided into two stages, with the loss functions for the two stages being:

[0083]

[0084] in, Indicates reconstruction error; The associated difference loss term is calculated using symmetric KL divergence. This is the prior correlation matrix; It is a sequence incidence matrix; It is a weighting coefficient; This indicates that gradient calculation has stopped; This represents the parameters that influence the prior association matrix; This represents the parameters that affect the sequence correlation matrix; In the first stage of the training process, the relevant parameters of sequence association are fixed, and only the network parameters that generate prior associations are updated to minimize the association difference between prior associations and sequence associations, so that the prior association matrix fits the current attention pattern. In the second stage of the training process, the relevant parameters of the prior association are fixed, and only the network parameters of the generated sequence association are updated to maximize the association difference between the prior association and the sequence association, while minimizing the reconstruction error, thereby encouraging the model to capture data-specific association structures. The definition of the reconstruction error is as follows:

[0085] in, Output the feature matrix of the sequence to the model; Input the feature matrix of the sequence into the model; W This represents the length of the sliding window.

[0086] The following is the calculation process for the correlation difference: First, define the first l Layer h The sequence association and prior association matrices of each attention head are as follows:

[0087] Before calculating the correlation variance, take the average for the long positions:

[0088]

[0089] Based on this, the definition is... To find the location i Key position j Attention weights; This represents the corresponding weight given prior.

[0090] To ensure numerical stability, we first perform probability normalization on each row:

[0091]

[0092] Definition of the first l Layer, query location i The differences in association are as follows:

[0093] The KL divergence calculation method is as follows:

[0094] The above correlation differences are averaged over the time and layer dimensions, and the correlation differences for the entire window are defined as follows:

[0095] During the testing phase, the abnormal score is compared with a detection threshold adaptively determined based on the statistical characteristics of the initial normal time period. If the abnormal score exceeds the detection threshold, it is marked as the beginning of a GNSS spoofing attack.

[0096] The anomaly score is composed of reconstruction error and correlation difference, and is defined as follows: During the testing phase, to obtain the time-by-time reconstruction loss of the samples:

[0097] Among them, coefficient This indicates that the average is taken across two acceleration feature dimensions; This means taking the diagonal elements of a matrix to form a vector; Meanwhile, in order to obtain the time-by-time correlation difference weights in the samples, the average is taken in the encoder layer dimension during the testing phase;

[0098] in, L This represents the number of encoder layers. Indicates the first l The correlation difference vector of each coding layer; During the testing phase, weights need to be assigned based on the correlation differences at each time step to reflect the degree of correlation deviation at each time step. The correlation differences are then negative and normalized using the Softmax function to obtain the correlation weight coefficients:

[0099] Multiplying the correlation weight coefficient by the value at the last moment of the reconstructed error sequence yields the GNSS deception anomaly score at the current moment:

[0100] in, Representing vectors The final time component, Representing vectors The final time component.

[0101] The detection threshold is determined by the maximum and minimum values ​​of the score benchmark segment, and the same threshold is applied across the entire test sequence. The specific setting method is as follows: Selecting the initial time before M The anomaly score corresponding to each time step is used as the baseline segment:

[0102] Record its maximum and minimum values ​​respectively:

[0103] The detection threshold is then set as follows:

[0104] in, This is the threshold adjustment coefficient.

[0105] Example 2: This embodiment further illustrates the solution and effects of the present invention through specific application examples.

[0106] This embodiment provides a ship GNSS spoofing detection method based on temporal-constrained prior association weights, and its algorithm architecture is as follows: Figure 2 As shown, this method consists of two stages: pose estimation and deception detection. The following sections will elaborate on the dataset construction and specific implementation steps of this method.

[0107] In this embodiment, the shipborne GNSS, IMU, and compass synchronization data collected during normal navigation are divided into a training set and a validation set in an 8:2 ratio. A deception attack is injected into another set of normal navigation data as a test set. The GNSS data includes receiver weekday seconds, longitude, latitude, and altitude data; the IMU data includes specific force acceleration and angular velocity in the carrier coordinate system; and the compass data includes the ship's heading. The GNSS sampling rate is 1Hz, the IMU sampling rate is 100Hz, and the compass sampling rate is 10Hz. The carrier coordinate system is selected as the forward right-downward direction, and the navigation coordinate system is selected as the northeast-northeast direction.

[0108] Step 1: Implementation of attitude estimation and coordinate transformation.

[0109] At the initial moment, the three-axis zero-bias prediction state of the gyroscope is set. Initialize to zero. Approximate the negative value of the specific force measured by the IMU accelerometer as the local gravity vector. The initial roll and pitch angles are then normalized and calculated; simultaneously, combined with the initial bow heading from the compass, they are transformed into initial state quaternions. The initial covariance matrix for the error state is then calculated. Set initial attitude error The standard deviation is 5°, and the initial gyroscope zero bias error is... The standard deviation is 1° / h.

[0110] Process noise includes white noise from gyroscope angular velocity measurements. and gyroscope zero-bias random walk drive noise The standard deviations of the two are respectively and Standard deviation of bow-head noise observed by compass The angle is 0.1°, and the local gravitational acceleration constant is... Set the acceleration update threshold. The value is set to 0.25, meaning that gravity direction measurements are only updated when the relative deviation between the acceleration modulus and the standard gravity is less than 25%; otherwise, the gravity direction updated at the previous moment is used. Simultaneously, a minimum effective standard deviation for gravity observation noise is set. The maximum effective standard deviation is 0.006. A value of 0.3 allows the filter to adaptively adjust the confidence level based on the intensity of the ship's motion acceleration within this range. In the time-progression process, state prediction and gravity direction updates are performed at the IMU frequency, while compass observation updates are executed at a frequency of 10 Hz.

[0111] Step 2: Extract the acceleration residual sequence and construct a multivariate time series sliding window.

[0112] like Figure 2As shown, the IMU acceleration is transformed from the carrier coordinate system to the navigation coordinate system. Then, to maintain temporal consistency with the GNSS acceleration obtained by second-order difference of the GNSS position sequence, the IMU acceleration in the navigation coordinate system is downsampled by 1Hz. The difference between the GNSS and IMU accelerations in the navigation coordinate system is then calculated, and the northward and eastward acceleration residuals are extracted to construct a two-dimensional feature sequence. To capture temporal anomalies, a sliding window length is set. W The value is 30, and the sliding step size is set to 1. This means that at any given time, the model will use the two-dimensional acceleration difference data from the past 30 time steps (i.e., 30 seconds) to form a 30×2 feature matrix X, which will be used as a single input to the deception detection model.

[0113] Step 3: Based on the characteristics of the ship's dynamic acceleration residual sequence, construct a deception detection model with temporal constraint prior association weights using Transformer.

[0114] The network architecture of the deception detection model is as follows Figure 3 As shown. The model input is the feature matrix X from step 2. The constructed deception detection model includes a data embedding layer, a multi-head self-attention encoder, and a feature reconstruction projection layer. The specific network structure and hyperparameter configuration are as follows: The data embedding layer consists of a one-dimensional convolutional layer and a positional encoding layer. The one-dimensional convolutional layer extracts local features from the input two-dimensional residual sequence, mapping the number of channels from 2 to 128. A kernel size of 3 is set, and circular padding is used to avoid loss of boundary information. The positional encoding employs standard sine and cosine positional encoding, fusing the extracted temporal and positional features before inputting them into the encoder.

[0115] Encoder layers L The hidden layer dimension is 3. d ff Set to 256; the network activation function is GELU, and in this embodiment, in order to retain all feature information, the random inactivation rate is set to 0.

[0116] In the multi-head attention mechanism of each coding layer, the number of heads H Set to 8. In addition to calculating sequence associations through standard scaled dot product, this embodiment constructs a temporally constrained prior association. Set the time step... i The varying reference bandwidth exhibits a non-linear attenuation trend, being wider at the beginning and narrower at the end. (Set shape parameters) p Setting the value to 2.2 makes the attenuation at the end steeper; setting the initial bandwidth The end bandwidth is 14.0. The value is 0.5. Latent variables are mapped to dynamic parameters with a dimension equal to the number of attention heads through linear projection, and bounded modulation is performed using an activation function. The modulation slope parameter is set. The value is 5, which limits the range of the modulation factor. It is 0.7. The value is set to 1.3, allowing for adaptive fine-tuning of the model without disrupting the overall convergence trend. The effective bandwidth at each point is taken. Then, through geometric mean Construct a symmetric diagonal bandwidth and calculate the elements in the Gaussian similarity matrix based on the adjusted distance. To prevent the prior correlation matrix from being overly sharp, a proportionality is introduced. Smooth mixing is performed with a globally uniform distribution of 0.05.

[0117] To ensure that the prior association matrix maintains symmetry while satisfying the condition that the sum of probabilities in each row is 1, this embodiment modifies the mixed similarity matrix. A Perform diagonal scaling iterations and row normalization. The scaling vector elements are initialized to 1, i.e. The number of iterations was set to 30, and a stable prior correlation matrix was finally output.

[0118] After processing by a three-layer encoder, the model extracts deep temporal dependency features. The 128-dimensional feature matrix output by the encoder is then input into a fully connected linear projection layer to re-map it back to two-dimensional physical space, resulting in the reconstructed northward and eastward acceleration residual sequences. This is used for subsequent calculation of the mean square reconstruction error.

[0119] Step 4: Experimental verification of the deception detection algorithm.

[0120] The model was trained using the Adam optimization algorithm with a learning rate of 10. -4 To prevent overfitting, a weight decay mechanism is introduced, with a decay coefficient set to 10. -6 The training batch size is set to 8, and the maximum number of training rounds is 50. An early stopping mechanism is introduced: if the validation set loss does not decrease for 5 consecutive rounds, training is terminated early and the current best model is saved.

[0121] During model training, a two-stage alternating adversarial training strategy is employed for each batch of data. Based on the calculated reconstruction loss, the first stage aims to minimize association discrepancies. In this stage, gradient propagation of sequence associations is blocked, and only the prior association network parameters are updated to fit the current sequence attention pattern to the prior association matrix. The second stage aims to maximize association discrepancies. In this stage, gradient propagation of prior associations is blocked, and the remaining network parameters are updated to expand the symmetric KL divergence between prior associations and sequence associations, while simultaneously minimizing the reconstruction error. To ensure that the reconstruction loss term and the association discrepancy loss term are on the same order of magnitude, this embodiment sets a correlation weight coefficient when adding or subtracting the two. 10 -4 .

[0122] During the testing phase, data injected with GNSS spoofing attacks was used as the test set. The test set was processed through the aforementioned attitude estimation, coordinate transformation, and sliding window truncation of residual sequences before being input into the trained spoofing detection model. To meet the requirements of real-time ship warnings, this embodiment does not calculate the mean score for the entire sliding window, but instead extracts the data features of the latest moment. After the model outputs, the mean square error of each time step is first calculated on the feature dimension to obtain the reconstruction error sequence within the window; simultaneously, the symmetric KL divergence between prior association and sequence association in the multi-head attention mechanism is calculated and averaged to obtain the time-by-time association difference sequence. Subsequently, the association difference sequence is negatively charged and normalized on the time dimension using the Softmax function, transforming it into time-by-time association weight coefficients. The last time-time element of this set of weight coefficients is multiplied by the last time-time element of the reconstruction error sequence to obtain the GNSS spoofing anomaly score for the current moment.

[0123] For anomaly detection, this embodiment employs a dynamic proportional threshold mechanism based on the initial normal state. In the initial testing phase, anomaly scores are extracted from the first 100 baseline time steps, and the maximum and minimum scores within this period are calculated to characterize the system's fluctuation baseline under normal navigation conditions. A globally fixed threshold is set using the product of the baseline maximum value, the fluctuation range, and a threshold adjustment coefficient. In this example, the threshold adjustment coefficient is set... The threshold is 30. In subsequent online detection, when the latest GNSS deception anomaly score exceeds this threshold, the ship is determined to be under deception attack and an alarm is triggered.

[0124] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, and not to limit them; although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some or all of the technical features; and these modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the scope of the technical solutions of the embodiments of the present invention.

Claims

1. A method for detecting ship GNSS deception based on temporal-constrained prior association weights, characterized in that, Includes the following steps: S1: The angular velocity measurement value of the gyroscope based on the IMU, the acceleration measurement value of the accelerometer and the bow observation value of the compass are used to perform attitude estimation using the error state extended Kalman filter to obtain the ship attitude estimation result; S2: Using the ship attitude estimation results, calculate the rotation matrix from the navigation coordinate system to the vehicle coordinate system; Based on the rotation matrix, the three-axis acceleration measured by the IMU in the carrier coordinate system is projected into the navigation coordinate system to obtain the IMU acceleration observation in the navigation coordinate system; Gravity compensation is applied to the IMU acceleration observations to obtain the IMU linear acceleration in the navigation coordinate system; S3: Perform second-order difference calculation on the GNSS position sequence output by the receiver to obtain the GNSS acceleration in the navigation coordinate system; then, subtract the GNSS acceleration at the corresponding time from the IMU linear acceleration in the navigation coordinate system, extract the northward acceleration residual and the eastward acceleration residual to form a two-dimensional feature, and truncate it through a multivariate time series sliding window to form the input sequence; S4: For the input sequence, a deception detection model with temporal constraint prior correlation weights is constructed based on Transformer. The deception detection model includes a data encoding module, a multi-head attention mechanism module, and a decoder module. The data encoding module uses one-dimensional convolution and positional encoding to map the input sequence to the latent space, enabling the deception detection model to perceive and capture the temporal information of each data point in the input sequence. The multi-head attention mechanism module extracts the correlation weights between each time step in the input sequence and outputs the sequence correlation matrix and a high-dimensional latent space representation, while constructing a temporal constraint prior correlation matrix that approximates GNSS deception features. The decoder module reconstructs the high-dimensional latent space representation output by the multi-head attention mechanism module back to the original physical feature space through a linear projection layer. S5: Train the deception detection model using the input sequence, then process the navigation data to be detected according to the process described in S1 to S3 to obtain the input sequence to be detected; input the input sequence to be detected into the trained deception detection model; perform deception detection by evaluating anomaly scores, wherein the anomaly scores include reconstruction error and association difference, wherein the reconstruction error is determined by the mean square error between the reconstructed sequence output by the deception detection model and the input sequence, and the association difference is determined by the symmetric KL divergence between sequence association and prior association under the multi-head attention mechanism, and the anomaly score is obtained by combining the reconstruction error and the association difference; The abnormal score is compared with a detection threshold adaptively determined based on the statistical characteristics of the initial normal period. If the abnormal score exceeds the detection threshold, it is marked as the beginning of a GNSS spoofing attack.

2. The ship GNSS spoofing detection method based on temporal constraint prior association weights according to claim 1, characterized in that, In S1, during the state prediction phase, the gyroscope is zero-bias compensation is performed and the attitude is updated based on the result of the zero-bias compensation. During the measurement update phase, the accelerometer observation vector is normalized to obtain the gravity direction estimate, and an adaptive noise covariance based on the gravity amplitude deviation is introduced to update the gravity direction.

3. The ship GNSS spoofing detection method based on temporal constraint prior association weights according to claim 2, characterized in that, The implementation process of "introducing adaptive noise covariance based on gravity amplitude deviation for gravity direction update" is as follows: Set minimum standard deviation, maximum standard deviation and threshold values; the minimum standard deviation corresponds to the background noise level of the normalized directional observation of the accelerometer under stationary or stable sailing conditions; the maximum standard deviation corresponds to the upper limit of the observation uncertainty under dynamic conditions, which is used to significantly reduce the observation weight under normal maneuvering conditions. The threshold value is preset with reference to the upper quantile of the normalized amplitude deviation of the vessel under typical sea conditions; To effectively suppress the observation weights at high dynamic moments, the adaptive noise standard deviation is dynamically adjusted between the minimum and maximum standard deviations based on the gravity amplitude deviation with a square growth law, and the adaptive standard deviation is calculated, thereby obtaining the adaptive noise covariance matrix. Determine whether the gravity amplitude deviation is less than or equal to the threshold value; If the gravity amplitude deviation is less than or equal to the threshold value, then the gravity direction is updated using the adaptive noise covariance matrix; If the gravity amplitude deviation is greater than the threshold value, then skip this gravity direction update.

4. The ship GNSS spoofing detection method based on temporal constraint prior association weights according to claim 1, characterized in that, In S2, the mathematical expression for "forming the input sequence by truncating a multivariate time series sliding window" is: in, This represents the input sequence, with a length of... ; for The two-dimensional acceleration residual vector at time step, whose elements are the northward acceleration residuals. and eastward acceleration residual .

5. The ship GNSS spoofing detection method based on temporal constraint prior association weights according to claim 1, characterized in that, In S4, the process of "constructing a temporal constraint prior correlation matrix for approximate GNSS deception features" is as follows: The multi-head attention mechanism module contains multiple parallel attention heads, each of which corresponds to an effective bandwidth at each time step, used to control the diffusion degree of the weights of the prior correlation matrix generated by the attention head at the corresponding time step. Since the input sequence is constructed based on a multivariate time series sliding window, the arrival of a GNSS spoofing attack will cause changes in the residual data at the end of the multivariate time series sliding window, resulting in a significant reduction in the cross-time step correlation between the residual data at the end of the multivariate time series sliding window and the residual data at the beginning of the window. To reflect the cross-time step correlation, the prior correlation matrix is ​​designed as a diagonal trend structure that is wide at the beginning and narrow at the end; the diagonal bandwidth at the beginning of the window is larger to accommodate the time correlation of normal acceleration residuals, and the bandwidth at the end gradually shrinks to show the local correlation characteristics of the end samples. The expression for the reference bandwidth, which monotonically changes with time and exhibits a trend of being wider at the beginning and narrower at the end, is defined as follows: in, For time step index, The length of the sliding window for the multivariate time series; p Shape parameters that represent how bandwidth changes over time; and This refers to bandwidth hyperparameters. For the first i The baseline bandwidth for each location; To allow the deception detection model to adaptively fine-tune without disrupting the trend, a bandwidth modulation parameter predicted by the linear mapping layer for encoded features is introduced. This parameter, combined with an activation function, constructs a bounded modulation factor. Within a set upper and lower limit, the baseline bandwidth is adaptively fine-tuned to obtain the effective bandwidth for each attention point at different locations. The bounded modulation factor is defined as follows: in, For attention head index; These are the modulation parameters predicted by the network; and Modulation range; The slope; For the first Height in position Effective bandwidth; To ensure the prior correlation matrix is ​​symmetric, a geometric mean is applied to the effective bandwidth; the paired bandwidths are obtained using the geometric mean, expressed as: in, Indicates position Symmetric bandwidth; After obtaining the paired bandwidths, a Gaussian similarity is constructed based on the distance between time points, expressed as: in, These are the elements of the prior similarity matrix; To avoid an overly sharp distribution in the prior correlation matrix, a uniform mixture term is introduced, expressed as: in, For elements The prior similarity matrix formed by; To achieve a uniform bottom ratio; It is a matrix of all 1s; Diagonal scaling and row normalization are performed on the prior similarity matrix to generate a probabilistic prior association matrix: in, For the first h Prior association matrix of the size; To scale vectors A diagonal matrix composed of main diagonal elements; for The Middle i Row vectors; for The Middle i Line number j Column elements; To The result of normalization is easier to correlate with the sequence. Perform distributional variation measurement; Middle elements Iterative updates are performed as follows: in, t This represents the number of iterations.

6. The ship GNSS spoofing detection method based on temporal constraint prior association weights according to claim 1, characterized in that, In step S5, the "training the deception detection model using the input sequence" employs a two-stage alternating adversarial training strategy, the process of which is as follows: For each batch of data, the training process is divided into two stages, with the loss functions for the two stages being: in, Indicates reconstruction error; The associated difference loss term is calculated using symmetric KL divergence. This is the prior correlation matrix; It is a sequence incidence matrix; It is a weighting coefficient; This indicates that gradient calculation has been stopped; The parameters that influence the prior association matrix; The parameters that affect the sequence correlation matrix; In the first stage of the training process, the relevant parameters of sequence association are fixed, and only the network parameters that generate prior associations are updated to minimize the association difference between prior associations and sequence associations, so that the prior association matrix fits the current attention pattern. In the second stage of the training process, the relevant parameters of the prior association are fixed, and only the network parameters of the generated sequence association are updated to maximize the association difference between the prior association and the sequence association, while minimizing the reconstruction error, thereby encouraging the model to capture data-specific association structures. The definition of the reconstruction error is as follows: in, Output the feature matrix of the sequence to the model; Input the feature matrix of the sequence into the model; W This represents the length of the sliding window.

7. The ship GNSS spoofing detection method based on temporal constraint prior association weights according to claim 1, characterized in that, The anomaly score is composed of reconstruction error and correlation difference, and is defined as follows: During the testing phase, to obtain the time-by-time reconstruction loss of the samples: Among them, coefficient This indicates that the average is taken across two acceleration feature dimensions; This means taking the diagonal elements of a matrix to form a vector; Meanwhile, in order to obtain the time-by-time correlation difference weights in the samples, the average is taken in the encoder layer dimension during the testing phase; in, L This represents the number of encoder layers. Indicates the first l The correlation difference vector of each coding layer; During the testing phase, weights need to be assigned based on the correlation differences at each time step to reflect the degree of correlation deviation at each time step. The correlation differences are then negative and normalized using the Softmax function to obtain the correlation weight coefficients: Multiplying the correlation weight coefficient by the value at the last moment of the reconstructed error sequence yields the GNSS deception anomaly score at the current moment: in, Representing vectors The final time component, Representing vectors The final time component.

Images