A Method for Detecting Abnormal Behavior of Ship Platforms Based on ADS-B Big Data

Abstract

This invention discloses an abnormal behavior detection method for ships based on an ADS-B big data platform, specifically relating to the field of ship navigation safety monitoring technology. The method involves acquiring automatic correlation surveillance broadcast data of ships in a target sea area, constructing a ship trajectory evolution matrix, and generating a ship micro-disturbance motion feature tensor. It also involves constructing a ship neighborhood interaction topology network and calculating the trajectory disturbance propagation weights between nodes to form an abnormal propagation probability matrix. Furthermore, it establishes a coupling response relationship between the local trajectory disturbance energy sequence and the neighborhood interaction intensity through temporal folding mapping, generating an abnormal behavior evolution prediction curve and calculating an abnormal energy aggregation index. Combined with the abnormal propagation probability matrix, it generates a risk distribution map of abnormal ship behavior, extracts implicit behavior indicators, and dynamically matches them with preset abnormal behavior patterns. This invention can effectively identify micro-scale abnormal navigation behavior in complex sea environments, improving the accuracy of abnormal ship behavior detection and early warning capabilities.

Description

Technical Field

[0001] This invention relates to the field of ship navigation safety monitoring technology, specifically to a method for detecting abnormal behavior of ships based on an ADS-B big data platform. Background Technology

[0002] With the widespread application of Automatic Identification Systems (AIS / ADS-B type broadcast navigation data) in ocean shipping, port scheduling, and maritime safety supervision, regulatory agencies can continuously acquire high-frequency broadcast data such as vessel position, heading, speed, and timestamps. However, in complex sea conditions or high-density waterways, some vessels may evade routine supervision by exhibiting slight heading fluctuations, short-period speed drifts, or abnormal track overlaps. Examples include short-term berthing and unberthing outside anchorages, covert transshipment, or probing approaches. These abnormal behaviors often manifest as track disturbances within extremely small scales and short time windows, making them difficult to identify using traditional methods based on rule-based thresholds or single track statistical features. In particular, they struggle to characterize the diffusion process of local track anomalies caused by multi-vehicle interactions, resulting in insufficient ability to identify early-stage anomalies. Summary of the Invention

[0003] The purpose of this invention is to provide an abnormal behavior detection method for ships based on ADS-B big data platform, so as to solve the shortcomings of the background technology.

[0004] To achieve the above objectives, the present invention provides the following technical solution: a method for detecting abnormal behavior of a ship platform based on ADS-B big data, comprising:

[0005] Acquire ADS-B broadcast data streams of ships in the target sea area, extract multi-dimensional track sequences, and construct the ship track evolution matrix through a sliding time window;

[0006] Based on the trajectory evolution matrix, the heading perturbation angular velocity, velocity fluctuation gradient, and trajectory curvature density are calculated and fused to generate the ship's perturbation motion feature tensor.

[0007] Based on the spatial distribution relationship of the perturbation motion feature tensor, a ship neighborhood interaction topology network is constructed, and the trajectory perturbation propagation weights between nodes are calculated to obtain the anomaly propagation probability matrix.

[0008] Based on the anomaly propagation probability matrix, the perturbation motion feature tensor is temporally folded and mapped to establish a coupled response model between the local track perturbation energy sequence and the neighborhood interaction intensity, thereby generating anomaly behavior evolution prediction curves.

[0009] The abnormal energy accumulation index of the target ship is calculated based on the abnormal behavior evolution prediction curve, and the abnormal energy accumulation index is nonlinearly coupled with the abnormal propagation probability matrix to generate a risk distribution map of abnormal ship behavior.

[0010] Multi-scale feature extraction is performed on the risk distribution map to obtain implicit behavioral indicators, including trajectory drift density, interactive perturbation spectrum and abnormal energy convergence, and dynamic matching is performed with preset abnormal behavior patterns.

[0011] When any implicit behavior indicator exceeds the set threshold, the abnormal behavior detection result of the corresponding ship and its associated ship interaction path are output.

[0012] Preferably, the step of constructing the ship's trajectory evolution matrix through a sliding time window includes:

[0013] Collect ADS-B broadcast data of ships in the target sea area over a continuous time period, merge and sort the data according to the ship's unique identification code, and extract latitude and longitude coordinates, heading angle, speed and timestamp information to form the original multidimensional track sequence;

[0014] The original multidimensional track sequence is subjected to time consistency correction and spatial coordinate unification processing. A standardized track sequence with uniform time intervals is generated by linear interpolation and abnormal jump point removal methods.

[0015] The standardized track sequence is segmented and sampled using a preset sliding time window. The changes in latitude and longitude, heading and speed within the same time window are arranged in a matrix according to time order to obtain a local track feature submatrix of the ship.

[0016] Multiple local track feature sub-matrices corresponding to continuous sliding time windows are superimposed and integrated according to the time dimension to construct a track evolution matrix that reflects the dynamic change law of short-term ship track.

[0017] Preferably, the steps of calculating the heading perturbation angular velocity, velocity fluctuation gradient, and trajectory curvature density based on the trajectory evolution matrix and generating the ship perturbation motion characteristic tensor include:

[0018] Extract the heading angle and speed change data of the continuous time series from the track evolution matrix, calculate the heading perturbation angular velocity sequence based on the heading angle difference and time interval between adjacent time sampling points, and calculate the speed fluctuation gradient sequence based on the speed change and time change rate.

[0019] The spatial displacement vector between adjacent track points is calculated based on the longitude and latitude changes in the track evolution matrix, and the trajectory curvature density sequence is calculated using a three-point trajectory fitting method.

[0020] The heading perturbation angular velocity sequence, velocity fluctuation gradient sequence, and trajectory curvature density sequence are normalized, and corresponding multidimensional feature matrices are constructed.

[0021] The multidimensional feature matrix is ​​reconstructed by tensor quantization according to the time dimension and the feature dimension to generate a perturbation motion feature tensor.

[0022] Preferably, the steps for constructing a ship neighborhood interaction topology network based on the spatial distribution relationship of the perturbation motion feature tensor include:

[0023] Based on the time index corresponding to the perturbation motion feature tensor, the longitude and latitude coordinates of each ship within the same time window are extracted, and the spatial distance between any two ships is calculated. Ship pairs with a spatial distance less than a preset neighborhood distance threshold are selected as candidate neighborhood ship sets. The similarity of the perturbation motion feature tensors of the ships in the candidate neighborhood ship sets is calculated. By calculating the feature distance between the three feature sequences of heading perturbation angular velocity, velocity fluctuation gradient, and trajectory curvature density, the perturbation motion association weights between ships are generated.

[0024] Using the unique identification code of a ship as a network node, the spatial distance relationship is used as the initial connection relationship, and the connection relationship is weighted according to the perturbation motion association weight to construct a ship neighborhood interaction relationship matrix; a ship neighborhood interaction topology network is generated according to the ship neighborhood interaction relationship matrix, where network nodes represent individual ships and network edge weights represent the perturbation motion interaction intensity between ships.

[0025] Preferably, the steps of calculating the propagation weights of trajectory perturbations between nodes and generating the anomaly propagation probability matrix include:

[0026] Based on the ship neighborhood interaction topology network, the connection edges between any two adjacent nodes are extracted, and the heading perturbation angular velocity, velocity fluctuation gradient, and trajectory curvature density feature values ​​of the corresponding ship in the perturbation motion feature tensor are read to construct the perturbation feature difference vector between node pairs; the trajectory perturbation similarity index between nodes is calculated based on the perturbation feature difference vector, and the trajectory perturbation similarity index is multiplied by the spatial adjacency weight in the ship neighborhood interaction topology network to obtain the trajectory perturbation propagation weight between nodes;

[0027] The trajectory disturbance propagation weights between all nodes are arranged in a matrix according to the node connection relationship to construct the trajectory disturbance propagation weight matrix between ship nodes; the trajectory disturbance propagation weight matrix is ​​subjected to normalized probability transformation to obtain the abnormal propagation probability matrix representing the diffusion trend of abnormal ship behavior.

[0028] Preferably, the steps of performing temporal folding mapping on the perturbation motion feature tensor based on the anomaly propagation probability matrix and generating anomaly behavior evolution prediction curves include:

[0029] Based on the anomaly propagation probability matrix, the propagation probability weights of each ship node in adjacent time windows are extracted, and the perturbation motion feature tensors of the corresponding ships are weighted and aggregated according to time order to obtain the track perturbation energy sequence characterizing the intensity of local motion perturbation of the ship.

[0030] The track disturbance energy sequence is synchronized and aligned with the neighborhood interaction intensity within the corresponding time window, and the disturbance energy values ​​of consecutive time windows are superimposed and mapped using a sliding time folding method to form a time evolution feature sequence.

[0031] Based on the time evolution characteristic sequence, a coupled response model between the local track disturbance energy sequence and the neighborhood interaction intensity is constructed, and the abnormal behavior evolution trend parameters are obtained by calculating the dynamic response coefficient between the two.

[0032] Based on the abnormal behavior evolution trend parameters, the disturbance energy changes in subsequent time windows are recursively calculated to generate a ship abnormal behavior evolution prediction curve.

[0033] Preferably, the steps for calculating the anomalous energy accumulation index of the target ship based on the anomalous behavior evolution prediction curve include:

[0034] Based on the abnormal behavior evolution prediction curve, the coupled response prediction values ​​within the continuous prediction time window are extracted, and the predicted perturbation energy sequence is constructed in chronological order.

[0035] Local peak identification is performed on the predicted perturbation energy sequence. Abnormal energy peak points are determined by calculating the energy change amplitude between adjacent time windows, and the distribution density of peak points in the predicted perturbation energy sequence is statistically analyzed to obtain candidate intervals for abnormal energy aggregation.

[0036] The time-weighted cumulative energy is calculated based on the disturbance energy value within the candidate interval of the abnormal energy accumulation, and the energy concentration coefficient is calculated in combination with the peak distribution density.

[0037] The time-weighted cumulative energy is coupled with the energy concentration coefficient to obtain the abnormal energy accumulation index, which represents the intensity of the abnormal behavior of the target ship.

[0038] Preferably, the step of performing nonlinear coupling calculations between the abnormal energy accumulation index and the abnormal propagation probability matrix to generate a risk distribution map of abnormal ship behavior includes:

[0039] Read the abnormal energy accumulation index of the target ship and extract the propagation probability sequence corresponding to the target ship from the abnormal propagation probability matrix. Use the propagation probability sequence as the spatial diffusion weight to perform neighborhood diffusion mapping on the abnormal energy accumulation index to obtain the initial abnormal energy diffusion vector.

[0040] The initial abnormal energy diffusion vector is subjected to nonlinear transformation, and an exponential amplification function is constructed to enhance the calculation of the high propagation probability region, thereby obtaining the risk intensity sequence of abnormal ship behavior.

[0041] The abnormal behavior risk intensity sequence of the ship is mapped to the grid coordinates of the target sea area according to the spatial location of the ship, and the cumulative value of risk intensity in each grid cell is calculated to generate a spatial risk intensity matrix;

[0042] Continuous spatial interpolation is performed on the spatial risk intensity matrix to construct a risk distribution map of abnormal ship behavior in the target sea area.

[0043] Preferably, the step of multi-scale feature extraction of the risk distribution map includes: constructing a multi-scale spatial analysis window based on the ship abnormal behavior risk distribution map, and calculating the risk intensity gradient change rate and spatial distribution density within different scale windows to obtain trajectory drift density features, which are used to characterize the degree of ship track deviation; extracting the risk intensity change sequence over time within the multi-scale spatial analysis window, and performing frequency domain transformation on the sequence to obtain interactive disturbance spectrum features; calculating the abnormal energy convergence index based on the spatial clustering degree of high-risk areas in the risk distribution map, and generating an abnormal energy convergence feature vector by statistically analyzing the spatial contraction rate of high-risk units in a continuous time window.

[0044] Preferably, the step of dynamically matching with a preset abnormal behavior pattern includes: fusing the trajectory drift density features, interactive disturbance spectrum features, and abnormal energy convergence feature vectors, and calculating the similarity with the behavioral features in the preset abnormal behavior pattern to achieve dynamic matching and identification of abnormal ship behavior patterns.

[0045] The technical effects and advantages provided by the present invention in the above technical solution are as follows:

[0046] 1. This invention constructs a track evolution matrix and a perturbation motion feature tensor based on ship automatic dependent surveillance broadcast data. This enables joint analysis of fine-grained motion characteristics such as ship heading perturbation angular velocity, velocity fluctuation gradient, and trajectory curvature density. Furthermore, it constructs a ship neighborhood interaction topology network and anomaly propagation probability matrix, thereby characterizing the propagation relationship of anomalous disturbances between ships in the spatial neighborhood. Compared to traditional methods based solely on single-ship trajectory thresholds, this invention can identify micro-scale track anomalies from the perspective of multi-ship interactions, significantly improving the ability to identify short-term, concealed anomalies. It is particularly suitable for anomaly detection in complex marine environments such as port peripheries, anchorage areas, and high-density waterways.

[0047] 2. This invention further analyzes the development trend of abnormal ship behavior by constructing abnormal behavior evolution prediction curves, abnormal energy accumulation indices, and risk distribution maps. It also obtains implicit behavioral indicators such as trajectory drift density, interaction disturbance spectrum, and abnormal energy convergence degree through multi-scale feature extraction methods. Furthermore, it achieves refined identification and spatial correlation analysis of abnormal behavior through abnormal behavior pattern matching and related ship interaction path extraction. This method can not only accurately identify abnormal behavior of individual ships but also reveal the diffusion path of abnormal behavior among neighboring ships, thereby providing a more forward-looking risk warning capability for maritime safety supervision and improving the intelligence level of maritime traffic management and safety supervision. Attached Figure Description

[0048] To more clearly illustrate the technical solutions in the embodiments of this application or the prior art, the drawings used in the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments recorded in this invention. For those skilled in the art, other drawings can be obtained based on these drawings.

[0049] Figure 1 This is a flowchart of an abnormal behavior detection method for a ship platform based on ADS-B big data, according to the present invention.

[0050] Figure 2 This is a flowchart of the method for generating a risk distribution map of abnormal ship behavior according to the present invention. Detailed Implementation

[0051] 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, not all embodiments. 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.

[0052] For examples, please refer to Figure 1 , 2 As shown in this embodiment, an abnormal behavior detection method for ships based on an ADS-B big data platform includes:

[0053] The system acquires ADS-B broadcast data streams of ships in the target sea area, extracts multi-dimensional track sequences, and constructs a track evolution matrix of ships through a sliding time window.

[0054] In this embodiment, firstly, the Automatic Dependent Surveillance (ADS) broadcast data of ships within a continuous time period in the target sea area is acquired. This ADS broadcast data is collected by a shore-based receiving station or satellite receiving terminal. Each data record contains at least a unique ship identifier, longitude coordinates, latitude coordinates, heading angle, speed, and timestamp information. Subsequently, all data records are categorized according to the unique ship identifier and sorted in ascending order based on the timestamp, thus obtaining a single ship's continuous navigation data sequence. After sorting, five fields—longitude coordinates, latitude coordinates, heading angle, speed, and timestamp—are extracted from each data record and combined in chronological order to form a multidimensional track sequence. The data vector at each moment is represented as F(t) = [LON(t), LAT(t), HDG(t), SPD(t), T(t)]; where LON(t) represents the longitude coordinate at time t, LAT(t) represents the latitude coordinate, HDG(t) represents the heading angle, SPD(t) represents the speed, and T(t) represents the timestamp. The above processing yields the original multidimensional track sequence arranged continuously in time, providing basic data for subsequent data correction processing.

[0055] After obtaining the original multidimensional track sequence, time consistency correction and spatial coordinate unification processing are performed on the original multidimensional track sequence. First, a unified time sampling interval Δt is set, which is set to 30 seconds according to the ship broadcast frequency. Then, the time difference ΔT between two adjacent data records is calculated. When ΔT is greater than the time sampling interval Δt, supplementary track points are generated using a linear interpolation algorithm. The linear interpolation calculation method is as follows: Where X represents longitude coordinates, latitude coordinates, heading angle, or speed parameters; ti and tj are the original data time points; and tk is the interpolation time point. After interpolation, abnormal jump points are removed from the spatial positions. Specifically, the spatial distance D between adjacent waypoints is calculated, and the theoretical maximum movement distance is calculated based on the speed SPD. When D > 2 × Dmax, the waypoint is identified as an outlier and deleted. After time interpolation and outlier removal, a standardized waypoint sequence with consistent time intervals and spatial continuity is obtained.

[0056] After obtaining the standardized track sequence, the standardized track sequence is segmented and sampled using a preset sliding time window. First, the sliding time window length W and sliding step size S are set, where the sliding time window length W is set to 10 minutes and the sliding step size S is set to 1 minute. For the track data within each sliding time window, the changes in longitude ΔLON, latitude ΔLAT, heading ΔHDG, and speed ΔSPD between adjacent track points are calculated sequentially. The ΔLON, ΔLAT, ΔHDG, and ΔSPD values ​​arranged in chronological order within the same sliding time window are combined to form a feature sequence, which is then matrixed according to the time dimension to form a ship local track feature submatrix M(k), where the rows of the matrix represent the time sequence position and the columns represent different track change characteristics.

[0057] After obtaining the ship local track feature sub-matrices corresponding to multiple sliding time windows, the ship local track feature sub-matrices are superimposed and integrated in the time dimension. Specifically, the multiple ship local track feature sub-matrices generated by consecutive sliding time windows are stacked in chronological order, and a time index is established in the third dimension of the matrix, thereby forming a three-dimensional track evolution matrix E. The structure of the track evolution matrix E is defined as: E={M(1),M(2),...,M(n)}; where M(k) represents the ship local track feature sub-matrix corresponding to the k-th sliding time window, and n represents the number of sliding time windows. The track evolution matrix constructed in the above manner can simultaneously reflect the dynamic evolution process of the ship's longitude, latitude, heading, and speed changes within multiple consecutive time windows, thus providing a continuous time-series data foundation for subsequent calculation of abnormal ship behavior characteristics.

[0058] Based on the trajectory evolution matrix, the heading perturbation angular velocity, velocity fluctuation gradient, and trajectory curvature density are calculated and fused to generate the ship's perturbation motion feature tensor.

[0059] In this embodiment, the ship's heading angle time series is extracted based on the track evolution matrix to calculate the heading perturbation angular velocity. The track evolution matrix contains heading angle data HDG(i) arranged in chronological order and corresponding timestamps T(i). First, the heading angle change ΔHDG(i) between adjacent time sampling points is calculated, and then the time interval ΔT(i) between adjacent time sampling points is calculated; the heading perturbation angular velocity ω(i) is calculated based on the above two parameters. Where ω(i) represents the rate of change of the heading angle per unit time, used to describe the degree of disturbance of the ship's heading over a short time scale. To avoid calculation errors caused by the heading angle spanning 360 degrees, when When the angle is greater than 180 degrees, the change in heading angle is normalized, i.e., ΔHDG(i) = ΔHDG(i) ± 360 degrees, so that the change in heading angle is within the range of -180 degrees to 180 degrees, thereby obtaining a continuous and stable heading perturbation angular velocity sequence.

[0060] After obtaining the heading perturbation angular velocity sequence, the velocity fluctuation gradient is calculated based on the velocity time sequence in the track evolution matrix. First, the velocity sequence SPD(i) is extracted, and the velocity change ΔSPD(i) between adjacent time sampling points is calculated. Then, the velocity change rate G(i) is calculated. Where ΔT(i) is the calculated time interval. To reflect the changing trend of velocity fluctuations in the time series, a local gradient is calculated for the velocity change rate sequence. Specifically, the gradient change value is calculated within a time window of length k, where the time window length k is set to 5 time sampling points. The velocity fluctuation gradient Grad(i) is calculated as follows: The velocity fluctuation gradient sequence obtained through the above calculations is used to characterize the intensity of acceleration and deceleration fluctuations of a ship over a short timescale.

[0061] After obtaining the heading perturbation angular velocity sequence and velocity fluctuation gradient sequence, the trajectory curvature density is further calculated based on the longitude and latitude coordinates in the trajectory evolution matrix. First, the longitude and latitude coordinates are converted into position coordinates X(i) and Y(i) in a plane coordinate system, calculated as follows: X(i) = R × cos(LAT(i)) × LON(i); Y(i) = R × LAT(i); where R represents the Earth's radius, taken as 6,371,000 meters. Then, three consecutive trajectory points P(i-1), P(i), and P(i+1) are selected, and the trajectory vectors V1 and V2 for two adjacent segments are calculated, expressed as: , The trajectory curvature K(i) is calculated using the cross product of the vectors and the magnitude of the vector: ;in This represents the absolute value of the cross product of two-dimensional vectors. The trajectory curvature density is defined as the change in curvature per unit track length, and it is calculated as follows: The trajectory curvature density sequence obtained through the above calculations is used to describe the degree of curvature of the ship's trajectory in local space.

[0062] After obtaining the heading perturbation angular velocity sequence, velocity fluctuation gradient sequence, and trajectory curvature density sequence, these three feature sequences are aligned using a unified time index and then numerically normalized. The normalization method employs a minimum-maximum value normalization algorithm, the calculation method of which is as follows: Where X(i) represents the original feature value, Xmin represents the minimum value in the feature sequence, and Xmax represents the maximum value in the feature sequence. This represents the normalized feature values. After normalization, the heading perturbation angular velocity sequence, velocity fluctuation gradient sequence, and trajectory curvature density sequence are matrix-arranged according to the time dimension, forming three sets of feature matrices. Subsequently, these three sets of feature matrices are stacked and combined in the third dimension to construct a three-dimensional structured data set, represented as T(i,j,k); where i represents the time series index, j represents the sliding time window position index, and k represents the feature type index, corresponding to the heading perturbation angular velocity, velocity fluctuation gradient, and trajectory curvature density, respectively. This forms a perturbation motion feature tensor used to describe the micro-scale motion disturbance state of the ship, providing basic feature data for subsequent ship abnormal behavior identification calculations.

[0063] Based on the spatial distribution relationship of the perturbation motion feature tensor, a ship neighborhood interaction topology network is constructed, and the trajectory perturbation propagation weights between nodes are calculated to obtain the anomaly propagation probability matrix.

[0064] In this embodiment, the position data of all ships within the same time window are extracted based on the time index corresponding to the perturbation motion feature tensor. Specifically, the longitude coordinate LON(i) and latitude coordinate LAT(i) of each ship within the corresponding time window are read from the perturbation motion feature tensor. Then, the longitude and latitude coordinates are converted into the geographic coordinate form required for planar distance calculation, and the spatial distance D(i,j) between any two ships is calculated using the spherical distance calculation method. Subsequently, a neighborhood distance threshold D0 is set. This threshold is obtained by statistically analyzing historical navigation data of the target sea area, specifically by calculating the distribution of distances between all ships and taking the 15th percentile as the neighborhood distance threshold. When D(i,j) ≤ D0, ship i and ship j are defined as candidate neighbor ship pairs, and all ship pairs that meet the condition constitute a candidate neighbor ship set.

[0065] After obtaining the candidate neighboring ship set, similarity is calculated for the perturbation motion feature tensors of the ships in the set to determine the perturbation motion association weights between ships. First, the heading perturbation angular velocity sequence, velocity fluctuation gradient sequence, and trajectory curvature density sequence of two ships in the candidate neighboring ship pair within the same time window are read. Then, the Euclidean distance between the three feature sequences is calculated as D1(i,j), D2(i,j), and D3(i,j), where ω represents the heading perturbation angular velocity sequence, Grad represents the velocity fluctuation gradient sequence, C represents the trajectory curvature density sequence, and k represents the time series index. Finally, the three distance values ​​are weighted and combined to obtain the comprehensive feature distance. Where α, β, and γ are weighting coefficients, set to 0.4, 0.3, and 0.3, respectively. Finally, the perturbation motion correlation weights are calculated using an exponential decay function. ; The perturbation motion correlation weight is used to represent the similarity degree of the perturbation motion states of two ships.

[0066] After obtaining the perturbation motion correlation weight, a node set is constructed with the ship unique identification code as the node identifier. For each pair of ships in the candidate neighborhood ship set, first establish an initial connection relationship based on the spatial proximity relationship, and convert the corresponding spatial distance into a spatial adjacency weight. The calculation method of the spatial adjacency weight is: ; where D(i,j) is the calculated spatial distance. Subsequently, the spatial adjacency weight and the perturbation motion correlation weight are multiplied to obtain the comprehensive interaction weight: A(i,j)=S(i,j)×W(i,j); Arrange the comprehensive interaction weights between all ships in the order of the ship unique identification code to construct a two-dimensional matrix A, where the matrix element A(i,j) represents the interaction intensity between ship i and ship j, thus forming a ship neighborhood interaction relationship matrix.

[0067] After obtaining the ship neighborhood interaction relationship matrix, a ship neighborhood interaction topology network is constructed according to the matrix structure. Specifically, map each ship unique identification code to a network node. When A(i,j) is greater than 0, establish a connection edge between node i and node j, and use A(i,j) as the weight of this connection edge. In order to eliminate the interference of weak association edges, a screening threshold T0 is set for the comprehensive interaction weight. The screening threshold is calculated by statistically averaging all non-zero values in the comprehensive interaction weight matrix. When A(i,j)≥T0, the connection relationship is retained, and when A(i,j)<T0, the connection relationship is deleted. Through the above processing, a ship neighborhood interaction topology network composed of a node set and weighted connection edges is obtained, where the nodes represent individual ships and the connection edge weights represent the perturbation motion interaction intensity between ships.

[0068] After constructing the ship neighborhood interaction topology network, extract the perturbation motion characteristics of the corresponding ships for any two adjacent nodes in the network. Specifically, read the heading perturbation angular velocity, speed fluctuation gradient, and trajectory curvature density eigenvalue in the perturbation motion feature tensor of the two ships, and calculate the difference between the corresponding features. The perturbation feature difference vector is expressed as: ; where this vector is used to characterize the difference degree of the perturbation motion states of two ships.

[0069] After obtaining the perturbation feature difference vector, calculate the trajectory perturbation similarity index between nodes according to the difference vector. First, calculate the modulus length L(i,j) of the perturbation feature difference vector, and then convert the difference value into a similarity index through a similarity function: Sim(i,j)=1 / (1+L(i,j)); Then perform a multiplication operation on the similarity index and the comprehensive interaction weight A(i,j) in the ship neighborhood interaction topology network: Where Pw(i,j) represents the propagation weight of trajectory disturbance between nodes. After obtaining the propagation weights of trajectory disturbance between all nodes, the propagation weights between nodes are arranged in order according to the unique identification code of the ship, constructing a two-dimensional matrix P. The element P(i,j) in matrix P represents the weight of the trajectory disturbance propagated from ship i to ship j. When there is no connection between nodes, the corresponding matrix element is set to 0. The complete trajectory disturbance propagation weight matrix is ​​constructed in the above manner.

[0070] After obtaining the trajectory perturbation propagation weight matrix, the matrix is ​​subjected to probability normalization. Specifically, each row of the matrix is ​​normalized so that the sum of the weights propagated from each node to its neighboring nodes is 1. The normalization calculation method is as follows: ;in The sum of all elements in the i-th row of the matrix is ​​represented by . The anomaly propagation probability matrix P' is obtained through the above calculation. The matrix element P'(i,j) represents the probability value of the trajectory disturbance of ship i propagating to ship j. This matrix is ​​used to describe the diffusion trend of abnormal ship behavior among neighboring ships.

[0071] Based on the anomaly propagation probability matrix, the perturbation motion feature tensor is temporally folded and mapped to establish a coupled response model between the local track perturbation energy sequence and the neighborhood interaction intensity, thereby generating anomaly behavior evolution prediction curves.

[0072] In this embodiment, the perturbation motion feature tensor is weighted and mapped according to the anomaly propagation probability matrix to construct the track perturbation energy sequence. First, the anomaly propagation probability P'(i,j) between the target ship node and its neighboring ship nodes is read from the anomaly propagation probability matrix. Then, three feature values ​​of the target ship in the perturbation motion feature tensor are extracted: heading perturbation angular velocity, velocity fluctuation gradient, and trajectory curvature density, denoted as ω(i,t), Grad(i,t), and C(i,t), respectively. Next, the local track perturbation energy E(i,t) at a single moment is calculated as follows: E(i,t) = a × |ω(i,t)| + b × |Grad(i,t)| + c × |C(i,t)|; where a, b, and c are feature weight coefficients, set to 0.4, 0.3, and 0.3, respectively. Finally, the neighborhood propagation weighted energy is calculated based on the anomaly propagation probability matrix. Where E(j,t) represents the track disturbance energy of neighboring ships at the same moment. The time-varying track disturbance energy sequence is obtained through the above calculation: S(i)={E'(i,1),E'(i,2),...,E'(i,n)}; this sequence is used to describe the change process of local motion disturbance intensity of the target ship under the influence of neighboring interaction.

[0073] After obtaining the track disturbance energy sequence, a time-series folding mapping process is performed on the sequence. First, a folding time window length W is set to 15 minutes, and the track disturbance energy sequence is divided into multiple time segments according to the time sampling interval. Then, the disturbance energy within each time segment is calculated using the folding accumulation method: F(i,k)=Σ(E'(i,t)); where t belongs to all time sampling points within the k-th time window. After completing the folding accumulation calculation, the folded energy values ​​are time-normalized using the following method: Where W represents the number of sampling points within the time window. The time-folded energy sequence is obtained through the above processing: This sequence is used to represent the evolution characteristics of the intensity of track disturbances of a ship within a continuous time window.

[0074] After obtaining the time-folded energy sequence, a coupled response model between the local track disturbance energy sequence and the neighborhood interaction intensity is further established. First, the neighborhood interaction intensity I(i,k) within each time window is calculated as follows: Where A(i,j) is the interaction weight in the ship neighborhood interaction matrix. Then, a coupled response function is constructed to characterize the dynamic relationship between disturbance energy and neighborhood interaction strength. The coupled response function is defined as: R(i,k)=F'(i,k)×I(i,k); where R(i,k) represents the coupled response value within the k-th time window. To characterize the response change trend, the response change rate between adjacent time windows is further calculated: The above calculations form a coupled response sequence. This sequence is used to describe the dynamic response relationship between local ship track disturbances and neighborhood interaction effects.

[0075] After obtaining the coupled response sequence, trend prediction calculations are performed on the coupled response sequence to generate a predicted curve for the evolution of anomalous behavior. First, a time series recursive prediction algorithm is used for prediction calculations. Specifically, the most recent L time windows are selected as prediction inputs, where L is set to 6. Then, the trend coefficient K(i) is calculated: Where ΔR(i,k) represents the rate of change of the coupled response. The coupled response value for future time windows is recursively predicted based on the trend coefficient. By repeating the above calculation over multiple future time windows, the predicted coupling response sequence can be obtained: Rpred(i)={Rpred(i,1),Rpred(i,2),...,Rpred(i,p)}; Finally, a continuously changing curve is plotted with the time window as the horizontal axis and the predicted coupling response value as the vertical axis, thereby generating a prediction curve for the evolution of abnormal ship behavior. This curve is used to characterize the evolution trend of abnormal ship behavior in future time periods.

[0076] The abnormal energy accumulation index of the target ship is calculated based on the abnormal behavior evolution prediction curve, and the abnormal energy accumulation index is nonlinearly coupled with the abnormal propagation probability matrix to generate a risk distribution map of abnormal ship behavior.

[0077] In this embodiment, the predicted values ​​of coupled responses within a continuous prediction time window are extracted based on the abnormal behavior evolution prediction curve. The abnormal behavior evolution prediction curve consists of a time series and corresponding predicted values ​​of coupled responses, where the time windows are arranged at fixed time intervals. First, the predicted value of the coupled response, Rpred(i,t), is read from the prediction curve, where i represents the target ship number and t represents the prediction time window index. Then, all the predicted values ​​of coupled responses are arranged in chronological order to form a predicted disturbance energy sequence: Spred(i) = {Rpred(i,1), Rpred(i,2), ..., Rpred(i,n)}; where n represents the number of prediction time windows. To ensure the comparability of energy values ​​between different ships, the predicted disturbance energy sequence is normalized using a minimum-maximum value normalization method. Where Smin represents the minimum value in the sequence, Smax represents the maximum value in the sequence, and Snorm(i,t) represents the normalized predicted perturbation energy value. The above processing yields a standardized predicted perturbation energy sequence, which is used for subsequent abnormal energy peak identification calculations.

[0078] After obtaining the standardized predicted perturbation energy sequence, local peaks are identified within the sequence. First, the energy variation amplitude between adjacent time windows is calculated: Then, an energy change threshold T1 is set. This energy change threshold is obtained through statistical analysis of historical track data. Specifically, the average value μ and standard deviation σ of historical disturbance energy changes are calculated, and then determined according to the formula T1 = μ + 1.5σ. When ΔS(i,t) > T1 and and When the time window t is reached, the energy value corresponding to the time window t is determined as the abnormal energy peak point.

[0079] After identifying all anomalous energy peaks, the density distribution of these peaks in the time series is statistically analyzed. A sliding statistical window of length K is defined, where K represents 5 time windows. The number of peaks, Pk, is calculated within each statistical window. When Pk ≥ 2, this time window is defined as a candidate interval for anomalous energy aggregation. The above steps yield a set of candidate intervals for anomalous energy aggregation.

[0080] Step 3 of the instruction manual: Calculate the time-weighted cumulative energy and the energy concentration coefficient.

[0081] After obtaining candidate intervals for abnormal energy accumulation, time-weighted cumulative calculation is performed on the disturbance energy within the candidate intervals. First, a time weighting function W(t) is defined as: W(t) = t / T; where t represents the time index within the interval, and T represents the length of the candidate interval. Then, the time-weighted cumulative energy within the candidate intervals is calculated: Where Snorm(i,t) is the normalized predicted perturbation energy value. After completing the time-weighted cumulative energy calculation, the energy concentration coefficient is further calculated. First, the average energy value Epeak of all peak points in the candidate interval and the average energy value Eavg of all energy values ​​in the candidate interval are calculated. Then, the energy concentration coefficient is calculated according to the following formula: Cenergy=Epeak / Eavg; the energy concentration coefficient is used to characterize the concentration of anomalous energy in the time series.

[0082] After obtaining the time-weighted cumulative energy and the energy concentration coefficient, the two are coupled to calculate the anomalous energy concentration index. The coupling calculation uses a product coupling function, and its calculation method is as follows: Where AEI(i) represents the anomalous energy accumulation index of the target vessel. To avoid extreme values ​​affecting subsequent calculations, the anomalous energy accumulation index is logarithmically compressed: The standardized abnormal energy accumulation index is obtained through the above calculations and is used to represent the intensity of abnormal behavior of the target ship.

[0083] After obtaining the abnormal energy accumulation index, it is coupled with the abnormal propagation probability matrix. First, the abnormal propagation probability matrix P'(i,j) is read, where i represents the target ship node and j represents neighboring ship nodes. Then, the data in the i-th row of the matrix is ​​extracted as the propagation probability sequence: Where m represents the number of neighboring ships. Then, the anomalous energy aggregation index AEInorm(i) is multiplied by the propagation probability sequence to obtain the initial anomalous energy diffusion vector: This vector represents the intensity of the diffusion of anomalous energy from the target ship to neighboring ships.

[0084] After obtaining the initial anomalous energy diffusion vector, a nonlinear transformation is applied to it to enhance the influence of the high propagation probability region. The nonlinear transformation uses an exponential amplification function, and its calculation formula is as follows: Where λ is the nonlinear amplification factor, and λ takes a value of 2. The risk intensity sequence of abnormal ship behavior is obtained through the above calculations: The risk intensity sequence represents the degree of risk that vessels in each neighboring region are affected by abnormal behavior.

[0085] After obtaining the risk intensity sequence, the risk intensity values ​​are mapped to the target sea area grid coordinates according to the ship's spatial location. First, the target sea area is divided into a regular grid, with each grid cell having a side length of 500 meters. Then, the grid cell containing the ship is determined based on its longitude and latitude coordinates, and the corresponding risk intensity value is accumulated into that grid cell. Let the grid cell number be g; then the risk intensity of that grid cell is calculated as: G(g) = ΣR(j); where R(j) is the risk intensity value of the ship falling into that grid cell. By accumulating the risk intensity values ​​for all ships, a two-dimensional spatial risk intensity matrix is ​​constructed.

[0086] After obtaining the spatial risk intensity matrix, continuous spatial interpolation is performed on the matrix to generate a risk distribution map. The interpolation method uses an inverse distance weighted interpolation algorithm, and its calculation method is as follows: Where Z(x) represents the risk value of the location to be calculated, Zi represents the known risk intensity value of the grid cell, di represents the distance between the location to be calculated and the center point of the grid cell, and p represents the distance weighting index, with a value of 2.

[0087] By performing continuous interpolation calculations over the entire target sea area, a continuous spatial risk intensity distribution is obtained. Subsequently, risk level intervals are divided according to the risk intensity values, and different risk levels are displayed with different colors, thus forming a risk distribution map of abnormal ship behavior, which is used to intuitively represent the spatial distribution of abnormal behavior within the target sea area.

[0088] Multi-scale feature extraction is performed on the risk distribution map to obtain implicit behavioral indicators, including trajectory drift density, interactive perturbation spectrum and abnormal energy convergence, and then dynamically matched with preset abnormal behavior patterns.

[0089] In this embodiment, a multi-scale spatial analysis is performed on the risk distribution map of abnormal ship behavior to extract trajectory drift density features. First, the target sea area corresponding to the risk distribution map is divided into spatial analysis windows of different scales, namely 500 meters, 1000 meters, and 2000 meters. Then, within each spatial analysis window, the risk intensity value G(g) of the grid cell is statistically analyzed, and the risk intensity gradient between adjacent grid cells is calculated. The risk intensity gradient is calculated as follows: Where G(g) represents the risk intensity of the current grid cell, G(gn) represents the risk intensity of the adjacent grid cells, and d represents the distance between the center points of two grid cells. Then, within each scale window, the number of grid cells with a risk intensity gradient greater than the threshold T2 is counted. The threshold T2 is obtained through statistical analysis of historical risk intensity gradient data, and its calculation method is as follows: Where μg represents the average historical risk intensity gradient, and σg represents the standard deviation. Finally, the ratio of the number of grid cells satisfying the conditions to the window area is calculated: Where Ng represents the number of grid cells that meet the conditions, and Aw represents the window area. The trajectory drift density feature obtained through the above calculation is used to represent the degree of deviation of the ship's track within a spatial range.

[0090] After obtaining the trajectory drift density characteristics, frequency domain analysis is performed on the risk intensity time series in the risk distribution map to extract the interactive perturbation spectrum characteristics. First, the risk intensity values ​​of each grid cell within a continuous time window are read in chronological order to construct the risk intensity time series: Subsequently, the Discrete Fourier Transform (DFT) algorithm was used to perform frequency domain transformation on the risk intensity time series. The DFT calculation method is as follows: Where t represents the time index, k represents the frequency index, and n represents the time series length. The amplitude of each frequency component is calculated as follows: Subsequently, the dominant frequency component fmax with the largest amplitude was identified, and the corresponding amplitude A(fmax) was used as the intensity index of the interactive disturbance spectrum. This index is used to represent the periodic variation characteristics of the interactive disturbance behavior between ships over time.

[0091] After obtaining the spectral characteristics of the interactive perturbation, spatial clustering analysis is performed on high-risk areas in the risk distribution map to calculate the convergence index of abnormal energy. First, a risk intensity threshold T3 is set. This threshold is determined by statistically analyzing all values ​​in the risk intensity matrix and taking their 90th percentile. When the risk intensity G(g) of a grid cell ≥ T3, the grid cell is defined as a high-risk cell. Subsequently, the spatial distribution area of ​​high-risk cells is statistically analyzed within a continuous time window. Let At represent the total area of ​​high-risk cells within time window t, then the spatial contraction rate between adjacent time windows is calculated as follows: Subsequently, the spatial contraction rate of the continuous time window was averaged and calculated: Where m represents the number of statistical time windows. The above calculation yields an anomaly energy convergence index, which is used to represent the concentration trend of anomaly energy within a spatial range.

[0092] After obtaining the trajectory drift density features, interactive perturbation spectrum features, and anomaly energy convergence index, feature fusion processing is performed on these three indices. First, a behavioral feature vector is constructed: V = [Ddrift, A(fmax), Cenergy]; then, a preset abnormal behavior pattern is established. This abnormal behavior pattern is constructed using historical abnormal trajectory data. Specifically, the same feature extraction process is performed on historical abnormal events to obtain multiple standard behavioral feature vectors, which are then used as pattern templates.

[0093] During dynamic matching, the cosine similarity calculation method is used to calculate the similarity between the current behavior feature vector and the pattern template: Where Vm represents the pattern template feature vector. If the calculated similarity value is greater than the matching threshold T4, the current ship behavior is determined to match the corresponding abnormal behavior pattern. The matching threshold T4 is set to 0.8 and determined through statistical analysis of historical abnormal behavior recognition experiments. The above method achieves dynamic matching and recognition of abnormal ship behavior patterns.

[0094] When any implicit behavior indicator exceeds the set threshold, the abnormal behavior detection result of the corresponding ship and its associated ship interaction path are output.

[0095] In this embodiment, threshold judgment processing is performed on three implicit behavioral indicators: trajectory drift density, interactive disturbance spectrum intensity, and anomalous energy convergence. First, long-term normal navigation data of the target sea area is collected, and statistical analysis is performed on the trajectory drift density, interactive disturbance spectrum intensity, and anomalous energy convergence under normal navigation conditions. Specifically, a historical data sample set is established for each indicator, and the corresponding indicator values ​​are recorded in chronological order. Then, the historical distribution characteristics of each indicator are calculated, and the corresponding anomaly judgment threshold is determined. The threshold is determined using the historical distribution quantile statistical method, that is, the historical samples are sorted according to their numerical values, and the high-order distribution interval is selected as the anomaly boundary. The trajectory drift density threshold is determined based on the higher proportion interval in the historical samples, the interactive disturbance spectrum intensity threshold is determined based on the high-frequency disturbance sample interval, and the anomalous energy convergence threshold is determined based on the energy concentration sample interval. After the thresholds are set, the three implicit behavioral indicators calculated at the current detection time are compared with their corresponding thresholds. When any indicator reaches or exceeds the corresponding threshold, the vessel is marked as a candidate for anomalous behavior.

[0096] After determining the thresholds for implicit behavioral indicators, vessels marked as candidates for anomalous behavior undergo behavioral confirmation processing. First, the vessel's trajectory drift density change records, interactive disturbance spectrum change records, and anomalous energy convergence change records are read within a continuous time window, and the consistency of these three indicators' trends is analyzed. Specifically, the stability of indicator changes within a continuous time window is statistically analyzed; when the same indicator remains within an anomalous range for multiple consecutive time windows, it is classified as a stable anomalous feature. Simultaneously, synchronous anomalous situations of multiple implicit behavioral indicators within the same time window are examined; when two or more indicators are simultaneously within an anomalous range within the same time window, the vessel's behavior is classified as confirmed anomalous behavior. This process avoids misjudgments caused by instantaneous fluctuations, thus obtaining reliable anomalous behavior detection results.

[0097] After confirming abnormal vessel behavior, the interaction paths of related vessels are extracted based on the previously constructed vessel neighborhood interaction topology network. First, using the abnormal vessel node as the starting node, neighboring vessel nodes with interactive connections to it are searched within the network. Then, the neighborhood is expanded layer by layer according to the connections, and the connection order and spatial position changes between nodes are recorded. For each connection, the interaction weights and spatial proximity relationships between the corresponding vessels are extracted, and the connected nodes are arranged sequentially according to the time window, forming a continuous vessel interaction trajectory path. This path reflects the interaction propagation relationship between the abnormal vessel and surrounding vessels and can demonstrate the propagation structure of abnormal behavior among neighboring vessels.

[0098] After obtaining the interaction paths of abnormal vessels and their associated vessels, the detection results are organized and output. First, abnormal behavior record data is generated, including the unique identification information of the abnormal vessel, the time window of the anomaly, the corresponding implicit behavioral index value, and the anomaly determination type. Then, the interaction paths of associated vessels are organized chronologically to form a complete interaction path record, including the positional changes and interaction relationships of the corresponding vessels at each node. Finally, the abnormal behavior records and interaction path records are combined to generate abnormal vessel behavior detection result data, which is used as the final output information to display the abnormal behavior of the target vessel and its interaction impact paths with surrounding vessels.

[0099] The above description is merely a specific embodiment of this application, but the scope of protection of this application 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 this application should be included within the scope of protection of this application.

Claims

1. A method for detecting abnormal behavior of ships based on an ADS-B big data platform, characterized in that: include: Acquire ADS-B broadcast data streams of ships in the target sea area, extract multi-dimensional track sequences, and construct the ship track evolution matrix through a sliding time window; Based on the trajectory evolution matrix, the heading perturbation angular velocity, velocity fluctuation gradient, and trajectory curvature density are calculated and fused to generate the ship's perturbation motion feature tensor. Based on the spatial distribution relationship of the perturbation motion feature tensor, a ship neighborhood interaction topology network is constructed, and the trajectory perturbation propagation weights between nodes are calculated to obtain the anomaly propagation probability matrix. The steps for constructing the ship neighborhood interaction topology network based on the spatial distribution relationship of the perturbation motion feature tensor include: Based on the time index corresponding to the perturbation motion feature tensor, the longitude and latitude coordinates of each ship within the same time window are extracted, and the spatial distance between any two ships is calculated. Ship pairs with a spatial distance less than a preset neighborhood distance threshold are selected as candidate neighborhood ship sets. The similarity of the perturbation motion feature tensors of the ships in the candidate neighborhood ship sets is calculated. By calculating the feature distance between the three feature sequences of heading perturbation angular velocity, velocity fluctuation gradient, and trajectory curvature density, the perturbation motion association weights between ships are generated. Using the unique identification code of a ship as a network node, the spatial distance relationship is used as the initial connection relationship, and the connection relationship is weighted according to the perturbation motion association weight to construct a ship neighborhood interaction relationship matrix; a ship neighborhood interaction topology network is generated according to the ship neighborhood interaction relationship matrix, where network nodes represent individual ships and network edge weights represent the intensity of perturbation motion interaction between ships; The steps for calculating the propagation weights of trajectory perturbations between nodes and generating the anomaly propagation probability matrix include: Based on the ship neighborhood interaction topology network, the connection edges between any two adjacent nodes are extracted, and the heading perturbation angular velocity, velocity fluctuation gradient, and trajectory curvature density feature values ​​of the corresponding ship in the perturbation motion feature tensor are read to construct the perturbation feature difference vector between node pairs; the trajectory perturbation similarity index between nodes is calculated based on the perturbation feature difference vector, and the trajectory perturbation similarity index is multiplied by the spatial adjacency weight in the ship neighborhood interaction topology network to obtain the trajectory perturbation propagation weight between nodes; The trajectory disturbance propagation weights between all nodes are arranged in a matrix according to the node connection relationship to construct the trajectory disturbance propagation weight matrix between ship nodes; the trajectory disturbance propagation weight matrix is ​​subjected to normalized probability transformation to obtain the abnormal propagation probability matrix representing the diffusion trend of abnormal ship behavior; Based on the anomaly propagation probability matrix, a temporal folding mapping is performed on the perturbation motion feature tensor to establish a coupled response model between the local track perturbation energy sequence and the neighborhood interaction intensity, generating anomaly behavior evolution prediction curves, including: Based on the anomaly propagation probability matrix, the propagation probability weights of each ship node in adjacent time windows are extracted, and the perturbation motion feature tensors of the corresponding ships are weighted and aggregated according to time order to obtain the track perturbation energy sequence characterizing the intensity of local motion perturbation of the ship. The track disturbance energy sequence is synchronized and aligned with the neighborhood interaction intensity within the corresponding time window, and the disturbance energy values ​​of consecutive time windows are superimposed and mapped using a sliding time folding method to form a time evolution feature sequence. Based on the time evolution characteristic sequence, a coupled response model between the local track disturbance energy sequence and the neighborhood interaction intensity is constructed, and the abnormal behavior evolution trend parameters are obtained by calculating the dynamic response coefficient between the two. Based on the abnormal behavior evolution trend parameters, the disturbance energy changes in subsequent time windows are recursively calculated to generate a ship abnormal behavior evolution prediction curve. The abnormal energy accumulation index of the target ship is calculated based on the abnormal behavior evolution prediction curve, and the abnormal energy accumulation index is nonlinearly coupled with the abnormal propagation probability matrix to generate a risk distribution map of abnormal ship behavior. Multi-scale feature extraction is performed on the risk distribution map to obtain implicit behavioral indicators, including trajectory drift density, interactive perturbation spectrum and abnormal energy convergence, and dynamic matching is performed with preset abnormal behavior patterns. When any implicit behavior indicator exceeds the set threshold, the abnormal behavior detection result of the corresponding ship and its associated ship interaction path are output.

2. The method for detecting abnormal behavior of ships based on ADS-B big data platform according to claim 1, characterized in that: The steps for constructing the ship's trajectory evolution matrix using a sliding time window include: Collect ADS-B broadcast data of ships in the target sea area over a continuous time period, merge and sort the data according to the ship's unique identification code, and extract latitude and longitude coordinates, heading angle, speed and timestamp information to form the original multidimensional track sequence; The original multidimensional track sequence is subjected to time consistency correction and spatial coordinate unification processing. A standardized track sequence with uniform time intervals is generated by linear interpolation and abnormal jump point removal methods. The standardized track sequence is segmented and sampled using a preset sliding time window. The changes in latitude and longitude, heading and speed within the same time window are arranged in a matrix according to time order to obtain a local track feature submatrix of the ship. Multiple local track feature sub-matrices corresponding to continuous sliding time windows are superimposed and integrated according to the time dimension to construct a track evolution matrix that reflects the dynamic change law of short-term ship track.

3. The method for detecting abnormal behavior of a ship based on an ADS-B big data platform according to claim 1, characterized in that: The steps for calculating the heading perturbation angular velocity, velocity fluctuation gradient, and trajectory curvature density based on the trajectory evolution matrix, and generating the ship perturbation motion characteristic tensor, include: Extract the heading angle and speed change data of the continuous time series from the track evolution matrix, calculate the heading perturbation angular velocity sequence based on the heading angle difference and time interval between adjacent time sampling points, and calculate the speed fluctuation gradient sequence based on the speed change and time change rate. The spatial displacement vector between adjacent track points is calculated based on the longitude and latitude changes in the track evolution matrix, and the trajectory curvature density sequence is calculated using a three-point trajectory fitting method. The heading perturbation angular velocity sequence, velocity fluctuation gradient sequence, and trajectory curvature density sequence are normalized, and corresponding multidimensional feature matrices are constructed. The multidimensional feature matrix is ​​reconstructed by tensor quantization according to the time dimension and the feature dimension to generate a perturbation motion feature tensor.

4. The method for detecting abnormal behavior of ships based on ADS-B big data platform according to claim 1, characterized in that: The steps for calculating the anomalous energy accumulation index of a target ship based on the anomalous behavior evolution prediction curve include: Based on the abnormal behavior evolution prediction curve, the coupled response prediction values ​​within the continuous prediction time window are extracted, and the predicted perturbation energy sequence is constructed in chronological order. Local peak identification is performed on the predicted perturbation energy sequence. Abnormal energy peak points are determined by calculating the energy change amplitude between adjacent time windows, and the distribution density of peak points in the predicted perturbation energy sequence is statistically analyzed to obtain candidate intervals for abnormal energy aggregation. The time-weighted cumulative energy is calculated based on the disturbance energy value within the candidate interval of the abnormal energy accumulation, and the energy concentration coefficient is calculated in combination with the peak distribution density. The time-weighted cumulative energy is coupled with the energy concentration coefficient to obtain the abnormal energy accumulation index, which represents the intensity of the abnormal behavior of the target ship.

5. The method for detecting abnormal behavior of ships based on ADS-B big data platform according to claim 4, characterized in that: The steps for performing nonlinear coupling calculations between the abnormal energy accumulation index and the abnormal propagation probability matrix to generate a risk distribution map of abnormal ship behavior include: Read the abnormal energy accumulation index of the target ship and extract the propagation probability sequence corresponding to the target ship from the abnormal propagation probability matrix. Use the propagation probability sequence as the spatial diffusion weight to perform neighborhood diffusion mapping on the abnormal energy accumulation index to obtain the initial abnormal energy diffusion vector. The initial abnormal energy diffusion vector is subjected to nonlinear transformation, and an exponential amplification function is constructed to enhance the calculation of the high propagation probability region, thereby obtaining the risk intensity sequence of abnormal ship behavior. The abnormal behavior risk intensity sequence of the ship is mapped to the grid coordinates of the target sea area according to the spatial location of the ship, and the cumulative value of risk intensity in each grid cell is calculated to generate a spatial risk intensity matrix; Continuous spatial interpolation is performed on the spatial risk intensity matrix to construct a risk distribution map of abnormal ship behavior in the target sea area.

6. The method for detecting abnormal behavior of a ship based on an ADS-B big data platform according to claim 1, characterized in that: The steps for multi-scale feature extraction from the risk distribution map include: constructing a multi-scale spatial analysis window based on the ship abnormal behavior risk distribution map, and statistically analyzing the risk intensity gradient change rate and spatial distribution density within different scale windows to calculate the trajectory drift density feature, which is used to characterize the degree of ship track deviation; extracting the risk intensity change sequence over time within the multi-scale spatial analysis window, and performing frequency domain transformation on the sequence to obtain the interactive disturbance spectrum feature; calculating the abnormal energy convergence index based on the spatial clustering degree of high-risk areas in the risk distribution map, and generating an abnormal energy convergence feature vector by statistically analyzing the spatial contraction rate of high-risk units within a continuous time window.

7. The method for detecting abnormal behavior of a ship based on an ADS-B big data platform according to claim 6, characterized in that: The steps for dynamically matching with preset abnormal behavior patterns include: fusing the trajectory drift density features, interactive disturbance spectrum features, and abnormal energy convergence feature vectors, and calculating the similarity with the behavioral features in the preset abnormal behavior patterns to achieve dynamic matching and identification of abnormal ship behavior patterns.

Images