Time Series Analysis Using Composite Multiscale Entropy

Abstract

Multiscale entropy (MSE) was recently developed to evaluate the complexity of time series over different time scales. Although the MSE algorithm has been successfully applied in a number of different fields, it encounters a problem in that the statistical reliability of the sample entropy (SampEn) of a coarse-grained series is reduced as a time scale factor is increased. Therefore, in this paper, the concept of a composite multiscale entropy (CMSE) is introduced to overcome this difficulty. Simulation results on both white noise and 1/f noise show that the CMSE provides higher entropy reliablity than the MSE approach for large time scale factors. On real data analysis, both the MSE and CMSE are applied to extract features from fault bearing vibration signals. Experimental results demonstrate that the proposed CMSE-based feature extractor provides higher separability than the MSE-based feature extractor.

1. Introduction

Quantifying the amount of regularity for a time series is an essential task in understanding the behavior of a system. One of the most popular regularity measurements for a time series is the sample entropy (SampEn) [1] which is an unbiased estimator of the conditional probability that two similar sequences of m consecutive data points (m is the embedded dimension) will remain similar when one more consecutive point is included [2]. The SampEn characterizes complexity strictly on a time scale defined by the sampling procedure which is used to obtain the time series under evaluation. However, the long-term structures in the time series cannot be captured by SampEn. In regard to this disadvantage, Costa proposed the multiscale entropy (MSE) algorithm [3], which uses sample entropies (SampEns) of a time series at multiple scales to tackle this problem. The MSE has been successfully applied to different research fields in the past decades. These applications include the analyses of the human gait dynamics [2], heart rate variability [3,4], electroencephalogram [5], postural control [6], vibration of rotary machine [7,8], rainfall time series [9], time series of river flow [10], electroseismic time series [11], time series of traffic flow [12], social dynamics [13], chatter in the milling process [14], and vibrations of a vehicle [15], etc.. These works demonstrate the effectiveness of the MSE algorithm for the analysis of the complex time series.

The conventional MSE algorithm consists of two steps: (1) a coarse-graining procedure is used to derive the representations of a system’s dynamics at different time scales; (2) the SampEn algorithm is used to quantify the regularity of a coarse-grained time series at each time scale factor. To obtain a reasonable entropy value by using SampEn, the time series length is suggested to be in the range of 10^m to 30^m [16]. As reported in [2,5], in case of m = 2, the SampEn is significantly independent of the time series length when the number of data points is larger than 750. However, for a shorter time series, the variance of the entropy estimator grows very fast as the number of data points is reduced. In the MSE algorithm, for an N points time series, the length of the coarse-grained time series at a scale factor τ is equal to N /τ. The larger the scale factor is, the shorter the coarse-grained time series is. Therefore, the variance of the entropy of the coarse-grained series estimated by SampEn increases as a time scale factor increases. In many practical applications, the data length is often very short and the variance of estimated entropy values at large scale factors would become large. Large variance of estimated entropy values leads to the reduction of reliability in distinguishing time series generated by different systems. In order to reduce the variance of estimated entropy values at large scales, a composite multiscale entropy (CMSE) algorithm is proposed in this paper. The effectiveness of the CMSE algorithm is evaluated through two synthetic noise signals and a real vibration data set provided by Case Western Reserve University (CWRU) [17].

2. Methods

2.1. Multiscale Entropy

Essentially, the MSE is used to compute the corresponding SampEn over a sequence of scale factors. For an one-dimensional time series, $x=\{x_1,x_2,\mathrm{...},x_N\}$, the coarse-grained time series, $y^{(\tau )}$, can be constructed at a scale factor of τ, according to the following equation [3]:

$$y_j^{(\tau )}=\frac{1}{\tau }{\displaystyle \sum _{i=(j-1)\tau +1}^{j\tau }{x}_i},1\le j\le \frac{N}{\tau }$$

(1)

As shown in Figure 1, the coarse-grained time series is divided into non-overlapping windows of length τ, and the data points inside each window are averaged. We then define the entropy measurement of each coarse-grained time series as the MSE value. In this paper, the SampEn [1] is used as the entropy measurement. The algorithm proposed in [18] is repeated here, and we refer to the algorithm as shown in Figure 2. In the whole study of this paper, we calculate MSE values from scale 1 to scale 20 (τ = 1 to 20), and the sample entropy of each coarse grained time series is calculated with m = 2 and r = 0.15σ [2], where σ denotes the standard deviation (SD) of the original time series.

Figure 1. Schematic illustration of the coarse-grained procedure. Modified from reference [3].

Figure 1. Schematic illustration of the coarse-grained procedure. Modified from reference [3].

Figure 2. Brute force method. Modified from reference [18].

Figure 2. Brute force method. Modified from reference [18].

Most of the entropy measurements are dependent on the length of time-series. Since the length of each coarse-grained time series is equal to that of the original time series divided by the scale factor, τ, the variance of entropy measurements grows as the length of coarse-grained time series is reduced. The estimation error of a conventional MSE algorithm would be very large at large scale factors. In the following section, the modified MSE algorithm, named composite multiscale entropy (CMSE), is proposed to overcome this drawback.

2.2. Composite Multiscale Entropy

As shown in Figure 3, there are two and three coarse-grained time series divided from the original time series for scale factors of 2 and 3 respectively. The kth coarse-grained time series for a scale factor of τ, $y_k^{(\tau )}=\left\{\begin{array}{cccc}{y}_{k,1}^{(\tau )}& y_{k,2}^{(\tau )}& \cdots & y_{k,p}^{(\tau )}\end{array}\right\}$ is defined as

$$y_{k,j}^{(\tau )}=\frac{1}{\tau }{\displaystyle \sum _{i=(j-1)\tau +k}^{j\tau +k-1}{x}_i},1\le j\le \frac{N}{\tau },1\le k\le \tau$$

(2)

Figure 3. Schematic illustration of the CMSE procedure.

Figure 3. Schematic illustration of the CMSE procedure.

In the conventional MSE algorithm, for each scale, the MSE is computed by only using the first coarse-grained time series, $y_1^{(\tau )}$:

$$MSE(x,\tau ,m,r)=SampEn\left(y_1^{(\tau )},m,r\right)$$

(3)

In the CMSE algorithm, at a scale factor of τ, the sample entropies of all coarse-grained time series are calculated and the CMSE value is defined as the means of τ entropy values. That is:

$$CMSE(x,\tau ,m,r)=\frac{1}{\tau }{\displaystyle \sum _{k=1}^{\tau }SampEn\left(y_k^{(\tau )},m,r\right)}$$

(4)

Figure 4 shows the flow charts of the MSE and CMSE algorithms for comparison. The Matlab code of CMSE is shown in Appendix A.

Figure 4. Flow charts of MSE and CMSE algorithms.

Figure 4. Flow charts of MSE and CMSE algorithms.

3. Comparative Study of MSE and CMSE

To evaluate the effectiveness of the CMSE, two synthetic noise signals, white and 1/f noises, and a real vibration data set were applied in comparison with that of the MSE in this section.

3.1. White Noise and 1/f Noise

The SampEn for coarse-grained white noise time series is computed by [4]:

$$SampEn=-\ln \left\{{\displaystyle {\int }_{-\infty }^{\infty }\frac{1}{2\sigma }\sqrt{\frac{\tau }{2\pi }}\left[erf\left(\frac{x+r}{\sigma \sqrt{2/\tau }}\right)-erf\left(\frac{x-r}{\sigma \sqrt{2/\tau }}\right)\right]\exp \left(\frac{-x^2\tau }{2{\sigma }^2}\right)dx}\right\}$$

(5)

where erf refers to the error function [4]. For 1/f noise, the analytic value of SampEn is 1.8 with N = 30,000 in all scales [4]. In order to further investigate the effect of different data lengths on the MSE and CMSE, we first test the MSE on simulated white noises with different data lengths. As shown in Figure 5a, for short time series, the estimated MSE values are significantly different from the analytic solutions. This significant error may reduce the reliability in distinguishing time series generated by different systems. We then applied the MSE on simulated 1/f noises with different data lengths. As shown in Figure 5b, the variance of the entropy estimator increases with the reduced data lengths and difference between analytic solutions and numerical solutions existing in all scales. Figure 6a,b show the entropies of white noise and 1/f noise by applying CMSE, respectively. In comparison with the estimation by the MSE, the variance of the entropy estimator can be improved by the CMSE evidently. However, the over estimation due to the shortage of the data length still exists when CMSE is applied to the 1/f noise.

Figure 5. MSE results of (a) white noise and (b) 1/f noise with different data lengths.

Figure 5. MSE results of (a) white noise and (b) 1/f noise with different data lengths.

Figure 6. CMSE results of (a) white noise and (b) 1/f noise with different data lengths.

Figure 6. CMSE results of (a) white noise and (b) 1/f noise with different data lengths.

The numerical results of white noise with two different data lengths (N = 2,000 and 10,000) are shown in Figure 7. The error bar at each scale indicates the SD of an entropy value which calculated 100 independent noise signals. For a scale factor of one, the MSE value is equal to the CMSE value because the coarse-grained time series is the same as the original time series. In all cases, the means of the entropy values have no significant difference between the MSE and CMSE. However, the SDs of the entropy values between the MSE and CMSE are different. For a longer length of white noise (N = 10,000, Figure 7b), the SD of the CMSE is slightly less than that of the MSE. For a shorter length of white noise (N = 2,000, Figure 7a), at a large scale, the SD of the CMSE can be reduced greatly. For instance, in the case of white noise with 2,000 data points, the SD of the MSE at a scale factor of 20 is 0.1033 while the SD of the CMSE is only 0.0658. Figure 8a,b show the results of the MSE and CMSE applied to the 1/f noise with 2,000 and 10,000 data points, respectively. The result of 1/f noise is similar to that of white noise; the CMSE can reduce the SDs of estimations.

Figure 7. MSE and CMSE results of white noise with data lengths (a) N = 2,000 and (b) N = 10,000.

Figure 7. MSE and CMSE results of white noise with data lengths (a) N = 2,000 and (b) N = 10,000.

Figure 8. MSE and CMSE results of 1/f noise with data lengths (a) N = 2,000 and (b) N = 10,000

Figure 8. MSE and CMSE results of 1/f noise with data lengths (a) N = 2,000 and (b) N = 10,000

Table 1 summarizes the SDs of the MSE and CMSE at different time scales. These results indicate that the entropy values calculated by the conventional MSE and CMSE algorithms are almost the same, but the CMSE can estimate entropy values more accurate than the MSE. This improvement is significant when the CMSE is utilized to analyze the time series with short data length.

Table 1. Standard deviations of the MSE and CMSE at different time scales.

Data Length	Signals	Methods	Scales
			1	3	5	7	9	11	13	15	17	19	20
2,000	white noise	MSE	0.026	0.049	0.059	0.064	0.072	0.067	0.084	0.076	0.091	0.091	0.103
		CMSE	0.026	0.030	0.037	0.043	0.046	0.051	0.054	0.054	0.057	0.063	0.066
	1/f noise	MSE	0.084	0.097	0.121	0.139	0.162	0.227	0.236	0.284	0.265	0.282	0.310
		CMSE	0.084	0.088	0.092	0.099	0.108	0.126	0.122	0.148	0.153	0.155	0.163
10,000	white noise	MSE	0.007	0.014	0.019	0.025	0.028	0.030	0.032	0.035	0.038	0.035	0.035
		CMSE	0.007	0.010	0.015	0.018	0.020	0.021	0.023	0.025	0.025	0.027	0.028
	1/f noise	MSE	0.069	0.069	0.070	0.073	0.072	0.071	0.078	0.074	0.086	0.085	0.080
		CMSE	0.069	0.068	0.068	0.067	0.069	0.069	0.068	0.069	0.072	0.072	0.070

Table 1. Standard deviations of the MSE and CMSE at different time scales.

Data Length	Signals	Methods	Scales
			1	3	5	7	9	11	13	15	17	19	20
2,000	white noise	MSE	0.026	0.049	0.059	0.064	0.072	0.067	0.084	0.076	0.091	0.091	0.103
		CMSE	0.026	0.030	0.037	0.043	0.046	0.051	0.054	0.054	0.057	0.063	0.066
	1/f noise	MSE	0.084	0.097	0.121	0.139	0.162	0.227	0.236	0.284	0.265	0.282	0.310
		CMSE	0.084	0.088	0.092	0.099	0.108	0.126	0.122	0.148	0.153	0.155	0.163
10,000	white noise	MSE	0.007	0.014	0.019	0.025	0.028	0.030	0.032	0.035	0.038	0.035	0.035
		CMSE	0.007	0.010	0.015	0.018	0.020	0.021	0.023	0.025	0.025	0.027	0.028
	1/f noise	MSE	0.069	0.069	0.070	0.073	0.072	0.071	0.078	0.074	0.086	0.085	0.080
		CMSE	0.069	0.068	0.068	0.067	0.069	0.069	0.068	0.069	0.072	0.072	0.070

3.2. Real Vibration Data

In order to validate the utility of the CMSE algorithm for real data, experimental analysis on bearing faults is carried out. All the bearing fault data used in this paper are obtained from the Case Western Reserve University (CWRU) Bearing Data Center [17]. The test stand is composed of a 2-horsepower motor and a dynamometer, which are connected by a torque transducer. The test bearings using electro-discharge machining with fault diameters of 7, 14, and 21 mils (1 mil is one thousandth of an inch) are used to detect single point faults. Bearing conditions of the experiments include normal states, ball faults, inner race faults and outer race faults located at 3, 6, and 12 o’clock positions which are at 0°, 270°, and 90° on the front section diagram of the bearing, respectively. In other words, from the cross section diagram of the bearing, the 3 o’clock position is parallel to the direction of the load zone; 6 and 12 o’clock positions are perpendicular to the load zone. Vibration data are collected by accelerometers which are placed at the 12 o’clock position at both the drive end and fan end of the motor housing. Digital data are collected at a sampling rate of 48,000 samples per second for drive end bearing experiments. The motor speeds controlled by motor load are set to be 1,730, 1,750, and 1,772 rpm.

In the experiments, the vibration signals were divided into several non-overlapping segments with a specified data length, N = 2,000. Each non-overlapping segment was regarded as one sample in the validation process. The numbers of samples for each bearing condition are listed in Table 2. Each sample is a time series with 2,000 data points. We then calculated the MSE and CMSE values up to scale 20 for each sample. Therefore, the dimension of sample in the feature space is 20 in the following experiments. Partial measured acceleration signals of vibrations at the six different conditions are shown in Figure 9. The MSE and CMSE of bearing data in specific condition is shown in Figure 10. For each condition, the means of the entropy estimator obtained by the CMSE are very similar to those obtained by the MSE, while less SDs are achieved by the CMSE. This consists with the analysis results of synthetic noise signals. The collected data in Figure 10 are dependent on the neighbor states sampled in the experiments.

Table 2. Numbers of data sets are corresponding to different faulted classes, defective levels and rotation speeds.

Shaft Speed / Defect Level	Rotation Speed (rpm)
	1730			1750			1772
	Fault diameter (mils)
Fault Class	7	14	21	7	14	21	7	14	21
Normal state	243			242			242
Ball	244	243	243	243	243	243	243	243	243
Inner race fault	243	242	244	243	244	245	243	191	242
Outer race fault (3)	243		242	243		243	242		245
Outer race fault (6)	244	244	244	243	243	244	243	242	244
Outer race fault (12)	242		243	241		243	241		243

Table 2. Numbers of data sets are corresponding to different faulted classes, defective levels and rotation speeds.

Shaft Speed / Defect Level	Rotation Speed (rpm)
	1730			1750			1772
	Fault diameter (mils)
Fault Class	7	14	21	7	14	21	7	14	21
Normal state	243			242			242
Ball	244	243	243	243	243	243	243	243	243
Inner race fault	243	242	244	243	244	245	243	191	242
Outer race fault (3)	243		242	243		243	242		245
Outer race fault (6)	244	244	244	243	243	244	243	242	244
Outer race fault (12)	242		243	241		243	241		243

Figure 9. Measured acceleration signals of vibrations in the time domain of six different bearing conditions (a) normal state, ball fault and inner race fault; (b) outer race faults at 3, 6, and 12 o’clock positions.

Figure 9. Measured acceleration signals of vibrations in the time domain of six different bearing conditions (a) normal state, ball fault and inner race fault; (b) outer race faults at 3, 6, and 12 o’clock positions.

Figure 10. MSE and CMSE results on bearing vibration data (1,730 rpm, 7 mils). (a) Normal state. (b) Ball fault. (c) Inner race fault. (d) Outer race fault (3 o’clock position). (e) Outer race fault (6 o’clock position). (f) Outer race fault (12 o’clock position)

Figure 10. MSE and CMSE results on bearing vibration data (1,730 rpm, 7 mils). (a) Normal state. (b) Ball fault. (c) Inner race fault. (d) Outer race fault (3 o’clock position). (e) Outer race fault (6 o’clock position). (f) Outer race fault (12 o’clock position)

3.3. Performance Assessment

In this subsection, we first use Mahalanobis distance to assess the effectiveness of the MSE and CMSE methods. Mahalanobis distance [19] is a popular method to measure the separation of two groups of samples. Assuming two groups with mean vectors ${\overline{s}}_1$ and ${\overline{s}}_2$, Mahalanobis distance is shown as the following equation [19]:

$$d=\sqrt{{({\overline{s}}_1-{\overline{s}}_2)}^T{C}^{-1}({\overline{s}}_1-{\overline{s}}_2)}$$

(6)

The pooled variance-covariance matrix C in equation (6) is shown below [19]:

$$C=\frac{1}{n_1+n_2-2}(n_1{C}_1+n_2{C}_2)$$

(7)

where $n_i$ is the number of samples of group i and $C_i$ is the covariance matrix of group i.

Table 3 shows Mahalanobis distances for six different distinguishing conditions. The motor speed was set at 1,730 rpm and the fault diameter was 7 mils. The faults are normal state (N), ball fault (B), inner race fault (I), outer race defects at 3, 6, and 12 o’clock positions (O3, O6, and O12). Larger Mahalanobis distance in the table represents the higher level of linear separability for two different groups [19]. Comparing N with B (N in fault class 1 and B in fault class 2 or N in fault class 2 and B in fault class 1), both values of the CMSE and MSE are high and the value of the CMSE is higher than that of MSE, indicating that the normal state can easily be distinguished from ball fault and the CMSE has the higher distinguishability. Therefore, it is obvious that the normal state can easily be distinguished from ball fault, inner race fault, and outer race faults located at 3 and 6 o’clock positions, but not easily distinguished from the outer race fault located at 12 o’clock position due to the smaller value in Table 3. The most indistinguishable conditions are (1) between ball fault and inner race fault and (2) between outer race faults located at 3 and 6 o’clock positions. In addition, in all cases, Mahalanobis distance of two different groups of features extracted by the CMSE algorithm is larger than that extracted by the MSE algorithm. Therefore, compared with the MSE, the CMSE as a feature extractor can have the higher distinguishability.

3.4. Fault Diagnosis Using an Artificial Neural Network

We built a fault diagnosis system based on a neural network and used the CMSE as a feature extractor contrasting with MSE. The aforementioned quantities in twenty scales were selected as the features for bearing fault diagnosis. The training of a neural network with bearing vibration data was performed by the MATLAB Neural Networks Toolbox V6.0.2. A three-layer backpropagation neural network was trained by the Levenberg–Marquardt algorithm [20]. The network had 20 nodes in the input and visible layers (each node corresponding to a scale of MSE and CMSE), 30 nodes in the hidden layer, and 4 or 6 nodes in the output layer dependence on how many fault conditions were considered. For training, a target mean square error of 0, a learning rate of 0.001, a minimum gradient of 10⁻¹⁰ and maximum iteration number of 1,000 were used. To improve generalization, the data sets were randomly divided by three parts: (1) training (50%), (2) validation (15%), and (3) testing (35%). The average accuracy of prediction for each experiment was quantified over 200 tests.

Table 3. Mahalanobis distances for six different distinguishing conditions.

Fault Class 1	Feature Extractor	Fault class 2
		N	B	I	O3	O6	O12
N	MSE		25.200	25.200	20.044	25.484	5.538
	CMSE		25.898	28.069	21.614	26.115	6.992
B	MSE	25.200		3.675	5.312	5.515	9.657
	CMSE	25.898		5.852	8.105	7.534	11.387
I	MSE	25.200	3.675		7.128	7.560	12.883
	CMSE	28.069	5.852		8.672	9.652	16.965
O3	MSE	20.044	5.312	7.128		5.934	9.303
	CMSE	21.614	8.105	8.672		6.605	13.239
O6	MSE	25.484	5.515	7.560	5.934		9.624
	CMSE	26.115	7.534	9.652	6.605		12.139
O12	MSE	5.538	9.657	12.883	9.303	9.624
	CMSE	6.992	11.387	16.965	13.239	12.139

Table 3. Mahalanobis distances for six different distinguishing conditions.

Fault Class 1	Feature Extractor	Fault class 2
		N	B	I	O3	O6	O12
N	MSE		25.200	25.200	20.044	25.484	5.538
	CMSE		25.898	28.069	21.614	26.115	6.992
B	MSE	25.200		3.675	5.312	5.515	9.657
	CMSE	25.898		5.852	8.105	7.534	11.387
I	MSE	25.200	3.675		7.128	7.560	12.883
	CMSE	28.069	5.852		8.672	9.652	16.965
O3	MSE	20.044	5.312	7.128		5.934	9.303
	CMSE	21.614	8.105	8.672		6.605	13.239
O6	MSE	25.484	5.515	7.560	5.934		9.624
	CMSE	26.115	7.534	9.652	6.605		12.139
O12	MSE	5.538	9.657	12.883	9.303	9.624
	CMSE	6.992	11.387	16.965	13.239	12.139

In this paper, we conducted nine experiments with a single operation speed and a single fault diameter. Table 4 lists the diagnostic accuracy results using the MSE and CMSE of each experiment. The experiments with the fault diameter of 21 mils have lower diagnosis accuracy than others. However, in these indistinguishable cases, the improvement by using the CMSE as a feature extractor is more obvious. Therefore, it can be inferred that the accuracy of bearing fault diagnosis can be enhanced by the CMSE. Furthermore, although CMSE is only applied to the univariate time series, it also can be applied to multivariate time series [21]. The proposed CMSE algorithm is for SampEn in this research. It is also for permutation entropy while the multiscale analysis is performed [22,23].

Table 4. Diagnostic accuracy results using MSE and CMSE.

Fault Diameter / Feature Extractor	Rotation Speed (rpm)
	1730			1750			1772
	7	14	21	7	14	21	7	14	21
MSE	97.33%	98.81%	96.77%	99.05%	98.01%	95.65%	99.34%	96.67%	95.89%
CMSE	99.29%	99.75%	98.26%	99.58%	99.86%	98.50%	99.91%	99.65%	98.42%

Table 4. Diagnostic accuracy results using MSE and CMSE.

Fault Diameter / Feature Extractor	Rotation Speed (rpm)
	1730			1750			1772
	7	14	21	7	14	21	7	14	21
MSE	97.33%	98.81%	96.77%	99.05%	98.01%	95.65%	99.34%	96.67%	95.89%
CMSE	99.29%	99.75%	98.26%	99.58%	99.86%	98.50%	99.91%	99.65%	98.42%

5. Conclusions

In this paper, the concept of CMSE is introduced for the analysis of the complexity of a time series. The proposed method presents better performance on short time series than the MSE. For the analysis of white noise and 1/f noise, simulation results show that the CMSE provides a more reliable estimation of entropy than the MSE. In addition, for the CMSE as the feature extractor of the bearing fault diagnosis system and the Mahalanobis distance used as a performance assessment, the simulation results show that the CMSE can enhance the linear distinguishability in comparison with the MSE. Experimental results also demonstrate that the proposed CMSE provides a higher accuracy of bearing fault diagnosis.