Fault Detection for CNC Machine Tools Using Auto-Associative Kernel Regression Based on Empirical Mode Decomposition

: In manufacturing processes using computerized numerical control (CNC) machines, machine tools are operated repeatedly for a long period for machining hard and difﬁcult-to-machine materials, such as stainless steel. These operating conditions frequently result in tool breakage. The failure of machine tools signiﬁcantly degrades the product quality and efﬁciency of the target process. To solve these problems, various studies have been conducted for detecting faults in machine tools. However, the most related studies used only the univariate signal obtained from CNC machines. The fault-detection methods using univariate signals have a limitation in that multivariate models cannot be applied. This can restrict in performance improvement of the fault detection. To address this problem, we employed empirical mode decomposition to construct a multivariate dataset from the univariate signal. Subsequently, auto-associative kernel regression was used to detect faults in the machine tool. To verify the proposed method, we obtained a univariate current signal measured from the machining center in an actual industrial plant. The experimental results demonstrate that the proposed method successfully detects faults in the actual machine tools.


Introduction
Computerized numerical control (CNC) refers to a method of automating the control of machine tools by inputting programmed milling information into a microcomputer without a manual operator.The programmed milling specification is stored in the memory of the computer to process large amounts of work efficiently.Furthermore, by flexibly controlling various milling conditions (rotating speed, cutting force, etc.), high-quality material products can be produced at low costs.However, when the milling process is running, machine tools frequently suffer from faults, because these are operated in extreme environments to cut materials [1].
A fault is defined as an unpermitted deviation of at least one characteristic property of a variable from an acceptable behavior [2,3].Although CNC machines have been well developed, they lack a function that diagnoses the condition of the tools and replaces faulty tools with new ones.The failure of machine tools requires maintenance time (also known as downtime) to replace the machine tool by stopping the machine.This may reduce the quality of the products and efficiency of the target process and increases the maintenance cost.In addition, if the damaged tool is used continuously to process the material, it may cause severe physical damage owing to the failure of the CNC machine.In fact, 79.6% of the maintenance time for CNC machines in modern industry is observed to be caused by damage to machine tools [4].To solve these problems in the industrial field, when the life Multivariate statistics for detecting abnormal scenarios have been studied intensively in the research on multivariate quality control [23].Principal component analysis (PCA), independent component analysis (ICA), and AAKR have commonly been used for multivariate statistics.PCA and ICA can efficiently manage multivariate process data through a dimensionality reduction.However, these algorithms have a limitation in terms of the distribution of the latent variables.For example, PCA assumes that the hidden variables follow a Gaussian distribution.These statistical assumptions can degrade model performance because of the difficulty in satisfying actual industrial process data.In contrast, AAKR can effectively detect failures in various industrial processes due to no restriction on the abovementioned assumptions.
AAKR is a nonparametric multivariate regression method (also known as lazy, memorybased, and instance-based learning) that compares the similarity between the training data stored in memory and query vector.Subsequently, high weights are assigned to the training vectors with high similarity to compute the estimated vector.AAKR is updated online using a local model; therefore, extensive studies were conducted in [24][25][26][27][28][29] online monitoring of multivariate processes because of the advantage that it can be applied to time-varying processes operating in different process modes.
Huang proposed EMD in 1998 for decomposing the unique features contained in a univariate signal into amplitude modulation (AM)/frequency modulation (FM) components [30].This method is not based on any assumption regarding linearity or stationary signals.Thus, many studies have been conducted to analyze nonlinear and non-stationary signals, such as biomedical signal [31,32] and industrial process data [33,34].EMD is compared with the short-time Fourier transform (STFT) and wavelet transform (WT), a classical time-frequency analysis method.STFT provides a constant resolution for the entire signal through a fixed window.However, the resolution has limitations (also known as the uncertainty principle) owing to the trade-off between frequency and time.Therefore, this method is suitable for quasi-stationary signals.WT can extract the time-frequency characteristics of a signal through expansion and transformation.When using WT, multiscale signal analysis is possible, and it is more suitable for non-stationary signals than STFT.Nevertheless, WT has a critical disadvantage in that its analysis results depend on the selection of the wavelet basis function.In other words, only the signal characteristics correlated with the shape of the selected wavelet function generate high-value coefficients, and other characteristics may be masked or completely ignored.In contrast to WT, EMD is a self-adaptive signal processing method [35].It can decompose the original signal into IMFs that reflect the vibrational mode features of the original signal without requiring a predefined basis function for signal analysis.Therefore, EMD is more suitable for analyzing industrial process data, which are nonlinear and non-stationary signals, than STFT and WT.

Preliminary
Most of the milling processes using a CNC machine have a smaller scale than chemical processes and power plants, and obtaining multivariate process data by installing a large number of sensors is difficult.EMD can extract some IMFs that can obtain multivariate data from univariate signals in the same time domain.This can address the inapplicability of multivariate fault detection models in the environment where only one signal is measured from the target system.In other words, multivariate statistical methods can be applied to univariate signals using EMD.In practice, PCA [36,37] and ICA [38,39] have been employed with EMD for process monitoring.In this study, EMD was used to extract IMFs from univariate signals measured using an actual CNC machine.Subsequently, we conducted an AAKR-based fault detection by constructing a multivariate dataset.
To demonstrate the efficiency of the proposed method, we used the current signal obtained from the long-term operation of an actual machine tool.Several tools of the same type were repeatedly utilized for straight and spiral cutting processes to obtain normal and fault data.The straight-cutting process provided tool breakage datasets.Artificial fault data were generated by applying disturbances (bias and drift) because the tools used in the spiral cutting process did not break.The proposed method was applied to the following two types of process datasets mentioned above: straight cutting process (actual tool breakage data) and spiral cutting process (artificial fault data through bias and drift).The experimental results demonstrated that the proposed method can effectively detect faults in a machine tool for two types of cutting processes.The contributions of this study are summarized as follows: 1.
Existing signal processing and deep learning methods used for detecting faults in machine tools require the characteristics of the machining process and massive data, respectively.Furthermore, when deep learning is used, the features of the original signal can be lost during the process of converting a signal into an image.In contrast, the proposed method can efficiently perform tool-state diagnosis by avoiding these problems.

2.
By using EMD, the proposed method has the advantage of applying a multivariate fault-detection model whose performance has been verified through prior related studies, even in a limited environment where only univariate signals can be acquired.

3.
AAKR was employed to detect the fault in a machine tool for the first time.In addition, it has never been used in combination with EMD to detect faults in various industrial processes.4.
To obtain the actual machine tool data, we repeatedly conducted some experiments through straight parallel and spiral circular cutting, and then the proposed method was validated using massive data obtained from an actual operating CNC machine.
The remainder of this paper is organized as follows.Section 2 describes the proposed method for fault detection of a machine tool.The target systems and datasets used in this study are described in Section 3. The experimental results and discussion are presented in Section 4. Finally, Section 5 presents our conclusions.

Fault Detection of Machine Tools Using EMD and AAKR
Figure 1 shows the overall procedure for the fault detection of machine tools using EMD and AAKR.First, IMFs are extracted from the univariate current signal measured in a CNC machine using EMD.In this study, the number of IMFs was determined by repeatedly conducting experiments (trial and error) with the best AAKR performance.Subsequently, AAKR is trained using a dataset consisting of two extracted IMFs and original current signal.The first step in AAKR is to calculate the similarity between training and query vectors.Among the various distance functions for computing similarity, the most used Euclidean distance function was selected in this study.Next, a Gaussian weight function is used to assign weights to the training data vector.In this step, the appropriate value of the bandwidth parameter h in the Gaussian weighting function is determined using k-fold cross-validation.Subsequently, we calculate the estimated and residual vectors, and then, to obtain the detection indices, SPE, the multivariate statistic, is computed as the square of the residual vector.The confidence limit (also called the threshold and control limit) for declaring the fault was determined using kernel density estimation (KDE) in this study.More details of the AAKR algorithm for fault detection are presented in Figure 2 and Section 3.2.

Empirical Mode Decomposition
The EMD technique was proposed by Huang [30] in 1998 and can decompose a signal into an IMF.As a nonlinear, multiresolution, self-adaptive decomposition technique, EMD has been an effective signal processing tool for signals emanating from nonlinear and nonstationary systems [39].The IMFs represent the natural oscillatory mode embedded in the signal and function as the basis functions, which are determined by the signal itself, rather than the pre-determined kernels [35].IMFs must satisfy the following two conditions [30]: (i) in the entire dataset, the number of extrema and number of zero crossings must either be equal or differ at most by one, and (ii) at any point, the mean value of the envelope defined by the local maxima and envelope defined by the local minima is zero.EMD is based on the assumption that any signal contains different simple IMF components.With this definition, a signal x(t) can be decomposed as follows [30]:

•
Step 1. Identify all local extrema of x(t).

•
Step 3. Verify that a 1 is correct for the IMF conditions.Ideally, if h 1 is an IMF, h 1 is the first component of x(t).

•
Step 4. If h 1 is not an IMF, h 1 is considered the input signal, and repeat Steps 1-3, then After repeated this shifting procedure k times, until h 1k is an IMF, that is then h 1k is designated as the first IMF component from the original signal.

•
Step 5. Separate c 1 from the remainder of the data by where r 1 is considered new data; the above processes are resumed, and the second IMF component c 2 can be obtained.We can repeat the above steps to achieve n-IMFs of signal x(t), as follows: The shifting process can be stopped by any of the following predefined stopping criteria: (i) the component, c n , or residue, r n , becomes so small that it is less than the predetermined value, and (ii) the residue r n becomes a monotonic function when no more IMF can be extracted.By summing Equations ( 5) and ( 6), we obtain Finally, EMD is completed, and we can decompose the original signal into n empirical modes.A residue, r n , can be either a mean trend or a constant of signal x(t).A more detailed explanation of the EMD is available in [30,35].

Fault-Detection Method Based on AAKR
Figure 2 shows the AAKR procedure used in this study.The fault-detection procedure using AAKR can be roughly divided into offline and online processes for training and testing the algorithm, respectively.In the offline process, we applied z-score normalization to the multivariate training dataset constructed by EMD to set the average and variance of each variable to 0 and 1, respectively.The bandwidth parameter of the weighting function was selected using k-fold cross-validation, where we used the Gaussian weighting function as the kernel function.The residual vector for the training data was then calculated using the leave-one-out method.The confidence limit (also known as the threshold value) for identifying faults was determined in advance.In this study, the threshold value for the fault declaration was obtained using the KDE.In the online process, similarity and weight vectors were computed for the query vector (testing data), and the estimated and residual vectors were computed based on them.Fault detection was then conducted by calculating multivariate statistics and comparing them with a pre-defined threshold in an offline process.A detailed description of fault detection is as following subsections.

AAKR Algorithm
AAKR is a multivariate regression method that stores training data in memory and assigns high weights to vectors having high similarity with the test vector to compute the score [40].When multivariate data are collected from the target system, data matrix X is constructed as follows: where j and p are the number of samples and process variables, respectively.In AAKR, if the training data matrix X is large, the algorithm requires longer training time.However, after the completion of the training phase, the appropriate weights were assigned to training vectors and AAKR can distinguish between normal and abnormal for query vectors.Generally, each variable of the multivariate data measured in the industrial process has a different mean and standard deviation.Therefore, if this dataset is trained on the model without data normalization, it can result in biased training results.For example, any variable with a large range of values is not important but can be considered more important than the other variables.To avoid this problem, data normalization is essential for setting equal mean and variance of all variables.In this study, we employed the z-score normalization technique to set the mean and variance of each variable to 0 and 1, respectively.The z-score normalization is defined as follows: where µ m and σ m are the average and variance corresponding to the mth variable of the training dataset, respectively.In the AAKR algorithm, the similarity between the training and query vectors is calculated using the distance function to compute the estimated vector for the query vector.Several distance functions (such as the Euclidean distance and the Mahalanobis distance) are available for measuring the distance.To select the proper distance function, we compared the performance of the AAKR by applying Euclidean distance and Mahalanobis distance.Mahalanobis distance function required longer than Euclidean distance function.In addition, application of Euclidean distance improved performance compared with Mahalanobis distance.In this study, therefore, we employed the Euclidean distance defined in (10).
where x(t) is the query vector for time t, and x l is the training vector.If the training and query vector are highly similar, the value of the distance function reduces.The magnitude of the weight for each data vector is determined according to the degree of similarity measured using the distance function.In this study, we employed a Gaussian weighting function to assign weights as follows: As the value of similarity increases, K h assigns less weight to training data x i .In (11), h is a bandwidth parameter of the weighting function.In AAKR, multivariate regression performance depends on the bandwidth parameter value [41].If the value of parameter h is too small, then we obtain crude estimation results because a small amount of data is used to estimate the query vector.However, if the value of parameter h is too large, we obtain smooth estimation results.The parameter value is typically set by following a trial-and-error procedure for validation data [42].For example, the objective function was used to select h with minimal error by varying the values h [43,44].In this study, the optimal bandwidth parameter was determined using k-fold cross-validation, which is a method of repeatedly performing training and validation to prevent bias in the training data.This method is widely used to tune the hyperparameters of the models with inherent uncertainties.A more detailed explanation of parameter determination is available in [40].In this study, the mean squared prediction error (MSPE) was calculated to determine the optimal parameter value as follows: where q is the number of repetitions up to 1 + (h max − h min )/∆h, and X s is a subset of k training data randomly divided by k without duplication.h max , h min , and ∆h are the maximum, minimum, and incremental values of h, respectively.Thus, we determine the value of bandwidth parameter h with the minimum MSPE as the optimal parameter value.When the proper parameter value is determined, the estimated vector x(t), and the residual vector, e(t), are calculated as follows: where K h (d j ) is the weight computed by the weighting function, and x(t) is the query vector at time t.The estimated vector, x(t), is obtained using the weighted average.The residual vector, e(t), can be obtained by calculating the difference between the query vector, x(t), and estimated vector, x(t).In this study, as presented in Algorithm 1, the residual vector for training data was calculated using the leave-one-out method [45].More details of the leave-one-out method are presented in [40,45].

Detection Indices and Confidence Limit
Hotelling's T 2 statistic and SPE are typically used as multivariate statistics for fault detection.In this study, SPE was employ as a multivariate statistic and is calculated using the residual vector e(t) as follows: If the target system is normal, the magnitude of SPE(t) is small.In contrast, the values of SPE(t) begin to increase rapidly when a fault occurs in the system.Thus, SPE(t) represents the condition of the target system for multivariate monitoring.In fault-detection steps using multivariate statistics, the threshold should be defined in advance.In this study, we employed KDE because it is widely used to define the threshold value for declaring the fault [29,46,47].In KDE, a univariate kernel estimator with kernel function K is defined as where K(•) is the kernel function; j is the number of samples; h is the smoothing parameter; and The kernel estimator, fh (•), is the sum of the bumps located in each sample.The shape of the bumps is determined using the kernel function K(•) [29].Various kernel functions are available, such as Epanechnikov, uniform, and Gaussian.Among these kernel functions, we employed the most commonly used Gaussian kernel function.In practice, the confidence limit determined in KDE is influenced by the value of the smoothing parameter, h.To estimate cumulative kernel function, we used the 'ksdensity' MATLAB function built into the Statistic and Machine Learning Toolbox.More details of KDE and its smoothing parameter are available in [48,49].

Performance Indices
In this study, the false alarm rate (FDR) and miss detection rate (MDR) (also known as type I and type II errors, respectively) were used as performance indices to quantify the fault detection performance of the proposed method, as summarized in Table 1.They are based on statistical hypothesis tests and are widely used to quantify the performance of fault detection.In Table 1, H 0 is the null hypothesis, which means that the target system is normal, and H 1 is an alternative hypothesis that means that the target system is abnormal.That is, FAR indicates the rates of false alarms, in which the model has detected a fault, but the fault has not occurred in the actual target system.MDR refers to the rates of miss detection, in which a fault has occurred in the actual target system but the model has not detected the fault.If the MDR value is high, the model does not detect faults occurring in the real target system.Therefore, MDR is considered to be much more important than FAR in this research field.

Data Acquisition
In this section, we introduce the target system and milling process used to measure the current data.Current data are collected from the Internet of Things (IoT) sensor in an actual CNC machining center.The sampling frequency of the current data was 6 Hz.The models of the CNC machining center and cutting tool were Fanuc Oi, Mynx 5400, and YG-1 SUS-CUT, respectively.We used a CNC machining center operated in an actual manufacturing plant and stainless steel to obtain milling process data.Figure 3 shows the CNC machining center used to obtain the milling process data.To obtain two types of milling process data, we conducted spiral circular and straight parallel cuttings, as illustrated in Figure 4.These cutting processes are primarily used to machine rigid materials during machining.In the spiral circular cutting process, we did not obtain tool breakage data despite repeatedly performing spiral circular cutting processes.Therefore, in this study, we generated artificial fault data by applying two types of disturbances (bias and drift) to normal data measured through a spiral cutting process.For the straight parallel cutting process, machine tools were broken during long-term operation.Thus, we first validated the proposed method by using two types of artificial fault data, and then conducted a real-world application by applying the proposed method to actual tool breakage data.The details of each machining method and the acquired dataset are described below.

Case 1. Spiral Circular Cutting Process
Spiral cutting was performed clockwise until the length of the center of the processed stainless steels materials was 16 mm.We used three tools with a diameter of 10 phi, and 73 stainless materials were processed during spiral circular cutting.Notwithstanding the long-term operation of an actual machine tool, we could not obtain actual fault data.Hence, to obtain the fault data in this study, three types of artificial faults for each bias and drift disturbance (i.e., a total of six cases) were applied to the normal data measured by spiral circular cutting, respectively.After sampling approximately 23,500 normal data points among the datasets obtained through the spiral cutting process, 16,500 samples were used for learning, and the remaining 7000 samples were used for testing.The normal dataset was sampled by considering the period immediately after replacing the tool with the new tool as the normal condition.To generate the artificial fault data, the bias and drift faults were applied to the test samples from the 4001st sample to the end.

•
Bias fault: a step change of the current variable by N was introduced from the 4001st sample to the end (N = 0.02, 0.025, and 0.03, respectively).

•
Drift fault: the original current variable was linearly increased from the 4001st sample by adding M(k − 4000) to the current variable of each sample in this range, where k is the sample number (M = 0.0005, 0.00075, and 0.001, respectively).
The bias fault data have the form of a step function.These fault data describe a scenario in which fluctuations in the measured data are significant because of a fatal fault in the actual industrial process.In contrast, the drift fault data are in the form of a ramp function.This describes a phenomenon in which the size of a fault rapidly increases owing to the small faults in the facility.The generated artificial fault data in this manner were used as the test data to verify the proposed method.

Case 2. Straight Parallel Cutting Process
Figure 5 shows the status of straight cutting of stainless-steel materials before and after milling.As shown in Figure 5a, after the straight parallel cutting of processed stainless steel material 86 times, straight cutting was performed under the same conditions on the opposite side of the processed material.In this way, we used five tools with a diameter of 10 phi, and 68 stainless materials were machined during the straight cutting process.Four of the five tools used in the straight cutting process broke during the process.Among the tool breakage cases, we used three fault datasets of tool breakage to verify the proposed method.Similar to the spiral circular process, the training dataset was sampled by considering the period immediately after replacing the used tool with a new tool as a normal state.
A summary of the datasets measured using the five tools is presented in Table 2.During the use of the first tool, several values were missing owing to communication errors between the machining center and IoT sensor.In addition, the first tool used in the straight cutting process did not break despite processing 18 materials being processed.The second tool cut a total of 16 materials, and broke when cutting the 34th material.The third tool suddenly broke when processing the first material.Because cause of breakage cannot be clearly explained, and it differed from the general failure of the running tool, this dataset was not used as test data.In these experiments, 14 and 15 materials were processed using the fourth and fifth tools, respectively, and each tool broke when the last material was processed.In this manner, we constituted the dataset used in this study for detecting faults in machine tools.

Experimental Results and Discussion
In this section, we explain the results of the fault detection of the proposed method for the actual tool data introduced in Section 3. To perform comparative studies, we employed the local outlier factor (LOF), PCA, and ICA as comparison methods.LOF is a distancebased method frequently used to detect faults in multivariate processes, and can calculate local outlier factors for multivariate data by k-nearest neighbors (k-NN) using a predefined distance function.In this study, we set the value of k to 10 and employed the Euclidean distance to compute the distance between the nearest neighbors.PCA and ICA are the most well-known multivariate statistical techniques used primarily in large-scale processes that require the analysis of many process variables.They are suitable models for multivariate process data, owing to dimensionality reduction for the selection of useful process variables.In this study, Hotelling's T 2 and SPE statistics were used as detection indices for PCAbased process monitoring.The T 2 statistic is calculated based on the principal components (PCs) selected for dimensionality reduction.In the case of ICA, we employed the SPE and I 2 d statistic that have commonly been used.The I 2 d statistic is computed to consider the dominant part of independent components (ICs).For more details regarding the I 2 d statistic, readers are invited to read the paper [46,47].In ICA, we used the FastICA algorithm, which is a very simple and highly efficient fixed-point algorithm proposed in [50].As mentioned above, the comparison methods are effective for a multivariate dataset.To consider the property, we combined EMD with comparative methods (i.e., EMD-LOF, EMD-PCA, and EMD-ICA) to verify the proposed method.The experimental results are discussed below.

Artificial Fault Cases (Spiral Circular Cutting Process)
In this subsection, we introduce the experimental results, in which artificial fault datasets generated by two types of disturbances were used to verify the performance of the proposed method.Artificial fault data were generated by applying drift and bias to normal data that were not used for training.Figure 6 shows the behavior of the faulty data generated by three types of artificial faults for each bias and drift disturbance.In these figures, the blue line describes the normal behaviors of the faulty variable, and the red line represents the abnormal behaviors of the faulty variable owing to disturbances.As shown in Figure 6a,c,e, the magnitude of current value has risen steeply at the 4001th sample by 0.02, 0.025, and 0.03, respectively.In Figure 6b,d,f, the current value gradually increases from time t = 4001 drift faults.From now on, only the results of the best cases for each disturbance will be explained to save space.The best cases for each bias and drift fault are N = 0.03 and M = 0.001, respectively.7a) and test datasets (Figure 7b,c) of AAKR.An appropriate number of IMFs was selected while increasing the number of IMFs from one to five in increments of one.When IMF 1 and IMF 2 were added to the dataset, AAKR exhibited the best fault-detection performance.That is, the final constituted dataset X = [bias faulty variable, IMF 1 , IMF 2 ].Thus, we constructed a multivariate dataset by extracting the IMFs from the univariate current signal.
Processes 2022, 10, x FOR PEER REVIEW 14 of Figure 8 shows the results of the bandwidth parameter selection using k-fold crossvalidation.In this study, we set the value of k to 10, and the values of hmax, hmin, and ∆h were set at 0.01, 1.3, and 0.01, respectively.As shown in the figure, the MSPE values decreased when h = 0.01, and then increased rapidly for the smallest value of h.From the kfold cross-validation, the appropriate bandwidth parameter value was determined to be h = 0.06. Figure 9 shows the empirical and estimated cumulative distribution functions (CDFs) for the detection indices.The blue and dashed green lines indicate the empirical and estimated CDFs, respectively.The inset of Figure 8 shows that the estimated CDF is adequately estimated for the empirical CDF.Therefore, we determined the confidence limits

Bias fault data
Drift fault data h = 0.06 Figure 8 shows the results of the bandwidth parameter selection using k-fold crossvalidation.In this study, we set the value of k to 10, and the values of h max , h min , and ∆h were set at 0.01, 1.3, and 0.01, respectively.As shown in the figure, the MSPE values decreased when h = 0.01, and then increased rapidly for the smallest value of h.From the k-fold cross-validation, the appropriate bandwidth parameter value was determined to be h = 0.06. Figure 9 shows the empirical and estimated cumulative distribution functions (CDFs) for the detection indices.The blue and dashed green lines indicate the empirical and estimated CDFs, respectively.The inset of Figure 8 shows that the estimated CDF is adequately estimated for the empirical CDF.Therefore, we determined the confidence limits with a significance level of α = 0.01 from the estimated CDFs.Figures 10 and 11 show monitoring charts that are the results of fault detection using the comparison and proposed method for bias and drift data, respectively.The black line denotes the detection indices, and the dotted purple line indicates the confidence limits determined using KDE.In Figures 10f and 11, each monitoring chart at time t = 1 to 5000 is enlarged because the detection indices for the normal region are very low compared with that for the abnormal area.First, as shown in Figure 10, the LOF values exceeded the confidence limit before the bias fault were occurred.Furthermore, many LOF values after failure were often lower than the threshold.In the case of PCA and ICA, T 2 (PCA) and SPE (ICA) statistics did not entirely detect faults in the machine tool.SPE (PCA) and I d (ICA) statistics can be distinguished between normal and abnormal areas, but each threshold is set to be high.This can lead to a high value of MDR.In contrast, the proposed method (EMD-AAKR) completely detected the faults in all areas except when a small degree of false alarm was declared at approximately time t = 2000.Next, in Figure 11, although the faults occurred from time 4001, the LOF values continuously exceeded the threshold from roughly t = 4500.In the case of PCA and ICA, T 2 (PCA) and SPE (ICA) statistics consecutively exceeded the threshold at approximately t = 4500.This increased the MDR value.SPE (PCA) and I d (ICA) were expected to exhibit a lower MDR than other statistics (T 2 and SPE statistics).As shown in Figure 11f, detection indices for proposed method quickly exceeded the threshold at t = 4000 compared to other comparative methods, but some false alarms were declared in normal area.The false alarms were negligible considering the improvement of MDR for the proposed method.The performances of the proposed and comparison methods for the two types of artificial fault data are summarized in Table 3.The lowest values in each performance index are expressed in bold.As shown in Table 3, if the magnitude of each bias and drift is larger, then the performance of the proposed method and comparison methods are better; in case of bias 3 (N = 0.03), SPE (PCA) and I d (ICA) have significantly improved performance compared to the other bias cases.LOF has a higher FAR than the other methods for both the bias and drift fault datasets.PCA and ICA have significantly lower FAR values for all datasets; however, the MDRs for T 2 (PCA) and SPE (ICA) were relatively high.In particular, the MDRs for all bias fault dataset were higher than 90%.This implies that the detection indices did not detect faults in the target system.In the case of drift fault dataset, SPE (PCA) and I d (ICA) achieved lower FARs and MDRs than other comparative methods, but the MDRs for bias fault data were still high.The proposed method outperformed the comparison methods; the FAR for proposed method were slightly higher than that of the other comparative methods except LOF, and the MDRs were the lowest than other methods.In particular, the MDR of the proposed method for the bias 3 (N = 0.03) was calculated as 0%.As the results of process monitoring, Table 3 shows that the proposed method successfully detects the artificial faults.In this subsection, we introduce the fault-detection results of the proposed method for the tool breakage dataset.As mentioned in Section 4.2., straight parallel cutting was performed with five machine tools, and four of them broke.Among the broken tools, the dataset for fourth tool was excluded from the test dataset because of a sudden failure, which was different from other tool failures and difficult to explain.Thus, we used three types of tool breakage data measured using the second, fourth, and fifth tools to demonstrate the performance of the proposed method.
Figure 12 shows the original and IMFs signals extracted from the normal and tool breakage data (second, fourth, and fifth tool) using EMD.In Figure 12, each first-row signal is the original signal, and the others are the IMF extracted using EMD from the original current signal.By comparing with each first-row signal, we can confirm that the fault signal gradually increased owing to tool breakage compared with the normal signal.In Figure 12, the signals contained in the red box are composed of a set of training and test data of the AAKR for fault detection.,In this case, the number of IMFs was determined using the method introduced in Section 5.1.Finally, we constructed the dataset X = [actual faulty variable, IMF1, IMF2] to verify the proposed method.
Figure 13 shows the results of the bandwidth parameter selection for straight parallel cutting data using k-fold cross-validation.We selected the value of k as 10 for k-fold cross validation.h max , h min , and ∆h were set to 0.01, 1.3, and 0.01, respectively.As shown in this figure, when the value of h increased from the predefined initial value h min , the corresponding MSPE decreased exponentially.Subsequently, MSPE attained minimal at h = 0.081.Therefore, we determined the value of the bandwidth parameter h to be 0.081.current signal.By comparing with each first-row signal, we can confirm that the fault signal gradually increased owing to tool breakage compared with the normal signal.In Figure 12, the signals contained in the red box are composed of a set of training and test data of the AAKR for fault detection.,In this case, the number of IMFs was determined using the method introduced in Section 5.1.Finally, we constructed the dataset X = [actual faulty variable, IMF1, IMF2] to verify the proposed method.Figure 13 shows the results of the bandwidth parameter selection for straight parallel cutting data using k-fold cross-validation.We selected the value of k as 10 for k-fold cross validation.hmax, hmin, and ∆h were set to 0.01, 1.3, and 0.01, respectively.As shown in this figure, when the value of h increased from the predefined initial value hmin, the corresponding MSPE decreased exponentially.Subsequently, MSPE attained minimal at h = 0.081.Therefore, we determined the value of the bandwidth parameter h to be 0.081.Next, we provide the results of fault detection for three tool-breakage datasets gathered from an actual CNC machine.As mentioned in Section 4.2., the datasets using second, fourth, and fifth tools were used to verify the proposed method.Figures 15-17 show the monitoring charts applying the comparison and proposed methods to the three types of tool breakage datasets, as mentioned above.The black line denotes the detection indices, and the dotted purple line represents the thresholds determined using KDE.In Figures 16f and 17c,d,f, the monitoring charts at t = 1 to 5000 are enlarged because the detection indices for the normal region were very small compared with those for the abnormal area.First, for the results of second tool data, the values of SPE (PCA), I d (ICA), and proposed method statistics began to increase at t = 400.From this, we can determine the occurrence of the first fault.However, the others did not identify the time points; the majority of LOF, T 2 (PCA), and SPE (ICA) statistics were below each threshold.Compared with the proposed method, since the confidence limits of SPE (PCA) and I d (ICA) statistics were set too low, FAR was expected to be high.In the case of the fourth tool dataset, we identified the start time (t = 280) of consecutive faults through the monitoring chart of SPE (PCA), I d (ICA), and SPE (AAKR) statistics.In addition, the proposed method was expected to have smaller false alarms.In contrast, LOF, T 2 (PCA), and SPE (ICA) statistics still did not detect the faults.As shown in Figure 17, the comparison methods outperformed each result for the second and fourth tool.Nevertheless, because LOF, T 2 (PCA), and SPE (ICA) could not detect the start time of continuous faults, the MDRs of these methods were expected to be high.On the other hand, SPE (PCA), I d (ICA), and proposed method identified the outset (t = 420) of the faults.In particular, the proposed method achieves a satisfactory performance for every tool dataset.
The performance indices of the proposed and comparison methods for the actual tool breakage datasets are summarized in Table 4.The lowest values for each performance index are written in bold.For the results of every tool dataset, the FARs and MDRs of the LOF, T 2 (PCA), and SPE (ICA) were calculated as very low and high, respectively; in the case of second and fourth tool, the MDR of these methods were approximately 100%.As can be seen from Figures 15 and 16, these methods could not detect the faults.SPE (PCA) and I d (ICA) achieved better performance than the other comparison methods; for the fifth tool dataset, these methods outperformed the proposed method, but FARs were the highest in the case of second and fourth tools.The proposed method for every tool dataset achieved satisfactory performance in comparison to SPE (PCA) and I d (ICA); although FAR for the second tool was high, except this case, both FAR and MDR were relatively low.

Figure 1 .
Figure 1.Framework of fault-detection method for machine tool based on EMD and AAKR.

Figure 2 .
Figure 2. Schematic block diagram of fault-detection method based on AAKR.

Algorithm 1 :
Leave-one-out method for calculating the residual vector Input: Training data X = [x 1 , . . ., x j ] h ← bandwidth parameter determined by k-fold cross validation for l from 1 to j X ← X/{x l } Calculate the distance function values d l (x l , x i ), l = {1, . . ., j}\{l} between data vectors in X and x i Generate weights K h (d l ), l = {1, . . ., j}\{l} of each data vector Obtain lth estimated vector xl Calculate residual vectors e l = x l − xl = [e 1 l , . . ., e p l ] T end return {e 1 , . . ., e j }

Figure 3 .
Figure 3. CNC machining center used in this study for data acquisition.

Figure 4 .
Figure 4. Blueprints of process shape for spiral circular cutting and straight parallel cutting: (a) for spiral circular cutting process; (b) for straight parallel cutting process.

Figure 5 .
Figure 5. Condition of straight parallel cutting and images of materials before and after the straight cutting process: (a) condition of straight cutting process; (b) before the straight milling process; (c) after the milling process.

Figure 7
Figure 7 depicts the IMF signals extracted from the normal and artificial fault data (bias and drift) using EMD.Each first-row signal in Figure 7 is the input signal of EMD, and the other row signals are the outputs of EMD.Some signals included in the red box in Figure 7 are the training (Figure7a) and test datasets (Figure7b,c) of AAKR.An appropriate number of IMFs was selected while increasing the number of IMFs from one to five in increments of one.When IMF 1 and IMF 2 were added to the dataset, AAKR exhibited the best fault-detection performance.That is, the final constituted dataset X = [bias faulty variable, IMF 1 , IMF 2 ].Thus, we constructed a multivariate dataset by extracting the IMFs from the univariate current signal.

Figure 7 .
Figure 7. (Spiral cutting data) Results of EMD for constructing the training and test datasets: (a) for normal training data, (b) for bias fault data, (c) for drift fault data.

Figure 8 .
Figure 8. (Spiral cutting data) Result of bandwidth parameter selection using k-fold cross validation.

Figure 7 .
Figure 7. (Spiral cutting data) Results of EMD for constructing the training and test datasets: (a) for normal training data, (b) for bias fault data, (c) for drift fault data.

Figure 8 .
Figure 8. (Spiral cutting data) Result of bandwidth parameter selection using k-fold cross validation.

Figure 9 .
Figure 9. (Spiral cutting data) Empirical and estimated CDFs of SPE for training samples.

Figure 12 .
Figure 12. (Straight cutting data) Results of EMD for constructing the training and test dataset: (a) for normal training data, (b) for second tool fault data, (c) for forth tool fault data, (d) for fifth tool fault data.

Figure 13 .
Figure 13.(Straight cutting data) Result of bandwidth parameter selection using k-fold cross validation.

Figure 14
Figure 14 depicts the empirical and estimated CDFs for detection indices.The blue

Figure 12 .
Figure 12. (Straight cutting data) Results of EMD for constructing the training and test dataset: (a) for normal training data, (b) for second tool fault data, (c) for forth tool fault data, (d) for fifth tool fault data.

Figure 13 .
Figure 13.(Straight cutting data) Result of bandwidth parameter selection using k-fold cross validation.

Figure 14
Figure 14 depicts the empirical and estimated CDFs for detection indices.The blue and dashed green lines indicate the empirical and estimated CDFs, respectively.In Figure 8, the enlarged figure shows that the estimated CDF was reasonably estimated for the empirical CDF.Therefore, we determined the confidence limits with significance level of α = 0.01 from the estimated CDFs.

Figure 14 .
Figure 14.(Straight cutting data) Empirical and estimated CDFs of SPE for training samples.
2. Extract the ith IMF candidate c i .

Table 1 .
The performance indices based on statistical hypothesis test.

Table 2 .
Summary of the datasets measured using five tools in straight cutting process.

Table 3 .
(Spiral cutting data.)Performance indices of comparison and proposed methods.

Table 4 .
(Straight cutting data) Performance indices of comparison and proposed methods.