1. Introduction
In the fault diagnosis, it is impossible to establish an accurate system model because of the complexity of the system and the interaction of the working characteristics of each component, so the errors of model deviation and the parameter deviation are inevitable. Moreover, it is greatly reduced the accuracy of fault diagnosis by the influence of environment and human disturbance factors. With the increasing requirements for the fault diagnosis precision of large-scale complex equipment, the fault diagnosis method of complex system in strong noise environment becomes one of the difficult problems in this field [1,2].
Particle filter (PF) is an approximate estimation of state based on Monte Carlo sampling theory, which can be theoretically applied to any nonlinear non-Gaussian system, so it has been widely researched and popularized in the field of fault diagnosis [3-5]. But in the application, the present PF method and its improved algorithms assumed that the statistical characteristic of noise is the known Gaussian distri¬bution. We cannot obtain the noise accurately in the real environment, and its statistical characteristics are unknown, which is not in conformity with the actual situation [6,7], and the accumulation of errors in the complex system can lead to lower or even divergence of the accuracy of the state estimation. To solve this problem, this paper proposed the cost evaluation of particle filter (CRPF) algorithm [8]. CRPF algorithm is based on the cost or risk of particles for the reference to achieve resampling and state estimation, the cost and risk calculation does not require known noise statistic characteristics. It has been widely used in the field of fault diagnosis and target tracking.
Particle filter (PF) is an approximate estimation of state based on Monte Carlo sampling theory, which can be theoretically applied to any nonlinear non-Gaussian system, so it has been widely researched and popularized in the field of fault diagnosis [3-5]. But in the application, the present PF method and its improved algorithms assumed that the statistical characteristic of noise is the known Gaussian distri¬bution. We cannot obtain the noise accurately in the real environment, and its statistical characteristics are unknown, which is not in conformity with the actual situation [6,7], and the accumulation of errors in the complex system can lead to lower or even divergence of the accuracy of the state estimation. To solve this problem, this paper proposed the cost evaluation of particle filter (CRPF) algorithm [8]. CRPF algorithm is based on the cost or risk of particles for the reference to achieve resampling and state estimation, the cost and risk calculation does not require known noise statistic characteristics. It has been widely used in the field of fault diagnosis and target tracking.
It can be seen that CRPF has been well used in solving problems of the nonlinear and unknown noise statistic characteristics, but there are still two deficiencies in practical application: (1) when the back¬ground noise is large, and the dimension of observation is high, it can lead to the sample risk and the enormous cost, and sample impoverishment will appear in the resampling process [10,13]; (2) CRPF takes the state predictive value of the current moment as the mean value and randomly selects the particles to complete the state updating within a certain range of Gaussian distribution [11,14]. However, as time goes on, the covariance of state updating density gradually tends to be constant, therefore the ability to modify time-varying noise diminishes, greatly reducing the ability of the algorithm to adapt to noise environment. The research of fault diagnosis method in complex noise environment is a basic scientific problem that must be solved in the research of fault diagnosis and fault prediction of modern complex systems. Therefore, this study is of great significance to improve the accuracy of diagnostic system and further improving the reliability of modern complex equipment fault diagnosis.
In this paper, the problem of sample impoverishment in resampling process is overcome by the opti¬mization of small-value sample particles to replace resampling process. At the same time, the variance of state transfer density is adaptively adjusted according to the correlation between current measured value and prior state. The reasonable region of state update is determined to improve the adaptability to time-varying noise. On this basis, the problem of nonlinear non-Gaussian system fault diagnosis is simulated and analyzed in the condition of strong noise.
2. The Cost Reference Particle Filter
The nonlinear and non-Gaussian stochastic state space model can be described as:
where [TeX:] $$x_{k}$$ represents the state vector at time [TeX:] $$k$$, [TeX:] $$y_{k}$$ represents the measurement vector at time [TeX:] $$k$$, [TeX:] $$g(\bullet)$$ and [TeX:] $$h(\bullet)$$ are the system transition function and measurement function which can be nonlinear functions. [TeX:] $$u_{k}$$ and [TeX:] $$w_{k}$$ are system noise and measurement noise.
In the CRPF algorithm, it used the cost function and risk function to represent the performance and quality of particles. Meanwhile, the forgetting factor is introduced. Based on the principle of cost minimization, a particle weight evaluation method with unknown statistical characteristics are given. The cost function is defined as:
The equation abbreviated as:[TeX:] $$\boldsymbol{C}_{k}^{i_{p}}=\lambda \boldsymbol{C}_{k-1}^{i_{p}}+\Delta \boldsymbol{C}_{k}^{i_{p}}$$. Where [TeX:] $$0 \leq \lambda \leq 1$$ is the forgetting factor. [TeX:] $$\Delta \mathcal{C}_{k}^{i_{p}}$$ is the “incremental cost function” and computed by [TeX:] $$\left\|y_{k}-h\left(x_{k}^{i_{p}}\right)\right\|^{q}$$ with [TeX:] $$q \geq 1$$.
the risk function is a prediction of the cost increment [TeX:] $$\Delta \mathcal{C}_{k}^{i_{p}}$$ , and the predictive cost function is defined as follows:
The probability mass function (PMF), also known as class weight, is computed by using the risk function as follows:
where [TeX:] $$\delta, \beta>0$$. [TeX:] $$\delta$$ is to ensure that the denominator is not zero. According to the above definition of parameters, the CRPF algorithm obtains the state estimation by recursive calculation through the steps of risk estimation, selection, particle transfer and cost update. The algorithm steps are as follows:
(1) At time k = 0, N samples were obtained from the prior distribution [TeX:] $$x_{0}^{i_{p}} \sim P_{0}\left(x_{0}\right)$$ , Set the cost [TeX:] $$C_{0}^{i_{p}}$$ and the transfer density covariance[TeX:] $$\sigma_{0}^{2, i_{p}}$$ of the particles initial time, then sample and the cost set are [TeX:] $$\left\{x_{0}^{i_{p}}, \mathcal{C}_{0}^{i_{p}}\right\}_{i_{p}=1}^{N}$$ , [TeX:] $$i_{p}=1,2, \cdots, N$$ , at time k = 0.
(2) recursive update:
① Calculate the risk function [TeX:] $$R_{k_{k}}^{i_{p}}$$ and PMF [TeX:] $$\tilde{\pi}_{k}^{i_{p}}$$ according to (4) and (5), respectively.
② Resampling. According to the value of [TeX:] $$\tilde{\pi}_{k}^{i_{p}}$$ , select N particles randomly, generate particle and cost set [TeX:] $$\left\{\begin{array}{l} \left.\bar{x}_{k-1}^{i_{p}}, \widetilde{\mathcal{C}}_{k-1}^{i_{p}}\right\}_{i_{p}=1}^{N} \\ \end{array}\right.$$.
③ particle update.
I is the unit matrix with the same dimension as x, and the [TeX:] $$\sigma_{k}^{2, i_{p}}$$ is given by
From the formula (7) we can see that the covariance of the transfer density [TeX:] $$\sigma_{k}^{2, i_{p}}$$ in the CRPF algorithm is adjusted online.
④ Calculate the cost function [TeX:] $$C_{k}^{i_{p}}$$ according to Eq. (2). Calculate [TeX:] $$\pi_{k}^{i_{p}}$$ according to Eq. (8) and normalize it.
(3) State estimation:
It can be seen from the above algorithm steps that no calculation of process noise and measurement noise is involved in the recursive estimation of the CRPF algorithm. Therefore, it is not necessary to know the statistical characteristics of noise in the calculation, which improves the adverse effects of random external disturbance on the basis of the evaluation of measurement likelihood in PF. However, when the system noise and measurement noise are larger, the corresponding risk and cost of each sample will be large, which can be seen from Eqs. (5) and (8), the small [TeX:] $$R_{k}^{i_{p}}$$ and [TeX:] $$C_{k}^{i_{p}}$$ samples of the class weights are much larger than the other sample class weights, resampled particles with large weights will be copied many times, and particles with small weights are discarded, resulting in the lack of sample diversity. In addition, from the particle updating process of Eqs. (6) and (7), we can see that the time update of the particles is to select new particles randomly in the determined probability region, the variance [TeX:] $$\sigma_{k}^{2, i_{p}}$$ of this region is recursive update based on the state estimates of the current moment and the variance of the previous moment. When [TeX:] $$k \rightarrow \infty$$, [TeX:] $$\frac{k-1}{k} \rightarrow 1$$, [TeX:] $$\frac{1}{k} \rightarrow 0$$, then [TeX:] $$\sigma_{k}^{2, i_{p}} \rightarrow \sigma_{k-1}^{2, i_{p}}$$. It can be seen that the current value of the system state loses the adjustment effect on the variance over time, the variance tends to be constant, which greatly reduces the accuracy of the algorithm when the time-variant noise or the noise fluctuates greatly.
3. Using Particle Optimization Instead of Resampling Process
Aiming at the problem of particle depletion in the process of resampling, we improve the crossover and mutation operation by considering the real-time performance and the effectiveness of the particles, and the credibility of the whole particle is improved by optimizing the selection of the larger risk particles. First, we divide the sample into large weight particles set and small weight particles set (according to the size of the class weight), and only the operation of the small weight particles are performed to reduce the calculation. For the crossover operation, the weighting particle is adopted to modify the small-value particle, which can enhance the adaptability of the small-value particle without changing the weighting particle. Therefore, the superiority of the weighting particle is maintained. Then, the modified small-value particle is mutated according to the given mutation probability. Weights of small particles evolve into large-valued particles so as to improve the problem of particles depletion.
Sort the particles in ascending order of [TeX:] $$\tilde{\pi}_{k}^{i_{p}}$$ 's value to get [TeX:] $$\left\{x_{k}^{\pi^{1}}, x_{k}^{\pi^{2}}, \ldots, x_{k}^{\pi^{N}}\right\}$$ , and then follow the steps below to optimize them.
(1) According to the value of [TeX:] $$\tilde{\pi}_{k}^{i_{p}}$$ , the particles are divided into large weight particles set [TeX:] $$C_{H}$$ and small weight particles set [TeX:] $$C_{L}$$.
The number of effective particles is defined as [TeX:] $$N_{e f f}=\left[1 / \sum_{i_{p}=1}^{N}\left(\pi_{k}^{i_{p}}\right)^{2}\right]$$ , and the PMF of the particle with number [TeX:] $$N_{e f f}$$ is set to the threshold [TeX:] $$\pi_{T}$$ , then compute as [TeX:] $$\pi_{T}=\tilde{\pi}_{k}^{i_{N e f f}} \propto \mu_{1}\left(R_{k}^{i_{N e f f}}\right)$$ and make [TeX:] $$x_{k L}^{l} \in C_{L}\left(l=1, \ldots, N_{L}\right)$$ , [TeX:] $$x_{k H}^{j} \in C_{H}\left(\mathrm{j}=1, \ldots, N_{H}\right)$$ , where [TeX:] $$N_{L}$$ and [TeX:] $$N_{H}$$ represent the number of the small weight particles and the large weight particles, respectively.
(2) The small weight particles are crossed,
[TeX:] $$x_{k S}^{l}$$ are the crossed particles. The larger the parameter [TeX:] $$\alpha \in[0,1]$$, indicates that more information is passed from [TeX:] $$x_{k L}^{l}$$ to the descendant particles [TeX:] $$x_{k S}^{l}$$ . If [TeX:] $$\alpha=1$$ , there is no cross operation, and [TeX:] $$x_{k H}^{j}$$ will be selected from [TeX:] $$C_{H}$$ randomly.
(3) Mutate the small weight particles after crossover operation,
where [TeX:] $$r_{l}$$ is a random variable about [TeX:] $$x_{k S}^{l}$$ , which can be obtained from the uniform distribution of [0,1], and [TeX:] $$p_{M}$$is the mutation probability.
The class weights of the particles optimized by the above steps are recalculated according to Eq. (5) to generate a new sample set of particle and cost [TeX:] $$\left\{\check{x}_{k-1}^{i_{p}}, \check{C}_{k-1}^{i_{p}}\right\}_{i_{p}=1}^{N}$$ . By optimizing the particles to replace the resampling process, the particles’ diversity is improved, and the computational cost is saved at the same time.
4. Adaptive Update the State Transition Density Covariance
From the analysis of part 1, the variance of state transition density directly affects the accuracy of state prediction. But in the CRPF algorithm, the recursive update of variance tends to be constant over time. In the case of random fluctuation of time-varying noise, it is likely that accurate state predictions will not be available in the determined small range. Therefore, this paper considers the correlation coefficient between the currently observed value and the prior state, adaptively adjusts the variance of state transition density, and enhances the ability of the algorithm to adapt to time-varying noise.
The correlation between the observed value and the prior state is expressed as:
The covariance of the state transition density is corrected according to the following formula:
where ω_1 and ω_2 are adaptive adjustment factors of the variance. [TeX:] $$\omega_{1}=\gamma, \omega_{2}=1-\gamma$$. If the correlation coefficient is large, indicating that the state noise is small, then [TeX:] $$\omega_{1}>\omega_{2}$$, and the variance is mainly determined by the variance of the previous moment, which can enhance the particle disturbance appropriately. Otherwise, if the correlation coefficient is small, indicating that the noise disturbance is larger, then [TeX:] $$\omega_{1}<\omega_{2}$$, the variance is mainly determined by the difference between the state estimate and the predicted value at the current time. That is to say, the noise disturbance at the current moment is mainly used to correct the variance, which improves the accuracy of the algorithm for the state estimation under the background of time-varying noise.
5. Steps of the Proposed Method
On the basis of the previous introduction, intelligent optimization resampling and state transition density covariance adaptive control method are introduced into the traditional CRPF algorithm, and the improved method is obtained. The detailed steps for the sequential estimation x are given in Table 1.
The algorithm steps proposed in this paper
6. Simulation and Analysis
The dynamic model of a 160-MW fuel unit [15] is used to validate the algorithm. The discrete equations of the model are as follows:
where, the relevant parameters are defined as,
where [TeX:] $$x_{1}$$ is the drum pressure, [TeX:] $$x_{2}$$ is the drum liquid density, [TeX:] $$u_{1}$$ is the fuel control valve; [TeX:] $$u_{2}$$ is the steam turbine regulating valve, [TeX:] $$u_{3}$$ is the feed water control valve, [TeX:] $$y_{2}$$ is the drum water level, [TeX:] $$v_{k}$$ is process noise, [TeX:] $$n_{k}$$ is measurement noise, they are all non-Gaussian noise with unknown statistical properties. The experiments of drum water level sensor faults are carried out and suppose there are 3 fault modes of sensor:
Fault mode 1: Constant deviation fault,
Fault mode 2: Constant gain fault,
Fault mode 3: Stuck fault,
Assuming that, when [TeX:] $$0 \leq t \leq 100 s$$ , the system is in normal mode, when [TeX:] $$100 s \leq t \leq 200 s$$ , the system is in fault mode 1, when [TeX:] $$200 s \leq t \leq 300 s$$ , the system is in fault mode 2, and when [TeX:] $$300 s \leq t \leq 400 s$$ , the system is in fault mode 3.
Set initial state [TeX:] $$x_{0}=\left[\begin{array}{ll} 108 & 428 \end{array}\right]^{\mathrm{T}}$$ , discrete step [TeX:] $$\Delta t=0.1 s$$ , particle number N=500, and sampling frequency of sensor is 1 Hz. Other parameters:[TeX:] $$a_{11}=0.0018, a_{12}=0.9, a_{13}=0.15, a_{21}=141, a_{22}=1.1, a_{23}=0.19, b_{1}=0.131, b_{2}=0.068, b_{3}=0.00154, b_{4}=0.8, b_{5}=25.6, b_{6}=1.0394, b_{7}=0.00123, b_{8}=0.854, b_{9}=0.147, b_{10}=45.59, b_{11}=2.514, b_{12}=2.096, u_{t}=[0.3 \quad 0.4 \quad 0.5] .$$. The initial prior distribution of the state is [TeX:] $$x_{0}^{i_{p}} \sim N\left(x_{0}, \Sigma_{0}\right), \quad \Sigma_{0}=\operatorname{diag}(0.01,0.01), \quad C_{0}^{i_{p}}=0, \sigma_{0}^{2, i_{p}}=\left[\begin{array}{ll} 0.01 & 0.01 \end{array}\right]^{T}$$ , residual is [TeX:] $$r=\left|\hat{y}_{2}-y_{2, \text { 正常 }}\right|$$ .
The mean absolute error and the effective particle number are chosen as the evaluation index of the algorithm, the mean absolute error is defined as,
[TeX:] $$x_{s k}$$ and [TeX:] $$\hat{x}_{s k}$$ are the actual value and estimated value of the step [TeX:] $$k$$ status of the [TeX:] $$s^{t h}$$ simulation, respectively, the total simulation time is [TeX:] $$N_{s}=50$$, and the time step is T=4000 in a simulation.
Experiment 1. Set the noise as gamma distribution, [TeX:] $$\left\{\begin{array}{l} v_{k} \sim 0.09 \times 10^{-4} \times \Gamma(0.25,0.5) \\ n_{k} \sim 0.16 \times 10^{-4} \times \Gamma(0.25,0.5) \end{array}\right.$$ , The experimental results of CRPF and improved CRPF are compared and analyzed as follows.
Experiment 2. Set the noise as Gaussian mixture distribution,
The experimental results of CRPF and improved CRPF are compared and analyzed as follows.
Experiments have been done based on the traditional CRPF method and the proposed method. Experiments 1 and 2, respectively use two different kinds of non-Gaussian random noise to simulate the traditional CRPF method and the proposed method. Each simulation result is the average value obtained by 4000 steps iteration calculation 50 times. Figs. 1 and 2 show the result of tracing the fault condition of the drum water level sensor with the traditional CRPF under two different kinds of noise, respectively. As shown in Figs. 1 and 2, the deviation between the estimated and actual values is large. Figs. 3 and 4 show the tracking results of the improved CRPF method under two kinds of noise. Due to the introduction of the difference between the current state value and the predicted value, the variance of the state transition density can be corrected in real time so that it can achieve fast tracking of random fault states. The simulation results show that the proposed method can achieve accurate tracking to different fault modes under two different kinds of noise. The residuals can also be seen clearly, in the failure time of 100s, 200s, and 300s, the residual of the proposed method has changed dramatically, which can provide an effective basis for the identification of different faults.
(a) The tracking results and (b) the residuals for fault status [TeX:] $$y_{2}$$ using CRPF in Experiment 1.
(a) The tracking results and (b) the residuals for fault status [TeX:] $$y_{2}$$ using CRPF in Experiment 2.
(a) The tracking results and (b) the residuals for fault status [TeX:] $$y_{2}$$ using the proposed method in Experiment 1.
(a) The tracking results and (b) the residuals for fault status [TeX:] $$y_{2}$$ using the proposed method in Experiment 2.
As shown in Tables 2 and 3, it can be seen that the tracking accuracy of [TeX:] $$x_{1}$$ and [TeX:] $$y_{2}$$ of the improved CRPF method is significantly higher than that of CRPF method in the strong noise environment. The tracking error of [TeX:] $$x_{2}$$ is great using the improved method, because the fault value of [TeX:] $$y_{2}$$ is used as the error measure to correct [TeX:] $$x_{2}$$, so [TeX:] $$x_{2}$$ deviates from the actual value and gradually approaches the fault state. However, as the number of iterations increases, the traditional CRPF method gradually loses the effect of the current measurement on the state. Thus, the tracking status does not change much in the event of a failure. This also shows that the improved method has greatly enhanced the adaptability to external disturbance and time-varying noise. At the same time, particle optimization can replace the resampling process to improve the particle leanness of the CRPF algorithm. From Tables 1 and 2, the number of effective particles is obviously increased.
Through the above simulation analysis, we obtain a completely consistent conclusion that the accuracy of the improved method has been significantly improved in the strong noise environment.
Performance of the two algorithms in Gamma noise environment
Performance of two algorithms in Gaussian mixture noise environment
7. Conclusions
The CRPF algorithm can recursively estimate the state of nonlinear non-Gaussian system under the condition of unknown noise statistical characteristics, and provide a reference method for the fault diagnosis of a complex nonlinear system. In this paper, we analyze the problem of particle impoverishment in the resampling process of CRPF algorithm and introduce the intelligent optimization operation. Replacing the resampling process by particle optimization improves the reliability of the whole particle, thereby improving the problem of particle depletion and saving computational cost. By analyzing the performance of CRPF algorithm for state estimation under strong noise interference, the correlation discriminant function is designed to adaptively adjust the variance of state transition density of CRPF, which overcomes the problem of large error or divergence of CRPF with time, because of the correction ability for strong noise and time-varying noise is gradually disappear. Finally, the fault detection for nonlinear non-Gaussian system is carried out in different noise environments. The results show that the proposed method improves the accuracy of fault diagnosis and has important theoretical significance in practical application.
Acknowledgement
This paper is funded by the National Natural Science Foundation of China (No. 61763028) and the Natural Science Foundation of Gansu (No. 1506RJZA105) and the open project of Key Laboratory of Gansu Advanced Control for Industrial Processes (No. XJK201805).