## 3 Variational Bayesian Approximation-Based Adaptive Kalman Filter Algorithm  
As we know, when the measurement noise covariance is unknown or time varying, the classical approach of solving the problem is to use adaptive filters. The Bayesian approach is a typical choice. But it is usually hard to get the analytical form due to complex probability density function or high dimension of integration. In [27], variational Bayesian (VB) approximation was firstly used for the Kalman filter to estimate the joint posterior distribution of the states and noise covariances.  
Given a discrete-time linear state space model as:  
\[X_{k+1} = F_{k}X_{k}+w_{k} \tag{17}\] \[y_{k} = H_{k}X_{k}+v_{k} \tag{18}\]  
where \(w_{k}\sim\mathcal{N}(0,Q_{k})\) is a Gaussian distributed process noise, \(v_{k}\sim\mathcal{N}(0,\Sigma_{k})\) is a Gaussian measurement noise with diagonal covariance \(\Sigma_{k}\). Note \(Q_{k}\) is assumed known but \(\Sigma_{k}\) is unknown.  
Assuming the state vector and the measurement noise covariance are independent, the joint posterior distribution of the state and covariance \(p(X_{k},\Sigma_{k}|y_{1:k})\) can be solved by VB approximation as follows:  
\[p(X_{k},\Sigma_{k}|y_{1:k})\approx Q_{x}(X_{k})Q_{\Sigma}(\Sigma_{k}) \tag{19}\]  
The VB approximation can now be formed by minimizing the Kullback-Leibler (KL) divergence between the separable approximation and the true posterior:  
\[KL\left[Q_{x}(X_{k})Q_{\Sigma}(\Sigma_{k})||p(X_{k},\Sigma_{k}|y _{1:k})\right]\] \[= \int Q_{x}(X_{k})Q_{\Sigma}(\Sigma_{k})\times\log\left(\frac{Q_{ x}(X_{k})Q_{\Sigma}(\Sigma_{k})}{p(X_{k},\Sigma_{k}|y_{1:k})}\right)dX_{k}d \Sigma_{k} \tag{20}\]  
Minimizing the KL divergence with respect to the probability densities \(Q_{x}(X_{k})\) and \(Q_{\Sigma}(\Sigma_{k})\) in turn, while keeping the other fixed, we can get the following equations:  
\[Q_{x}(X_{k})\propto exp\left(\int\log p(y_{k},X_{k},\Sigma_{k}|y_{1:k-1})Q_{ \Sigma}(\Sigma_{k})d\Sigma_{k}\right) \tag{21}\]  
\[Q_{\Sigma}(\Sigma_{k})\propto exp\left(\int\log p(y_{k},X_{k},\Sigma_{k}|y_{1:k -1})Q_{x}(X_{k})dX_{k}\right) \tag{22}\]  
Computing the above equations, we can get the following densities [27]:  
\[Q_{x}(X_{k})=\mathcal{N}(X_{k}|m_{k},P_{k}) \tag{23}\]  
\[Q_{\Sigma}(\Sigma_{k})=\prod\text{Inv-Gamma}(\sigma_{k,i}^{2}|\alpha_{k,i}, \beta_{k,i}) \tag{24}\]where the parameters \(m_{k}\), \(P_{k}\), \(\alpha_{k,i}\), \(\beta_{k,i}\) can be calculated as:  
\[\begin{split} m_{k}&=m_{k}^{-}+P_{k}^{-}H_{k}^{T}(H_{k} P_{k}^{-}H_{k}^{T}+\hat{\Sigma}_{k})^{-1}(y_{k}-H_{k}m_{k}^{-})\\ P_{k}&=P_{k}^{-}-P_{k}^{-}H_{k}^{T}(H_{k}P_{k}^{-}H_{k}^{ T}+\hat{\Sigma}_{k})^{-1}H_{k}P_{k}^{-}\\ \alpha_{k,i}&=1/2+\alpha_{k-1,i}\\ \beta_{k,i}&=\beta_{k-1,i}+\frac{1}{2}[(y_{k}-H_{k} m_{k})_{i}^{2}+(H_{k}P_{k}H_{k}^{T})_{ii}]\end{split} \tag{25}\]  
where \(i=1,...,d\) and \(d\) denote the dimensionality of the measurement. \(m_{k}^{-}\) and \(P_{k}^{-}\) are the predicted state estimate and its covariance, respectively. Here we assume each component of the measurement noise variance is mutually independent, and then the covariance matrix is diagonal, estimated as:  
\[\hat{\Sigma}_{k}=diag(\beta_{k,1}/\alpha_{k,1},\ldots,\beta_{k,d}/\alpha_{k,d}) \tag{26}\]  
Furthermore, in order to describe the dynamics of the measurement noise variance, the inverse Gamma distribution parameters are assumed to change by a scale factor \(\rho\in(0,1]\). The formulas are given as:  
\[\begin{split}\alpha_{k,i}^{-}&=\rho_{i}\alpha_{k-1,i}\\ \beta_{k,i}^{-}&=\rho_{i}\beta_{k-1,i}\end{split} \tag{27}\]  
Note that the value \(\rho=1\) corresponds to stationary variances and lower values represent larger time-fluctuations.  
The above are the basic equations of adaptive Kalman filter based on VB approximation. Moreover, when the state equation and the measurement equation are nonlinear, the VB method can be rewritten in the EKF framework.