J. Cent. South Univ. Technol. (2009) 16: 0815-0820

DOI: 10.1007/s11771-009-0135-y

An improved particle filtering algorithm based on observation inversion optimal sampling

HU Zhen-tao(胡振涛), PAN Quan(潘 泉), YANG Feng(杨 峰), CHENG Yong-mei(程咏梅)

(College of Automation, Northwestern Polytechnical University, Xi’an 710072, China)

Abstract: According to the effective sampling of particles and the particles impoverishment caused by re-sampling in particle filter, an improved particle filtering algorithm based on observation inversion optimal sampling was proposed. Firstly, virtual observations were generated from the latest observation, and two sampling strategies were presented. Then, the previous time particles were sampled by utilizing the function inversion relationship between observation and system state. Finally, the current time particles were generated on the basis of the previous time particles and the system one-step state transition model. By the above method, sampling particles can make full use of the latest observation information and the priori modeling information, so that they further approximate the true state. The theoretical analysis and experimental results show that the new algorithm filtering accuracy and real-time outperform obviously the standard particle filter, the extended Kalman particle filter and the unscented particle filter.

Key words: particle filter; proposal distribution; re-sampling; observation inversion

1 Introduction

In recent years, Monte Carlo simulation technique has been a growing research field due to improving computer performances. The main breakthrough comes from particle filter (PF) proposed by GORDON et al [1-2]. Based on the concept of sequential importance sampling and the recursive Bayesian principle, particle filter is particularly useful in dealing with nonlinear and non-Gaussian problems [3-5]. The underlying principle of the methodology is to approximate the posterior probability density distribution with random measures composed of particles and their associated weights. It is known that the PF suffers from the degeneracy phenomenon, where all but one of the normalized importance weights is very close to zero after a few iterations [6-7]. The degeneracy means that a large computational effort is wasted to update trajectories whose contribution to the final estimate is almost zero.

In order to resolve particle degeneracy, the following two methods are applied: the re-sampling and the optimization of proposal distribution. While re-sampling process restrains the degeneracy phenomenon to some extent, it also introduces a new problem: the particle impoverishment [8-9]. A simple way for the problem is to increase sampling particle numbers. But it makes rapid expansion of computational complexity. Markov Chain Monte Carlo (MCMC) process is applied after re-sampling, which makes particles set tend to stationary distribution so as to decrease correlation among particles and improve particles diversity [9]. Nevertheless, the probability transition number is large in order to ensure a better convergence, which weakens real-time. Some experts and scholars attempt to introduce swarm intelligence technologies into PF algorithm, such as genetic algorithm and evolutionary strategy [10-11]. Particles diversity after re-sampling is improved by the particles crossover step or the particles mutation step. The precondition of better estimation result obtained is that some parameters are selected reasonably, which requires more prior information. The kernel function is introduced to approximate the posterior probability density distribution in the regularized particle filter (RPF), which can make the distribution of sampling particles much smoother because particles are sampled in continuous distribution. The disadvantage is that the regularized process is difficult for high dimensional system. In addition, particle set after re-sampling is not the asymptotially unbiased estimation for posterior probability density distribution [12]. The purpose of designing better proposal distribution is to improve the effectiveness of sampling particles by means of the effective utilization of the latest observation information, which is implemented by combining PF with some suboptimal filters (EKF, IEKF, UKF and GHF). Some existing algorithms, for example, EPF, IEPF, UPF and GHPF [13-16], are their implementations. Nevertheless, EKF and IEKF can diverge in strong nonlinearity cases, while UKF and GHF are unsuitable to non-Gaussian noise model and have more computational complexity. These disadvantages are also introduced into the improved algorithms. In auxiliary particle filter (APF), these most promising particles whose prediction likelihood degrees are large at the previous time, are augmented to the current time based on particle weights after re-sampling. This realizes the improvement of effectiveness of sampling particles to some degree. When process noise is small, APF is superior to PF. Conversely, the advantage is guaranteed difficultly [17].

In this work, a novel proposal distribution optimizing design principle is presented firstly, and then an improved particle filtering algorithm based on observation inversion optimal sampling (IPF) is given.

2 Particle filter

Consider the following nonlinear state space model:

x_k+1=f(x_k)+w_k                                 (1)

z_k=h(x_k)+v_k                                  (2)

where  x_k and z_k are the state variable and observation at time k, respectively; f and h are known nonlinear evolution functions for state in system model and observation model, respectively; w_k and v_k are system noise and observation noise sequences, respectively, which independently meet identical distribution. The problem of state estimation can be solved by calculating the posterior probability density function p(x_k|z_1:k) of state variable x_k at time k based on all the available data of observation sequence z_1:k={z₁, z₂, …, z_k}. Because the complete information of sequential estimation is in p(x_k|z_1:k), some parameters of which system state estimation need can be obtained, such as mean and variance. The main idea is to approximate p(x_k|z_1:k) with sampling particles in PF, whose mathematical description can be written as

                   (3)

where  δ(?) is Dirac’s delta function (δ(x)=1 for x=0 and δ(x)=0 otherwise), and { i=1, 2, …, N, N→∞} represents sampling particle, that is sampled directly from p(x_k|z_1:k). Unfortunately, p(x_k|z_1:k) is unknown generally, and the above sampling process is often impossible to implement. The difficulty can be circumvented by sampling particles   i=1, 2,

…, N} with associated importance weights, =1,

from a known and easy-to-sample proposal distribution q(x_k|z_1:k). The above process is described as the importance sampling. Here, the associated importance weight of particle is defined by

∝                    (4)

To further depict the generation of  the proposal distribution q(x_k|z_1:k) is factorized as follows:

q(x_k|z_1?k)=q(x_k|x_k-1, z_1:k)q(x_k-1|z_1:k-1)                (5)

It is known that  is generated by augmenting each of the existing sampling  sampled from the proposal distribution q(x_k-1|z_1:k-1) with the new state sampled from q(x_k|x_k-1, z_1:k). In order to obtain the updating equation of particle weight  p(x_k|z_1:k) is expressed in terms of p(z_k|x_k), p(x_k|x_k-1) and p(x_k-1|z_1:k-1). Noting that

p(x_k|z_1:k)=p(z_k|x_k, z_1:k-1)p(x_k|z_1:k-1)/p(z_k|z_1:k-1)∝

p(z_k|x_k)p(x_k|x_k-1)p(x_k-1|z_1:k-1)             (6)

Under assumptions that states are subjected to a Markov process and the observations are conditionally independent, according to Eqs.(4)-(6), the particle weight is given by

     (7)

In the practical application, the proposal distribution is commonly selected as

                   (8)

Hereby, weight expression is written as

                           (9)

After assessing particle weights, re-sampling step is applied to solving particle degeneracy. At present, re-sampling methods are mainly in the following categories [1, 6]: the residual re-sampling, the system re-sampling, and the polynomial re-sampling. That is standard particle filter and also known as bootstrap filter.

3 Particle filtering algorithm based on observa- tion inversion optimal sampling

The research shows that the estimation precision of the standard particle filter is restricted under the high observation accuracy conditions. Particles sampled from the above proposal distribution cannot consider the utilization of the latest observation information. There is a large deviation between particles sampled from the proposal distribution and those sampled from true posterior distribution, so only a few particles have significant importance weights evaluated by their likelihood functions, especially when the likelihood function locates in the tail of state transition probability density function or is too narrow compared with the state transition probability density. That is, the deviation is particularly evident when observation model has a quite high precision. In this case, the majority of particles become invalid samples for small weights, and low sampling efficiency even leads to filtering divergence. With the development of sensor manufacturing technology, high-precision observation occasions can be frequently encountered in practical engineering. From the above analysis, a novel particle filtering based on observation inversion optimal sampling is proposed as follows.

First, N virtual observations  are generated based on the latest observation z_k. Here, two kinds of sampling strategies are proposed. The first strategy is that  is sampled from normal distribution function N₁(z_k, σ), where σ is a small value, so  can centralize on the circumference of true observation. Then,  is generated by the inversion of  based on the function relationship between z_k and x_k-1. Finally,  is sampled in terms of the system one-step state transition model. The second strategy is that  is generated by increasing the sampling speed of sensor, and then  is sampled based on the anterior principle. Hence, the expression of proposal distribution is modified as follows:

 (10)

where  reflects the prior modeling information;  and reflect redundant information and complementary information coming from the latest observation, respectively. When particles are sampled from the proposal distribution expressed in Eq.(10), they may make fully use of the latest observation information and the prior modeling information. Under the high observation accuracy condition for non-linear system, the less the information lost in the generation process of  the higher the approximation of  for true state. Combining Eq.(10) with Eq.(7), weightis evaluated by

       (11)

In accordance with the above construction process of IPF, the main differences between IPF and PF lie in the following two aspects: the generation of virtual observations and the application of observation inversion strategy. In essence, the generation of virtual observations is a kind of sampling process for the latest observation. It realizes the preservation of redundant information and the extension of complementary information from current observation. In addition, it provides condition of the improvement for sampling particle diversity. According to the function relationship between observation and state, an observation inversion strategy is used to implement the absorption and the utilization for the above two kinds of information, and  is generated. On this basis,  is sampled by combining  and model prior information. In the IPF, the observation information and the model prior information are effectively utilized so as to implement further improvement of the effectiveness of sampling particles. The IPF algorithm flow is expressed by the pseudo-code, and the concrete algorithm description is as follows.

Initialization step: at time k=0

(1) For i=1, 2, …, N

Draw state ～ p(x₀)

Evaluate weight ～

Normalize weight ～

End

(2) Importance sampling step: for k≥1

For i=1, 2, …, N

Generate～N₁(z_k, σ) though increasing obser- vation speed

 is sampled according to function relation- ship between  and

 is sampled from

Evaluate the importance weight

Normalize weight ～

End

(3) Apply re-sampling algorithm

Set weight 1/N

(4) Output step:

4 Simulation result and analysis

In order to illustrate the advantages of IPF, two examples are given, which respectively aim at the state estimation of nonlinear/non-Gaussian single-dimensional system and nonlinear high-dimensional system.

4.1 Example 1

The system model is taken from Ref.[6] as follows.

                 (12)

＞

≤

                        (13)

where  w_k is a GA(3, 2) random variable modeling the process noise. The observation noise is drawn from a Gaussian distribution N₁(0, 0.000 1). In this case, the number of Monte Carlo run M is 50. The true initialization state of the target is x₀=1. In every run, the total simulation time t_total is 60 s, sampling interval t is  1 s and the number of particles N is 200, and the sample

average form  is considered as the

algorithm output. The root mean square error σ_RMSE is used as the performance evaluation index of algorithm precision, which is defined as

                 (14)

Figs.1 and 2 show the estimation of the system over single simulation and σ_RMSE over 50 times Monte Carlo simulations based on four algorithms, respectively. It is shown that the estimation precision of IPF is evidently superior to that of PF, EPF and UPF. Especially, the state model and observation model are all at nonlinear stage.

Fig.1 Comparison of state estimation with true result

Fig.2 Comparison of σ_RMSE for different filters

The quantitative comparison of the average time, the mean and variance of σ_RMSE of state estimation is shown in Table 1. It further verifies that the performance of the IPF is superior to that of PF, EPF and UPF. The mostly reason is that new sampling method can make full use of advantages of the latest observation information and the prior modeling information for system, so that the current time particles more effectively approximate the posterior probability density distribution. In addition, the process of observation inversion sampling avoids the Gaussian noise distribution restriction, and overcomes the disadvantage that the EPF and the UPF are sensitive to non-Gaussian noise.

Table 1 Comparison for average time, mean and variance of dimensionless σ_RMSE over 50 independent runs

4.2 Example 2

The observations of two-coordinate radar are used to realize target tracking. The target moves within the horizontal-vertical plane. The target state evolution model and the observation evolution model are as follows:

                          (15)

          (16)

where  x_k,  y_k and  denote the position component and velocity component of target state based on horizontal and vertical directions, respectively; system state transition matrix

 and  u_k is the

system noise vector, which meets zero mean Gaussian white noise process with standard deviation Q_k=       I; sampling interval is t=1 s, and system noise matrix

 v_k is observation noise vector

which meets zero mean Gaussian white noise process

with standard deviation  the noise

standard deviation of radial distance component is R_r=  10 m, and  the noise standard deviation of azimuth angle component is R_θ=0.2?. The true initialization state of the target is X₀=[2 000 80 4 000 50]^T, and the variance

initialization is  where P_k,x=  P_k,y= The number of Monte Carlo

run is M=50. In every run, the total simulation time is t_total=50 s, and the number of particles is N=200. The simulation results are as follows.

Fig.3 shows the σ_RMSE of position estimation in horizontal and vertical directions over 50 times Monte Carlo. Table 2 shows the quantitative comparison of the average time, the mean of σ_RMSE of position estimation in horizontal and vertical directions. The simulation significance illustrates the effectiveness of IPF for high-dimensional nonlinear system state estimation. In accordance with the results given in Fig.3, the estimation precision of IPF is superior to that of PF, EPF and UPF, and the estimation result of PF shows filter divergence phenomenon. Though IPF increases less than PF in computational complexity, it is far superior to EPF and UPF. In addition, it also shows that the real time of IPF is better than that of EPF and UPF when system state dimensional is high.

Fig.3 σ_RMSE of position estimation in horizontal (a) and vertical (b) directions

Table 2 Comparison for average time, mean of RMSE over 50 independent runs

5 Conclusions

(1) An improved particle filtering algorithm based on observation inversion optimal sampling is proposed.

(2) Compared with the classical PF and its variants, IPF has advantages as follows. Firstly, IPF adopts the basic framework of PF, so it preserves the advantage to solve nonlinear and non-Gaussian problems. Secondly, the model prior information and the latest observation information are utilized effectively in the generation process of due to the application of virtual observations and the observation inversion strategy, which makes particles move to high likelihood area and smooth the distribution of sampling particles, so as to improve the effectiveness and diversity of sampling particles. In addition, IPF also overcomes the rapid expansion of computational complexity caused by MCMC process and swarm intelligence method which are used to promote particle diversity, and the adverse effect caused by suboptimal filter which is used to optimize proposal distribution. Finally, the construction of IPF algorithm is based on existing information. Meanwhile, it also possesses explicit physical meaning and is easy to realize in engineering.

(3) A novel approach for the improvement of effectiveness and diversity of particles in IPF is presented, but it should meet the precondition that the analytical solution ofexists based on observation inversion strategy. The condition, which is weakened or overcome, will be the research emphasis at the next step.

References

[1] GORDON N J, SALMOND D J, SMITH A F M. Novel approach to non-linear/non-Gaussian Bayesian state estimation [J]. IEE Proceedings on Radar, Sonar and Navigation, 1993, 140(3): 107-113.

[2] ARULAMPALAM S, MASKELL S, GORDON N, ?LAPP T. A tutorial on particle filters for online non-linear/non-Gaussian Bayesian tracking [J]. IEEE Trans on Signal Processing, 2002, 50(2): 174-188.

[3] BOLIC M, DJURIC P M, HONG S J. Re-sampling algorithms and architectures for distributed particle filters [J]. IEEE Trans on Signal Processing, 2005, 53(7): 2442-2450.

[4] GIREMUS A, TOURNERET J Y, CALMETTES V. A particle filtering approach for joint detection/estimation of multi-path effects on GPS measurements [J]. IEEE Trans on Signal Processing, 2007, 55(4): 1275-1285.

[5] DUAN Zhuo-hua, CAI Zi-xing, YU Jin-xia, ZOU Xiao-bing. Particle filter based fault diagnosis for inertial navigation system of mobile robot [J]. Journal of Central South University: Science and Technology, 2005, 36(4): 642-647. (in Chinese)

[6] RISTIC B, ARULAMPALAM S, GORDON N. Beyond the Kalman filter: Particle filters for tracking applications [M]. London: Artech House Publishers, 2004.

[7] IOANNIS K, DARRYL M, ANTONIA P S. Sequential Monte Carlo methods for tracking multiple targets with deterministic and stochastic constraints [J]. IEEE Trans on Signal Processing, 2008, 56(3): 937-948.

[8] CAPPE O, GODSILL S J, MOULINES E. An overview of existing methods and recent advances in sequential Monte Carlo [J]. Proceedings of the IEEE, 2007, 95(5): 899-924.

[9] ZHAI Y, YEARY M. Implementing particle filters with Metropolis- Hastings algorithms [C]// Proceedings of Region 5 Conference: Annual Technical and Leadership Workshop. Tokyo, 2004: 149-152.

[10] UOSAKI K, HATANAKA T. Evolution strategies based particle filters for fault detection [C]// Proceedings of the IEEE Symposium on Computational Intelligence in Image and Signal Processing. Honolulu, 2007: 58-65.

[11] LI Cui-yun, JI Hong-bing. A new particle filter with GA-MCMC re-sampling [C]// Proceedings of International Conference on Wavelet Analysis and Pattern Recognition. Beijing, 2007: 146-150.

[12] GIREMUS A, TOURNERET J Y, DJURIC P M. An improved regularized particle filter for GPS/INS integration [C]// Proceedings of IEEE 6th Workshop on Signal Processing Advances in Wireless Communications. New York, 2005: 1013-1017.

[13] GUO Wen-yan, HAN Chong-zhao, LEI Ming. Improved unscented particle filter for nonlinear Bayesian estimation [C]// Proceedings of 10th International Conference on Information Fusion. Quebec, 2007: 1-6.

[14] LI Qian, FENG Jin-fu, PENG Zhi-zhuang, LU Qing, LIANG Xiao-long. An iterated extend Kalman particle filter for multi-sensor based on pseudo sequential fusion [C]// Proceedings of IEEE International Conference on Robotics and Biomimetics. Sanya, 2007: 1534-1539.

[15] LI Liang-qun, JI Hong-bing, LUO Jun-hui. The iterated extended Kalman particle filter [C]// Proceedings of IEEE International Symposium on Communications and Information Technology. Beijing, 2005: 1213-1216.

[16] KARLSSON R. Particle filter for positioning and tracking applications [D]. Linkoping: Linkoping University, 2005.

[17] SMITH L, AITKEN V. Analysis and comparison of the generic and auxiliary particle filtering frameworks [C]// Proceedings of IEEE Canadian Conference on Electrical and Computer Engineering. Ottawa, 2006: 2124-2127.

(Edited by CHEN Wei-ping)

Foundation item: Project(60634030) supported by the Key Project of the National Natural Science Foundation of China; Project(60702066) supported by the National Natural Science Foundation of China; Project (2007ZC53037) supported by Aviation Science Foundation of China; Project(CASC0214) supported by the Space-Flight Innovation Foundation of China

Received date: 2008-11-15; Accepted date: 2009-03-22

Corresponding author: HU Zhen-tao, Doctoral candidate; Tel: +86-13669277886; E-mail: guchenshou@yahoo.com.cn