Prediction of remaining useful life for lithium‐ion battery based on particle filter with residual resampling

Accurate prediction of the remaining useful life for lithium‐ion battery is beneficial to prolong the life of the battery and increase safety. With the capacity degradation curve obtained from the data of the battery charge and discharge experiment, the remaining useful life of the battery was predicted by using particle filter. In order to improve the prediction accuracy, the particle filter with residual resampling method is used to overcome the lack of particle diversity which has an important effect on the accuracy of state estimation. Compared with the prediction result of the extended Kalman filter, it was found that the precision and stability of particle filter are better than those of extended Kalman filter. The research results presented in this paper provide some suggestions for the health monitoring of power battery for electric vehicles.

are mass transfer, conduction, and electrochemical reaction, respectively. The processes above are organically unified by the governing equation in. 7 The equivalent circuit model represents the external characteristics of the battery by simulating the variation rules of the battery voltage, current, resistance, etc. It does not require much knowledge of electrochemistry for the researchers. 8 Empirical model is the mathematical model about system parameters and system variables, and it is the empirical function conforming to the laws obtained from the statistical analysis of system parameters and variables. 9 Due to the small amount of calculation, the empirical model can be used for online prediction and it is mostly used in the prediction of remaining useful life for the batteries, such as the exponential increase model based on the impedance, exponential decrease model based on the capacity, and polynomial model. 10 On the basis of these models, many filtering algorithms have been developed, such as unscented Kalman filtering and particle filtering. 11 Data-driven methods are focused on making short-time prediction for the battery based on their historical data. 12 Chao Hu et al put forward a kind of forecasting method which integrates data-driven. 13 It combines multiple algorithms by introducing the weight. The k-fold cross-validation is used to calculate the prediction error and three weighting schemes are proposed: they are based on accuracy, diversity and optimization, respectively. The results show that no matter which weighting scheme is used, the prediction accuracy of RUL is higher than that of single prediction algorithm. In addition, the prediction robustness can be improved. Patil et al proposed a multi-node support vector machine (SVM) method, which integrates classical models and regression models to improve the efficiency and achieve sufficient prediction accuracy under different working conditions. 14 The data-driven method requires a large amount of historical sample data and uses statistical as well as machine learning techniques to extract features from performance data such as current, voltage, time, and impedance to estimate the RUL of the battery. 15 The vector machine, neural network, and other intelligent algorithms used in the data-driven approach require a lot of computation, so how to reduce the computation complexity and improve the computational efficiency is a research focus. 16 The method based on fusion can be the cooperation of multiple forecasting methods or the weighted average of multiple forecasting methods. Dong Wang et al A prediction method with associated vector machine and three-hypothetical parameter capacity degradation model is proposed to predict RUL. 17 The associated vector machine is used to obtain the associated vector which can help to find the representative training vector, including the cycle of the associated vector and the predicted value in the cycle. Three-parameter capacity degradation model is used to fit the predicted value of the associated vector cycle. Then, the RUL value is obtained by extrapolation. Yang Chang et al proposed a new fusion method based on the idea of error correction, including unscented Kalman filter, complete empirical mode decomposition, and associated vector machine. 18 The general idea is to use unscented Kalman filter to forecast RUL and prediction error first. The prediction error is decomposed by complete empirical mode. Then, the reconstructed error sequence is obtained. At last, the associated vector machine is used to learn from the error sequence to get the predicted new error which is fed back to the unscented Kalman filter for correction.
For the prediction algorithm, more historical data are, more accurate the prediction is. However, in the early stage of the use of the new battery, the amount of data available is limited, so the accuracy of the prediction algorithm is not good and sometimes with large deviation from the actual value. It is a common phenomenon that the prediction becomes more and more accurate after the battery is used for a period of time. The degradation mechanism of the battery is very complex. Ignoring the operating conditions, the other factors like temperature, charging rate, and overcharge rate all will affect the battery life. If the actual operating conditions are taken into account, then there are more factors. Therefore, it is not practical to build a universally applicable observation model of battery degradation. Even if this complex observation model is built, it is unlikely to be used in cars due to the limited memory and computing speed of the battery management system. Also, the timeliness cannot be guaranteed. After the above analysis, it is obtained that the filtering methods can play an important role in the research of battery remaining useful life prediction, because they do not need too much data and the dynamic equation of the battery system can be used as the degradation observation model. Compared with Kalman filter, particle filter has better nonlinear adaptability and is not limited by system noise. Therefore, this paper studies the prediction of remaining useful life of lithium-ion battery by particle filter.

| Principle of particle filter
For the state space model of the unknown system, the essence of particle filter is to approximate the posterior probability (P( x k |z k ))of the unknown system. It uses Monte Carlo sampling method. With the prior distribution of the system and Bayesian estimation, distribution of the samples approximates the posterior distribution of the system. If we use particle filter to predict the posterior probability density (P( x k |z k )) of the system at time k, it is assumed that the posterior probability density (P( x k−1 |z k−1 )) of the system at time k-1 is known. If the prior distribution (P( x k )) of the system is known, then samples can be randomly selected according to the prior distribution. After obtaining the marginal | 1117 PAN et Al. probability (P( z k )) of the system at time k, the distribution represented by the samples can be modified. The corrected sample distribution can approximate the posterior distribution of the system at time k. The more samples there are, the more the predicted posterior distribution approaches the true posterior distribution.
Assume that the state equation and observation equation of the dynamic system are as follows: where f k is the state transition function, k−1 is the process noise at time k-1 and x k is the system state at time k.
where h k is the observation function of system state x k . k is the observed noise at time k. z k is the measurement at time k.
The algorithm flow of particle filter is as follows and shown in Figure 1: 1. Algorithm parameters setting: The number of particles N, process noise, observed noise, driving matrix, initial state value, etc 2. Initialization of particle set: The initial state x i 0 of each particle is obtained by sampling from the initial distribution P( x 0 ). Each particle has a weight i 0 . {x i 0 , i 0 } refers to the state and weight for the ith particle at the initial moment. The state and weight of the ith particle are equal in this step, namely {x i 0 , 1∕N};

Importance sampling
If the posterior distribution (P( x 0:k |z 1:k )) is known, the sampled particles obtained from it is most reliable. However, this distribution is multivariable and nonstandard, so it is difficult to directly sample from this distribution. Sampling from the prior distribution or introducing a similar distribution easy to sample whose probability density function is known is a choice. This distribution is called importance distribution.

Weight calculation
The calculation of the weight needs the observation value of the current moment Without this value, the state update for particle filter cannot be carried out, so the observation value cannot be optimized. The basic formula for weight calculation is as follows: where k (x 0:k ) is unnormalized importance weights. This formula is not a recursive formula. When the observation is updated, the importance weight needs to be calculated again and the calculation amount will increase with iterations, which will cause trouble in practical application.
Assume that the current state does not depend on future observations. It means only filtering is considered and smoothness is not taken into consideration. The importance function can be written as follows: Assume that the state conforms to the Markov process, it can be obtained: q(x 0:k |z 1:k ) , when formula 4 ~ 6 is substituted into formula 3, the recursive weight formula is obtained as follows: If the prior distribution p(x k | x k−1 ) is select as the importance distribution, it can be obtained: Substitute the formula 8 into formula 7, it can be obtained as follows: The weights need to be normalized:

Resampling
Resampling is to copy or discard particles ) according their weight. Then, the weight of all particles is reset to 1/N. The purpose of resampling is to avoid the problem of particle degradation.

Output:
Choose the mathematical expectation of the sampled particles as output:

Cyclic iteration
Return to step (3) and repeat steps (3) to (6) until the last measurement.

| Resampling algorithm
In order to solve the problem of particle degradation, it is necessary to screen the particles to eliminate some particles with large deviation and realize the "survival of the better" to suppress the particle degradation. The specific implementation method is using resampling algorithm. The reasons for particle degradation are as follows: firstly, the number of sampled particles is insufficient; secondly, because the accurate prior distribution of the system is difficult to get, the sampling method itself is not reasonable; thirdly, due to the reason that particles with large sampling error may have larger weight, these particles are duplicated in large quantities in the iterative process. Make k is the corresponding weight for each particle. At the first sampling, the particle set of particle filter is {x ( 1 ) k , 1∕N}. N is the number of particles. The addition of new observed values helps modify the distribution represented by the particle set. Before the next importance sampling, the particle set is resampled according to the revised results. The weight of each particle becomes 1/N again after resampling. The following are some common resampling algorithms:

Random resampling
The principle of this algorithm is shown in Figure 2. Assume that the number of particles is N. The interval of the cumulative distribution function cdf is divided into subintervals with number of N. A random number obeying uniform distribution over [0,1] is generated. Depending on the range of the random number, the output will activate different responses. In other words, the greater the length of the interval in which the random number falls, the greater the weight of the corresponding particle is, the more times the particle is indexed. In turn, particles with small weights may be discarded.

Polynomial resampling
The distribution function of the discrete random variable X is expressed in the form of probability accumulation: Random number with uniform distribution over [0,1] is generated according to the above equation:  Copy or discard the particles according to their weight and the set after resampling is a. System resampling It is similar to random resampling algorithm. The algorithm steps are as follows: a. Divide (0, 1] into N intervals with the law:

a. Residual resampling
Residual resampling is proposed on the basis of random resampling. Its implementation is shown as the following pseudocode:

Residual resampling algorithm
Where N is the number of particles, i is the weight of each particle before resampling, "Index" stores the sequence number of particles, [x] refers to the operation of rounding x, N_res is the number of particles that need to be copied, N-N_res is the number of particles directly retained, W i is the cumulative distribution by calculating the sum of the elements: 1 , 2 , . . . , i , u is the random number over [0,1], "OutIndex" stores the sequence number of the particles that need to be copied or retained.
Random resampling, polynomial resampling and system resampling all determine whether the particle is retained or not according to the random number that conforms to the uniform distribution over [0,1]. However, the basic idea of residual resampling is to retain the particles with large weight directly, then the insufficient particles are supplemented with randomly selected re-weighted particles. In addition, residual resampling can maintain a certain precision while reducing the amount of calculation. Therefore, in this paper, particle filter with residual resampling is used to predict the 1∕N , remaining useful life of lithium-ion battery. The predicted results are compared with those obtained by particle filter with random resampling and extended Kalman filter to prove the validity and superiority of the particle filter with residual resampling.

| Experiment purpose
The charge and discharge cycle test of lithium-ion battery refers to that the battery is continuously charged and discharged after the corresponding condition is set. Data about capacity, resistance, current, voltage, and other parameters can be collected in each cycle. Test purpose and test equipment are as follows: The propose are as follows: 1. The data used for estimating the model state parameter are called training set. 2. The data compared with the predicted results to verify the feasibility of the prognostic approach are called the verification set.
The test equipment includes an incubator (see it in Figure 3) and a battery performance test cabinet (see it in Figure 4). The tested batteries are shown in Figure 5.
The four batteries used in the test are 18 650 batteries produced by Shenhuo Co. LTD in Shenzhen. The main parameters of the battery are shown in Table 1.

| Experiment settings and results
At first, charge the battery with constant current. Then, choose constant voltage to charge the battery. The specific steps are as follows:  (2) to (4) and stop the test when the capacity of the batteries drops to at least 75% of the rated capacity. Figures 6 and 7 show the capacity degradation curves of the tested batteries. It can be seen from the curves that the four batteries show two different degradation mechanisms. The degradation mechanism of battery 1 and 2 is the same. Battery 3 and 4 have another degradation mechanism. This is due to the different discharge cut-off voltage of four batteries. Lower cutoff voltage results in faster battery capacity degradation, such as battery 1 and battery 2. With the capacity degradation data, we can move on.

| Capacity degradation observation model of the Lithium-ion battery
There are many degradation observation models for lithiumion batteries. For example, RB Wright et al proposed an The more parameters the model has, the higher the model accuracy is. However, the model will be more complex. In addition to the conventional parameters, the model above also considers the influence of temperature and SOC, so the accuracy is high and the calculation is complex.
In the continuous charge and discharge, the maximum capacity that can be released or stored by the battery will decrease continuously. Saha et al found that the capacity of lithium-ion batteries decays exponentially. 20 He et al proposed an exponential decay observation model to fit the degradation curve of the battery. 21 Where Q k is the capacity of the battery at kth cycle. a, b, c, and d are the parameters of the model. k is the variable of the model, and it represents the number of the cycle. In practice, the parameters of the model are variable. Choose these four parameters as the state of the battery degradation dynamic model: Then, the state equation and observation equation of the model can be expressed as follows: Where X(k) is the state vector, a(k), b(k), c(k), d( k) are state parameters of the observe model and they are dimensionless. a ( k), b (k), c (k), d (k) are process noise.
( k) is the observation noise. Q k is the observed capacity. k is the number of cycles. This model is simpler than the impedance increase model, but it still maintains high accuracy.
In order to investigate whether the capacity degradation data of the battery used in this experiment can be approximated by the double exponential model, the actual capacity degradation data of battery 1, battery 2, battery 3, and battery 4 are imported into MATLAB and fitted with the double   Tables 2 and 3. From the figures above, it can be found that a good regression with a small fitting error for the experimental data can be gained using the double exponential model. Therefore, the degradation trend of battery capacity can be well predicted as long as the parameters of the double exponential model are accurately estimated. Therefore, the double exponential model is chosen as the degradation model of the battery in this paper.

| Parameter settings
The failure threshold of the capacity for the batteries is set as 75% of the rated capacity. The calculation of battery remaining useful life can be expressed by the following formula: Where RUL represents the remaining useful life. k EOL is the number of cycles when the battery failure threshold is reached. k EOM is the number of cycles which is at the end point of filter tracking and also the starting point of the prediction.
The method to determine the initial values of state parameter a, b, c, and d is as follows: import the experimental data about capacity degradation of battery 1 and battery 2 into the curve fitting tool in MATLAB and obtain the fitting values. Their averages are taken as the initial values for the parameters of the particle filter, which is used to do the remaining capacity prediction and remaining useful life prediction of battery 2. The fitting values are shown in the Table 2. The initial values for the parameters of the particle filter used to predict the remaining capacity and remaining useful life of battery 4 are determined by the same method, and the corresponding fitting values are shown in the Table 3.
It should be noted that the initial value is an approximate estimate which has deviation from the actual one, but it will gradually approach the actual value in the filtering process. It means the parameters of the model are constantly changing during the tracking iteration.
The process noise covariance is 1 × 10 − 6 × diag [0.01, 0.01, 0.01, 0.01]. The observation noise covariance R = 0.001. The state transition matrix Φ = diag [1, 1, 1, 1 ]. The prediction is classified into three groups, and the starting point k EOM is set as the 40th, 50th, and 60th cycle, respectively. The prediction results are compared with that based on extended Kalman filter under the same conditions.

| Prediction results
Kalman filtering algorithm is the best filter for discrete linear system but the noise is required to conform to gaussian distribution. Improved Kalman filter, such as extended Kalman filter, can be used to deal with nonlinear problems. Extended Kalman filter is a locally linearized Kalman filter. Taylor formula is used to expand the nonlinear equation by retaining the first two order terms and omitting the terms above the second order. The essence of extended Kalman filter is a linear filter, so the precision will decrease when dealing with complex nonlinear problems. Particle filter approximates the posterior probability distribution of the system by randomly scattering points. The initial particle distribution is completely random. With the emergence of new measurements, the value of the particles is corrected slowly. The particles with larger errors are eliminated and the particles with smaller error are replicated. The distribution represented by the particles will gradually be close to the true posterior distribution. Finally, the optimal estimation will be obtained. So even if the dynamic system is nonlinear, particle filter can still have good tracking performance. Figure 10 is the tracking curve of extended Kalman filter (EKF) and residual resampling particle filter (PF). In the figure, the black curve is the actual degradation curve of battery 2. The pink curve is the tracking curve of EKF, and the blue curve is the tracking curve of PF. It can be seen that the blue curve is closer to the black curve which means PF has better tracking performance than EKF.  Figures 11-13 show the prediction results of EKF and PF for the battery 2, respectively, when the starting point is selected as 40, 50, and 60. Figure 14 shows the prediction error comparison curves of two different algorithms at each prediction starting point. Figures 15-17 show the prediction results of EKF and PF for the battery 4, respectively, when the starting point is selected as 40, 50, and 60. Figure 18 shows the prediction error comparison curves of two different algorithms at each prediction starting point.

| Validity of residual resampling algorithm
With the knowledge of hierarchical statistics, the random variable U generates a random number obeying uniform distribution over [0, 1]. The interval in which the random number falls is called the response interval. According to this, it is determined to retain which particle. The essence of random resampling is to index particles with larger weight more to retain them and discard ones with smaller weight. The Figures 19-21 show the comparison of the prediction results for battery 2 using the particle filter with random resampling and residual resampling respectively. It can be seen that the curves obtained from the particle filter with residual resampling are closer to the true capacity degradation curves, which means the residual resampling is better than the random resampling in remaining capacity prediction in this paper. From the Table 4, it can be found that when the prediction starting point is the 40th and 60th cycle, the particle filter with residual resampling has a higher prediction accuracy of RUL than the particle filter with random resampling which are 1.79% and 2.78% percent higher, respectively. The results above show that the residual resampling algorithm is slightly better than random resampling in this study.

| Sensitivity analysis of prediction error
The effect of parameter initialization on the prediction error is researched in this part. When the prediction starting point is the 40th cycle and the remaining capacity of the battery 2 is predicted, by keeping b, c, and d unchanged, the change rule about the prediction error is studied by changing a with the step size of 0.00001. The results are shown as the picture at the top in the Figure 22: when a is selected as the fitting value (0.0002084), the root mean square prediction error is in the trough of the prediction error curve. Slight adjustment to a (near −0.0002084) can improve the accuracy in the process of capacity prediction and reduce the RMSE to a certain extent. However, due to the influence of each parameter initialization on the prediction error is interrelated, namely changing the initial value of a may make change the influence trend of other parameters on the prediction error so that all the parameters may not be relatively optimal at the same time. In addition, changing the initial value of a at some prediction starting point does not necessarily improve the prediction accuracy at another prediction starting point. Therefore, it is suitable that the initial value of a is set to −0.0002084 which is the fitting value because the prediction error is in the trough of the error curve and it is within the acceptable range. In addition, the influences of parameter a on the prediction error at different starting point are shown as the other pictures in the Figure 22. The influences of parameters b, c, and d on the prediction error when the starting point is 40 are shown in the Figure 23. From the figures above, it can be found that when the fitting values are used as the initial values for the parameters of the particle filter, the prediction errors are mostly in the trough of the error curves, which further indicates the validity of selecting the fitting values as all parameters' initial values of the particle filter. At this time, the prediction error can meet the demand of the prediction accuracy.

| Evaluation index
The following three evaluation indexes are used to comprehensively evaluate the accuracy of the prediction results:

Absolute error (AE)
The effect of absolute error is to reflect the true magnitude of the error between the predicted remaining useful life and the actual remaining useful life.

Relative error (RE)
The relative error reflects the reliability of the predictive results. The smaller the relative error is, the more reliable the result is.

Stability error (SE)
The stability error is calculated by the formula of rootmean-square error. The effect of the stability error is to reflect the stability of the prediction. The smaller the value is, the better the stability is. where Q k is the actual capacity at the kth cycle. Q k is the predicted capacity at the kth cycle. k EOL is the number of cycles when the capacity reaches the failure threshold. k EOM is the number of cycle where the prediction starts.

| Results analysis
Based on the capacity prediction curves and prediction error curves in Figures 11-18, qualitative comparison about the prediction accuracy of PF and EKF algorithms can be obtained. Take battery 2 for example, in the Figures 11-14, it can be seen that the prediction error curves of the PF are mostly below those of EKF and closer to the X axis in Figure 14. Also, the capacity prediction curves of PF are closer to the true capacity degradation curves than those of EKF in Figures 11-13. The conclusions above indicate remaining capacity prediction accuracy of PF is higher. The same conclusions can be drawn for battery 4 from the Figures 15-18. In addition, quantitative comparison about the prediction accuracy between PF and EKF can be gained with the evaluation indexes which are calculated as shown in Tables 5 and 6. It can be found that when the prediction starting point is selected as the 40th cycle, the 50th cycle and the 60th cycle, respectively, for battery 2, the remaining useful life prediction of PF is 3.58%, 2.17%, and 2.78% higher than that of EKF respectively; for battery 4, the F I G U R E 1 6 The comparison of prediction result between particle filter (PF) and extended Kalman filter (EKF) at the 50th cycle for battery 4 F I G U R E 1 7 The comparison of prediction result between particle filter (PF) and extended Kalman filter (EKF) at the 60th cycle for battery 4 | 1129 PAN et Al.
accuracy of PF is 0.93%, 2.07%, and 1.15% higher than EKF, respectively. In general, PF has higher prediction accuracy of remaining capacity and RUL than EKF.
Furthermore, from group 1 to group 3, as the starting points of the prediction move back, the absolute error, relative error, and stability error all decrease which indicates that the accuracy and stability of the prediction are improved. There are two reasons: 1. The backward movement of the starting point means that more actual capacity data can be used and then the final parameter estimates are more reliable; F I G U R E 1 8 The comparison of prediction error curves at different prediction starting points for battery 4 F I G U R E 1 9 Comparison of particle filter with residual resampling and random resampling respectively at the 40th cycle 2. The backward movement of the starting point also means that the range needed to be predict is reduced and then the corresponding uncertainty of prediction is decreased. So the accuracy and stability of the prediction are improved.
Whichever the group is, when the starting point of the prediction is the same, the absolute error, relative error and stability error of PF are smaller than those of EKF. F I G U R E 2 0 Comparison of particle filter with residual resampling and random resampling respectively at the 50th cycle F I G U R E 2 1 Comparison of particle filter with residual resampling and random resampling respectively at the 60th cycle