Self-Tuning Process Noise in Variational Bayesian Adaptive Kalman Filter for Target Tracking

: Many practical systems, such as target tracking, navigation systems, autonomous vehicles, and other applications, are usually applied in dynamic conditions. Thus, the actual noise statistics characteristics of these systems are generally time varying and unknown, which will deteriorate the state estimation accuracy of the Kalman ﬁ lter (KF) and even cause ﬁ lter diverging. To address this issue, this paper proposes an adaptive process noise covariance ( Q k )-based variational Bayesian adaptive Kalman ﬁ lter (AQ-VBAKF) algorithm. Firstly, the adaptive factor is introduced to self-tune the process noise covariance; the adaptive factor is obtained based on the innovation sequences, which can adapt to the input measurement values. Then, the VB solution is applied to approximate the time variant and unknown measurement noise covariance. Therefore, this proposed algorithm can adjust the process noise covariance and the measurement noise covariance simultaneously based on the variable input signals, which can improve the self-adaptive ability of the state estimation ﬁ lter in dynamic conditions. According to the dynamic target tracking test results, the proposed AQ-VBAKF outperforms several other existing ﬁ ltering methods in estimation accuracy, robustness, and computational e ﬃ ciency.


Introduction
In many practical applications, such as target tracking, positioning, autonomous vehicles, and so on, state estimation has been widely deployed [1][2][3][4].The KF (Kalman filter) can obtain the optimal estimation for the minimum mean square error (MMSE) if the noise covariances are accurate under the linear Gaussian condition [5].In the KF iteration, noise covariances are preselected and kept constant throughout the whole filtering process.However, when the practical systems work in dynamic environments, the constant noise covariances cannot reflect the real dynamic circumstances.For example, autonomous vehicle systems [6] inherently work in dynamic environments, where the sensor measurement, environment conditions, and vehicle dynamics are time varying and difficult to obtain.In this condition, the accurate noise statistics of practical systems are usually unknown.The actual noise covariances vary with the changes in environments.Thus, constant noise covariances for KF may degrade the filtering performance in dynamic conditions.Therefore, advanced estimation methods are necessary for autonomous vehicle systems, target tracking systems, and other applications.
Accordingly, to address this issue, several adaptive KF (AKF) methods have been developed.The AKF algorithms include the multiple mode AKF [7], the Sage-Husa AKF [8], the fading AKF [9], and the variational Bayesian-based AKF (VBAKF) [10].The multiple mode AKF requires a bank of KFs operating simultaneously, resulting in a large computational burden.The Sage-Husa AKF uses the covariance matching method to estimate the maximum posterior noise, which cannot ensure that it will converge to the accurate noise covariances.The fading AKF introduces the fading factor to decrease the weight of the current measurement value.However, the fading AKF usually needs multiple fading factors, which are not easy to calculate [11].
The variational Bayesian-based AKF (VBAKF) has been widely utilized, as it can use a factorized free-form distribution to approximate the unknown measurement noise covariance matrix (MNCM) and the state vector's joint posterior distribution [12,13].However, in the VBAKF algorithm, the process noise covariance matrix (PNCM) is predefined and kept constant, which cannot reflect the actual noise environment and, thus, will deteriorate the filtering performance in time-varying situations.Huang et al. [14] proposed the VBAKF-PR algorithm to estimate the online MNCM and the predicted state error covariance matrix together, and more accurate results were obtained.However, the VBAKF-PR approach needs a predefined PNCM and has a large computational burden.In reference [15], the predicted state error covariance matrix was estimated, but the PNCM stayed constant.In reference [16], to overcome the issue that the PNCM's distribution and the system state vector's distribution are non-conjugate, black-box variational inference was used to solve this problem.In [17], a sliding window variational AKF was proposed, which relies on approximating the smoothing posterior distribution.However, these algorithms either have a high computational burden or the PNCM stays constant, which will decrease the filtering performance in dynamic conditions.Since the inaccurate predefined PNCM is used for each iteration, this may cause a decline in estimation accuracy.
To enhance the VBAKF's performance in dynamic systems, an adaptive PNCM ( k

Q
)-based variational Bayesian adaptive KF (AQ-VBAKF) method is proposed in this paper.
The main contributions of this paper are as follows: 1. To further improve the VBAKF algorithm in dynamic systems and to improve the adaptability of PNCM, the proposed AQ-VBAKF algorithm can not only estimate the MNCM online, but can also self-tune the PNCM in real time.Therefore, the AQ-VBAKF can adapt to the input variable signals in the linear dynamic system.2. The VBAKF algorithm can only estimate the MNCM, but the PNCM stays constant during the whole filtering process, which will decrease the filtering performance.In the proposed algorithm, the adaptive factor is introduced to tune the PNCM in real time, meaning the fixed PNCM does not cause a decline in the VBAKF estimation accuracy at each iteration.3. The proposed algorithm is evaluated using the target tracking simulation.Based on a series of experiments, the proposed AQ-VBAKF can improve the filtering accuracy and robustness in comparison with other filtering methods under dynamic conditions.Moreover, the computational burden is not heavy, thus it can be applied in practical systems.
The structure of the paper is summarized as follows.The state-space model and problem statement are described in Section 2. The proposed AQ-VBAKF approach is illustrated in Section 3. The AQ-VBAKF method is verified by target tracking tests in Section 4. The conclusions are drawn in Section 5.

State-Space Model and Problem Statement
For the linear Gaussian discrete-time system, the state-space model is demonstrated as follows: where Equation ( The KF algorithm includes two procedures, which are the prediction process and the correction process.The prediction procedure is expressed as [18]: Then, the correction process is defined as [18]: where In dynamic conditions, such as integration navigation systems, target tracking systems, and autonomous vehicles, the actual noise covariances vary with the changes in environments and the accurate noise statistics of practical systems are usually unknown.The inaccurate and constant noise covariances k Q and k R in the KF may cause subop- timal estimations and even divergence.Therefore, an adaptive KF algorithm for the inaccurate noise covariance matrices is necessary under practical dynamic conditions.Accordingly, a novel AQ-VBAKF approach is proposed to self-tune the variable and unknown MNCM and PNCM under dynamic conditions.

Variational Bayesian Approximation for MNCM ( k R )
According to Equations ( 1) and (2), the probability form for the state-space model is [19]: where N( | , )  μ Σ represents the multivariate Gaussian distribution with a mean vector μ and covariance matrix Σ .Because k R is unknown, the VB algorithm is applied to acquire the posterior distribution However, it is usually not analytical to handle [12].The approximation scheme VB is used to approximate the intractable posterior inference in the Bayesian model [12].
In linear and Gaussian systems, k x is considered to be independent of the NCM.
The joint distribution at the (k−1)-th epoch is decomposed as: x z is assumed as the Gaussian distribution, and can be considered as an inverse Wishart (IW) distribution, which is expressed as [16]: x z x x P (10) where IW( | , ) w  W represents the IW distribution.w denotes the degrees of freedom's number.W represents the inverse scale matrix.
Employing the VB inference, the joint posterior distribution can be approximated as: where ( ) q  is the approximation for ( ) p  .( ) k q x and ( ) k q R are obtained based on the minimization of the Kullback-Leibler (KL) divergence between the real posterior distribution and the approximation ( ) ( ) where KL( m( )||n ( )) x x represents the KL divergence between m(x) and n(x).
The optimal solution of Equation ( 13) can be obtained as follows [20]: where E[•] represents the expectation operator.Because of the Bayesian rule, it can be given that: Due to x , then: Based on Equations ( 8), (10), and (11), x R z z is decomposed as [21]: Accordingly, we can have: where tr[•] denotes the matrix's trace.
x n and z n denote the dimensions of k x and k z , respectively.C represents a constant independent of k x and k R .
Substituting Equations ( 17) and (19) into Equation ( 14), we can obtain the updated where The superscript s denotes the s-th VB iteration.The dimensions of Substituting Equations ( 17) and (19) into Equation ( 15), ( ) k q x is updated as: where In the KF estimation, when the state estimations are accurate, the following equation is held [22]: where I represents the identity matrix and k d denotes innovation sequences, that is: where k+ j k d d should be satisfied [22].Therefore, according to Equation (29), k K should be selected to establish the fol- lowing equation: By substituting Equation (5) into Equation (31), it can be given that: To satisfy Equation (32), the following equation should be held: Because of Equation (28), then it can then be determined that: Therefore, the adaptive factor k  can be obtained as: To keep the whole filtering process stable, the adaptive factor k  should be larger than 1; that is: When the system experiences sudden changes, and the adaptive factor is larger than 1, then the corresponding Kalman gain will become larger based on Equation ( 5).Therefore, according to Equation ( 6), more weight will be added to the current measurement data and less weight will be given to the past state information, and then the state estimation will be more accurate.Therefore, k  is obtained with: The innovation covariance k C is defined as: In general, it is difficult to obtain the actual innovation covariance.Thus, we use the window function method to estimate the innovation covariance [23]: where ˆk C denotes the estimated innovation covariance, and N > 1 represents the window size and is a scalar.

Proposed Algorithm
To summarize, the structure diagram of the proposed algorithm is shown in Figure 1.The detailed algorithm is demonstrated in Algorithm 1, where   represents the forgetting factor and S denotes the VB iteration number.

Time update for state values:
1: ˆk k x can be obtained via Equation (3).

The adaptive factor update and the time update for
| 1 k k P : 2: The adaptive factor is computed using Equations (30), (37), and (39).

P
is obtained with Equation (28).

Simulation Verification and Analysis
In this section, we adopt a target tracking test to assess the availability and superiority of the proposed AQ-VBAKF algorithm.The target tracking example's state equation and measurement equation are: where x y x y    x is the system state vector.Then, ( , ) 0.5 [0.5cos( ) 6.5] 0.5 where 2 I is a two-dimensional identity matrix.T denotes the simulation time and T = 1000 s.Additionally, q is set as 1 m 2 /s 3 and r is set as 100 m 2 .Moreover, the initial val-

Estimation Performance Analysis
To verify the filtering accuracy, we compare the proposed AQ-VBAKF with five other filtering estimation methods: 1. KFTCM algorithm: the KF is iterated with true covariance matrices (KFTCM), which act as a reference.

VBKF algorithm:
k R is updated using the VB method and k Q is constant through- out the whole process.
In the simulation, the initial values are set to 0|0 diag([100, 100, 100, 100])  P and T 0 =[100, 100, 10, 10] x . Table 1 shows the different algorithms' parameter settings.S is the VB iteration times.The iteration times should be proper to ensure the VB approximate converges to the optimal results.Generally, VB iteration times should be more than 6 to make the filter converge, thus S is chosen as 15 [24].N represents the window size for the AQKF and AQ-VBAKF algorithms, which is analyzed in the simulation later.  to ensure that the posterior and prior PDFs have the same forms.τ is used for the VBAKF-PR algorithm.x y and ˆ( , ) x y represent the true and the estimated positions at the l-th MC test, respectively.T is set to 1000 s.For velocity, the definitions of RMSE and ARMSE are similar.Figure 2 demonstrates the RMSEs for different algorithms, where the KFTCM has the smallest RMSEs, and acts as a reference.In comparison with the KFNCM, AQKF, VBAKF, and VBAKF-PR algorithms, the proposed AQ-VBAKF has the smallest RMSEs for position and velocity.The VBAKF has the worst performance, which indicates that the VBAKF is even worse than the KFNCM when the initial values are far away from the true values.The VBAKF cannot guarantee an optimal outcome in any case.Due to the adaptive predicted state error covariance matrix, the VBAKF-PR shows higher accuracy than the VBAKF, but is also worse than the KFNCM due to the initial value deviating from the true values.The AQKF is slighter worse than the proposed AQ-VBAKF since the AQKF lacks the ability to update the MNCM ( k R ).

Robustness Analysis for PNCM and MNCM Initial Value Settings
To  M is set to 20.Moreover, Table 1 demonstrates other parameter settings, including ξ, S, N, and τ.In Figure 4, as we can see, the ARMSEs for position and velocity stay small and stable in most cases.However, when the parameters are in the range of ( , ) [100, 1000] [0.1, 1] , the ARMSEs for position and velocity increase, which indicates that the filtering estimations are not so accurate.Figure 5 shows the VBAKF-PR's results.In the range of ( , ) [0.1, 1000] [0.1, 1] , as the parameter β decreases, the ARMSEs increase.Therefore, the AQ-VBAKF is more robust than the VBAKF-PR in the range of ( , ) [0.1, 100] [0.1, 1] , in which the initial values of the PNCM ( 0 Q ) and MNCM ( 0 R ) deviate from the true value.Therefore, when the initial values of β are smaller than 1, which represents that the initial MNCM ( 0 R ) deviates from the true value, the proposed AQ-VBAKF algorithm is more robust than the VBAKF-PR.In many practical systems, the initial state values are usually unknown.In this condition, when the initial settings of the MNCM and PNCM values differ from the true values, the proposed AQ-VBAKF performs better in self-adaptability and robustness than the VBAKF-PR algorithm.

Analysis of Window Size Setting
To analyze the effect of the window size on the proposed AQ-VBAKF algorithm, some further simulations are performed.The window size is set as N = 1:1:25, and other parameters are the same as those in Table 1.The ARMSEs of different window sizes are shown in Figure 6.The KFTCM has the smallest ARMSEs of position and velocity, which acts as a reference.In Figure 6a, the ARMSEs of the position for the proposed AQ-VBAKF are smaller than that for the VBAKF-PR and KFNCM when N = 1:1:25.In Figure 6b, the ARMSEs of the position for the proposed AQ-VBAKF are smaller than that for the VBAKF-PR and KFNCM when N ≥ 4. As we can see, the ARMSEs for both position and velocity cannot decrease as the window size N increases when N = 1:1:25.When the window size N is larger than 10, the ARMSEs for both the position and velocity estimations tend to be stable.
In dynamic conditions, the window size cannot be too large.A small window size can help the innovation covariance adapt to the time-varying input signals [25].On the other hand, the window size should be large enough to enhance the innovation covariance estimation stability.Moreover, if the window size is too small, the computational burden will be heavier.Accordingly, taking into account the innovation covariance estimation accuracy and computational efficiency simultaneously in the dynamic environment, the window size N should be set as 4 ≤ N ≤ 10.

Computational Efficiency Comparison
To verify the computational efficiency, Table 2 shows the total simulation time for the different algorithms.The parameter settings are shown in Table 1, and c M = 1000.The proposed method adds only 9% computational time compared with the VBAKF method and is only about half of the VBAKF-PR simulation time.Therefore, the proposed method has almost the same computational burden compared with the VBAKF algorithm and the computational complexity is not high.To further evaluate the computational complexity, the total floating-point operations of each step for the proposed AQ-VBAKF were analyzed and are shown in Table 3, where x n and z n are the dimensions of the state vector and measurement vector, respectively.S is the VB iteration time.
Based on [14,17], the total floating-point operations of VBAKF-PR are: According to Equation (49), with the increase in the VB iteration time S, the computational burden of VBAKF-PR is larger than the proposed AQ-VBAKF, which is consistent with the computational time results.

Conclusions
This paper proposes a novel AQ-VBAKF algorithm that can adjust the PNCM and MNCM simultaneously.This algorithm introduces the adaptive factor to self-tune the PNCM in real time based on the input signals.Furthermore, the MNCM is estimated online using the VB approximation.According to the dynamic target tracking test results, the proposed AQ-VBAKF shows more accurate estimation and better robustness in comparison with other filtering methods.Moreover, for a general case in practical systems, when the initial PNCM and MNCM are not appropriately set and deviate from the true values, the proposed algorithm has better self-adaptability and robustness than other algorithms.The proposed algorithm's computational burden is not heavy; thus, it can be applied in practical engineering applications.
For future prospects, the proposed filtering estimation algorithm can be applied in dynamic conditions, such as autonomous vehicles, target tracking, positioning, and so on.For example, autonomous vehicles inherently work under dynamic environments, and the environment noise statistics are unknown.The adaptability of the proposed AQ-VBAKF algorithm can adjust process and measurement noise covariances dynamically, which can significantly enhance the vehicle's perception accuracy and response.Our research provides a theoretical basis in terms of dynamic conditions to enhance the filtering estimation accuracy.In future work, the proposed algorithm will be employed in positioning systems to enhance positioning accuracy.

k C R and k C
x denote the constants independent of k R and k x , respectively.

R
fact, the prior noise information of the sys- tem is usually unknown; thus, the setting of 0 Q and 0 R are usually not accurate.More- over, α and β are both variables to test the robustness of the algorithms with different 0 Q and 0 values.ξ is used for VB inference initialization, as shown in Algorithm 1

M
represents the Monte Carlo (MC) test number, and c M = 1000.( , )

Figure 2 .
Figure 2. RMSEs for different filtering algorithms: (a) the RMSEs of position for different filtering algorithms; (b) the RMSEs of velocity for different filtering algorithms.

Figure 3a shows the
Figure 3a shows the mean values of k  in the 1000 Monte Carlo (MC) tests based on

Figure 3 .
Figure 3. Values of k  in the proposed AQ-VBAKF algorithm: (a) the mean values of k  ; (b) the 100th MC test values of k  ; (c) the 500th MC test values of k  .
assess the proposed AQ-VBAKF's robustness, the ARMSEs of the AQ-VBAKF algorithm for different initial PNCM ( 0 Q ) and MNCM ( 0 R ) values were calculated and are presented in Figure 4.For comparison, Figure 5 demonstrates the ARMSEs for the VBAKF-PR.α and β are selected in the range of ( , ) [0.1, 1000] [0.1, 1000]     .To reduce the time consumption and keep the results reliable, c

Figure 6 .
Figure 6.ARMSEs of position and velocity under different window sizes N: (a) ARMSEs of position under different window sizes; (b) ARMSEs of velocity under different window sizes.
1) denotes the system model, and Equation (2) denotes the measurement model.k   represents the sampling time point.

Table 1 .
Parameters settings for different algorithms.
The RMSE (root mean square error) and ARMSE (average RMSE) are used to validate the state estimation accuracy based on the following formulas:

Table 2 .
Total computational time for different algorithms.