Intelligent expert system fault diagnosis based on PCA–SPE–CNN classifier

PCA-SPE-based fault monitoring

PCA-SPE is a data processing method that combines principal component analysis (PCA) and square prediction error (SPE). PCA is used to downscale high-dimensional data by extracting the main features in the data and transforming complex multidimensional data into a few major components, thus reducing the dimensionality and noise of the data. SPE, on the other hand, is used to measure the prediction error of the downscaled data, i.e., the difference between the original data and the reconstructed data of the PCA model. By analyzing the SPE values, it is possible to identify anomalous samples or points of failure, as these data points are often quite different from the predictions of the PCA model. Thus, the PCA-SPE approach not only effectively simplifies the data, but also maintains sensitivity to anomalous or faulty data while downgrading the dimensions.

The compressor data $X$ ∈ $R_{n\times m}$ (where n and m denote the number and dimension of samples, respectively) were collected following the procedure of Wang & Hu (2025). Each sample $x_k$ was measured at $t_k=t_0+kd$, with $t_0$ and $d$ representing the start time and sampling interval.

The monitoring samples are continuously obtained from $X$ using a sliding window of size 100, i.e., data fragments of the form { $x_k$,…, $x_{k+100}$}. Each frame subseries is denoted as $X$ [ $t_k$: $t_{k+100}$]. The data within each window are normalized using the mean and variance computed from the training data. Specifically, each column of $X$ is scaled to have a mean of 0 and a root mean square of 1, resulting in the normalized matrix $Y$.

Normal running data $H=\{x\left(t_1\right)...x\left(t_{1500}\right)$}, having the same dimension as the monitoring samples, are selected from the offline dataset. After error-tolerant filtering (Shaolin, Na & Wenming, 2016), each column of $H$ is normalized to have a mean of 0 and a root mean square of 1, resulting in the matrix $\hat{X}$.

• (1)

  The covariance matrices of the normalized training $U_a$ and monitoring data $U_s$ are computed:

• (1) $$\begin{array}{c}{U}_a=\frac{1}{n-1}{\hat{X}}^T\end{array}\hat{X}$$ (2) $$U_s=\frac{1}{h-1}{Y}^{T}Y$$

• where n and h are the numbers of rows in the training and monitoring samples, respectively, and the resulting covariance matrix is a 12 × 12 characteristic matrix.

• (2)

  The eigenvalues are sorted from large to small, and the top k features with a sum of more than 91% eigenvalues are selected for PCA dimensionality reduction. (3) $$\frac{{\sum }_{i=1}^k{\lambda }_i}{{\sum }_{i=1}^m{\lambda }_i}\ge \mathrm{0.91.}$$

The first k eigenvalues from large to small form a diagonal matrix $S_{k\ast k}$, and k corresponding eigenvectors form a matrix $P_{m\ast k}$.

The formulations in Eqs. (4)–(11) are directly based on the approach proposed by Wang & Hu (2025). The training set diagonal matrix and the corresponding eigenvectors form the matrix:

(4) $$S_{k_1\ast k_1}=diag\left({\lambda }_1,{\lambda }_2,\cdots ,{\lambda }_{k_1}\right).$$

(5) $$P_{m\ast k_1}=\left[p_1,p_2,\cdots ,p_{k_1}\right].$$

The diagonal matrix of the monitoring set and the corresponding feature vector form a matrix:

(6) $${\hat{S}}_{k_2\ast k_2}=diag\left({\hat{\lambda }}_1,{\hat{\lambda }}_2,\cdots ,{\hat{\lambda }}_{k_2}\right).$$

(7) $${\hat{P}}_{m\ast k_2}=\left[{\hat{p}}_1,{\hat{p}}_2,\cdots ,{\hat{p}}_{k_2}\right].$$

For $\hat{X}$, $Y$ after dimensionality reduction, the number of samples remains unchanged, but the number of features becomes K. After dimensionality reduction, it is:

(8) $${\tilde{X}}_{n\ast k1}={\hat{X}}_{n\ast m}\ast P_{m\ast k_1}.$$

(9) $${\tilde{Y}}_{h\ast k_2}=Y\ast {\hat{P}}_{m\ast k_2}.$$

The formula for calculating the matrix of $\stackrel{⌢}{X},\stackrel{⌢}{Y}$ obtained by reconstructing $\hat{X}$ and $Y$ is as follows:

(10) $$\begin{array}{c}\stackrel{⌢}{X}={\tilde{x}}_{n\ast k_1}{P}_{m\ast k_1}^T=X_{n\ast m}{P}_{m\ast k_1}{P}_{m\ast k_1}^T\end{array}$$

(11) $$\begin{array}{c}\stackrel{⌢}{Y}={\tilde{Y}}_{h\ast k_2}{\hat{P}}_{m\ast k_2}^T=Y_{\mathrm{h}\ast \mathrm{m}}{\hat{P}}_{m\ast k_2}{\hat{P}}_{m\ast k_2}^T\end{array}$$where n is the number of rows in the training set, h is the number of rows in the monitoring set, m is the number of corresponding key variables, and k is the number of corresponding variables after dimensionality reduction.

The training set in this study consists of 1,500 sampled values, while the monitoring sample adopts a sliding window containing 100 data points. These two samples are statistically independent and both follow a normal distribution. Additionally, the variable $\stackrel{⌢}{X}$ is assumed to follow a zero-mean normal distribution. We used a PCA-based process monitoring methodology with the SPE statistic (Bakdi & Kouadri, 2017; Tamura & Tsujita, 2007) to make the determination.

The control limits for the SPE statistic were proposed by Jackson & Mudholkar (1979). $\stackrel{⌢}{Y}$ conforms to the multivariate normal distribution, and the SPE statistic (Zhou, Park & Liu, 2016) is established for $\stackrel{⌢}{Y}$, that is, for ${\stackrel{⌢}{Y}}_{i,:}$ (i = 1, 2 … h). Construct statistic Q:

(12) $$\begin{array}{c}Q={\stackrel{⌢}{Y}}_{i,:}\ast \left(I_{m\ast m}-{\hat{P}}_{m\ast k_2}{\hat{P}}_{m\ast k_2}^T\right)\ast {{\stackrel{⌢}{Y}}_{i,:}}^T\end{array}$$where ${\stackrel{⌢}{Y}}_{i,:}$ is the row vector of row i of $\stackrel{⌢}{Y}$, $I_{m\ast m}$ is the identity matrix, and ${\hat{P}}_{m\ast k_2}$ is the matrix ${\hat{P}}_{m\ast k_2}$ of k eigenvectors corresponding to k eigenvalues from large to small from (1) and (2).

Since $\stackrel{⌢}{X}$ and $\stackrel{⌢}{Y}$ are assumed to follow a normal distribution with zero mean, and the corresponding sampled data points are mutually independent, it has been established in prior studies (Jackson & Mudholkar, 1979) that Eq. (13) can be approximated by a normal distribution.

(13) $$(Q/{\theta }_1{)}^{h_0}\sim N\left[1+\frac{{\theta }_2{h}_0(h_0-1)}{{{\theta }_1}^2},\frac{2{\theta }_2{h_0}^2}{{{\theta }_1}^2}\right].$$

(14) $$\mathrm{W}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{e}{\theta }_r=\sum _{j=k1+1}^m{\lambda }_j^r,r=1,2,3.$$

(15) $$h_0=1-\frac{2{\theta }_1{\theta }_3}{3{\theta }_2^2}.$$

Under these assumptions, the variable $c$ can be approximated as normally distributed with a mean of zero. Its expected value is calculated as:

(16) $$c={\displaystyle \frac{{\theta }_1[{(Q/{\theta }_1)}^{h_0}-1-{\theta }_2{h}_0(h_0-1)/{{\theta }_1}^2]}{\sqrt{2{\theta }_2{h_0}^2}}.}$$

Accordingly, the control limit for the Q statistic is given by:

(17) $$Q_{\alpha }={\theta }_1{\left[{\displaystyle \frac{c_{\alpha }{h}_0\sqrt{2{\theta }_2}}{{\theta }_1}+1+{\displaystyle \frac{{\theta }_2{h}_0\left(h_0-1\right)}{{\theta }_1^2}}}\right]}^{{\displaystyle \frac{1}{h_0}}}.$$In this expression, $c_{\alpha }$ denotes the critical value of the standard normal distribution at a 99.6% confidence level. The parameter $\lambda$ corresponds to the eigenvalue of the covariance matrix $R_{ta}$, and $h_0$ is an empirical correction factor introduced to fine-tune the control limit.

According to the methodology outlined by Wang & Hu (2025), data were collected and processed using a sliding window strategy combined with standard normalization. However, unlike the original approach, the construction of the monitoring matrix differs in this study. The subsequent steps, including PCA transformation and SPE-based monitoring, follow the original framework, with necessary adjustments made to parameter settings and evaluation criteria.

• (1)

  According to Eqs. (12) and (17), the SPE statistic $Q$ and the statistical limit $Q_{\alpha }$.

• (2)

  Fault determination

• If the system is functioning normally, the $Q$-value of the sample should meet $Q$ < $Q_{\alpha }$, otherwise it is considered to be faulty.

• If $Q$ < $Q_{\alpha }$, the corresponding fault data $\hat{Y}$. (Please refer to the testx1. csv file in the Supplemental Document).

CNN module

The CNN module consists of a CNN layer, a pooling layer, and a fully connected layer, as shown in Fig. 1, which are described in detail below (Wang & Hu, 2025).

The CNN layer involves two primary operations: convolution and activation. The input to the convolution operation is a 3D tensor Z ∈ $R^{w\times h\times c}$, where w, h, and c denote the width, height, and number of channels, respectively. In the special case where c = 1, the input simplifies to a 2D form $\hat{Y}$ ∈ $R^{w\times h\times c}$, which is used as input to the proposed model. Specifically, in the two-dimensional input data, the horizontal axis corresponds to the feature dimension, while the vertical axis represents time, with their scales determined by the number of attributes and the sliding window width, respectively. To clearly describe the convolution operation, four hyperparameters are introduced: the number of convolutional kernels K, kernel size F, stride S, and padding P. To ensure that the output maintains the same shape as the input, the input data $\hat{Y}$ ∈ $R^{w\times h\times c}$ is zero-padded accordingly. As a result, the transformed data ${\hat{Y}}^{\ast }$ ∈ $R^{(w+\mathrm{p})\times (h+p)\times c}$ retains the original dimensionality. Following the padding, the convolutional slicing operation is performed as follows:

(18) $${\hat{Y}}_{t_w,t_h}={\hat{Y}}^{\ast }\left[t_w:t_w+F,t_h:t_h+F,:\right]$$where $t_w$: $t_w$ + F represents the operation of extracting the subset located between row $t_w$ and row $t_w$ + F, $t_h$ : $t_h$ + F represents the operation of extracting the subset located between column $t_h$ and column $t_h$ + F, and : represents all the data in the direction of the cutting channel.

The convolution and activation methods presented in this section are based entirely on the framework proposed by Wang & Hu (2025), $t_w$ and $t_h$ start at 1 and change with the number of steps S. According to ( $t_w$, $t_h$), data ${\hat{Y}}_{t_w,t_h}$ can be extracted from ${\hat{Y}}^{\ast }$. Let the convolution kernel be denoted as kernel ∈ $R^{F\times c}$, let $T_{t_w,t_h}\in R^{F\times c}$ denote the convolution kernel and ${\hat{Y}}_{t_w,t_h}$ the corresponding input subregion. The element-wise multiplication between $T_{t_w,t_h}\in R^{F\times c}$ and ${\hat{Y}}_{t_w,t_h}$ is performed, and the result is calculated as follows (Wang & Hu, 2025):

(19) $$T_{t_h}=Kernel\times {\hat{Y}}_{t_h}.$$

Accordingly, the output $d_{t_h}$ is obtained through the following operation (Wang & Hu, 2025):

(20) $${\mathrm{d}}_{t_h}=\sum _{j=0}^F\sum _{k=0}^C{T}_{t_w,t_h}\left[j,k\right].$$In this expression, $T_{t_h}\left[j,k\right]$ represents an individual element from the set $T_{t_w,t_h}$, and the resulting feature map $\tilde{D}\in R^{w\ast \times h\ast }$ is constructed from $d_{t_w,t_h}$, ordered according to the dimensions $t_w$ and $t_h$. It is worth noting that

(21) $$h^{\ast }={\displaystyle \frac{h-F+2\ast P}{S}+1.}$$

For each convolutional kernel, a corresponding feature map D is generated. Repeating this operation across all K kernels produces a new three-dimensional output tensor $\hat{Z}\in R^{w\ast \times h\ast \times K}$. Following the convolution operation, an activation function is applied to introduce nonlinearity into the network. This nonlinearity enhances the model’s capacity to learn complex representations and improves the separability of features. In this study, the Rectified Linear Unit (ReLU), a widely adopted activation function, is used to activate each element of the output tensor $\hat{Z}$. The activation is defined as Wang & Hu (2025):

(22) $$A\left[i,j,k\right]=f\left(\hat{Z}\left[i,j,k\right]\right)=max\{0,\hat{Z}\left[i,j,k\right]\}.$$Here, $\hat{Z}[i,j,k]$ denotes an element of $\hat{Z}$, and $A[i,j,k]$ denotes an element of $A\in R^{w\ast \times h\ast \times K}$, representing the values after the activation operation.

Dropout Regularization:

The Dropout layer helps prevent overfitting by randomly setting a subset of input elements to zero during training. This disrupts co-adaptations between neurons and promotes model generalization. The operation of Dropout is mathematically defined as:

(23) $$\mathrm{D}\mathrm{r}\mathrm{o}\mathrm{p}\mathrm{o}\mathrm{u}\mathrm{t}\left(\mathrm{L},p\right)=L\cdot M$$where L denotes the input data, p is the dropout rate, and M is a binary mask matrix sampled from a Bernoulli distribution.

Fully Connected Layer (Dense):

The fully connected layer transforms the output from the global average pooling layer into class probabilities using a weight matrix and a bias vector. Its operation is defined by the following formula:

(24) $$\mathrm{D}\mathrm{e}\mathrm{n}\mathrm{s}\mathrm{e}\left(G,W,b\right)=\mathrm{s}\mathrm{o}\mathrm{f}\mathrm{t}\mathrm{m}\mathrm{a}\mathrm{x}\left(G\cdot W+b\right)$$where G represents the result from the global average pooling layer, W is the associated weight matrix, b denotes the bias vector, and $\mathrm{s}\mathrm{o}\mathrm{f}\mathrm{t}\mathrm{m}\mathrm{a}\mathrm{x}$ ( $\cdot$) is the activation function used to convert raw outputs into normalized class probabilities.

Establishment of knowledge base framework

The causes of faults, troubleshooting methods and their corresponding rule patterns are shown in Table 1 below. The four rules correspond to four faults, namely 1 (compressor vibration or noise), 2 (dry gas seal failure), and 3 (oil pressure too low), 4 (the driver does not start). According to the corresponding inference rules, faults are identified according to the rule with the largest number of matched causes.

Rule mode and fault diagnosis

In order to facilitate the storage of expert knowledge and experience and the corresponding rule relationship in the later maintenance system in the form of documents, the automatic diagnosis mode is mainly based on the improved CNN classifier proposed in this article, and the corresponding rules of the expert system are used to confirm the corresponding faults, and then the follow-up analysis and diagnosis is carried out according to the expert knowledge. Implemented data-driven expert system automatic diagnosis, reducing manual maintenance. Figure 2 shows the automatic confirmation failure mode.

If you want to further reason and confirm the detection results, this article adds a manual auxiliary fault detection mode based on expert experience and knowledge. Manual rule reasoning uses the number of corresponding input condition matches to determine the most probable faults based on prior expert knowledge experience and other rule reasoning to determine the most probable faults with the highest number of corresponding input conditions, and Fig. 3 shows the manual confirmation of fault patterns.

The expert system automatically matches the expert experience in the knowledge base based on early automatically identified fault modes to achieve automatic fault diagnosis. The diagnosis results can be used as training data for classifier training, so that it can learn from actual operations, adjust according to the changing environment, and realize automatic updating of the expert system knowledge base. The diagnostic results are shown in Fig. 4 below.

Introduction to data settings

For each fault and normal state, the compressor status is represented in both a training data set and a test data set. Each training data set contains 1,500 samples, while each test data set includes 217 samples and corresponding tag data. Each sample consists of 12 features. Therefore, the total number of data points in the original training set is 1,500 × 12 = 18,000, and in the original test set, it is 217 × 12 = 2,604. The data sampling frequency is 1 sample per minute. To ensure effective learning outcomes, the training samples, which include both faulty and normal samples, are randomly shuffled and normalized prior to processing.

Network architecture setup

The input shape of the CNN network is set to a batch size of 50 along the vertical axis, with 12 features along the horizontal axis. The convolution kernels have a size of 1 × 3, and the network employs a series of 16, 16, 64, 64, 256, and 256 convolution kernels. The activation function used is the ReLU function, and the pooling layers include both maximum pooling and average pooling. The network structure consists of three consecutive convolutional layers followed by maximum pooling, with an average pooling layer at the end. To improve performance and prevent overfitting, regularization techniques are applied, and fully connected layers are used for fault classification. The model is trained for 25 epochs, with a total training time of 6 s, which is sufficient to meet the monitoring needs. The hardware reference graphics card is 3,060 or above, the CPU is i5 eighth generation or above, and the hard disk is 50 GB or above.

Parameter tuning

In the data preprocessing stage, we use a sliding window approach to process the monitoring time series. The selection of window length was empirically evaluated by implementing a grid search strategy on the training set, taking into account the sensitivity to transient faults and robustness to noise interference. Classification accuracy and SPE statistical stability were used as measurement criteria to ultimately determine the optimal window length. In addition, during the fault detection phase, residual statistics are obtained via PCA-SPE, with thresholds derived from the chi-square distribution of normal data SPE values. These thresholds are fine-tuned through confidence level adjustments and empirical testing to ensure both consistency and adaptability to varying data distributions.

The network structure is tailored to the nature of the 12-dimensional input, which corresponds to two fault types. A compact convolutional design with small kernels is adopted to emphasize local patterns and hierarchical feature extraction (Wang & Hu, 2025). The architecture, iteratively optimized through experiments, comprises two convolutional layers followed by pooling, supporting better generalization and resilience against noise. To reduce overfitting and strengthen feature learning, regularization techniques and a triple-layer stacking scheme are applied. When stacked across multiple layers, these kernels enable the extraction of broader, high-level features, thereby improving the model’s ability to capture global contextual information (Wang & Hu, 2025). Although larger kernels can directly capture a broader range of features, they require more computational resources, so small kernels are preferred, and the network is deepened to achieve the same results efficiently. This strategy, widely used in successful networks, ensures high computational efficiency. A lower initial learning rate of 0.0002 is chosen, with the Adam optimizer adjusting it during training. A moderate dropout rate of 0.3 is used to prevent overfitting while maintaining model performance, and the number of filters is given an initial value space. Finally, grid search techniques are employed to optimize the model’s parameters, ensuring strong and efficient performance.

This article uses independently divided training and testing sets, and implements 5-fold cross validation on the training set. Specifically, the train0.csv dataset is divided into five parts, with one part selected as the validation set and the remaining four parts as the training set, while maintaining consistent category distribution. After five rounds of iteration, five sets of training/validation results can be obtained, and their mean (± standard deviation) can be taken as the generalization performance evaluation indicator of the model on the training set. The test set testx1.csv is retained until the final unified evaluation to ensure subject independence, that is, data from the same device/subject will not appear in both the training and testing sets at the same time. The monitoring sample is constructed by dividing the original signal through a sliding window (window length 100, step size 1) and filtering it using the PCA-SPE method. To control the randomness in the experiment, all experiments were fixed with a random seed (seed = 42 in NumPy and TensorFlow), ensuring the reproducibility of the results.

Each experiment was run independently 10 times. The evaluation indicators for model performance include classification accuracy, precision, recall, F1 score, and AUC, and the classification effect is visualized through a confusion matrix.

Classification accuracy and performance analysis of PCA-SPE-CNN classifier

The performance of the PCA-SPE-CNN classifiers has been validated using data, and their performance metrics were comprehensively evaluated. These metrics include accuracy, precision, recall, F1 score, AUC-Receiver Operating Characteristic (ROC), and the ROC curve, all of which were assessed using k-fold cross-validation. The ROC curve, derived from signal detection theory, illustrates the relationship between the true positive rate (TPR) and the false positive rate (FPR) of the classifier, with TPR plotted on the x-axis and FPR on the y-axis. The validation and evaluation results, based on these metrics, are presented and discussed below.

To evaluate classification performance, four key terms are used: true positive (TP), true negative (TN), false positive (FP), and false negative (NP). TP refers to correctly identified positive samples, while TN denotes correctly identified negative samples. FP represents negative samples misclassified as positive, and NP refers to positive samples misclassified as negative. Based on these definitions, the True Positive Rate (TPR) is calculated as TP/(TP + NP), reflecting the proportion of actual positives correctly identified. Similarly, the False Positive Rate (FPR) is FP/(FP + TN), indicating the proportion of actual negatives incorrectly classified as positive.

First, offline normal data is used to train the model. A sliding window is then applied to perform fault monitoring on real-time data. According to Eqs. (12) and (17), the SPE statistic and its corresponding control limit are calculated. Fault detection is conducted by comparing these two values: data points with SPE values exceeding the control limit are identified as faults. The monitoring results are illustrated in Fig. 5. The corresponding normal data filtering part is shown in Fig. 6.

After training the PCA-SPE-CNN algorithm, a confusion matrix was constructed, and the improved algorithm achieved 100% classification accuracy. The classification accuracy for both the training and verification sets, as well as the model’s classification results, are presented in Fig. 7. Additionally, the model underwent an evaluation based on the ROC curve, which is derived from the theory of signal detection. The ROC curve illustrates the relationship between the true positive rate (TPR) and the false positive rate (FPR) of the classifier (with TPR on the x-axis and FPR on the y-axis). This curve reflects the model’s performance across different thresholds, as shown in Fig. 8. The model also underwent k-fold cross-validation, as demonstrated in Fig. 9. To comprehensively assess the model’s performance, an analysis was conducted on various metrics, including accuracy, precision, recall, F1 score, and AUC-ROC. The results of this analysis are presented in Table 2 (Wang & Hu, 2025). Furthermore, the Bootstrap method was applied to estimate the performance metrics, enhancing the robustness of the analysis. By resampling 1,000 times, this method ensured that the model performed consistently under different validation set samples. Finally, the mean and 95% confidence intervals for the performance indicators were calculated to assess the variability in model performance under different sampling conditions. The calculated results are as follows: mean AUC: 1.000, 95% confidence interval (CI) [1.000–1.000]. mean accuracy: 1.00, 95% CI [1.00–1.000].

Mean AUC: this refers to the average AUC value of the model across all Bootstrap resamples. A value closer to 1 indicates better discriminative ability of the model.

95% CI for AUC: The 95% confidence interval for the AUC represents the range in which the model’s AUC value is expected to fall with 95% probability under different resampling conditions. A smaller interval suggests that the model’s performance is relatively stable, while a larger interval indicates significant fluctuations in performance across different data. Mean accuracy: this is the average accuracy of the model’s predictions across all resamples. 95% CI for accuracy: The 95% confidence interval for accuracy shows the fluctuation range of the model’s accuracy under various resampling conditions.

Classification accuracy and performance analysis of three classifiers

Due to the significant differences in the selected data, it is difficult to distinguish between known edge normal data and abnormal data. After training CNN, LSTM, and Transformer algorithms, a confusion matrix was constructed, and complex model comparisons were made through accuracy, precision, recall, and F1 scores. According to the confusion matrix, the classification performance is far inferior to the improved PCA-SPE-CNN model, and the accuracy, recall, and F1 score of the three compared models are not as good as the improved model. Table 3 shows their metric evaluations, while Figs. 10, 11, and 12 present their classification results. By comparing the above, the best solution can be achieved.

Expert system automatic diagnosis and analysis

After the new data is input, the system will perform feature matching and pattern recognition to determine whether existing rules are applicable. Apply relevant rules from the rule library to new data through a rule engine to generate new knowledge. Integrate the generated new knowledge with the content in the existing knowledge base to avoid redundancy and maintain consistency. Once the new knowledge is successfully integrated with the existing knowledge base, the system automatically updates the knowledge base without manual intervention. The expert system’s diagnostic results are shown in Fig. 13 below.