A DPCA-based online fault indicator for gear faults using three-direction vibration signals

2.1. Theory of PCA

PCA is a well-known and multi-variable feature extraction method based on 2-order statistics. It determines a set of orthogonal vectors, called loading vectors, which can be ordered by the amount of variance explained in the loading vector directions. By PCA, the high-dimensional correlated measurements are transformed into the low-dimensional uncorrelated variables.

Consider an $n$-by-$m$ history data matrix $\mathbf{X}$ which has been scaled to zero mean and unit variance to avoid problems associated with some measurement having large values and others with small ones, where $n$ is the sample size and $m$ is the number of the measured variables. By applying PCA, the matrix $\mathbf{X}$ is decomposed into the summation of cross products of $m$ score vectors ${\mathbf{t}}_i\in {\mathfrak{R}}^n$ and $m$ loading vectors ${\mathbf{p}}_i\in {\mathfrak{R}}^m$, that is:

where ${\mathbf{t}}_i$ and ${\mathbf{p}}_i$ are the $i$th score vector and loading vector respectively, and each ${\mathbf{p}}_i$ is a unit vector; $\mathbf{T}={\left[\begin{array}{llll}{\mathbf{t}}_1& {\mathbf{t}}_2& \cdots & {\mathbf{t}}_m\end{array}\right]}_{n\times m}$ is called the score matrix and $\mathbf{P}={\left[\begin{array}{llll}{\mathbf{p}}_1& {\mathbf{p}}_2& \cdots & {\mathbf{p}}_m\end{array}\right]}_{m\times m}$, an orthogonal transformation matrix, is called the loading matrix. Both score vectors and loading vectors are orthognals. If ${\mathbf{t}}_i$ is sorted in decrease according to its length, namely $‖{\mathbf{t}}_1‖>‖{\mathbf{t}}_2‖>\cdots ‖{\mathbf{t}}_m‖$, the first principal component explains the greatest amount of variability. The first loading vector ${\mathbf{p}}_1$ defines the direction of the greatest variance, ${\mathbf{p}}_2$ does the second greatest variance, and so on. Therefore, the full principal components decomposition of $\mathbf{X}$ can be given as:

If only the first $l$ ($l\ll m$) principal componets are kept, the Eq. (1) can be rewritten as:

where ${\mathbf{E}}_l$ is a matrix of residuals which represents the variability corresponding to process noise. ${\mathbf{P}}_l^{}$ composed of the first $l$ loading vectors is usually called the PCA model. Using ${\mathbf{P}}_l^{}$, an original $m$-dimensional space is transposed into two orthogonal subspaces: the $l$-dimensional principal subspace $\widehat{S}$ and the $(m-l)$-dimensional residual subspace $\tilde{S}$. Their relationship can be expressed as:

Once the PCA model is built via the singular value decomposition (SVD), a new measured vector ${\mathbf{x}}_{new,k}^{}$ at $k$th sampling time can be mapped into the principal subspace, and corresponding $1\times l$ score vector ${\mathbf{t}}_{new,k}^{}$ is given as following:

So, the estimated value of ${\mathbf{x}}_{new,k}^{}$, ${\widehat{\mathbf{x}}}_{new,k}$, is:

Therefore, the estimation error ${\mathbf{e}}_{new,k}$ is obtained by:

When monitoring if a system is normal by PCA based method, Hotelling $T^2$ statistic in principal subspaces and SPE statistic (also known as $Q$ statistic) in residual subspace are used. The $T^2$ statistic and the SPE statistic of the $k$th new sample are respectively represented as:

where $\lambda$ is an $l$-by-$l$ the diagonal matrix composed of the first $l$ eigenvalues of the matrix ${\mathbf{X}}^T\mathbf{X}/n$. $T^2$ statistic reflects changes of variables through fluctuation of principal component vector module of the PCA model. SPE is the sum of errors between the data and its estimates, so it describes the deviation degree of the measured value not to be captured by the PCA model at a certain time. In general, as long as one or both of the two statistics exceeds their own control limits, the system is out control. In other words, some fault occurs.

2.2. Dynamic principal component analysis

In general, static monitoring, such as PCA-based method, doesn’t function well for auto-correlated data. So, the dynamic PCA (DPCA) is developed to handle the process measurements with dynamic relationship [25]. The difference between DPCA and PCA is mainly reflected in the composition of modeling data. In DPCA, when the data are stacked with the current observation vector and the previous $d$ observations, the modeling data matrix is formed as:

where $\mathbf{x}\left(k\right)={\left[\begin{array}{llll}{x}_{1,k}& x_{2,k}& \cdots & x_{m,k}\end{array}\right]}^T$ is the $m$-dimensional observation vector at time point $k$; $\mathbf{x}(k-d)$ is the observation vector at time point $k$-$d$. The rest of the process of building the DPCA model is the same as that of PCA. The residuals of the DPCA model are much more uncorrelated than those of the traditional statistic PCA model. This means that DPCA is much better than the traditional static PCA in detecting the fault occurrence from serially correlated data.

3.1. Fault indicator

As the previous section described, SPE statistic reflects how well each sample conforms to the PCA model, and $T^2$ statistic is a measure of the variation of each sample within the PCA model. Although, preliminary research results indicated that SPE statistic, rather than $T^2$ statistic, could be used to monitor the rolling bear by vibration analysis, the types of faults cannot be identified by only these two statistics. Hence, aiming at directly monitoring and diagnosing faults of the gear, a new fault indicator is proposed based on two statistics. In addition, the vibration signals in each direction may contain information that reflects the fault characteristics. Therefore, in order to improve the accuracy of fault detection and diagnosis, a tri-axis acceleration sensor is used to collect radial-direction, axial-direction and rotational vibration signals denoted as $\left\{\mathbf{x}\left(t\right),\mathbf{y}\left(t\right),\mathbf{z}\left(t\right)\right\}$.

As we all know, the vibration signals of rotating parts, such as gear and bearing, are periodic and non-stationary, which leads to values of $T^2$ and SPE is not as smooth as ones of the stationary process, such as Tennessee Eastman (TE) process in chemical industry. According to the characteristics of gear vibration signals and the properties of these two statistics, the fault indicator of the $i$th sample, $F{I}_i$, can be defined as:

where, $n$ is the number of subgroups. That is, each sample is divided into $n$ equal parts. $T_i^2\left(j\right)$ and $SP{E}_i\left(j\right)$ are respectively the values of $T^2$and SPE of the $j$th subgroup in $i$th sample. It can be seen that the fault indicator describes the distribution characteristics of two statistics in the two-dimensional feature space and overcomes the fault cannot be effectively diagnosed by $T^2$and SPE because of overlapping of their range. For different faults, there are obvious differences in their distributions. Hence, a fault can be detected and diagnosed by monitoring FI.

The overall strategy of the proposed approach is summarized in Fig. 1. It mainly includes two phases: offline building and online applying.

3.2. Offline building DPCA model

As shown Fig.1, two major tasks need to be completed in offline modeling phase. One is to build the DPCA model which is developed to obtain two statistics as the features of the gear. The other is to determine the fault diagnosis thresholds which can be used to distinguish the fault types of the gear.

The main steps of offline building the DPCA model are summarized as follows:

Step 1: Acquire 3-D vibration signals $\left\{\mathbf{x}\left(t\right),\mathbf{y}\left(t\right),\mathbf{z}\left(t\right)\right\}$ when the gear without any faults is operated under the unload condition.

Step 2: Determine the size of a sample $m$ and the number of subgroups $n$. That is to say, each sample with $m$ points is equal divided $n$ subgroups without overlap. Each subgroup has $d=m/n$ points. Hence, a $j$th vector belonging to the $i$th sample is organized in the following form:

If there are $N$ samples, according to Eq. (10), an $nN$-by-3$d$ building data matrix ${\mathbf{X}}_d$ is given by:

Fig. 1The strategy of the proposed approach

Obviously, each row of ${\mathbf{X}}_d$ is the measurement variables with a time-lagged window $d-1$.

Step 3: Compute the sample mean vector ${\mathbf{x}}_d^{}$ and the sample standard deviation vector ${\mathbf{s}}_d^{}$. Then, normalize ${\mathbf{X}}_d$ to become a matrix ${\mathbf{X}}_d^{\mathrm{\text{'}}}$ with zero mean and unit variance using ${\mathbf{x}}_d^{}$ and ${\mathbf{s}}_d^{}$.

Step 4: Apply PCA to ${\mathbf{X}}_d^{\mathrm{\text{'}}}$. A DFCA model of the gear system, denoted as ${\mathbf{P}}_l$, is built after determining the number of principal componets $l$.

Step 5: Save the information of the DPCA model, including the loading matix ${\mathbf{P}}_l$, ${\mathbf{x}}_d^{}$ and ${\mathbf{s}}_d^{}$.

When there are clear distinctions among the different faults, the threshold based diagnosis method can be used to certain the fault types. So, the second work that needs to be done offline is to set the fault diagnosis thresholds. Three-sigma (3$\sigma$) quality level translates to a process yield of 99.73 % under the assumption of normality. That it is:

So, $\mu +3\sigma$ can be set as the fault warning threshold (FWT) to detect whether any failure occurs. The FWT can be computed based on the FIs of the healthy state:

where, $\mu$ and ${\sigma }^2$ are the mean and variance of the fault indicators of health state. If FIs of the vibration signals are greater than a fault warning threshold, it can say there is a fault.

3.3. Online fault monitoring and diagnosing via FI

Once the DPCA model and the thresholds are determined, the gear can be online monitored and diagnosed. The central idea in online fault monitoring and diagnosing using the FI is summarized as follows:

Step 1: Download the information of the DPCA model and the fault diagnosis thresholds.

Step 2: Acquire 3-D vibration signals in real time. Each sample has $m$ points.

Step 3: Obtain $n$ new subgrups according to Eq. (12) and scale respectively them with ${\mathbf{x}}_d^{}$ and ${\mathbf{s}}_d^{}$ which were obtained at offline modeling phase.

Step 4: Apply DPCA model on each scaled subgroups to obtain $T_i^2\left(j\right)$ and $SP{E}_i\left(j\right)$ ($j=\mathrm{1,2},\dots ,n$) according to Eqs. (8-9).

Step 5: Compute FI value of the current sample according to Eq. (11).

Step 6: Determine the gear state via the thresholds. Namely, detect and diagnose for the running gear.

Step 7: Return Step 2 and monitor the gear at the next new samples.

4.1. Fault simulation rig and 3D vibration data collection

In this section, to evaluate the performance of the proposed method, the experiments were carried out at a QPZZ-II rotating machinery fault simulation rig as shown in Fig. 2, which manufactured by Jiangsu Qianpeng Diagnosis Engineering Co., Ltd. The gear test part of the QPZZ-II rig is made of a gearbox, a 5 Nm (Rated Torque) and 0.5A (current) magnetic powder brake (MPB), a 0.55 kW power and 1450 rpm (Rated Speed) AC motor and a speed and load controller. Here, the controller is not shown in Fig. 2. The gearbox has a 55-tooth pinion (driving gear) and 75-tooth wheel (driven gear), which modules both are 2 mm. As shown in Fig. 3, the root-crack fault (RCF) and the one tooth-broken fault (TBF) were introduced into driven gear using the electron discharge machining (EDM) method. So, the driven gear has three conditions including healthy gear (HG).

Fig. 2QPZZ-II Rotating machinery fault simulation rig

The hardware of the data acquisition system mainly includes an 8-slot NI CompactDAQ Ethernet chassis (cDAQ-9188XT), a NI 9234 dynamic signal acquisition module, a B&K 4056 tri-axis acceleration sensor and a laptop. The data acquisition software was programmed by the National Instruments Labview.

The data sampling frequency was set as 5.12 kHz. The radial load for the gear box was respectively set as unload, 50 % load and full-load. The tri-axis acceleration sensor was placed on the top of the driven gear of the gearbox to acquire the gear’s radial-direction, axial-direction and rotational vibration signals.

Fig. 3Test gears with faults

a) One tooth broken

b) Tooth root crack

4.2. The results and analysis using fault indicator based on 3D vibration signals

Here, the time-domain waveforms of the raw three-direction vibration signals of all gear conditions under unload condition are respectively plotted in Figs. 4-6 when the rotation speed of AC motor was 1450 rpm.

Vibration signals of the healthy gear under unload condition at 1450 rpm are selected to build the DPCA model. Each 1024 sampling points are as a group. At the same time, a group is divided into 16 sub-groups in order to construct the modeling data set again. Obviously, each sub-group has 64 sampling points. There are 60 groups. Namely, $N=$60, $m=$1024, $n=$16 and $d=$64. So, a 960-by-192 building data matrix is obtained. After applying PCA, the original space with 192 dimensions is decomposed into a 52-dimensional principal subspace and a 140-dimensional residual subspace. After modeling, the fault warning threshold (FWT) is determined to be 1.1558 according to Eq. (15). If the fault indictor is greater than 1.1558, a fault occurs in the gear.

Fig. 4Three-direction vibration signals from a good gear

Fig. 53D vibration signals from a root crack gear

Fig. 63D vibration signals from a tooth broken gear

Fig. 7Fault indicator of testing samples under unloading at 1450 rpm

In order to test the proposed method affected by load changes, the testing samples of three gear states operated three different loading conditions at 1450 rpm are used. Fault indicators under unloading condition are calculated and plotted in Fig. 7. In this figure, lines with ‘o’, ‘*’ and ‘·’ indicate that the values of fault indicators come from a healthy gear, a gear with RCF and a gear with TBF, respectively. As shown in Fig. 7, faults are easily detected by the fault warning threshold (see the long dashed line in Fig. 7), and the accuracy reaches 100 %. To identify the types of faults, another threshold is determined based on the fault indicators of TBF since there is clear boundary between two faults. In this work, it is 2.5290. Table 1 lists the fault detection rates and fault recognition rates of testing data under all the loading conditions. As seen in Table 1, overall accuracy of fault detection is up to 100 %, and fault diagnosis rate of two faults is greater than 96 %. So, the proposed method is effective and robust to load variations.

Further, the vibration data at 1305 rpm, 1160 rpm, 1015 rpm, 870 rpm and 725 rpm are used to test the proposed method. Figs. 8-12 give the results at different driven speeds using the model built at 1450 rpm. Fig. 8 shows that the accuracy of gear health monitoring is up to 100 %, and TBF’s and RCF’s diagnose accuracies are respectively 92 % and 82 %. The result suggests that the fault of the gear can be diagnosed effectively by the fault indicator when the speed is reduced by 10 %. When the speed is reduced by 20 % and 30 %, although the failure types cannot be effectively diagnosed, as shown in Figs. 9-10, only 1 % of testing samples is not correctly monitored by FWT. The failure of gear can still be detected by resetting FWT (0.5498) when the speed drops to a lower level (see Figs. 11-12). So, the DPCA model at 1450 rpm can be used to monitor the health of gears at different speeds.

Fig. 8Fault indicator of testing samples at 1305 rpm

Fig. 9Fault indicator of testing samples at 1160 rpm

Fig. 10Fault indicator of testing samples at 1015 rpm

Fig. 11Fault indicator of testing samples at 870 rpm

Fig. 12Fault indicator of testing samples at 725 rpm

The proposed method is of less computation, and its average elapsed time for making fault detection only is 0.0047 second. In contrast, EMD spends 0.3402 second in processing a sample with 3×1024 points. And VMD needs 0.6387 second to get 6 modes. The proposed method greatly reduces the computation time comparing with EMD-based and VMD-based so that it can be applied to real-time fault detection.

4.3. Compared with traditional DPCA based method

For the purposed of comparison, the traditional control limit monitoring method using $T^2$ statistic and SPE statistic is used to determine whether faults are present. There are two cases to consider: one is to monitor the gear based on $T^2$ and SPE statistics of each subgroup with 64 points, the other is based on the statistics of each group with 1024 points. Detailed determining method of control limits of two statistics can be found in Ref. [22].

For the first case, SPE statistic and $T^2$ statistic of 1600 subgroups for all gear conditions are shown in Fig. 9. Here, 95 % and 99 % control limits of $T^2$statistic are respectively 82.7132 and 72.4265. 95 % and 99 % control limits of SPE statistic are respectively 44.0892 and 39.0683. As the two charts show, it can be observed that the $T^2$statistic cannot detect whether the gear is failure using the method of control limit because almost all $T^2$ values of the tooth broken gear and most $T^2$values of the crack fault gear both are within its control limits. At the same time, the SPE statistic is also not suitable to monitor the running state of gear because it will lead to a higher false alarm rate. As a consequence, monitoring gear states by mean of each subgroup’ two statistics doesn’t work well and it is invalid for gear fault detection and diagnosis.

In the second case, the SPE’s mean and the $T^2$’s mean value of each group including 16 subgroups replace its subgroup’s SPE and $T^2$ values, and the corresponding results are shown in Fig. 13.

Fig. 13Two statistics of the subgroups under unload condition at 1450 rpm

a) SPE statistic

b) Hotelling $T^2$statistic

Fig. 14 shows that the performance of fault detection can be improved significantly using the mean of SPE, and the $T^2$ can also monitor the gear running state. However, by the mean values of two statistics, fault reasons cannot be directly identified due to the similar distribution characteristics of the tooth broken fault and the root crack fault. In other words, the kinds of gear faults cannot be distinguished using SPE statistic. Therefore, the mean of SPE can effectively monitor the health of the gear, but cannot be used to identify the types of faults.

Fig. 14The mean values of two statistics of each group under unload condition at 1450 rpm

a) SPE statistic

b) Hotelling $T^2$statistic

4.4. Compared with the proposed method using single direction signals

Furthermore, in order to verify the superiority of FI based method using tri-axis vibration signals, the same building and testing data are used to construct 3 DPCA models only using each direction signals. All the FIs of testing samples are computed and the corresponding results are displayed in Fig. 15.

Fig. 15(a) shows that HG and TBF are completely separated, but there are overlaps between RCF and HG, between RCF and TBF. Under unload condition, the overlaps between RCF and HG is the most serious so that RCF doesn’t discern from HG. So, the monitoring and diagnosis of the gear cannot be realized only using radial-direction signals. In Fig. 15(b) and Fig. 15(c), the discrimination between HG and the other faults is obviously improved, but there are still serious overlaps between RCF and TBF. The above experimental results show the FIs based on only one direction vibration signals can’t be used to diagnose the gear faults.

Fig. 15The fault indicator of testing samples using one direction vibration signal at 1450 rpm

a) The results using radial-direction vibration signals

b) The results using axial-direction vibration signals

c) The results using rotational vibration signals