Multi-fading factor and update monitoring strategy adaptive Kalman filter based variational Bayesian with inaccurate time-varying process and measure- ment noise covariance matrices

Aiming at the problem that the performance of Adaptive Kalman filter estimation will be affected when the statistical characteristics of the process and measurement noise matrix are inaccurate and time-varying in the linear Gaussian state-space model, an algorithm of Multi-fading factor and update monitoring strategy adaptive Kalman filter based variational Bayesian is proposed. Inverse Wishart distribution is selected as the measurement noise model, the system state vector and measurement noise covariance matrix are estimated with the variational Bayesian method. The process noise covariance matrix is estimated by the maximum a posteriori principle, and the update monitoring strategy with adjustment factors is used to maintain the positive semidefinite of the updated matrix. The above optimal estimation results are introduced as time-varying parameters into the multiple fading factors to improve the estimation accuracy of the one-step state predicted covariance matrix. The application of the proposed algorithm in target tracking is simulated. The results show that compared with the current filters, the proposed filtering algorithm has better accuracy and convergence performance, and realizes the simultaneous estimation of inaccurate time-varying process and measurement noise covariance matrices.


Introduction
In many practical engineering applications, the actual values of the required state variables are often not directly available. For example, when radar detects an air target, it can calculate the target distance based on information such as reflected waves. Still, there is random interference in the radar detection process, resulting in random noise in the observation signal. In this case, it is impossible to obtain the required state variables accurately, and these state variables can only be estimated or predicted based on the observed signal. In linear systems, the Kalman filter is the optimal filter[1]. With the development of computer technology, the calculation requirements and complexity of Kalman filtering no longer become obstacles to its application [2]. At present, the Kalman filtering theory has been widely used in tracking, navigation, guidance, and other areas [3][4][5][6][7][8][9].
The application of the Kalman filter requires prior knowledge of the mathematical model of the system and the statistical characteristics of noise. Still, in many practical application problems, they are unknown or only partially known [10][11][12]. If inaccurate mathematical models or statistical noise characteristics are used to design the Kalman filter, the performance of the filter will be degraded, resulting in larger estimation errors, and even filter divergence. To solve this problem, various adaptive Kalman filters have been produced [13,14].
Sage-Husa filter (SH-KF) is widely used because of its simple algorithm, which can estimate the first and second moments of noise online [15]. However, the Sage-Husa adaptive noise estimator has problems such as a large amount of calculation and easy divergence of state estimation [16]. And some literature pointed out that the covariance matrix of process noise and observation noise cannot be estimated dynamically in real-time by the Sage-Husa adaptive estimator at the same time, and can only estimate the other noise covariance matrix when one noise covariance matrix is known [17,18]. The maximum likelihood-based adaptive filtering method (ML-KF) can evaluate and correct the second-order moments of the noise statistics online, but it needs to rely on accurate innovation covariance estimation, and ML-KF requires a large sliding window to obtain a precise estimate of the noise covariance matrix, which makes it only theoretically for the time-varying noise covariance matrix estimation [19,20]. Strong tracking Kalman filter (ST-KF) is an adaptive filtering algorithm based on the principle of residual orthogonality. It adjusts the weight of new measurement data by adding an estimate of the one-step predictive covariance matrix. It has strong robustness regarding model parameter mismatch, lower sensitivity to the statistical characteristics of noise and initial values, and a strong ability to track sudden changes. But its adjustment ability for each filtering channel is the same, and PNCM and MNCM are not estimated [21,22] In recent years, many scholars have introduced the variational Bayesian machine learning method into the KF algorithm and proposed adaptive Kalman filter (VB-KF) algorithm based on the variational Bayesian approach, which is an approximation of the Bayesian method. By choosing a suitable conjugate prior distribution, the slowly time-varying measurement noise covariance can be estimated [23][24][25]. Literature [26] proposed a variational adaptive Kalman filter,(R-VBKF) but only estimated the measurement noise covariance matrix; the accuracy is not satisfactory enough. In the algorithm presented in the literature [27], the state predicted error covariance matrix PECM and the measurement noise covariance matrix MNCM are estimated, but the process noise covariance matrix PNCM is not directly assessed.
Aiming at the linear Gaussian state-space model with slow time-varying covariance of process and measurement noise, taking into account the estimation accuracy, convergence performance, robustness, and the realization of simultaneous estimation of noise covariance matrices, the multifading factor and update monitoring strategy adaptive Kalman filter based variational Bayesian (MFMS-VBAKF) is proposed. Its feasibility is proved by simulation experiments.
The main structure is as follows. In Section 2, we modelled the problem mathematically. In Section 3, the multi-fading factor and update monitoring strategy adaptive Kalman filter based variational Bayesian is derived. In Section 4, the proposed algorithm is compared with the existing algorithm through the simulation of the target tracking application, and the excellent performance of the proposed algorithm is proved. Finally, Section 5 summarizes our conclusions.

Problem modelling
Consider the following discrete linear stochastic system of the state-space model where (1) and (2) are process and measurement equations, respectively. is discrete-time, ∈ × is the state vector of the system at time , ∈ × is the measurement signal vector of the corresponding state. ∈ × is the state-transition matrix, ∈ × is the measurement matrix.
∈ and ∈ are uncorrelated white Gaussian noise with zero mean vectors and covariance matrices and respectively. The initial state 0 is assumed to be a Gaussian distribution with mean vector ̂0 and the covariance matrix 0 . 0 is uncorrelated to and at any time [1].
For linear Gaussian state-space models, the Kalman filter (KF) algorithm is an optimal estimation filter algorithm. If the noise covariance matrices and are fully known, KF estimates the state vector through the measurement information of 1: , and the estimation accuracy is satisfactory. However, the performance of the KF algorithm overly depends on the prior knowledge of the noise statistics. If the time-varying noise covariance matrices and are unknown or inaccurate, the accuracy of the KF algorithm will decrease, and even cause the estimation to diverge. Besides, when most existing adaptive Kalman filter algorithms estimate PNCM and MNCM at the same time, the filtering will diverge. Therefore, A Multi-fading factor and update monitoring strategy adaptive Kalman filter Based variational Bayesian with inaccurate time-varying PNCM and MNCM is proposed.

3.The proposed Multi-fading factor and update monitoring strategy adaptive Kalman filter Based variational Bayesian
In the VBAKF algorithm, the independent state vector and the measurement noise covariance matrix are regarded as the parameters to be estimated.

Adaptive Kalman filter based variational Bayesian (VBAKF)
3.1.1 Prediction process and the selection of the prior distribution In the traditional Kalman filter framework. The Gaussian distributions are selected as the distributions of one-step predicted Probability Density Function (PDF) ( | 1:k−1 ) and likelihood PDF ( | k ): ( | 1: −1 ) = ( −1 ;̂: −1 , : −1 ), where ( ; , ) is the Gaussian distribution, and represent the mean and variance of the distribution, respectively. The PDF of the Gaussian distribution is: According to Equation (1), the predicted state vector ̂: −1 and the corresponding one-step predicted error covariance matrix (PECM) : −1 can be described as: where ̂: −1 and : −1 represent the state estimation at time − 1 and the corresponding estimation error covariance matrix, respectively. (. ) represents the transpose of the matrix. Among them, it is assumed that the real PECM is unknown due to the complex environmental factors of target tracking, which will lead to an inaccurate : −1 in equation (7). The estimation methods for and : −1 will be given in the next two sections. In this section, the aim is to infer with . For the purpose, the conjugate prior distributions need to be first selected for inaccurate MNCM since the conjugacy can ensure that the posterior distribution and the prior distribution maintain the same functional form.
According to Bayesian statistical theory, if the Gaussian distribution has a known mean, the conjugate prior distribution of the covariance matrix can be regarded as the inverse Wishart (IW) distribution [28] . −1 is the inverse matrix of a positive definite matrix . If −1 follows the Wishart-distribution ( −1 ; , −1 ), the matrix follows the IW distribution: In equation (8), is a symmetric positive definite random matrix, distribution parameter is a dof parameter, is a symmetric positive definite matrix; d is the dimension of , (. ) represents a multivariate gamma function,  [29]. Since is the covariance matrix of the Gaussian PDF, the prior distribution of ( | 1: −1 ) can be written as IW distribution: ( | 1: −1 ) = ( : −1 ;̂: −1 ,̂: −1 ), (9) wherê: −1 and ̂: −1 denote the dof parameter and scale matrix of ( | −1 ), respectively. Next, the value of ̂: −1 and ̂: −1 need to be assigned.
Owing to the Bayesian theorem, the prior distribution ( | −1 ) can be written as: where ( −1 | 1: −1 )) is the posterior probability density function (PDF) of the MNCM −1 . Utilizing (9), the distribution of posterior PDF ( −1 | 1: −2 ) can be replaced as inverse Wishart distribution, due to the prior distribution ( −1 | 1: −1 ) of MNCM −1 is selected as inverse Wishart distribution, ( −1 | 1: −1 ) can be written as: To guarantee ( | 1: −1 ) also obeys the inverse Wishart distribution, we introduce a changing factor to modify the one-step predicted values of the distribution parameters ̂: −1 and ̂: −1 . The formulas are as follows: among them, m is the dimension of , ∈ (0,1], the time-varying measurement noise covariance matrix can be changed with a certain probability distribution, and control the posterior and prior probability density functions to have the same distribution. In addition, the initial PDF distribution of MNCM 0 is also assumed as inverse Wishart distribution. ( 0 ) = ( 0 ;̂0 :0 ,̂0 :0 ).At the initial moment, to formulate the prior information of the measurement noise covariance matrix, the mean value of 0 is considered as the initial fixed measurement noise covariance matrix ̃0 , i.e., Assuming that the prior distribution of the joint probability density function of the state variable and the MNCM is the product of the Gaussian distribution and the inverse Wishart distribution, the prediction process can be defined as: ( , | 1: −1 ) = ( | 1: −1 ) ( | −1 ) = ( −1 ;̂: −1 , : −1 ) ( −1 ;̂− 1: −1 ,̂− 1: −1 ) (15)

Variational update process
Aiming at estimating the state and the MNCM , their joint posterior PDF ( , | 1: ) need to be calculated. However, the analytical solution of this joint posterior PDF cannot be obtained directly. The variational Bayesian method is utilized to find an approximate PDF of a free form as follows [30]: where (. ) means the approximate posterior PDF of (. ) .Minimizing the Kullback-Leibler Divergence (KLD) between the approximation posterior PDF ( ) , ( ) and the true joint posterior ( , | 1: −1 ) is used to form the VB-approximation [30]: The divergence function KLD(.) is defined as: Preprints (www.preprints.org) | NOT PEER-REVIEWED | Posted: 7 December 2020 doi:10.20944/preprints202012.0164.v1 Combine with (17) and (18), the optimal solution of (16) is derived as: log ( ) = [log ( , , 1: )] + , (20) where log(. ) stands for natural logarithm calculation, [. ] denotes the expectation calculation of the approximate posterior PDF of the variable , and represent the constants of variable , and MNCM , respectively. The solutions of equations (19) and (20) can not be solved directly since ( ) and ( ) are coupled. Therefore, the fixed-point iteration method is introduced to calculate the solution to these parameters.
The further form of equation (20) can be derived as (See Appendix A for details): where: ( ) (. ) represents the approximate probability distribution of (. ) at the -th iteration, (. ) is the calculation of the matrix trace, is a constant related to which independent of the distribution form; is the dimension of the real observation matrix. Define the expectation part of equation (21) as ( ) , and expand it as: it can be seen that ( +1) ( ) obeys a new inverse Wishart distribution form as follows, and the distribution parameters ̂( +1) and ̂( +1) are respectively as follows: Similarly, the logarithmic expression of the approximate distribution of the system state is as follows: where ( +1) [ −1 ] −1 is given by: The likelihood PDF ( | ) in equation (4) after updating the ( + 1)-th iteration can be derived as follows: The corrected measurement noise covariance matrix (MNCM) ̂( +1) can be written as: Since ( ) obeys the Gaussian distribution as ( ) = ( ;̂, ). Combining with the standard Kalman filter framework, the gain matrix ( +1) , system state ̂( +1) , and state covariance ( +1) in the variational measurement update are corrected as follows respectively: Analyzing the above derivation, it can be seen that the implicit solution of the variational update formula is constituted by the equations (22), (24), (25) and (29)-(32). The expected maximum approach is used to iteratively calculate ( ) and ( ) to update the parameters and to be estimated continuously. When ( ) and ( ) are closer to ( , −1 | 1: ), the KL divergence value of formula (17) is smaller, and the estimation results of the parameter to be estimated adaptively approaching to the true value until the iteration of the variational update process finished . At this time, the optimal estimation results of parameters ̂ and ̂ to be estimated at time can be calculated as follows ( is the number of fixed-point iterations):

Update monitoring strategy based on maximum a posterior (MAP) for estimating PNCM
Some existing adaptive filtering methods estimate process noise covariance matrix (PNCM) and measurement noise covariance matrix (MNCM) at the same time, it is easy to cause the accuracy of the estimated value of the state to decrease or even diverge. This is caused by the value of becoming negative definite matrix during the estimation process [17]. Aiming at realizing the simultaneous estimation of PNCM , MNCM and improving the estimation accuracy of the state vector . An Update monitoring strategy based on maximum a posterior (MAP) for estimating PNCM is proposed. According to the state-space model as equation (1) and (2), in paper [15], the maximum a posteriori suboptimal unbiased estimation method based on the noise statistics of measurement { 1 , 2 , 3 , ⋯ } for estimating PNCM is given as: further deduced as follows: where is the attenuation factor. In ∑ [. ] =0 of equation (35), each item is multiplied by the weight coefficient , instead of the original weight coefficient. The time-varying process noise covariance matrix (PNCM) estimation method is obtained, and the recursive algorithm is derived as: Equation (7), (30)-(32) and (38) constitutes the VBAKF algorithm that simultaneously estimates PNCM , MNCM , state vector and one-step predicted state error covariance matrix (PECM) : −1 . However, through simulation experiments, we found that the estimation of the PNCM by the above algorithm is prone to abnormality; that is, it loses semi-definiteness, which leads to filtering divergence.
To solve this problem, based on not losing the information of the original process noise estimation algorithm (38), an update monitoring strategy of process noise parameters is designed. Firstly, it is judged whether calculated by equation (38) is a positive semi-definite matrix. Then the adjustment factor is introduced to update the process noise estimation parameters to ensure that the corrected PNCM meets the requirements.
The right side of the equation (38) is shifted as follows: Generally speaking, the selection of the initial value of the state error covariance matrix is randomly, and the deviation from the ideal value is large in the initial stage of filtering, resulting in the theoretical process noise covariance matrix ̂ determined by equation (39) being much larger than . Therefore, it is necessary to introduce an adjustment factor to attenuate the effect of the state error covariance matrix at the initial moment, to avoid the indefiniteness of the estimated value of the PNCM, to prevent the filter from diverging. The update monitoring strategy of the PNCM is as follows: where ≥ 1 is a positive integer (the initial value is 1), is the adjustment factor, and the value of is related to the state error variance matrix as follows: The specific process for update monitoring strategy is as follows: monitoring the process noise covariance matrix calculated by equation (38), and judge whether ̂ is a positive semi-definite matrix to determine whether ̂ needs to be updated. If not, output ̂. Otherwise, turn to equation (41), set = 1, equation (41) is used to update the process noise parameters. Continue the monitoring of the updated estimated process noise covariance matrix to determine whether it is necessary to continue to update. If it is necessary, take = + 1, equation (41) is used to recalculate ̂. The loop is executed until ̂ is a positive semi-definite matrix. End the update of the process noise covariance matrix at the current moment. The flowchart of One-time step of the update monitoring strategy is shown in Figure (1).
So far, combining with traditional Kalman filter framework, the VBAKF algorithm with update monitoring strategy is derived to estimate system state vector , process noise covariance matrix (PNCM) , measurement noise covariance matrix (MNCM) , one-step predicted state error covariance matrix (PECM) : −1 , and state error covariance matrix at the same time.

Improved by introducing multiple fading factors
If the statistical characteristics of the process and measurement noise are time-varying, the convergence speed of VBAKF will slow down, and there will be a particular error in the estimation result, which will be reflected by the residual sequence [31]. In view of this, to improve the accuracy of estimation, the multiple fading factor is introduced to realize the correction of the one-step predicted error covariance matrix (PECM) : −1 . equation (7) can be rewritten as: adjust the gain ( ) in real-time to keep orthogonal, forcing the filter to keep track of the actual state of the system. Thereby the tracking ability is improved. The calculation method of multiple fading factors is as follows: where is the dimension of the state vector , 0 = • , the value of is determined by the system prior information. The formula of is as follows: where is the -th element of the main diagonal of , the calculation formulas of and are as follows: In equation (48), is the weakening factor, is unknown and can be estimated by the following formula: where ∈ (0,1] is the forgetting factor. Since is used to correct the one-step predicted error covariance matrix (PECM) in the prediction step of the filtering algorithm, the initial value ̂1 of the MNCM ̂ must be set in advance, assuming that ̂1 also obeys the inverse Wishart distribution. The estimated value of measurement noise covariance matrix ̂ in time update step is defined as: In equation (51), the variation range of the slowly time-varying measurement noise covariance matrix is small, and the estimated value ̂− 1 at the previous time still has a great reference value for the current time estimation. Therefore, ̂− 1 estimated by the variational update recursively of the previous time is used as ̂ at time > 1. ̂ is used as the time-varying parameter of to modify the one-step predicted error covariance matrix (PECM)

4.Simulations and results
The The simulation environment is set as follows: = 1000 s is the total simulation time, is a parameter related to process noise, is a parameter related to measurement noise. The fixed PNCM and MNCM are set as ̃= 4 and ̃= 2 , respectively, where and are the prior confidence parameters used to adjust the initial fixed noise covariance matrix.
This paper compares MFMS-VBAKF and true noise covariance matrix Kalman filter (TCMKF), fixed noise covariance matrix Kalman filter (FCMKF), ST-KF, ML-KF, SH-KF and VBKF-R algorithms. Table 1 lists the estimated parameters and parameter settings of the existing algorithms. All algorithms are programmed using MATLAB R2018a, and the simulation program runs on a computer with Intel® Core™ i5-6300HQ CPU at 2.30 GHz and 8GB of RAM.
Aim to evaluate the accuracy of system state estimation, the root mean square error and average root mean square error of position and velocity are regarded as performance indicators, which are defined as follows:  Table 2 lists the average root mean square error (ARMSE) of different KF filtering algorithms: according to the data in Table 2, it can be found that in comparison with other adaptive Kalman filter algorithms, MFMS-VBAKF algorithm has the smallest ARMSE and the highest accuracy in estimating target position and velocity.
To evaluate the accuracy of the estimation of one-step predicted state error covariance matrix PECM and the noise covariance matrices PNCM, MNCM, the square root of the normalized Frobenius norm (SRNFN) and the averaged SRNFN (ASRNFN) are used as the measure of error, which are defined as: where ̂, : −1 and ̂: −1 represents the true value and estimated value of the noise covariance matrix or one-step predicted state error covariance matrix in the -th Monte Carlo experiment, respectively. The SRNFN and ASRNFN of the estimation result of PECM are shown in Figure 3 and Table 3, respectively.  It can be clearly seen that, compared with the existing adaptive KF algorithm, if the noise covariance matrices are slowly time-varying, the SRNFN of the MFMS-VBAKF algorithm is smaller than the SRNFN of the current algorithm. Compared with R-VBKF with similar performance, the ASRNFN of MFMS-VBAKF is reduced by 5.45%.     Figure 5 and Table 4 show respectively the SRNFNs and the ASRNFNs of PNCM from the existing filters and MFMS-VBAKF algorithm. It can be seen that the SRNFN and ASRNFN of the proposed MFMS-VBAKF are both smaller than the current filters. Thus, the MFMS-VBAKF has better estimation accuracy and satisfactory convergence speed in PNCM estimation.
Next, we compare and analyze the influence of the values of two critical parameters (changing factor and forgetting factor ) in the MFMS-VBAKF algorithm on the estimated effect.   when = 0.95, the MFMS-VBAKF algorithm has the best estimation accuracy and convergence performance.
For the sake of testing the robustness of the adaptive correction capability of the MFMS-VBAKF algorithm when the fixed noise covariance matrices are set to different initial values, the priori confidence parameters and are set to change in combination within the grid area of ( , ) ∈ [0.1, 800] × [0.1, 800]. The ARMSEs estimated by the algorithm for position and velocity are displayed in Figure 7.  It can be analyzed from Figure 7 that the ARMSEs of position and velocity estimation are flat in a large area of the set grid, and the estimation results are close to the actual values. However, the initial setting values of the fixed noise covariance matrices in the extremely narrow area on the right edge of the grid are too different from the actual values, which lead to unsatisfactory performance of the estimation results. This is caused by the variational Bayesian method that can only guarantee local convergence. In general, the estimated effects of the MFMS-VBAKF algorithm can converge to near the actual values, with excellent robust performance.

Conclusions
This paper presents a Multi-fading factor and update monitoring strategy adaptive Kalman filter based variational Bayesian with the inaccurate time-varying process and measurement noise covariance matrices. The model of measurement error is defined as the Inverse Wishart distribution, and the variational Bayesian method is used to recursively estimate the measurement noise covariance matrix and the system state, the estimation results of the two can be approximated to the true value. The process noise covariance matrix is estimated by the maximum a posteriori principle, and the update monitoring strategy with adjustment factors is used to guarantee the positive semidefiniteness of the updated matrix. The estimated value of measurement noise covariance obtained by the variational iteration recursion and the estimated value of process noise covariance obtained by updating the monitoring strategy are used as time-varying parameters of multiple fading factors, which can be corrected to obtain more accurate state predicted error covariance. Variational Bayesian, update monitoring strategy and multi-fading factors complement each other, which not only enhances the responsiveness of target tracking, but also improves the estimation accuracy of variational iteration recursion. The simulation results show that the proposed MFMS-VBAKF algorithm realizes the simultaneous estimation of the process noise covariance matrix and the measurement noise covariance matrix, and has achieved satisfactory results in terms of estimation accuracy, convergence performance, and robustness.