Robust Interacting Multiple Model Filter Based on Student’s t-Distribution for Heavy-Tailed Measurement Noises

In maneuvering target tracking applications, the performance of the traditional interacting multiple model (IMM) filter deteriorates seriously under heavy-tailed measurement noises which are induced by outliers. A robust IMM filter utilizing Student’s t-distribution is proposed to handle the heavy-tailed measurement noises in this paper. The measurement noises are treated as Student’s t-distribution, whose degrees of freedom (dof) and scale matrix are assumed to be governed by gamma and inverse Wishart distributions, respectively. The mixing distributions of the target state, dof, and scale matrix are achieved through the interacting strategy of IMM filter. These mixing distributions are used for the initialization of time prediction. The posterior distributions of the target state, dof, and scale matrix conditioned on each mode are obtained by employing variational Bayesian approach. Then, the target state, dof, and scale matrix parameters are jointly estimated. A variational method is also given to estimate the mode probability. The unscented transform is utilized to solve the nonlinear estimation problem. Simulation results show that the proposed filter improves the estimation accuracy of target state and mode probability over existing filters under heavy-tailed measurement noises.


Introduction
The Kalman filter is widely used for target tracking due to the low computational complexity for real time processing and the statistic optimality under a linear state space model with Gaussian noises. Because of the model uncertainty of target maneuvering and the heavy-tailed measurement noises deduced by outliers, the performance of Kalman filter may break down dramatically, and even the convergence cannot be guaranteed. Addressing model uncertainty and investigating robust filter against heavy-tailed measurement noises are two important research topics in the target tracking field.
Model uncertainty in maneuvering target tracking exhibits that the motion model of target is time varying and that the switching time between different models is unknown. Due to the model uncertainty, the model used by Kalman filter may not match the actual dynamics so that the filter may gradually diverge. Multiple-model (MM) methods are considered as the mainstream approach to maneuvering target tracking under model uncertainty [1]. These methods establish a finite set of models to describe the different motion patterns for maneuvering target, and achieve the ultimate estimate of state by a certain combination of the estimates from each model. In the literature, the MM methods can be categorized into two different types: fixed-structure MM (FSMM) methods which use a fixed model set at all times and variable-structure MM (VSMM) methods in which the model set is not fixed. The generalized pseudo-Bayesian estimator of order n (GPBn) filter and the interacting multiple model (IMM) filter [2] are the typical FSMM methods. Among them, the IMM filter is considered to be the best compromise between complexity and performance, and has been successfully applied in a large number of tracking applications [3][4][5][6][7][8]. In IMM filter, several sub-filters operate in parallel and cooperate with each other through an interacting strategy, leading to improved performance of estimation. For the FSMM methods, if the model set used does not match the set of actual models of target movement, the performance can be degraded. To cope with this problem, several VSMM methods have been proposed, such as likely mode set (LMS) filter [9], expected-mode augmentation (EMA) filter [10], equivalent-model augmentation (EqMA) filter [11], and hybrid grid multiple model (HGMM) filter [12]. These filters use different model-set adaptation schemes to adjust the model set. For example, the EqMA filter augments the model set by a new variable model depending directly on the history information from the recent models, while the HGMM filter adapts the model set by the hybrid grid scheme which consists of a fixed coarse grid and an adaptive fine one. In recent years, some researchers have proposed Sequential Monte Carlo (SMC) methods [13,14] for dealing with the problem of model uncertainty in target tracking. These methods maintain several particle filters (PFs) in parallel and fuse the results of PFs to provide the posterior probability of model and the posterior distribution of state. In addition, the SMC method in [14] uses a scale parameter to detect the model change and does not require the mode transition probabilities. This strategy of handling model change is quite different from the above mentioned FSMM and VSMM methods.
Besides the model uncertainty, the outlier measurements existing in real tracking scenarios should also be handled. Outliers can induce heavy-tailed non-Gaussian noises, and the performance of Kalman filter may be seriously degraded since its underlying optimality criterion is minimum mean square error, which is very sensitive to heavy-tailed noises. A large number of robust filters have been proposed to improve the robustness against non-Gaussian noises. The Huber-based filter [15][16][17][18][19] is a classical robust method based on a combined minimum l 1 and l 2 norm estimation. The maximum correntropy Kalman filter (MCKF) [20][21][22][23][24] is another effective way to suppress the impact of non-Gaussian noises. It uses the robust maximum correntropy criterion (MCC) as the optimality criterion. Simulation experiments show that the MCKF outperforms Huber-based filter in estimation accuracy [20,[22][23][24]. However, neither Huber-based filter nor MCKF takes advantage of the inherent heavy-tailed feature of the noise distribution. Thus, their estimation precision is limited. The Student's t-distribution is robust to outliers since it has heavier tail compared with Gaussian distribution. Based on this idea, many robust filters utilizing Student's t-distribution has been proposed [25][26][27][28][29][30]. The Student's t filter [25,26] is derived by modeling the posterior distribution of the state, the heavy-tailed process and measurement noises as Student's t-distributions. The filter recursions are designed according to the framework of Bayesian filter. This filter may induce significant bias since an approximated method of moment matching is employed to prevent the growth of the degrees of freedom (dof) of Student's t-distribution. The variational Bayesian (VB) based Student's t filter [27][28][29][30] represents the Student's t-distribution of noises as an infinite mixture of Gaussians, and employs VB approach to jointly estimate the state and the unknown parameters of the Student's t-distribution. This filter can achieve higher estimation accuracy than the Student's t filter [29,30].
In maneuvering target tracking, the problems of model uncertainty and heavy-tailed measurement noises could exist simultaneously. Investigating robust IMM filter to cope with both model uncertainty and heavy-tailed noises is meaningful. As far as we know, Shen et al. [31] considered utilizing Student's t-distribution to improve the robustness of IMM filter for the first time. They modeled the heavy-tailed measurement noises as Student's t-distribution and employed IMM and VB approaches to estimate the target state, the probability of motion mode and the parameters of noises. However, the dof of Student's t-distribution, which is determined by the estimated parameters of an auxiliary random variable, becomes very large after several filtering steps. Then, the Student's t-distribution converges to Gaussian distribution and loses the heavy-tailed property. Therefore, the robustness of this filter against heavy-tailed measurement noises is limited. In addition, the scale matrix of Student's t-distribution Sensors 2019, 19, 4830 3 of 20 in this filter is assumed to be known exactly. However, this premise is usually unavailable since the scale matrix is unknown or time varying in many actual applications. Furthermore, this filter is based on linear systems and should be generalized to nonlinear systems which are more common in maneuvering target tracking.
Motivated by the above discussions, we propose a new robust IMM filter to better handle the heavy-tailed measurement noises and the system nonlinearity. The posterior distribution of target state is approximated by Gaussian distribution, and the measurement noises are modeled as Student's t-distribution. The unknown dof and scale matrix of the Student's t-distribution are assumed to obey gamma and inverse Wishart distributions, respectively. Then, the state and the parameters of gamma and inverse Wishart distributions for each mode are interacted via IMM mixing, and they are jointly estimated by using VB approach. An approximation method is given to derive the mode-conditioned predicted likelihood which is used for calculating the mode probability. The unscented transform (UT) is employed to tackle the system nonlinearity. The growth of dof estimates is prevented and the estimation accuracy of state is improved over existing filters under heavy-tailed measurement noises as shown in our simulation example.

System Model and Assumptions
Consider the following jump Markov system where k is the scan time index, x k ∈ R n and z k ∈ R m are the state and measurement vectors respectively, r k is the mode of target movement taking value in a finite set {1, 2, · · · , M}, M is the total number of dynamic models, F k−1 (x k−1 , r k ) and H k (x k , r k ) are the state transition function and measurement function based on the mode r k respectively, w k−1 (r k ) is the process noise vector based on the mode r k , and ε k is the measurement noise vector. The mode r k is assumed to be a homogeneous Markov chain with the transition probability π ij = p(r k = j|r k−1 = i), and w k−1 (r k ) is assumed to be uncorrelated zero-mean Gaussian noises with covariance matrix Q k−1 (r k ). Throughout this paper, we abbreviate respectively. Then, the transition probability density function (PDF) of state is given by , where N(·; m, P) denotes a Gaussian PDF with mean m and covariance matrix P.
When measurement outliers occur, the distribution of measurement noises has heavy-tailed non-Gaussian characteristic. The performance of traditional IMM filter based on Gaussian noise assumption may dramatically degrade since it is very sensitive to heavy-tailed noises. Since the Student's t-distribution has heavier tail and is more robust to outliers than Gaussian distribution, we model the measurement noises as Student's t-distribution, i.e., ε k ∼ Std(ε k ; 0, R k , v k ), where Std(·; µ, R, v) stands for a Student's t PDF with mean µ, scale matrix R and dof v. The Student's t-distribution can be represented by an infinite mixture of Gaussian distributions as [32] Std(x; µ, where λ is an auxiliary random variable, G(·; a, b) denotes a Gamma PDF with shape parameter a and rate parameter b. The expression for G(·; a, b) is where Γ(a) = +∞ 0 exp(−t)t a−1 dt is the Gamma function. Then, the likelihood PDF conditioned on the mode r k = i can be written as The tail thickness of Student's t-distribution is controlled by the dof. When the dof decreases, the tail becomes thicker. When the dof goes to infinity, the Student's t-distribution approaches a Gaussian distribution. Shen et al. [31] assumed that the posterior distribution of the auxiliary variable λ k is a Gamma distribution, i.e., p(λ k |r k = i, z 1:k ) = G(λ k ;α i k ,β i k ) and used the VB approach to estimate the parametersα i k andβ i k . However, the estimates of bothα i k andβ i k which determine the dof of Student's t-distribution become large after several filtering steps. Thus, the distribution of measurement noises loses the heavy-tailed property and the robustness of the filter cannot be guaranteed. To prevent the growth of the dof estimates, we assume that the posterior distribution of the dof v k is a Gamma distribution, i.e., p(v k |r k = i, z 1:k ) = G(v k ;â i k ,b i k ), and estimate the parametersâ i k andb i k instead. Furthermore, in order to deal with the uncertainty of scale matrix R k , we choose inverse Wishart distribution as the posterior distribution of R k , i.e., p(R k |r k = i, z 1:k ) = IW(R k ;t i k , T i k ), where IW(·; t, T) stands for an inverse Wishart PDF with degree t and inverse scale matrix T. The PDF of inverse Wishart distribution is defined by where Γ m (0.5t) = π 0.25m(m−1) m−1 i=0 Γ(0.5(t − i)) Then, the unknown statistic of the scale matrix R k can be determined by estimating the parameterst i k and T i k .

Model Interaction
We assume that the mode-conditioned posterior distributions of the state x k−1 , the scale matrix R k−1 and the dof v k−1 at time k − 1 are Gaussian, inverse Wishart, and Gamma distributions, respectively, and that they are mutually independent. Then, the joint mode-conditioned posterior PDF of x k−1 , R k−1 and v k−1 can be expressed as Let u i k = p(r k = i|z 1:k ) and u i k|k−1 = p(r k = i|z 1:k−1 ) denote the posterior and predicted mode probabilities respectively, and u ij k−1 = p(r k−1 = i|r k = j, z 1:k−1 ) denote the mixing probability, we have The mixing PDF of x k−1 , R k−1 and v k−1 is given by We approximate the above summed PDF in the right side of (6) by a single one as By matching the first and the second-order moments between (6) and (7) for the state distribution, we obtain the mixing statex 0i k−1 and the covariance matrix P 0i k−1 aŝ However, the mixing parameterst 0i k−1 and T 0i k−1 cannot be computed analytically by matching the first two moments of inverse Wishart distributions between (6) and (7) due to the complexity of second-order moment. We adopt the method of minimizing weighted Kullback-Leibler (KL) divergence in [33] to overcome this difficulty and derive thatt 0i For the Gamma distribution, we also use the method of matching the first two moments like the state. According to (6), the conditional mean and covariance of v k−1 are given by Since the mean and covariance of a Gamma PDF G(·; a, b) are a/b and a/b 2 respectively, we have On the other hand, according to (7), we have By combining (12)-(15), the mixing parametersâ 0i k−1 andb 0i k−1 for the Gamma distribution can be calculated asâ

Time Prediction
The time prediction step is to derive the mode-conditioned predicted PDF p(x k , R k , v k |r k = i, z 1:k−1 ) at time k by using the mixing PDF p(x k−1 , R k−1 , v k−1 |r k = i, z 1:k−1 ) and the dynamical models at time k − 1.
The mode-conditioned predicted PDF can be calculated by the following Chapman-Kolmogorov equation: We assume that the dynamical models of the state, the scale matrix and the dof are independent, i.e., According to (7), (18) and (19), we have We can see that the mode-conditioned predicted distributions of the state, the scale matrix and the dof are also independent. These distributions are assumed to be Gaussian, inverse Wishart and Gamma distributions, respectively, which have the same form of PDF as the posterior distributions. Then, the mode-conditioned predicted PDF can be represented by By matching the first and the second-order moments between (20) and (21) for the state, we obtain the predicted mean and the covariance matrix aŝ However, it is hard to compute the predicted parameters of inverse Wishart and Gamma distributions since the dynamical models of the scale matrix and the dof are usually unknown in practice. Here, a heuristic dynamical model is chosen as in [34] by introducing a forgetting factor ρ ∈ (0, 1], which indicates the extent of parameter fluctuation. Then, the predicted parameters are given by When ρ = 1, the parameters are kept stationary. When ρ is close to 0, the parameters are high time-fluctuation.

Measurement Update
The measurement update step is to compute the mode-conditioned posterior PDF p(Θ|z 1:k , r k = i) and the posterior mode probability u i k when a new measurement z k is collected at time k. However, deriving the mode-conditioned posterior PDF is analytically intractable. We employ the VB approach [35] to obtain an approximated solution.
According to the VB approach, the mode-conditioned posterior PDF is approximated by a free form factored PDF: The posterior PDFs q i (x k ), q i (R k ), q i (λ k ) and q i (v k ) are calculated by minimizing the KL divergence as where KL(q(·)||p(·)) q(x) ln(q(x)/p(x))dx is the KL divergence. The solutions of the above optimal problem satisfy ln where c x , c R , c λ and c v are constants independent of x k , R k , λ k and v k , respectively. According to (2) and (21), we have Substituting (1), (3) and the expression of Gaussian PDF into (32), we obtain the detailed expression of p(Θ, z 1: where c denotes a constant with respect to Θ. Substituting (33) into (28) yields where R k is given by where E q i (·) denotes the expectation with respect to the PDF q i (·). R k can be deemed as modified covariance matrix of measurements taking into account the uncertainty of both and R k and λ k . Then, According to the construction of nonlinear Gaussian filter in [36], q i (x k ) can be approximated by a Gaussian PDF, i.e., According to (3), we can deduce that q i (R k ) is an inverse Wishart PDF with the following expression: Substituting (33) into (30) yields According to (1), q i (λ k ) is a Gamma PDF as Substituting (33) into (31), and using the Stirling's approximation ln The unknown expectations in (34), (41), (44), (46), (47) and (50) are given by The posterior mode probability can be calculated by In (57), the predicted mode probability u i k|k−1 is given by (4), and the predicted likelihood p(z k |r k = i, z 1:k−1 ) can be derived by However, the above integration is computationally infeasible since x k , R k , λ k and v k are coupled in p(Θ, z k |r k = i, z 1:k−1 ) according to (32). We use the similar variational approach in [37] to derive an approximated predicted likelihood. The predicted likelihood can be rewritten as Since the KL divergence term in the right side of the above equation is minimized by VB approach, the predicted likelihood can be approximated by Substituting (35) The final estimate of the state is given by a probabilistically weighted average of all the mode-conditioned estimates, i.e.,x and the corresponding covariance matrix is calculated as

Approximated Gaussian Integrations Based on Unscented Transform
Due to the nonlinearity of the state transition function F i k−1 (x k−1 ) and the measurement function H i k (x k ), the Gaussian integrations in (22), (23), (37)-(39), (53) and (54) cannot be computed analytically. In this paper, the UT [38] is employed to calculate the Gaussian integrations approximately.
For the Gaussian integrations in (22) and (23), 2n + 1 sigma points are generated from the mixing statex 0i k−1 and the covariance matrix P 0i k−1 as , where ( √ P) p denotes the p-th column of the matrix square root of P, λ = α 2 (n + κ) − n is a scaling parameter, α controls the divergence of the sigma points and is usually set to a small positive value (e.g., 0.01), κ is a secondary scaling parameter which is usually set to zero. Then, the predicted statê x i k|k−1 and the covariance matrix P i k|k−1 are calculated aŝ where the weights of sigma points are given by , where β is used to incorporate the prior information of the distribution. β = 2 is optimal under Gaussian distributions. For the Gaussian integrations in (37)-(39), the sigma points generated from the predicted statê x i k|k−1 and the covariance matrix P i k|k−1 are given by Then, the predicted measurement meanẑ i k|k−1 , the covariance matrix P i zz,k|k−1 and the cross-covariance matrix P i xz,k|k−1 are calculated asẑ We generate the sigma points from the posterior statex i k and the covariance matrix P i k as Then, the Gaussian integrations in (53) and (54) are calculated aŝ With the above derivations, the proposed robust IMM filter is summarized as follows: Step 1: Choose initial estimatesx i 0 , P i 0 , u i 0 ,â i 0 ,b i 0 ,t i 0 and T i 0 for each mode i Set π ij and ρ. Let k = 1 Step 2: Calculate the predicted mode probability u j k|k−1 and the mixing mode probability u ij k−1 using (4) and (5). Then, calculate the mixing quantitiesx 0i k−1 , P 0i k−1 ,t 0i k−1 , T 0i k−1 ,â 0i k−1 andb 0i k−1 for each mode i using (8)-(11), (16) and (17).

Until converged
Step 5: Update the posterior mode probability u i k using (57) and (59).

Simulation Example
This section presents a two-dimensional maneuvering target tracking scenario with a period of 200 s to demonstrate the performance of the proposed filter. A maneuvering target moves following two models: constant velocity (CV) model and coordinated turn (CT) model. The modes representing the CV model and CT model at time k are denoted by r k = 1 and r k = 2, respectively. The state of the CV model is x k = (x k , y k , . x k , . y k ) T including the components of the position (x k , y k ) and the velocity ( . x k , . y k ). The dynamics of the CV model is given by where ∆t is the sampling period. The covariance matrix of the process noises w 1 k−1 is given by where q 1 x is the power spectral density. The state of the CT model is x k = (x k , y k , . x k , . y k , ω k ) T , where ω k is the turn rate. The dynamics of the CT model is given by and the covariance matrix of the process noises w 2 k−1 is given by where q 2 x and q ω are the power spectral densities corresponding to (x k , y k , . x k , . y k ) and ω k , respectively. A sensor located at (x s , y s ) collects noisy range and azimuth measurements of the target according to the equation To illustrate the performance of the proposed filter, we consider the following four cases under different conditions of measurement noises: Case A: Gaussian noises This case is used to simulate the situation that the actual covariance matrix is deviated from the nominal one R . Case C: Contaminated Gaussian noises (1 )N( ; , ) N( ; ,100 ) can be seen as a normal distribution, N( ; ,100 ) can be considered as a perturbing distribution with much larger covariance matrix due to outliers, and [0,1] δ ∈ is a perturbing parameter that represents the extent of the perturbation. This case is used to simulate the corrupted measurements by outliers. We set 0.1 δ = .  This case is used to simulate the situation that the actual covariance matrix is deviated from the nominal one R.
Both outliers and covariance mistuning are simultaneously simulated in this case. δ is also set to 0.1. Four existing filters including VB based Student's t filter (VBStdF) [30] utilizing only single CT model, IMM filter (IMMF) [1], VB based IMM filter (IMMVBF) estimating the unknown covariance matrix of measurement noises [33], and VB based IMM filter modeling measurement noises as Student's t-distribution (IMMVBStdF) [31] are compared with the proposed filter. 1000 Monte Carlo (MC) runs are carried out for each case of measurement noises above. All the filters are implemented in MATLAB on an Intel i7 3.6GHz processor.
The initial estimates of the degree and the inverse scale matrix with respect to the inverse Wishart distribution in VBStdF, IMMVBF and the proposed filter are chosen ast 0 = 7 and T 0 = (t 0 − m − 1)R. All the initial estimates of the Gamma distribution parameters in VBStdF, IMMVBStdF and the proposed filter are set to 0.5. For all the VB approach based filters, ρ = 1 − exp(−4) is used, and the VB iteration is terminated when the difference of position estimates between two adjacent iterations is less than 1e-6. The mode probabilities are initialized as u 1 0 = 0.1 and u 2 0 = 0.9, and the mode transition probability is set to π ij = 0.99 0.01 0.01 0.99 for all the IMM type filters. The initial covariance matrix of the state estimate is set to P 0 = diag(100 m 2 , 100 m 2 , 1 ( m/s) 2 , 1 ( m/s) 2 , 1 × 10 −5 (rad/s) 2 ), and the initial estimate of the statex 0 is chosen randomly from the Gaussian distribution N(·; x 0 , P 0 ) for all the filters. The root mean square errors (RMSEs) and average RMSEs (ARMSEs) of position, velocity and turn rate are used to evaluate the estimation accuracy of the filters. The RMSE and ARMSE of position are defined as where P is the total count of MC runs, (x                  In Case A, the estimation accuracy of the proposed filter is lower than IMMF and IMMVBF. The reason is that the assumed Student's t-distribution for measurement noises deviates from the actual Gaussian distribution. The performance of the proposed filter in estimation accuracy is also poorer than IMMVBStdF. That is because the scale matrix of Student's t-distribution is required to be estimated additionally for the proposed filter, while the scale matrix is known exactly for IMMVBStdF. However, the accuracy degradation of the proposed filter is not obvious compared with IMMF, IMMVBF, and IMMVBStdF.
In Case B, IMMF has larger RMSEs than IMMVBF, IMMVBStdF, and the proposed filter since IMMF utilizes mistuned covariance matrix of measurement noises, while the IMMVBF and the proposed filters can learn the covariance matrix adaptively. The IMMVBStdF seems to have the ability of capturing the unknown covariance matrix either. The performance of the proposed filter in estimation accuracy is also slightly worse than IMMVBF and IMMVBStdF as Case A.
In Case C, the proposed filter outperforms other filters obviously in estimation accuracy. IMMF and IMMVBF are based on the assumption of Gaussian noises so that they cannot cope with heavy-tailed measurement noises. Figure 6 shows the estimates of the dof parametersα 1 k ,β 1 k ,α 2 k andβ 2 k versus time for IMMVBStdF and the proposed filter in one MC run. We can see that the estimates of dof parameters for IMMVBStdF become large after several time steps. Therefore, the distribution of measurement noises converges to a Gaussian distribution and loses the heavy-tailed property. The estimates of dof parameters for the proposed filter maintain small values overall, except some high peaks forβ 1 k andβ 2 k curves caused by outliers. Therefore, the proposed filter is more robust against heavy-tailed measurement noises than IMMVBStdF. In addition, we can see that the IMMVBF and IMMVBStdF are more robust than IMMF.
The results under Case D are similar to Case C. The 1000 MC runs averaged estimates of mode probabilities for CV model are shown in Figure 7. In Case A, we can see that the estimation results are nearly the same for all the filters. In Case B, the accuracy of IMMF estimates is lower than other filters and the response to mode changes for IMMF is more lagging than other filters. The reason is also that the covariance matrix of the measurement noises utilized by IMMF is mistuned. In Case C and D, the proposed filter achieves more accuracy estimates of mode probabilities than other filters, indicating that the proposed filter can also improve the estimation accuracy of mode probability under heavy-tailed measurement noises.  The 1000 MC runs averaged estimates of mode probabilities for CV model are shown in Figure  7. In Case A, we can see that the estimation results are nearly the same for all the filters. In Case B, the accuracy of IMMF estimates is lower than other filters and the response to mode changes for IMMF is more lagging than other filters. The reason is also that the covariance matrix of the measurement noises utilized by IMMF is mistuned. In Case C and D, the proposed filter achieves more accuracy estimates of mode probabilities than other filters, indicating that the proposed filter can also improve the estimation accuracy of mode probability under heavy-tailed measurement noises. We use the number of floating points operations (flops) [39] to measure the computational complexity of filter. The computational complexity of VBStdF, IMMF, IMMVBF, IMMVBStdF and the proposed filter are listed in Table 2. We can see that the proposed filter has the highest computational complexity among the filters. According to Table 2, the computational complexity of We use the number of floating points operations (flops) [39] to measure the computational complexity of filter. The computational complexity of VBStdF, IMMF, IMMVBF, IMMVBStdF and the proposed filter are listed in Table 2. We can see that the proposed filter has the highest computational complexity among the filters. According to Table 2, the computational complexity of proposed filter is higher than VBStdF mainly because of multiple model operations, and is higher than IMMF mainly because of fixed point iterations. The proposed filter requires more flops compared with IMMVBF and IMMVBStdF since the proposed filter involves more parameters to be estimated in the filter recursions. Table 2. Computational complexity of filters, where m denotes the dimension of measurements, n denotes the dimension of target state, M denotes the number of models for interacting multiple model type filters, and N denotes the number of fixed point iterations for variational Bayesian approach based filters.

Conclusions
We propose a robust IMM filter against heavy-tailed measurement noises for maneuvering target tracking. The heavy-tailed measurement noises are treated as Student's t-distribution, and the unknown dof and scale matrix are assumed to be governed by Gamma and inverse Wishart distributions, respectively. Then, the filter recursions for the target state and the parameters of Gamma and inverse Wishart distributions are designed in the IMM framework. In the model interaction step, the mixing distributions of the state, dof, and scale matrix are achieved by using the methods of matching the first two moments and minimizing weighted KL divergence. In the measurement update step, the state and the unknown parameters are jointly estimated by employing the VB approach. The problem of system nonlinearity is solved by utilizing UT to compute the Gaussian integrations approximately. Simulation results show that the proposed filter outperforms other related filters in terms of the estimation accuracy for the state and the mode probability under heavy-tailed measurement noises.
The main contributions of this paper are summarized as follows. Firstly, a robust IMM filter is proposed to properly handle both model uncertainty and outlier measurements. Secondly, by treating the unknown dof of Student's t-distribution of heavy-tailed measurement noises as a Gamma distributed random variable, the heavy-tailed property of Student's t-distribution is maintained and the robustness is improved over the existing filter. Finally, the problem of system nonlinearity is addressed by using UT technology.