A maneuvering target tracking method based on hybrid reasoning

Abstract

The present invention provides a maneuvering target tracking method based on hybrid inference. This method, based on variational Bayesian and sequential Monte Carlo hybrid inference, first designs a proposed distribution based on initial variational parameters. This proposed distribution is sampled to generate samples. The generated samples are then used to calculate a surrogate evidence lower bound. Finally, the variational parameters are re-optimized using stochastic gradient descent, resulting in continuous iterative updating of the sampled particles and variational parameters. This method transforms nonlinear state estimation and maneuvering model parameter identification into an evidence lower bound optimization problem. Stochastic gradient descent is used to achieve joint optimization of target state estimation and model parameter identification, enabling simultaneous state estimation and model parameter identification. This effectively addresses nonlinear filtering issues and improves nonlinear maneuvering target tracking performance. The method is advantageous in that it can adapt to nonlinear dynamic systems with unknown model parameters.

Description

Technical Field

[0001] The invention belongs to the field of target tracking and relates to a maneuvering target tracking method based on hybrid reasoning. Background Art

[0002] Using remote sensing devices such as radar, sonar, and optoelectronics to continuously track target positions in a timely, accurate, and diverse manner has wide applications in areas such as long-range early warning and unmanned driving. The inherent nonlinear characteristics of tracking systems (such as nonlinear motion of orbiting targets and nonlinear coordinate transformations between target states and sensor measurements), as well as the mismatch between system modeling and actual target motion caused by non-cooperative target maneuvers, severely impact tracking system performance. Consequently, high-precision and stable tracking of nonlinear maneuvering targets has long been a key and challenging task for radar, sonar, and optoelectronic sensing and detection systems. The nonlinear maneuvering target tracking problem can be decomposed into two subproblems: nonlinear filtering and adaptive filtering. Most maneuvering target tracking methods employ a sequential processing architecture, combining nonlinear state estimation with model selection or parameter identification. These two subproblems are coupled, with nonlinear filtering accuracy affecting model identification accuracy and model identification irregularities affecting filtering effectiveness. A joint optimization framework, which simultaneously performs state estimation and model parameter identification, often achieves better performance than sequential processing. However, solving the joint posterior probability density function presents a complex integral problem, making closed-loop analytical solutions difficult to obtain, necessitating the development of approximate inference methods.

[0003] Existing approximate inference methods fall into two main categories: sampling-based stochastic approximation methods, such as Markov Monte Carlo (MCMC) and Sequential Monte Carlo (SMC), and optimization-based deterministic approximation methods, such as Variational Bayesian (VB) and Expectation Propagation (EP). Both types of methods can be used to solve nonlinear filtering and adaptive robust estimation problems.

[0004] MCMC and SMC based on sampling approximation are highly scalable and can effectively handle nonlinear problems, but are only suitable for small-scale problems. When dealing with nonlinear estimation problems of high-dimensional latent variables, the rapid increase in the number of particles leads to a large amount of computation. VB based on deterministic approximation uses analytical approximation of the posterior probability density function, which is computationally less complex and easier to handle large-scale problems than MCMC, such as the joint estimation of high-dimensional latent variables. However, it requires the construction of a conjugate prior distribution of the joint latent variables, and nonlinear problems make it impossible for standard VB to maximize the lower bound of evidence using the coordinate ascent method.

[0005] Hybrid inference methods based on sampling approximation and deterministic approximation can effectively utilize the advantages of both approximation methods and are becoming a research hotspot in the international field of machine learning. However, in current research, this hybrid inference method is mostly only applied in the field of statistics. Summary of the Invention

[0006] In order to solve the problems in existing maneuvering target tracking methods that nonlinear filtering accuracy affects the filtering effect and the solution process is complicated, the present invention proposes a maneuvering target tracking method based on hybrid reasoning. This method effectively handles nonlinear problems based on variational Bayesian and sequential Monte Carlo reasoning, improves the nonlinear maneuvering target tracking performance, and is suitable for nonlinear dynamic systems with unknown model parameters.

[0007] The technical solution of the present invention is:

[0008] A method for tracking a maneuvering target based on hybrid reasoning is particularly pointed out to include the following steps:

[0009] Step 1: System modeling and problem description:

[0010] Step 1.1: Establish a discrete-time nonlinear dynamic system model:

[0011] x k =f k (x k-1 )+G k v k

[0012] y k =h k (x k )+w k

[0013] in, represent the target state and sensor measurement respectively, f k (·),h k () represent the nonlinear state transfer function and measurement function, G k is the control matrix of process noise, process noise and sensor measurement noise

[0014] All obey zero-mean Gaussian white noise, n x and n y They represent the dimensions of the system state x and sensor observation y, respectively, and the subscript k represents the time;

[0015] Step 1.2: Select the parameter of the mismatch between the prior system and the actual system when the target maneuvers as the unknown model parameter θ k ,θ k Include: Q k and R k are the covariance matrices of process noise and sensor measurement noise, respectively, and P 0|0 are the mean and variance of the initial state x0 respectively;

[0016] Step 1.3, solve the joint probability density function p(x k ,θ k |y 1:k ), for the target state x k and unknown model parameters θ k Perform joint estimation;

[0017] Step 2: nonlinear filtering;

[0018] Step 2.1, according to the initial variation parameter λ1, the proposed distribution is designed as r(x k ,θ k ;λ k ), where λ k is x k ,θ k Hyperparameters of

[0019] Step 2.2: Generate samples from the proposed distribution and perform particle prediction;

[0020] At the initial time k = 1, sample N particles from the initial proposed distribution r(x1;λ1) and calculate the particle weight

[0021] When k>1 at the subsequent moment, the particles at the previous moment are first resampled according to the probability, and then the particles are resampled according to the recommended distribution r(x k ,θ k ;λ k ) generates one-step predictions for all particles, and the particle weights are updated as Step 2.3: Use the particle weights generated in step 2.2 to calculate the lower bound of the smoothed substitution evidence. and filtered surrogate evidence lower bound

[0022]

[0023]

[0024] Step 2.4: Use stochastic gradient descent to iteratively optimize the lower bound of the alternative evidence. The optimization objective function is:

[0025]

[0026]

[0027] Repeat steps 2.2 to 2.4, and terminate the iteration when the number of iterations is reached or the lower bound converges;

[0028] Step 3: Target state estimation;

[0029] According to the optimized target parameter λ 1:k or λ kRepeat step 2.2 to obtain the time series and the particle set at the current moment The target state estimation is obtained through particle weighted calculation, namely:

[0030]

[0031]

[0032] Furthermore, in step 1.1, let the initial state be x0, which obeys Gaussian distribution with mean and variance respectively. P 0|0 , x0, v k and w k Independent of each other.

[0033] Furthermore, in step 1.3, under the framework of Bayesian estimation, the joint probability density function p(x k ,θ k |y 1:k ) method is: given an initial value p(x0,θ0), recursively perform the following two steps:

[0034] p(x k ,θ k |y 1:k-1 )=∫p(x k |x k-1 ,θ k )p(θ k |θ k-1 )p(x k-1 ,θ k-1 |y 1:k-1 )dx k-1 dθ k-1 ;p(x k ,θ k |y k )∝p(y k |x k ,θ k )p(x k ,θ k |y k-1 ).

[0035] Furthermore, in step 2.2, at the subsequent time k>1, first the particle at the previous time By particle weight Resample the probability, and then make a one-step prediction on the resampled particles, that is, The particle weight is updated to

[0036] Furthermore, the smoothed surrogate evidence lower bound and filtered surrogate evidence lower bound Calculated by the following formula:

[0037]

[0038] Calculate the unbiased estimate of the sequence marginal likelihood function based on the particle weights generated in step 2.2 and a biased but consistent estimate of the likelihood function at the current moment,

[0039]

[0040]

[0041] Then the smoothed surrogate evidence lower bound and filtered alternative evidence lower bounds They are:

[0042]

[0043]

[0044] Furthermore, in step 2.4, the Monte Carlo MC is used to obtain and its unbiased estimate of the gradient, iteratively optimize λ k ,Right now

[0045]

[0046] Where n = 1,...,N SGD represents the number of iterations, represents the gradient operator, the learning rate Satisfies the Robbins-Monro condition and is obtained using the Adam method.

[0047] Furthermore, in step 2.4, in order to reduce the estimated covariance, the reparameterization technique is used to obtain the gradient The calculation formula is as follows:

[0048]

[0049]

[0050] Where L<<N is the number of sampling particles of the random gradient.

[0051] The beneficial effects of the present invention are:

[0052] 1. The maneuvering target tracking method provided by the present invention is based on variational Bayesian and sequential Monte Carlo hybrid reasoning to achieve simultaneous state estimation and model parameter identification.

[0053] 2. The method of the present invention transforms nonlinear state estimation and maneuver model parameter identification into an optimization problem of the evidence lower bound. It designs a parameterized variational proposal distribution under the framework of Sequential Monte Carlo, constructs a closed-loop iterative mechanism of sampling particle weights and proposal distribution hyperparameters, and uses stochastic gradient descent to achieve joint optimization of target state estimation and model parameter identification, effectively dealing with the nonlinear filtering problem and improving the nonlinear maneuvering target tracking performance.

[0054] 3. The method of the present invention is applicable to nonlinear dynamic systems with unknown model parameters.

[0055] 4. Through typical maneuvering target tracking simulation data and actual data tests on sea target tracking, it is shown that the method of the present invention has obvious advantages over the existing sequential Monte Carlo-based methods or variational Bayesian-based methods and extended Kalman filter methods, and has high tracking accuracy. BRIEF DESCRIPTION OF THE DRAWINGS

[0056] Figure 1 It is a flowchart of the method for tracking a maneuvering target based on hybrid reasoning of the present invention;

[0057] Figure 2 This is a target motion trajectory diagram in scene 1 according to an embodiment of the present invention;

[0058] Figure 3 This is a RMSE comparison curve diagram of the embodiment of the present invention in scenario 1;

[0059] Figure 4 This is a NPE comparison curve chart of an embodiment of the present invention in scenario 1;

[0060] Figure 5 This is a target motion trajectory diagram in scene 2 according to an embodiment of the present invention;

[0061] Figure 6 This is a RMSE comparison curve diagram of the embodiment of the present invention in scenario 2;

[0062] Figure 7 This is a NPE comparison curve chart of an embodiment of the present invention in scenario 2;

[0063] Figure 8 This is a diagram of the actual motion trajectory of the target in the actual measurement scene according to the embodiment of the present invention;

[0064] Figure 9 This is a RMSE comparison curve diagram of the embodiment of the present invention in the actual measurement scene;

[0065] Figure 10 The following is a comparison curve of the NPE of the embodiment of the present invention in the actual measurement scenario. DETAILED DESCRIPTION

[0066] The method of the present invention is described in detail below with reference to the accompanying drawings and specific embodiments.

[0067] See also Figure 1 — Figure 10 To address the issues in existing nonlinear maneuvering target tracking methods, where nonlinear filtering accuracy affects the filtering effect and the solution process is complex, the present invention proposes a maneuvering target tracking method based on hybrid inference. Nonlinear maneuvering targets include moving objects such as ships, vehicles, and aircraft that can be detected by sensors (e.g., radar). The method specifically includes the following steps:

[0068] Step 1: System modeling and problem description,

[0069] Step 1.1: State space modeling:

[0070] Discrete-time nonlinear dynamic system models:

[0071] x k =f k (x k-1 )+G k v k

[0072] y k =h k (x k )+w k

[0073] in, Represent the target state and sensor measurement respectively. k (·),·),h k (·)h k () represent the nonlinear state transfer function and measurement function respectively. Assume that the process noise and measurement noise All obey zero-mean Gaussian white noise, and their noise covariance matrices are Q k and R k , G k is the control matrix of process noise. x and n y They represent the dimensions of the system state x and sensor observation y, respectively, and the subscript k represents the time. Assume that the initial state x0 follows a Gaussian distribution with a mean and variance of P 0|0 Furthermore, assuming that x0, v k and w k Independent of each other. Step 1.2: Determine the unknown model parameters:

[0074] Since the target maneuver will cause the prior system model parameter θ k mismatch with the actual system, then θ k Include:

[0075]

[0076] Step 1.3: Problem Description:

[0077] In the framework of Bayesian estimation, solve the joint probability density function p(x k ,θ k |y 1:k ), for the target state x k and model parameters θ k The joint estimation is performed as follows: after giving the initial value p(x0,θ0), recursively perform the following two steps:

[0078] p(x k ,θ k |y 1:k-1 )=∫p(x k |x k-1 ,θ k )p(θ k |θ k-1 )p(x k-1 ,θ k-1 |y 1:k-1 )dx k-1 dθ k-1 p(x k ,θ k |y k )∝p(y k |x k ,θ k )p(x k ,θ k |y k-1 )

[0079] Step 2: Perform nonlinear filtering;

[0080] Step 2.1: Design the proposed distribution r(x k ,θ k ;λ k ), where λ k is x k ,θ k Hyperparameters of

[0081] According to the prior function f(x k ,θ k ) Design suggestion distribution r(x k ,θ k ;λ k ), where λ k is x k ,θ k Hyperparameters of

[0082] Step 2.2: Particle propagation according to the proposed distribution:

[0083] At the initial moment, k=1, N particles are sampled from r(x1;λ1), that is, The particle weight is At the subsequent moment, k>1, first the particle at the previous moment Resampling is performed according to probability, and in the embodiment, according to particle weight Resample and then make a one-step prediction on the resampled particles, i.e. The particle weight is updated to

[0084] Step 2.3: Calculate the lower bound of the alternative evidence;

[0085] Since the standard smooth evidence lower bound B(λ 1:k ) and the filtering evidence lower bound B(λ k ) is difficult to calculate, so and As the lower bounds of the smoothing surrogate evidence and the filtering surrogate evidence, respectively, the following inequalities hold:

[0086]

[0087]

[0088] Calculate an unbiased estimate of the sequence marginal likelihood function based on the particle weights generated in step 2.1

[0089]

[0090] And the biased but consistent estimate of the likelihood function at the current moment, that is

[0091]

[0092] Therefore, the lower bounds of the smoothing surrogate evidence and the filtering surrogate evidence are:

[0093]

[0094]

[0095] Step 2.4: Use the stochastic gradient descent method to re-optimize the next level of alternative evidence, so that the sampling particles and variational parameters are continuously iteratively updated.

[0096] The optimization objective function is:

[0097]

[0098]

[0099] Weight w k With the hyperparameter λ k The stochastic gradient descent method is used to optimize the lower bound of the alternative evidence and the Monte Carlo MC is used to obtain and its unbiased estimate of the gradient, iteratively optimize λ k ,Right now

[0100]

[0101] Where n = 1,...,N SGD represents the number of iterations, represents the gradient operator, the learning rate Satisfies the Robbins-Monro condition and is obtained using the Adam method.

[0102] In order to reduce the estimated covariance, the reparameterization technique is used to calculate the gradient The basic idea is to assume that the distribution s and parameter λ k It doesn’t matter. For any deterministic function f, there is x k =f(x k-1 ,y k ,ε k ;λ k ), we can first sample ε from s k ~s(ε k ), obtain a set of particles Substitute into the deterministic function f to obtain a new particle

[0103] The gradient calculation formula is as follows:

[0104]

[0105]

[0106] Where L<<N is the number of sampling particles of the random gradient.

[0107] Repeat steps 2.2 to 2.4 and terminate the iteration when the number of iterations or convergence condition is reached.

[0108] Step 3: Estimate the target state.

[0109] According to the optimized parameter λ 1:k or λ k Repeat step 2.2 to obtain the time series and the particle set at the current moment The target state estimation is obtained through particle weighted sum calculation.

[0110] The following takes the target performing a turning maneuver as an example to explain this method in detail.

[0111] Scenario 1 is a single target tracking scenario considering time-varying turning maneuvers. Its motion trajectory is shown in the figure below. Figure 2 The target is moving in a straight line and turning at a constant speed. The radar coordinates are (55000m, 55000m). The initial state x0 is [30000m, 100m / s, 30000m, 0m / s]. T The total simulation step is 400s, the radar sampling period is T = 2s, and the target makes different degrees of turning maneuvers at k = {60, 120, 150, 214, 240, 272, 300, 338, 360}s, and the corresponding turning angular rate ω k ={0,3,0,5.6,0,8.6,0,-7.25,0,7.25} o / s. System noise variance η 2 =10m 2 / s 4 , measurement noise covariance R=diag(10000m 2 / s 2 ,2.7e-5rad 2 ).

[0112] Scenario 2: The total simulation step length is 200s, the radar sampling period is T = 1s, the target makes different degrees of turning maneuvers at k = {30, 60, 75, 107, 120, 136, 150, 169, 180}s, and the system noise variance η 2 =7.5m 2 / s 4 , the rest of the parameters are the same as the previous scene, the motion trajectory is as follows Figure 5 shown.

[0113] Test scenario: A shore-based navigation radar performs omnidirectional tracking of a vessel within a 20km radius of the sea surface. The radar operates in the 9.4GHz-9.5GHz band, with a sampling period of T = 3s, a radial range accuracy of 4.88m, and an azimuth accuracy of 9.5mrad. The target's actual trajectory is provided by the vessel's GPS. The target's actual trajectory is shown in the figure below. Figure 8 shown.

[0114] The specific tracking steps are as follows:

[0115] Step 1: System modeling and problem description;

[0116] Step 1.1 Build a turning maneuver model.

[0117] Consider a non-cooperative target tracking problem with a turning maneuver, where the turning angular velocity ω k Unknown. Let the target state Represents the position x of the target in the xy plane coordinate system at time k k,y k and speed The target's motion model is as follows:

[0118] x k =F k x k +G k v k

[0119] Among them, the state transfer matrix F k and process noise control matrix G k as follows:

[0120]

[0121]

[0122] Where T represents the sensor sampling period. Assume that the system process noise v k It has zero mean and variance η k 2 Gaussian white noise of I2, that is When the target turning angular velocity ω k = 0, the model is transformed into a uniform linear motion model.

[0123] Two-dimensional radar is used to track non-cooperative targets, and radar measurement k =[r k ,α k ] T Represents the radial distance and azimuth of the target with the radar as the reference origin. Its nonlinear measurement model is as follows:

[0124]

[0125] Among them, s x ,s y Represents the position of the radar on the x and y axes respectively, and the sensor measurement noise w k It also obeys the zero-mean Gaussian white noise.

[0126] Step 1.2, determine the unknown model parameters;

[0127] When the target makes a turning maneuver, the nonlinear state transfer function f k The turning angular velocity ω in (·) k And the corresponding covariance matrix Q k is unknown and time-varying, while the rest of the parameters are known.

[0128] Consider two simulation cases. In scenario 1, the turning maneuver angular velocity ω is k Unknown, process noise covariance Q k Known, the unknown model parameter θk ={ω k}; In scenario 2, the angular velocity of the turning maneuver is ω k and process noise covariance Q k are unknown, and the unknown model parameters θ k ={ω k ,Q k}.

[0129] Step 1.3: Joint posterior distribution.

[0130] For the target state x k and unknown model parameters θ k Perform joint estimation. Under the framework of Bayesian estimation, solve the joint probability density function p(x k ,θ k |y 1:k ), after giving the initial value p(x0,θ0), recursively perform the following two steps:

[0131] p(x k ,θ k |y 1:k-1 )=∫p(x k |x k-1 ,θ k )p(θ k |θ k-1 )p(x k-1 ,θ k-1 |y 1:k-1 )dx k-1 dθ k-1 p(x k ,θ k |y k )∝p(y k |x k ,θ k )p(x k ,θ k |y k-1 )

[0132] Step 2: Nonlinear filtering.

[0133] Step 2.1 designs a proposed distribution based on the initial variational parameters.

[0134] Prior function f(x k ,θ k )Set the recommended distribution r(x k ,θ k ;λ k ), where λ k is x k ,θ k Hyperparameters of . Figure 2 and Figure 8In the experimental scenario 1 and the measured scenario in the experimental part, for the case where there are unknown parameters in the nonlinear state transfer function, the unknown model parameter θ k ={ω k}, let the turning maneuver angular velocity ω k Obeying Gaussian distribution, its mean and variance is unknown and time-varying, so the hyperparameter Accordingly, the parameterized variational proposal distribution is

[0135] Reference Figure 5 , for the case where both the state transfer function and the model noise have unknown parameters, we further assume that the process noise covariance Q k Unknown, model parameters θ k ={ω k ,Q k}, and assuming that the system noise variance η 2 The obedience parameter is (a k ,b k ) of the inverse gamma distribution IG(a k / 2,b k / 2), so the hyperparameter Accordingly, the parameterized variational proposal distribution is The parameter a k =a k-1 +n x ,parameter

[0136] Step 2.2: Sample the proposed distribution to generate samples and particle prediction.

[0137] At the initial moment, k=1, N particles are sampled from r(x1;λ1), that is, The particle weight is

[0138] At the subsequent moment, k>1, first the particle at the previous moment By particle weight Resample the probability, and then make a one-step prediction on the resampled particles, that is, The particle weight is updated to

[0139] Step 2.3: Use the particle weights generated in step 2.2 to calculate the lower bound of the smoothed substitution evidence. and filtered surrogate evidence lower bound

[0140] Since the standard smooth evidence lower bound B(λ 1:k ) and the filtering evidence lower bound B(λ k) is difficult to calculate, so we need to find an alternative lower bound that is easy to handle. The unbiased estimate of the marginal likelihood function of the sequence is calculated by the particle weights generated by SMC. Right now

[0141]

[0142] And the biased but consistent estimate of the likelihood function at the current moment, that is

[0143]

[0144] Taking the logarithm of the former and applying Jensen's inequality, we get

[0145]

[0146] make Expand the above formula to get

[0147]

[0148] in,

[0149] make The following inequality holds:

[0150]

[0151]

[0152] and are the lower bounds of the smoothing surrogate evidence and the filtering surrogate evidence, respectively, and both are functions of the particle weights.

[0153] Step 2.4 uses stochastic gradient descent to optimize the lower bound of the smoothing surrogate evidence and the lower bound of the filtering surrogate evidence.

[0154] The optimization objective functions are:

[0155]

[0156]

[0157] Weight w k With the hyperparameter λ k The stochastic gradient descent method is used to optimize the lower bound of the alternative evidence and the Monte Carlo MC is used to obtain and its unbiased estimate of the gradient, iteratively optimize λ k ,Right now

[0158]

[0159] Where n = 1,...,NSGD represents the number of iterations, represents the gradient operator, the learning rate Satisfies the Robbins-Monro condition and is obtained using the Adam method.

[0160] gradient The reparameterization technique (reducing the estimated covariance) is used to calculate the result. The basic idea is: assume that the distribution s and the parameter λ k It doesn’t matter. For any deterministic function f, there is x k =f(x k-1 ,y k ,ε k ;λ k ), then first sample ε from distribution s k ~s(ε k ), obtain a set of particles Substitute into the deterministic function f to obtain a new particle The gradient calculation formula is as follows:

[0161]

[0162]

[0163] Where L<<N is the number of sampling particles of the random gradient.

[0164] The iteration is terminated when the number of iterations is reached or the lower bound converges.

[0165] Step 2.5: Estimate the target state.

[0166] According to the optimized parameter λ 1:k or λ k Repeat step 2.2 for particle propagation to obtain the particle set of the time series and the current moment Then the target state estimation is obtained by particle weighted sum calculation, that is,

[0167]

[0168]

[0169] Maneuvering target tracking assessment:

[0170] The following comparative evaluation of the nonlinear maneuver filtering effects of three methods, namely adaptive parameter estimation filter APE, extended Kalman filter method EKF and variational adaptive filter VBAKF, is conducted.

[0171] Evaluation metrics include the root mean square error (RMSE) of the radial distance to the target and the Normalized Position Error (NPE) relative to the measurement accuracy. Smaller RMSE and NPE indicate better filtering effectiveness, while an NPE less than 1 indicates effective filtering relative to the measurement. In the proposed method, the batch method for optimizing the lower bound of the smoothing evidence is called BVSMC, while the recursive method for optimizing the lower bound of the filtering evidence is called RVSMC.

[0172] The simulation results and analysis are as follows:

[0173] In scenario 1, considering the unknown turning angular velocity and the known noise covariance, the recommended distribution is r1(·). In order to meet the requirements of real-time processing, the BVSMC method considers the sliding window processing structure and assumes that the variational hyperparameters within the window are constant and the sliding window length is 3. BVSMC, RVSMC and APE all use 10,000 particles. The number of optimized particles for BVSMC and RVSMC is L=4, the iteration step is 0.01, the number of BVSMC iterations is 100, the number of RVSMC iterations is 10, and the initial optimization parameter λ0={10 -8 ,0.0001. The probability of change point in APE β=0.05, the smoothing parameter h of kernel density estimation 2 =0.01, in order to better turn, the turning angular velocity variance in EKF is set to 0.01rad 2 / s 2 , system noise variance η 2 =22m 2 / s 4 The number of Monte Carlo simulations is 100.

[0174] Figure 3 The following figure compares the RMSE curves of different methods. The results show that when the target moves at a constant speed or turns with a small angular rate (time periods A, B, D, and E), all methods can track the target well. When the target turns with a large angular rate (time periods F, G, H, I, and J), the performance of the EKF method decreases, while the other three methods can maintain stable tracking results.

[0175] Figure 4NPE curve comparison chart showing different methods. The results indicate that when the target is moving at a constant speed or making a turning motion with a small angular rate (time periods A, B, D, E), all methods generally have a filtering effect; when the target makes a turning motion with a large angular rate (time periods F, G, H, I, J), the performance of the EKF method degrades, and the overall filtering accuracy is lower than the measurement accuracy, while other methods can maintain a good filtering effect. At the same time, when the target makes a turning maneuver, the NPE peak values of RVSMC and BVSMC are smaller than those of APE, indicating that the filtering robustness of the VB-SMC-based method is better than that of APE. Table 1 records the average RMSE and NPE of each method under different turning maneuver angular rates of the target.

[0176] Table 1 Performance comparison statistical table of each method under different maneuver modes (in the case of the proposed distribution r1(·))

[0177]

[0178] The results show that when the target is moving in a straight line at a constant speed or making a small turning maneuver (A, B, C, D, E), all methods have a good filtering effect; when the turning maneuver is large (F, G, H, I, J), the filtering effect of EKF drops sharply, and at the same time, the filtering error of EKF accumulates and grows, while other adaptive filtering methods APE, RVSMC, and BVSMC can still maintain a robust performance, and the performance of RVSMC and BVSMC is generally better than that of APE. Generally speaking, the RMSE tracking accuracy of the position from high to low is: BVSMC (146.23m) < RVSMC (154.21m) < APE (169.53m) < EKF (183.10m) < measurement (190.49m).

[0179] At the same time, the NPE of the method of the present invention in each time period is less than 1, indicating that the method has high filtering effectiveness for target maneuvers. The filtering failure of EKF during large maneuvers is mainly because the dimension expansion method considers the slow-changing dynamic characteristics of the turning angular velocity. When the target makes a large maneuver turning, there is a time delay in the estimation of the angular velocity, resulting in an increase in the filtering error; APE, RVSMC, and BVMSC directly optimize the turning angular rate, enabling them to obtain a better estimation effect during large turning maneuvers. Since BVSMC comprehensively utilizes the information of multiple frames of measurements, its filtering effect is better than the sequential filtering methods APE and RVSMC, but its computational complexity is also the highest.

[0180] In Scenario 2: Consider the situation where both the turning angular rate and the process noise covariance are unknown, that is, the proposed distribution is r2(). The performance of the VBAKF method depends on the setting of the initial covariance Q0. A smaller Q0 will cause the method to diverge, so the system noise variance is set to 100×η 2, The number of iterations is 50. The initial hyperparameters a0 and b0 of APE, RVSMC, and BVSMC are both 5 and 15. The number of iterations for RVSMC and BVSMC is 10, and the parameters of other methods are the same as those in r1(). Figure 6 The RMSE curve graphs representing different methods. The results show that the estimation performance of VBAKF varies, being sometimes good and sometimes bad, as the maneuvering model changes. The main reason is that the filtering performance of VBAKF highly depends on the process noise covariance in the initial iteration, and the noise covariance estimation is a gradual adjustment process, making the method less robust when the model mutates. While other adaptive filtering methods, APE, RVSMC, and BVSMC, can still maintain robust performance, and the performance of RVSMC and BVSMC is overall better than that of APE.

[0181] Figure 7 The NPE curve graphs representing different methods are shown. The results show that the NPE of VBAKF has large fluctuations. The errors of APE, BVSMC, and RVSMC are all below the measurement error throughout the entire segment, and the NPE of RVSMC and BVSMC is less than that of APE, indicating that the filtering effects of RVSMC and BVSMC are better than that of APE. The error curves of RVSMC and BVSMC are similar, but the NPE of BVSMC is slightly less than that of RVMSC. Table 2 records the average RMSE and NPE of each method for the target under different turning maneuver angular rates. Table 2 is a statistical comparison table of the performance of each method under different maneuvering modes considering the case where both the turning angular rate and the process noise covariance are unknown, that is, when the proposed distribution is r2().

[0182] Table 2 Statistical comparison table of the performance of each method under different maneuvering modes (in the case of proposed distribution r2())

[0183]

[0184] The results show that the VBAKF method has poor robustness, being sometimes good and sometimes bad. While other adaptive filtering methods, APE, RVSMC, and BVSMC, can still maintain good robustness, and the performance of RVSMC and BVSMC is overall better than that of APE, and the filtering effect of BVSMC is the best. Generally speaking, the RMSE tracking accuracy of the position from high to low is: BVSMC (142.81m) < RVSMC (147.86m) < APE (157.08m) < VBAKF (166.72m) < Measurement (189.24m). This tracking result is consistent with the result of Scenario 1.

[0185] In the actual measurement scenario, the parameter settings are as follows: The number of optimized particles L of RVSMC and BVSMC is 4, the iterative optimization step size is 0.01, the number of iterations is 50, and the initial optimization parameter λ0 = {10 -8 , 0.0001}. The η of RVSMC, BVSMC, and APE2 = 0.1 m 2 / s 4 , the η of EKF and VBAKF 2 = 2 m 2 / s 4 , the variance of the turning angular velocity of EKF is 0.0001 rad 2 / s 2 , the number of iterations of VBAKF is 20.

[0186] Figure 9 and Figure 10 are respectively the comparison diagrams of the position RMSE and NPE curves of different methods. The tracking results of VBAKF and EKF are similar, and the curves basically coincide in the middle and later segments. The tracking results of APE and BVSMC are similar, but BVSMC performs well at the initial moment. The tracking accuracy of RVSMC is the highest, especially the filtering accuracy is greatly improved in the middle and later segments. Table 3 records the average RMSE of the position estimation for the entire segment of the measured data.

[0187] Table 3 Average RMSE of the Estimation of Measured Data

[0188] Measurement APE VBAKF EKF RVSMC BVSMC 7.66 7.50 7.56 7.58 6.10 7.35

[0189] As can be seen from Table 3, all methods can track the target. The tracking accuracies of APE, VBAKF, EKF, and BVSMC are similar, but the filtering effect is weak, and the accuracy is only improved by about 1 - 4%. The filtering effect of RVSMC is the most obvious, and the accuracy is improved by 20.4%. Generally speaking, the RMSE tracking accuracies of the positions from high to low are: RVSMC (6.10 m) < BVSMC (7.35 m) < APE (7.50 m) < VBAKF (7.56 m) < EKF (7.58 m) < measurement (7.66 m). The reason why the tracking accuracy of BVSMC in this scenario is not as good as that in the simulation scenario is that the target has been making an approximate turning motion throughout the entire segment, and the turning speed at each moment is not the same, resulting in the invalidation of the assumption that the variational hyperparameters in the window of BVSMC are constant. Therefore, the performance of the sequential filtering method RVSMC is better than that of the batch processing method BVSMC.

[0190] In the present invention, a maneuvering target tracking method based on the hybrid inference of variational Bayesian and sequential Monte Carlo is adopted, which effectively deals with the nonlinear problems during the target maneuvering process, improves the tracking performance of nonlinear maneuvering targets, and the method of the present invention is applicable to nonlinear dynamic systems with unknown model parameters.

Claims

1. A maneuvering target tracking method based on hybrid reasoning, characterized by: The following steps are involved: Step 1: System modeling and problem description: Step 1.1: Establish a discrete-time nonlinear dynamic system model: x k =f k (x k-1 )+G k v k y k =h k (x k )+w k in, represent the target state and sensor measurement respectively, f k (·),h k (·) represent the nonlinear state transfer function and measurement function, G k is the control matrix of process noise, process noise and sensor measurement noise All obey zero-mean Gaussian white noise, n x and n y They represent the dimensions of the system state x and sensor observation y, respectively, and the subscript k represents the time; Step 1.2: Select the parameter of the mismatch between the prior system and the actual system during the target maneuver as the unknown model parameter θ k ,θ k Include: Q k and R k are the covariance matrices of process noise and sensor measurement noise, respectively, and P 0|0 are the mean and variance of the initial state x0 respectively; Step 1.3, solve the joint probability density function p(x k ,θ k |y 1:k ), for the target state x k and unknown model parameters θ k Perform joint estimation; Step 2: Nonlinear filtering: Step 2.1, according to the initial variation parameter λ1, the proposed distribution is designed as r(x k ,θ k ;λ k ), where λ k is x k ,θ k Hyperparameters of Step 2.2: Generate samples from the proposed distribution and perform particle prediction; At the initial moment, k=1: sample N particles from the initial proposal distribution r(x1;λ1) and calculate the particle weight At the subsequent moment, k>1: first resample the particles at the previous moment according to the probability, and then resample the particles according to the proposed distribution r(x k ,θ k ;λ k ) generates one-step predictions for all particles, and the particle weights are updated as Step 2.3: Use the particle weights generated in step 2.2 to calculate the lower bound of the smoothed substitution evidence. and filtered surrogate evidence lower bound Step 2.4: Use stochastic gradient descent to iteratively optimize the lower bound of the alternative evidence. The optimization objective function is: Repeat steps 2.2 to 2.4, and terminate the iteration when the number of iterations is reached or the lower bound converges; Step 3: Target state estimation; According to the optimized target parameter λ 1:k or λ k Repeat step 2.2 to obtain the time series and the particle set at the current moment The target state estimation is obtained through particle weighted calculation, namely:

2. The method for tracking a maneuvering target based on hybrid inference according to claim 1, wherein: In step 1.1, let the initial state be x0, which obeys Gaussian distribution with mean and variance respectively. P 0|0 , x0, v k and w k Independent of each other.

3. The method for tracking a maneuvering target based on hybrid inference according to claim 1, wherein: In step 1.3, under the framework of Bayesian estimation, solve the joint probability density function p(x k ,θ k |y 1:k ) method is: given an initial value p(x0,θ0), recursively perform the following two steps: p(x k ,i k |y 1:k-1 )=∫p(x k |x k-1 ,i k )p(θ k |θ k-1 )p(x k-1 ,i k-1 |y 1:k-1 )dx k-1 dth k-1 ; p(x k ,θ k |y k )∝p(y k |x k ,θ k )p(x k ,θ k |y k-1 )。 4. The method for tracking a maneuvering target based on hybrid reasoning according to claim 1, wherein: In step 2.2, at the subsequent time k>1, first the particle By particle weight Resample the probability, and then make a one-step prediction on the resampled particles, that is, The particle weight is updated to 5. The method for tracking a maneuvering target based on hybrid reasoning according to claim 1, wherein: In step 2.3, the smoothed surrogate evidence lower bound and filtered surrogate evidence lower bound Calculated by the following formula: Calculate the unbiased estimate of the sequence marginal likelihood function based on the particle weights generated in step 2.2 and a biased but consistent estimate of the likelihood function at the current moment, Then the smoothed surrogate evidence lower bound and filtered alternative evidence lower bounds They are:

6. The method for tracking a maneuvering target based on hybrid reasoning according to claim 1, wherein: In step 2.4, we obtain and its unbiased estimate of the gradient, iteratively optimize λ k ,Right now Where n = 1,...,N SGD represents the number of iterations, represents the gradient operator, the learning rate Satisfies the Robbins-Monro condition and is obtained using the Adam method.

7. The method for tracking a maneuvering target based on hybrid reasoning according to claim 6, wherein: gradient Using the reparameterization technique, it is calculated using the following formula: Where L<<N is the number of sampling particles of the random gradient.

Images