A ship waypoint prediction method and system based on bridging distribution and a storage medium

Abstract

A ship waypoint prediction method and system based on bridging distribution, and a storage medium, the method comprising: acquiring historical trajectory observation data of a ship and a preset future waypoint set; constructing a ship motion model based on a linear stochastic differential equation, and solving to obtain a closed-form solution of a probability distribution representing a ship motion state; using a segmented parameter estimation strategy, dividing the historical trajectory observation data into multiple continuous time segments, estimating the model parameters in each time segment respectively, and obtaining dynamically adapted segmented estimation parameters; based on a Bayesian framework, converting each waypoint in the future waypoint set into a pseudo-observation value, combining the segmented estimation parameters and the historical trajectory observation data, processing the time of arrival and the mean velocity vector through marginalization, and calculating the posterior probability distribution of each preset waypoint; and according to the posterior probability distribution, screening the waypoint with the highest probability as the optimal predicted waypoint.

Description

Technical Field

[0001] This application relates to the field of intelligent shipping technology, and in particular to a method, system and storage medium for predicting ship waypoints based on bridging distribution. Background Technology

[0002] In fields such as maritime traffic management, intelligent shipping, route planning, and collision avoidance, accurate prediction of a ship's future waypoint is a core technological support for ensuring navigation safety and improving shipping efficiency. As the global shipping industry transforms towards intelligence and automation, higher demands are being placed on the accuracy and real-time performance of ship trajectory prediction.

[0003] Current ship waypoint prediction technologies mainly revolve around kinematic modeling, machine learning fitting, and statistical inference. Among these, kinematic model-based methods have become one of the mainstream approaches due to their clear physical meaning and strong interpretability. Traditional ship kinematic models often use linear differential equations to describe state changes, but they neglect the stochastic nature of marine environmental noise, leading to the accumulation of prediction errors over time. In statistical inference methods, Kalman filtering and its variants are widely used in ship state prediction because they can effectively handle the state estimation problem of noisy linear systems. However, traditional Kalman filtering assumes fixed model parameters, which cannot adapt to the time-varying characteristics of parameters caused by environmental changes or maneuver adjustments during ship navigation; at the same time, its approximation of nonlinear systems introduces additional errors, limiting the accuracy of waypoint prediction.

[0004] Regarding the inference framework for waypoint prediction, existing technologies mostly employ a single-point estimation model, outputting only a single predicted waypoint without considering various random disturbances during navigation, thus failing to provide uncertainty references for decision-making. Furthermore, existing Bayesian inference methods suffer from high computational complexity when dealing with multivariate coupling and arrival time uncertainties, making it difficult to strike a balance between prediction accuracy and real-time performance. Some methods reduce computational load by simplifying model assumptions, but at the expense of prediction accuracy; while high-precision inference methods often rely on complex numerical calculations, resulting in low computational efficiency and failing to meet the real-time requirements of ship navigation systems.

[0005] Therefore, there is an urgent need for a ship waypoint prediction technology that can take into account both random disturbance characterization and efficient and accurate inference. Summary of the Invention

[0006] To address the shortcomings of existing technologies, this application aims to provide a method, system, and storage medium for predicting ship waypoints based on bridging distribution, which achieves accurate probability prediction and uncertainty quantification of future ship waypoints, striking a balance between prediction accuracy and real-time performance.

[0007] To achieve the above objectives, this application provides a ship waypoint prediction method based on bridging distribution, comprising the following steps: Acquire historical trajectory observation data of ships and a preset set of future waypoints; A ship motion model is constructed based on linear stochastic differential equations, and the closed-form solution of the linear stochastic differential equations is obtained. The closed-form solution is used to characterize the probability distribution of the ship's motion state. A segmented parameter estimation strategy is adopted, which divides the historical trajectory observation data into multiple continuous time segments, and estimates the model parameters in each time segment to obtain dynamically adapted segmented estimation parameters. Based on the Bayesian framework, each waypoint in the preset set of future waypoints is transformed into a pseudo-observation. Combined with the segmented estimation parameters and the historical trajectory observation data, the arrival time and mean velocity vector are marginalized to calculate the posterior probability distribution of each preset waypoint. Based on the posterior probability distribution, the waypoint with the highest probability is selected as the optimal predicted waypoint output.

[0008] Furthermore, the step of constructing a ship motion model based on linear stochastic differential equations and solving for the closed-form solution of the linear stochastic differential equations further includes: Define the ship state vector The state vector includes position components and their corresponding velocity components in an s-dimensional coordinate system, where s is a positive integer; Constructing model equations ,in , , , The mean velocity vector, The response coefficient matrix, The noise intensity matrix is... For the s-Vivener process; Solve the closed-form solution of the model equations using Ito integrals: ,in This is a Gaussian noise term. This is the state transition function. For the input response function, At the initial moment, t For any subsequent moment.

[0009] Furthermore, when the response coefficient matrix When reversible, the state transition function is simplified through eigenvalue decomposition. With the input response function The calculation is expressed as: , , in R is the recovery coefficient matrix. eigenvector matrix, For the eigenvalue matrix, The recovery coefficient matrix eigenvectors.

[0010] Furthermore, when the covariance matrix in the closed-form solution of the linear stochastic differential equation... When a closed-form solution cannot be obtained, matrix fractional decomposition is used for solving the problem. Then calculate the covariance matrix: .

[0011] Furthermore, the step of employing a segmented parameter estimation strategy to divide the historical trajectory observation data into multiple continuous time segments, and estimating the model parameters within each time segment to obtain dynamically adapted segmented estimation parameters further includes: The historical trajectory is divided into K consecutive time segments according to the preset time window length. Model parameters within each time segment Approximate time invariance, where The recovery coefficient matrix eigenvectors, The mean velocity vector, This is the noise intensity matrix; For each time segment Construct the log-likelihood function The parameter estimates for this time segment are obtained by maximizing the log-likelihood function. ,in The number of observations in the k-th time segment. Time interval The covariance matrix under the following conditions , express transpose, Time interval The corresponding state transition matrix, Time interval The corresponding input response matrix, where s is the coordinate system dimension.

[0012] Furthermore, the step of converting each waypoint in the preset future waypoint set into a pseudo-observation value based on the Bayesian framework, and calculating the posterior probability distribution of each preset waypoint by combining the segmented estimation parameters and the historical trajectory observation data and marginalizing the arrival time and mean velocity vector, further includes: Constructing pseudo-observation equations Set the preset waypoints Converted into pseudo-observations, where This is a pseudo-observation matrix. ; By integrating the arrival time T and the mean velocity vector v, the conditional probability of historical trajectory observation data for the preset waypoint is calculated. Calculate the posterior probability of each waypoint using Bayes' theorem: ,in, waypoints The posterior probability, Historical trajectory observation data For the preset waypoint The conditional probability, waypoints Prior probability.

[0013] Furthermore, when integrating the arrival time, Simpson's rule is used for numerical integration, and a fixed number of integration points are selected to cover the arrival time interval.

[0014] Furthermore, the step of selecting the waypoint with the highest probability as the optimal predicted waypoint based on the posterior probability distribution also includes the step of performing inference optimization using a fixed γ inference scheme, Monte Carlo approximate inference, or variational Bayesian inference, where γ is the eigenvalue vector of the recovery coefficient matrix.

[0015] To achieve the above objectives, this application also provides a ship waypoint prediction system based on bridging distribution, comprising: The ship motion modeling module is used to construct a dynamic model describing ship motion based on stochastic differential equations and solve it to obtain a closed-form solution of the ship's motion state. The parameter estimation module is used to divide the ship's historical trajectory into multiple continuous time segments using a segmented parameter estimation strategy, and to estimate the model parameters in each segment to obtain dynamically adapted segmented estimated parameters. The waypoint inference module is used to convert preset future waypoints into pseudo-observations based on a Bayesian framework, and to calculate the posterior probability distribution of each preset waypoint by combining the segmented estimation parameters and historical trajectory observation data. The inference optimization module is used to provide optimization schemes for fixed γ inference or Monte Carlo approximate inference to adapt to different accuracy and real-time requirements, and output the optimal predicted waypoints based on the posterior probability distribution.

[0016] To achieve the above objectives, this application also provides a computer-readable storage medium storing a computer program that is loaded and executed by a processor to implement the bridging distribution-based ship waypoint prediction method described above.

[0017] The ship waypoint prediction method based on bridging distribution provided in this application achieves accurate probability prediction and uncertainty quantification of future ship waypoints by integrating stochastic process modeling, piecewise dynamic parameter estimation and Bayesian probability inference.

[0018] Other features and advantages of this application will be set forth in the following description, and will be apparent in part from the description, or may be learned by practicing this application. Attached Figure Description

[0019] The accompanying drawings are provided to further illustrate the present application and form part of the specification. Together with the embodiments of the present application, they serve to explain the present application but do not constitute a limitation thereof. In the drawings: Figure 1 This is a flowchart of a ship waypoint prediction method based on bridging distribution according to an embodiment of this application; Figure 2 This is a schematic diagram of the structure of a ship waypoint prediction system based on bridging distribution according to an embodiment of this application. Detailed Implementation

[0020] The preferred embodiments of this application are described below with reference to the accompanying drawings. It should be understood that the preferred embodiments described herein are for illustration and explanation only and are not intended to limit this application.

[0021] Embodiments of this application will now be described in more detail with reference to the accompanying drawings. While some embodiments of this application are shown in the drawings, it should be understood that this application can be implemented in various forms and should not be construed as limited to the embodiments set forth herein. Rather, these embodiments are provided to provide a more thorough and complete understanding of this application. It should be understood that the drawings and embodiments of this application are for illustrative purposes only and are not intended to limit the scope of protection of this application.

[0022] The term "comprising" and its variations as used in this application are open-ended, meaning "including but not limited to". The term "based on" means "at least partially based on". The term "one embodiment" means "at least one embodiment"; the term "another embodiment" means "at least one additional embodiment"; and the term "some embodiments" means "at least some embodiments".

[0023] It should be noted that the terms "first" and "second" may be used in this application only to distinguish different devices, components or parts, and are not used to define the order of functions performed by these devices, components or parts or their interdependence.

[0024] It should be noted that the terms "one" and "more" used in this application are illustrative rather than restrictive, and those skilled in the art should understand that, unless explicitly stated otherwise in the context, they should be understood as "one or more". "More" should be understood as two or more.

[0025] In this application, "bridging distribution" refers to the conditional joint Gaussian probability distribution that connects the posterior distribution of the ship's historical observation data with the distribution of future waypoint constraints. It is the core mathematical carrier for achieving the probabilistic fusion of historical motion information, real-time state estimation, and future waypoint constraints. Its mathematical definition is: given the model's response coefficient matrix... Eigenvalue vector γ, waypoint arrival time T, future waypoint pseudo-observations Compared with historical trajectory observation data At time T, the ship's state vector Compared to the state vector at the nth time step The joint posterior distribution probability is calculated. This distribution uses the state transition matrix and input response matrix of the Ornstein-Uhlenbeck process as a probabilistic bridge, anchoring one end to the current state posterior obtained from historical observations through Kalman filtering, and the other end carrying the pseudo-observation constraints of future waypoint transformations. Based on this bridging distribution, the conditional prior distribution of the mean velocity vector v can be further derived. This mathematically establishes a closed-loop link between parameter estimation and waypoint inference.

[0026] To make the objectives, technical solutions, and advantages of this application clearer, the embodiments of this application will be described in further detail below with reference to the accompanying drawings.

[0027] Example 1 Figure 1 The flowchart below shows the ship waypoint prediction method based on bridging distribution according to an embodiment of this application. Figure 1 The embodiments of this application will be described in further detail.

[0028] First, in step 101, the historical trajectory observation data of the ship to be predicted, the initial model parameters, the preset set of future waypoints, and the time range of waypoint arrival are obtained.

[0029] In the embodiments of this application, historical trajectory observation data of the ship to be predicted are obtained. (Including time series) ), pre-trained initial motion model (including initial parameters) ), and a pre-set set of future waypoints .

[0030] In some exemplary implementations, historical trajectory observation data of the ship (such as latitude, longitude, speed, timestamp, etc.) is acquired from a GNSS receiver and preprocessed into a UTM plane rectangular coordinate system (dimension s=2, eastward). North ), to obtain the observation sequence (In this embodiment, the number of observations) (1Hz sampling for 24 hours). Receive initial model parameters, including: initial feature vector. Initial mean velocity vector Initial noise intensity matrix Initial recovery coefficient matrix Outliers (velocities > 10 m / s, position jumps > 100 m) were removed using the 3σ criterion (σ being the noise intensity matrix), and missing data were supplemented using linear interpolation to ensure the continuity of the time series. The core model parameters were determined as follows: coordinate system dimension s = 2, observation matrix... (4×4 identity matrix), observation noise covariance matrix Waypoint arrival time range False observation noise covariance matrix A small amount of historical trajectory data was selected as calibration data for distribution statistics in the parameter estimation and inference process.

[0031] Step 102: Construct a ship motion model based on linear stochastic differential equations.

[0032] In the embodiments of this application, a linear stochastic differential equation (SDE) describing the random motion of a ship is constructed based on the Ornstein-Uhlenbeck (OU) process. This SDE accurately characterizes the dynamic changes and random disturbances of the ship's position and velocity. The process involves constructing the SDE and solving it to obtain a closed-form solution, which characterizes the probability distribution of the ship's motion state. This step specifically includes: Step 1021, define the state vector in the s-dimensional coordinate system as: in, For positional components, For the corresponding velocity component, s∈Z + In this embodiment, s=2 (corresponding to a Cartesian coordinate system or a latitude and longitude coordinate system), and the state vector... , in , UTM coordinates (unit: m) , The corresponding velocity (unit: m / s). The position and velocity of the first valid observation point of the historical trajectory are taken to obtain the initialization. state of time .

[0033] In step 1022, the linear stochastic differential equation (SDE) is constructed as follows: in, ; It is an s×s identity matrix. The mean velocity vector, , The response coefficient matrix, The noise intensity matrix is... This is a standard S-Vivener process, used to characterize random disturbances such as those in the marine environment. In this embodiment: .

[0034] Wiener process Generate standard two-dimensional Gaussian noise using Python's numpy.random.normal method, ensuring... , Here, E represents the mathematical expectation, which means that random disturbances in the marine environment, under statistical average, will not cause continuous and systematic speed or position shifts to the ship. The disturbances will only fluctuate around the zero mean, which is consistent with the characteristics of random noise in actual navigation.

[0035] In step 1023, the closed-form solution of SDE is obtained by Itō integral: in, This is a Gaussian noise term. This is the state transition function. For the input response function, ( It follows a pattern with a mean of 0 and a covariance matrix of . (normal distribution) , , At the initial moment, t This refers to the current or any subsequent moment.

[0036] When the recovery coefficient matrix When it is reversible, eigenvalue decomposition is used to simplify. and The calculations are expressed as follows: in, R is the eigenvector matrix of ρ, Γ=diag(γ), and γ is the eigenvalue vector of ρ.

[0037] In this embodiment, for Perform eigenvalue decomposition: , , Substitute The calculation yields: .

[0038] Step 1024, calculate the covariance matrix. : ; If a closed-form solution cannot be obtained, matrix fractional factorization can be used for efficient solution: .

[0039] Specifically, construct the augmented matrix: ; Substituting A and B, we get Thus, the specific value of M can be obtained; Calculate the matrix exponent exp(M×( - ): Calculated using the scipy.linalg.expm function in SciPy. - =1s, we get exp(M×1); Solve : ,in (Covariance matrix at initial time t0); calculate Calculate using NumPy's matrix inversion function numpy.linalg.inv Then perform matrix multiplication to obtain .here It is the covariance matrix at a time interval Δt = 1 second, that is, the covariance matrix of random noise accumulated after 1 second from the initial time.

[0040] The final output model state vector Discretized sequence { , ,..., }, Covariance matrix sequence { , ,..., }, used for subsequent parameter estimation.

[0041] In step 103, piecewise parameter estimation is performed.

[0042] In the embodiments of this application, a piecewise maximum likelihood estimation (MLE) strategy is adopted, dividing the historical trajectory into multiple continuous time segments, and performing maximum likelihood estimation on the model parameters within each segment to obtain dynamically adapted piecewise estimation parameters. This step specifically includes: Step 1031: Segment the historical trajectory. Divide the historical trajectory into k consecutive time segments according to a preset time window length. Model parameters within each time segment When approximated, the trajectory remains unchanged and is approximately piecewise linear.

[0043] In the embodiments of this application, the determination of the time window length needs to strike an optimal balance between the time-varying adaptability of parameters and the accuracy of MLE estimation. In stable navigation scenarios, ship speed and heading fluctuations are small, environmental disturbances such as wind and current are stable, and the time-varying nature of motion parameters is weak. A longer window can be used to improve the robustness of MLE estimation and reduce noise interference with parameters. In dynamic navigation scenarios, ship control commands are frequently adjusted, environmental disturbances are severe, and the time-varying nature of motion parameters is strong. A shorter window must be used to ensure that the parameters within the window are approximately time-invariant and to avoid error accumulation. The selection of the time window length is also constrained by the sample size of the statistical estimation. For example, if AIS data with a sampling frequency of 1Hz is used, the lower limit of the window length is 30s. If low sampling frequency data (such as AIS data with 1 point per minute) is used, the window length needs to be extended accordingly to ensure that the sample size within the window meets the degree of freedom requirements for MLE estimation. The premise for simplifying the model closure solution through eigenvalue decomposition is that the reversion coefficient matrix ρ is invertible, that is, all eigenvalues ​​of ρ... >0. Window length directly affects the MLE estimation result of ρ: a window that is too short will lead to... The estimation results in a non-positive value, making ρ non-invertible and preventing the eigenvalue simplification. Therefore, the choice of window length must be made while simultaneously verifying that the estimated ρ satisfies the invertibility condition.

[0044] In this embodiment, the time window length is set to... (1 hour), dividing the 24-hour trajectory into Each time segment is divided into time segments. , (Include (Number of observation points). Segment boundary processing: Adjacent segments overlap by 100 data points (approximately 1.7 minutes) to avoid abrupt changes in boundary parameters. The overlapping part takes the last data point of the previous segment and the first data point of the next segment.

[0045] Step 1032, observation model adaptation. If the historical trajectory observation data directly corresponds to the ship's state (such as high-precision GPS data), then the state data is directly used for parameter estimation; otherwise, a linear Gaussian observation model is used. Process observation data.

[0046] In this embodiment, since GNSS data directly outputs position and velocity, the observation model is simplified to... : Observation matrix H= (4×4 identity matrix), meaning that the observed values ​​directly correspond to the state vector; Observation noise (In this embodiment) (is normally distributed) (with the initial time) covariance matrix Consistent (determined by sensor accuracy); Generate observation model output Divide each time segment ,calculate (j=1,...,3600) Using numpy.random.normal(0, )generate.

[0047] Step 1033, Likelihood Function Construction and Optimization. The likelihood function for the k-th time segment is the product of the conditional Gaussian distributions at each time step: in, The distribution is Gaussian, and the mean and covariance are determined by the closed-form solution of the model. Parameter estimates are obtained by maximizing the likelihood function. This enables dynamic adaptation of model parameters. The corresponding log-likelihood function is: in, , , for transpose, for The reverse, This represents the number of observations within the k-th time segment (i.e., the number of sampling points for ship historical trajectory observation data). Time interval The covariance matrix under the following conditions Time interval The corresponding state transition matrix, Time interval The corresponding input response matrix.

[0048] In this embodiment, the condition distribution is as follows: Where Δt = 1s, = , = , = (Calculated during the modeling phase). The log-likelihood function (to avoid numerical underflow) is: in , , These represent the state transition matrix and input response matrix with a time interval of 1 second, respectively. |for The determinant is calculated using numpy.linalg.det.

[0049] In the embodiments of this application, the Newton-Raphson method is used to maximize the log-likelihood function and estimate the parameters. In each iteration, the gradient The numerical differentiation calculation is performed as follows: Initialization parameters: ; Calculate gradient Numerical differentiation and gradient step size are used. For each parameter component (The i-th parameter), calculate , It is a unit vector; Calculate the Hessian matrix Similarly, numerical differentiation is used. ; The parameter update formula is: ; Convergence condition: The convergence condition of the Newton-Raphson method is set as the 2-norm of the parameter update being less than 1 / 2. That is, when Stop iteration when the time is reached, and output the piecewise estimated parameters of the k-th time segment. ,in, The eigenvector values ​​of the k-th piecewise recovery coefficient matrix ρ determine how quickly the ship's speed recovers to its average speed. Let be the mean velocity vector of the ships in the k-th segment, representing the steady-state average speed of the ships in that segment; The k-th segment noise intensity matrix quantifies the amplitude of random disturbances to navigation caused by external environmental factors such as sea breezes and ocean currents.

[0050] The iteration result of the first time segment in this embodiment is: .

[0051] Step 1034, parameter smoothing. A moving average smoothing process is performed on the estimated parameters for the 24 segments: Smoothing window size = 3, which is the piecewise estimated parameter after smoothing. ; First and last segment processing: In the first time segment In the 24th time segment ; Output smoothed piecewise estimated parameters and observation model output sequence.

[0052] In step 104, Bayesian inference is initialized and the preset waypoint probability is calculated.

[0053] In the embodiments of this application, a waypoint probability inference process is designed based on a Bayesian framework. Through pseudo-observation construction, variable marginalization, and numerical integration, observation data and prior information about future waypoints are fused to calculate the probability distribution of each preset waypoint. This step specifically includes: Step 1041, Preset Waypoints and Prior Settings. Define the set of preset future waypoints. Each waypoint For the target position at future time T, set the prior probability of waypoint D. Prior distribution probability with arrival time T (Uniform distribution). The uniform distribution parameters for arrival time T are: and (Estimated based on ship speed and waypoint distance).

[0054] In this embodiment, five candidate waypoints are set, namely, a preset set of future waypoints D={ , , , }, specific coordinates: (386500, 3569500) (388200, 3570800) (389800, 3572100) (391000, 3573500) (391800, 3574800) (unit: m). Prior probability of waypoint D. (Uniform prior), set using numpy.ones(5)*0.2. The prior probability distribution of arrival time T is: (1 to 3 hours), covering the time range required for the entire route.

[0055] Step 1042: Convert future waypoints into pseudo-observations and construct a pseudo-observation matrix. The pseudo-observation equation is as follows: in, It is an s×2s-dimensional invertible pseudo-observation matrix (only the location-related dimension is extracted). It is the state vector at time K. ( It follows a covariance matrix with a mean of 0 and a pseudo-observation noise of... (normal distribution) The pseudo-observation noise matrix at time K (set by waypoint accuracy requirements) is used to incorporate waypoint constraints into the inference process. , .

[0056] In this embodiment, The pseudo-observation output is: for each candidate waypoint ,generate , (Coordinates of the l-th waypoint); where For the l-th candidate waypoint The observation vector corresponds to the arrival time K. The coordinates are for the eastward position. These are the coordinates of the northward position. Let l be the eastward UTM coordinates of the l-th candidate waypoint. Here are the northward UTM coordinates of the l-th candidate waypoint.

[0057] Step 1043, Variable Marginalization: Historical trajectory observation data are obtained by integrally marginalizing the arrival time T and the mean velocity vector v. For candidate waypoints conditional probability : Where, p It is obtained by sequentially calculating prediction error decomposition (PEDs). Based on the model's state constraints, the distribution is derived as Gaussian. This indicates that the target waypoint is predetermined. Under the conditions of arrival time T, mean velocity v, and so on, the existing complete historical observation sequence was obtained. Conditional likelihood probability; Indicates that at a fixed target waypoint At arrival time T, the mean velocity v follows a conditional probability distribution.

[0058] Specifically, the joint probability decomposition yields: ,in Let be the measured observation vector of the ship at time j; The conditional distribution of v is: Wherein, the mean vector of the conditional distribution of mean velocity v The covariance matrix of the conditional distribution of mean velocity v Derivation using pseudo-observation constraints: in, For the nth time moment The time difference to the time T of arrival at the future waypoint, i.e. ; , It is an s×4 dimensional mapping matrix; Let be the mean vector of the bridging joint posterior distribution. Let be the covariance matrix of the bridging joint posterior distribution; This is the position change-mean velocity mapping coefficient matrix. The initial velocity-mean velocity correction coefficient matrix is ​​a diagonal matrix of s×s dimensions, which is calculated from the eigenvector γ of the recovery coefficient matrix and the time difference Δt.

[0059] Substitution (2×4 matrix), Δt=T- The specific forms of η and ζ are constructed as diagonal matrices; Sequential computation of prediction error decomposition (PEDs): ,in, To observe data from known historical trajectories Candidate waypoints Under the conditions of arrival time T and mean velocity v, the observation at time n The mean of the Gaussian distribution; To observe at time n under the same known conditions The prediction error covariance matrix; , , and Calculated using Woodbury's formula: in, The mean of the predicted state from the Kalman filter (KF) is... Let be the predicted state covariance matrix of the Kalman filter (KF). Preset waypoints coordinate vector, The pseudo-observation noise covariance matrix; To add the l-th candidate waypoint After the pseudo-observation constraint, the posterior predicted mean of the ship's state vector at time n; For binding After constraints, the ship state covariance matrix at time n represents the uncertainty of the state estimation; Its function is to convert the current ship's state into a positional dimension based on arrival time T, which is then used to introduce candidate waypoints. constraint; To introduce The constrained corrected gain matrix.

[0060] In this embodiment, the prior state of the Kalman filter is initialized. and Based on the observations at the initial time After completing the first state update, the following results were obtained: and .

[0061] In step 1044, numerical integration: Simpson's rule is used to numerically integrate the arrival time T, selecting a fixed number of quadrature (integration) points. Coverage arrival time interval [ (In this embodiment) , ), to balance computational accuracy and efficiency.

[0062] In this embodiment, the integral over arrival time T covers the interval [3600, 10800]: 10 quadrature points are selected, i.e., q = 10. =3600+800×(i-1) (i=1,...,10), step size h=800s; Applying the integral formula: in, ( (The prior probability distribution at arrival time T). Calculate each waypoint conditional probability And output it.

[0063] In the embodiments of this application, for each preset waypoint Iterate through the quadrature points of arrival time. (In this embodiment, q=10, step size h=800s), perform the following operations: Parameter selection and sampling: If a fixed γ inference scheme is used, select the latest segment. As a fixed value (in this embodiment) If the Monte Carlo (MC) approximation scheme is used, the kernel density estimation (KDE) approximation can be applied. Sample M Sample (M=1000 in this example).

[0064] Derivation of the conditional priors of v: Based on the joint distribution of model state constraints and pseudo-observation updates, the following derivation is performed. Gaussian distribution parameters and The mean and covariance of the joint distribution are: Sequential computation of prediction error decomposition (PEDs): calculating conditional distribution probabilities using the Woodbury formula. mean Covariance ; This leads to the decomposition of the prediction error: in Let be the observation noise covariance matrix at time n; Marginalization v and integration T: By integrating and marginalizing the mean velocity vector v, we obtain the velocity at a given candidate waypoint. Complete ship history observation sequence under the condition of arrival time T Marginal likelihood of occurrence : in, The model is derived to be a Gaussian distribution based on state constraints.

[0065] Numerical integration of arrival time T is performed using Simpson's rule: in, .

[0066] Calculate the posterior probability of each waypoint using Bayes' theorem: Posterior probability After normalization, we get: Finally, the posterior probability distribution of each preset waypoint is output, and the waypoint with the highest probability is selected as the optimal predicted waypoint. The probability distribution is also output to quantify the prediction uncertainty. If the prediction accuracy does not meet the requirements, parameters such as the time window length, the number of quadrature points, and the number of MC samples are adjusted, and the parameter estimation and waypoint inference process is re-executed until the application requirements are met.

[0067] Furthermore, the fixed γ inference scheme includes two strategies: fixed prior and variable prior. The fixed prior is set at time t1. Through sequential updates Optimize likelihood calculation: Combine Kalman filtering (KF) for prediction and updating. Prediction step: ; Update steps: in Kalman gain; Posterior of sequential update v: ;in, in, To introduce After constraints, observation The constant term in the expression is independent of the mean velocity v; and They are observed Given a preset waypoint Given the arrival time T, the posterior estimates of the mean vector and covariance matrix of the mean velocity vector v; The distribution used to describe the uncertainty of v under the current estimate is the complete prior distribution for the next round of parameter updates. When new observations are received... At that time, the above parameters will be updated to and .

[0068] The variable prior is at each time point Dynamically adjust prior distribution The revised posterior formula is as follows: .

[0069] Subsequent update steps are consistent with the fixed prior strategy; only the following needs to be added: and Replace with the adjusted posterior mean and covariance.

[0070] Monte Carlo (MC) approximation inference includes: when γ is the parameter to be estimated, For non-Gaussian distributions, kernel density estimation (KDE) is used to approximate the results based on parameter estimation. Through sampling ,calculate: .

[0071] In this embodiment, when γ is the parameter to be estimated, p(v|) (T) Non-Gaussian, using MC sampling; Kernel density estimation (KDE) approximates p(γ| ,T): estimated in pieces For the sample, a Gaussian kernel is selected as the kernel function, with a bandwidth h=0.001, implemented using scipy.stats.gaussian_kde; Sampling γ samples: from p(γ∣) M=1000 samples were sampled in T) (m=1,...,1000); Calculate per sample: For each Calculate p(y1:n|) using the method of inferring with a fixed γ. ,T, ); Average ; The Simpson integral arrives at time T, yielding p(y1:n|). ).

[0072] Regarding optimal waypoint selection and uncertainty quantification, the posterior probability is calculated using Bayes' theorem, including: Posterior probability: After normalization, we get ; Optimal waypoint: Filtering p(D= The waypoint with the largest value of |y1:n) is used as the optimal predicted waypoint. In this embodiment, The posterior probability is the highest (approximately 0.42). Uncertainty quantification: Calculating the standard deviation of the posterior probability for each waypoint Output waypoint probability distribution Sum of standard deviations ,..., }

[0073] Regarding parameter adjustment signal feedback, if the prediction accuracy does not meet the requirements (position error > 50m), a parameter adjustment signal is generated, including: Adjustment rules: If If the time window is >20m, shorten the time window from 3600s to 3000s (50 minutes) and re-perform parameter estimation; if the MC inference is still not satisfied, increase the number of samples M to 2000; Feedback path: Use the adjusted time window length, MC sampling number, and other parameters for ship motion modeling and parameter estimation, and restart the process.

[0074] The final output is "optimal predicted waypoint, waypoint probability distribution and standard deviation", and the parameter adjustment signal is used for ship motion modeling to complete the closed-loop data flow.

[0075] The method described in this application has the following beneficial effects: (1) Excellent modeling accuracy: By using the Ornstein-Uhlenbeck (OU) process to construct linear stochastic differential equations, it not only characterizes the random disturbances of the marine environment, but also fits the ship speed change law through mean recovery characteristics, solving the problems of traditional linear models ignoring noise and nonlinear models having poor real-time performance; the piecewise maximum likelihood estimation (MLE) dynamically adapts to the time-varying characteristics of parameters, avoids the accumulation of errors caused by global parameters, and improves the modeling adaptability for different navigation stages.

[0076] (2) Balancing reliability and practicality: By integrating pseudo-observation construction and numerical integration through a Bayesian framework, future waypoint constraints are incorporated into the inference process, while quantifying prediction uncertainty, thus overcoming the limitation that traditional single-point estimation cannot provide risk reference; at the same time, it provides two schemes: fixed γ inference (efficient) and Monte Carlo (MC) approximate inference (high accuracy), combined with Kalman filtering (KF) sequential updates, in 1-3 hour prediction scenarios, the position error is ≤50m and the arrival time error is ≤10 minutes, balancing real-time performance and accuracy requirements.

[0077] (3) Strong generalization and adaptability: This method does not depend on specific ship types or water environments. Through parameter calibration, it can be adapted to different ships such as bulk carriers and container ships, as well as various navigation scenarios such as nearshore and inland waterways. The algorithm process is modularly designed, and parameters such as time window and number of integration points can be flexibly adjusted to adapt to different accuracy and computing power requirements.

[0078] (4) Excellent engineering feasibility: The reasoning process does not involve complex calculations, and the orthogonal transformation and weight fusion do not increase the additional computing power overhead. It can be deployed on general-purpose hardware. The output waypoint probability distribution and uncertainty index can directly support actual decision-making such as maritime scheduling and route planning, and have outstanding practical value.

[0079] Example 2 In the embodiments of this application, a ship waypoint prediction system based on bridging distribution is also provided to implement the ship waypoint prediction method based on bridging distribution described in Embodiment 1.

[0080] like Figure 2 As shown in the embodiment of this application, the ship waypoint prediction system 200 based on bridging distribution includes the following modules: The ship motion modeling module 201 is used to construct linear stochastic differential equations based on the Ornstein-Uhlenbeck (OU) process and solve them to obtain closed-form solutions of the ship's motion state. Specifically, this module receives historical trajectory observation data of the ship and initial model parameters (γ, v, σ), and defines the state vector. And construct model equations The closed-form solution is obtained by using Itō integrals. The final output model includes its state vector (position / velocity) and covariance matrix. .

[0081] The parameter estimation module 202 employs a piecewise maximum likelihood estimation (MLE) strategy to divide the ship's historical trajectory into multiple continuous time segments, estimating the model parameters for each segment to obtain dynamically adapted piecewise estimated parameters. Specifically, this module receives the state vector and covariance matrix of the model output by the ship motion modeling module. The historical trajectory is divided into k segments according to a time window. For each segment, a log-likelihood function is constructed and maximized to estimate the parameters. Output piecewise estimated parameters and observation model output .

[0082] The waypoint inference module 203 is used to transform preset future waypoints into pseudo-observations based on a Bayesian framework, and calculate the posterior probability distribution of each preset waypoint by combining segmented estimation parameters and historical trajectory observation data. Specifically, this module receives the segmented estimation parameters and observation model output from the parameter estimation module 202, as well as the preset waypoint set D, the waypoint prior probability p(D), and the arrival time T distribution, and constructs pseudo-observation equations. Variable marginalization and Simpson numerical integration are used to output intermediate results related to the waypoint probability distribution, namely the conditional probability of each preset waypoint. And the corresponding piecewise estimation parameters.

[0083] The inference optimization module 204 provides optimization schemes for fixed γ inference or Monte Carlo (MC) approximate inference to adapt to different accuracy and real-time requirements, and outputs the optimal predicted waypoint based on the posterior probability distribution. Specifically, this module receives intermediate results output by the waypoint inference module 203 and selects the latest segment... As a fixed value (fixed γ inference) or using KDE sampling MC approximation inference, combined with Kalman filtering for state prediction and update, the final output is the optimal predicted waypoint and waypoint probability distribution (quantification uncertainty), and at the same time, parameter adjustment signals can be generated and fed back to the ship motion modeling module 201.

[0084] The above four modules are connected in sequence, and the data flow is as follows: ship motion modeling module 201 → parameter estimation module 202 → waypoint inference module 203 → inference optimization module 204. The inference optimization module 204 can send parameter adjustment signals to the ship motion modeling module 201 to form a closed-loop feedback.

[0085] Understandably, without departing from the core idea of ​​this application, "based on bridging distribution fusion dynamic parameter estimation and efficient Bayesian inference", the technical solution can be locally equivalently transformed. These alternative methods can all achieve the core functions of accurate waypoint probability prediction, uncertainty quantification and time-varying parameter adaptation, while maintaining the core advantages of the original solution.

[0086] Specifically, regarding the ship motion modeling module, in addition to the Ornstein-Uhlenbeck (OU) process, a fractional Ornstein-Uhlenbeck (fOU) process can be used to enhance the long-memory adaptation of the ship's navigation state through fractional differential operators, making it more suitable for long-term low-speed navigation scenarios near the shore; alternatively, the Cox-Ingersoll-Ross (CIR) process can be used to avoid the unreasonable negative velocity or noise intensity problems that may occur in the OU process, and to adapt to the modeling requirements of non-negative state vectors.

[0087] Furthermore, in terms of calculating the covariance matrix, the original matrix fractional decomposition can be replaced by the Kalman filter recursive formula, which eliminates the need for complex matrix exponentiation operations, improves real-time performance by more than 30%, and is more suitable for the deployment needs of edge devices with limited computing power.

[0088] Regarding the parameter estimation module, the original piecewise maximum likelihood estimation (MLE) can be replaced with Bayesian estimation. By setting the Gamma prior distribution of γ and incorporating prior information such as ship type, the estimation robustness is improved by 25% in scenarios with small amounts of data. Alternatively, the recursive least squares (RLS) algorithm can be used to update parameter weights in real time without the need for a fixed time window, making it more suitable for scenarios with rapid time-varying parameters caused by sudden water flow and wind direction.

[0089] Furthermore, in the trajectory segmentation strategy, the fixed time window can be replaced with an adaptive window, based on the rate of change of velocity (threshold 0.5 m / s). 2 The window length is dynamically adjusted using the model residual as an indicator. When the navigation status is stable, the window is extended to 2 hours, and when the status changes abruptly, it is shortened to 10 minutes, effectively balancing the time-varying nature of parameters and the estimation accuracy.

[0090] Regarding the waypoint inference module, the original pseudo-observations constrained only by position can be expanded to multi-constrained pseudo-observations, and the pseudo-observation matrix... It also constrains position and arrival speed, adapting to scenarios with multiple requirements for arrival status, such as anchorage.

[0091] In numerical integration, Simpson's rule can be replaced by Gauss-Legendre integration, which requires only 5-8 integration points to achieve the same accuracy and significantly reduces the amount of computation.

[0092] Regarding the inference optimization module, Monte Carlo (MC) approximate inference can be replaced by variational Bayesian (VB) inference. By optimizing the evidence lower bound (ELBO), an approximate analytical solution for the posterior distribution is obtained. The calculation speed is 3-5 times faster than MC inference, which is suitable for the real-time prediction requirements of the edge.

[0093] The standard Kalman filter (KF) can be replaced by the unscented Kalman filter (UKF), which processes nonlinear state transitions through unscented transformation. This improves the state estimation accuracy by 15-20% under strong disturbances and is more suitable for nonlinear observation scenarios such as radar and sonar.

[0094] The above alternative methods are all mainstream and efficient implementation paths in engineering implementation, and do not deviate from the core logic of this application. They can be flexibly selected according to specific navigation scenarios (such as long-term / short-term prediction, strong / weak disturbance environment, edge / cloud deployment).

[0095] Example 3 In the embodiments of this application, a computer-readable storage medium is also provided, which stores a computer program, wherein the computer program is configured to execute the steps in the embodiments of the ship waypoint prediction method based on bridging distribution as described above when it runs.

[0096] In this embodiment, the aforementioned computer-readable storage medium may include, but is not limited to, various media capable of storing computer programs, such as USB flash drives, read-only memory (ROM), random access memory (RAM), portable hard drives, magnetic disks, or optical disks.

[0097] It will be understood by those skilled in the art that the above descriptions are merely preferred embodiments of this application and are not intended to limit this application. Although this application has been described in detail with reference to the foregoing embodiments, those skilled in the art can still modify the technical solutions described in the foregoing embodiments or make equivalent substitutions for some of the technical features. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of this application should be included within the protection scope of this application.

Claims

1. A method for predicting ship waypoints based on bridging distribution, characterized in that, Includes the following steps: Acquire historical trajectory observation data of ships and a preset set of future waypoints; A ship motion model is constructed based on linear stochastic differential equations, and the closed-form solution of the linear stochastic differential equations is obtained. The closed-form solution is used to characterize the probability distribution of the ship's motion state. A segmented parameter estimation strategy is adopted, which divides the historical trajectory observation data into multiple continuous time segments, and estimates the model parameters in each time segment to obtain dynamically adapted segmented estimation parameters. Based on the Bayesian framework, each waypoint in the preset set of future waypoints is transformed into a pseudo-observation. Combined with the segmented estimation parameters and the historical trajectory observation data, the arrival time and mean velocity vector are marginalized to calculate the posterior probability distribution of each preset waypoint. Based on the posterior probability distribution, the waypoint with the highest probability is selected as the optimal predicted waypoint output.

2. The ship waypoint prediction method based on bridging distribution according to claim 1, characterized in that, The step of constructing a ship motion model based on linear stochastic differential equations and solving the closed-form solution of the linear stochastic differential equations further includes: Define the ship state vector The state vector includes position components and their corresponding velocity components in an s-dimensional coordinate system, where s is a positive integer; Constructing model equations ,in , , , It is an s×s identity matrix. The mean velocity vector, The response coefficient matrix, The noise intensity matrix is... For the s-Vivener process; Solve the closed-form solution of the model equations using Ito integrals: ,in This is a Gaussian noise term. This is the state transition function. For the input response function, At the initial moment, t For any subsequent moment.

3. The ship waypoint prediction method based on bridging distribution according to claim 2, characterized in that, When the recovery coefficient matrix When reversible, the state transition function is simplified through eigenvalue decomposition. With the input response function The calculation is expressed as: ， ， in R is the recovery coefficient matrix. eigenvector matrix, For the eigenvalue matrix, The recovery coefficient matrix eigenvectors.

4. The ship waypoint prediction method based on bridging distribution according to claim 2, characterized in that, When the covariance matrix in the closed-form solution of the linear stochastic differential equation When a closed-form solution cannot be obtained, matrix fractional decomposition or Kalman filtering recursive formulas are used to solve the problem. The expression for solving the matrix fractional decomposition is as follows: ； ； in Let represent the covariance matrix at the initial time t0.

5. The ship waypoint prediction method based on bridging distribution according to claim 1, characterized in that, The step of employing a segmented parameter estimation strategy, dividing the historical trajectory observation data into multiple continuous time segments, and estimating the model parameters within each time segment to obtain dynamically adapted segmented estimation parameters, further includes: The historical trajectory is divided into k consecutive time segments according to the preset time window length. Model parameters within each time segment Approximate time invariance, where The recovery coefficient matrix eigenvectors, The mean velocity vector, This is the noise intensity matrix; For each time segment Construct the log-likelihood function The parameter estimates for this time segment are obtained by maximizing the log-likelihood function. ,in The number of observations in the k-th time segment. Time interval The covariance matrix under the following conditions, Gaussian noise term , express transpose, Time interval The corresponding state transition matrix, Time interval The corresponding input response matrix, where s is the coordinate system dimension.

6. The ship waypoint prediction method based on bridging distribution according to claim 1, characterized in that, The step of converting each waypoint in the preset future waypoint set into a pseudo-observation value based on the Bayesian framework, and calculating the posterior probability distribution of each preset waypoint by combining the segmented estimation parameters and the historical trajectory observation data and marginalizing the arrival time and mean velocity vector, further includes: Constructing pseudo-observation equations Set the preset waypoints Converted into pseudo-observations, where This is a pseudo-observation matrix. Let K be the pseudo-observation noise matrix. Let K be the state vector at time K; By integrating the arrival time T and the mean velocity vector v, the conditional probability of historical trajectory observation data for the preset waypoint is calculated. Calculate the posterior probability of each waypoint using Bayes' theorem: ,in, waypoints The posterior probability, Historical trajectory observation data For the preset waypoint The conditional probability, waypoints The prior probability.

7. The ship waypoint prediction method based on bridging distribution according to claim 6, characterized in that, When integrating the arrival time, Simpson's rule is used for numerical integration, and a fixed number of integration points are selected to cover the arrival time interval.

8. The ship waypoint prediction method based on bridging distribution according to claim 1, characterized in that, The step of selecting the waypoint with the highest probability as the optimal predicted waypoint based on the posterior probability distribution further includes the step of performing inference optimization using a fixed γ inference scheme, Monte Carlo approximate inference, or variational Bayesian inference, where γ is the eigenvalue vector of the recovery coefficient matrix.

9. A ship waypoint prediction system based on bridging distribution, characterized in that, include: The ship motion modeling module is used to construct a dynamic model describing ship motion based on stochastic differential equations and solve it to obtain a closed-form solution of the ship's motion state. The parameter estimation module is used to divide the ship's historical trajectory into multiple continuous time segments using a segmented parameter estimation strategy, and to estimate the model parameters in each segment to obtain dynamically adapted segmented estimated parameters. The waypoint inference module is used to convert preset future waypoints into pseudo-observations based on a Bayesian framework, and to calculate the posterior probability distribution of each preset waypoint by combining the segmented estimation parameters and historical trajectory observation data. The inference optimization module is used to provide optimization schemes for fixed γ inference or Monte Carlo approximate inference to adapt to different accuracy and real-time requirements, and output the optimal predicted waypoints based on the posterior probability distribution.

10. A computer-readable storage medium, characterized in that, The storage medium stores a computer program, which is loaded and executed by a processor to implement the ship waypoint prediction method based on bridging distribution as described in any one of claims 1 to 8.

Images