Operating status prediction method for gas storage compressors based on BiLSTM

2.1. Singular spectrum analysis

SSA is a non-parametric time-series analysis technique that effectively extracts trend, periodic, and noise components from raw monitoring signals. The method first transforms the original sequence into a Hankel matrix. Singular Value Decomposition (SVD) is then applied to obtain a set of decomposed components. These components can be selectively reconstructed to retain meaningful information.

To begin, a sliding window and diagonal averaging are used to convert the time-series data into a Hankel matrix. SVD is then performed on this matrix, yielding the left singular vectors, singular values, and right singular vectors. By selecting the dominant singular values and their corresponding vectors, the key components of the original sequence can be reconstructed.

As a data-driven denoising and decomposition method, SSA separates the main trend components from low-energy, unstructured noise components. The trend components are typically represented by the first few high-energy singular values. This process enables different frequency elements within the signal to be isolated, while retaining only the principal components to reconstruct a cleaner signal with an improved signal-to-noise ratio.

2.2. Convolutional neural networks

CNNs are a type of feedforward neural network widely used for feature extraction in time-series analysis. The core principle of CNNs is to apply convolution kernels along the input sequence to capture features at different scales.

Fig. 1Flowchart of 1D CNN for time series processing

In this study, the raw time-series data are first passed through one-dimensional convolution layers to extract local temporal features. A ReLU activation function is then applied to introduce nonlinearity and enhance feature expressiveness. Subsequently, a max-pooling layer reduces the dimensionality of the feature maps while preserving the most informative components. The processed feature representations are then forwarded to subsequent network layers for further analysis and prediction. The flowchart of 1D-CNN for time series processing used in this study is shown in Fig. 1.

2.3. Bidirectional long short-term memory networks

LSTM networks are a specialized form of Recurrent Neural Networks (RNNs). They are designed to address the vanishing-gradient and exploding-gradient problems commonly encountered when modeling long-term temporal dependencies. LSTMs achieve this by introducing memory cells and gating mechanisms that regulate the flow of information. The memory cell stores long-term information, while the input, forget, and output gates determine how much information is added, retained, or released at each time step.

In a LSTM unit, the input gate controls the extent to which new information enters the memory cell, the forget gate determines how much historical information is preserved, and the output gate regulates how the internal state contributes to the next layer. This structure enables LSTMs to effectively capture long-range dependencies within sequential data. The structure of a LSTM unit is presented in Fig. 2.

Fig. 2The structure of a LSTM unit

The training of an LSTM model usually adopts the Back Propagation Through Time (BPTT) algorithm, which is similar to the classical Back Propagation (BP) algorithm. When dealing with long sequences, BPTT is prone to gradient vanishing or gradient explosion, and LSTMs alleviate this problem effectively through the gating mechanism.

Fig. 3Flowchart of the BiLSTM neural network

BiLSTM is an extension of the standard LSTM. A BiLSTM layer consists of two LSTM layers: one processes the sequence in the forward direction from beginning to end, and the other processes it in the backward direction from end to beginning. At each time step, the hidden states from the forward and backward LSTMs are concatenated or summed to obtain the final hidden state. While a standard LSTM can only use past information, a BiLSTM can exploit both past and future information, and therefore provides better performance for modeling and predicting time-series signals. A schematic diagram of the BiLSTM layer structure is shown in Fig. 3.

3.1. Signal decomposition and trend extraction

The raw data sequence is first processed by the trend-extraction module to obtain the trend component used for model training. The original time-series signal is expressed as:

where $x_i$ represents the data points in time step $i$. This time-series sequence serves as the fundamental dataset for training the model.

A suitable window length $L$ is selected to construct the trajectory matrix $X$, which is defined as follows:

where $K=N-L+1$. For a time-series of length $N$, $L$ is typically chosen to be proportional to the dominant periodicity but not exceeding $N/2$ to ensure a sufficient number of lagged vectors ($K=N-L+1$). SVD is then applied, yielding the following factorization:

where $U$ is the matrix of left singular vectors, $\mathrm{\Sigma }$ is the diagonal matrix containing the singular values, and $V^T$ is the matrix of right singular vectors. After SVD, $X$ can be expressed as the sum of several rank-one matrices:

where ${\sigma }_i$ denotes the singular value, $U_i$ and $V_i^T$ are the $i$th columns of the left and right singular vector matrices, respectively.

To extract the trend component, the first one or several components with the largest singular values are typically selected for reconstruction. The reconstructed trajectory matrix is then converted back into the form of the original time series through diagonal averaging.

To enable automatic identification and selection of the principal components and further enhance the intelligence of SSA analysis, the K-means clustering algorithm is introduced to cluster the reconstructed components obtained from decomposition. For each reconstructed component, features such as variance and the value of the autocorrelation function are extracted to form feature vectors. These vectors are then clustered using the K-means algorithm, grouping components with similar statistical characteristics into the same cluster. The number of clusters is set to $k=$ 2, corresponding to a “principal component cluster” and a “noise cluster.” Components representing dominant trends or significant periodic structures are assigned to the “principal component” cluster, whereas components with low energy and lacking periodic structure are categorized as “noise.” Finally, by reconstructing the components belonging to the principal cluster, the useful information in the time series is effectively preserved, achieving automatic extraction of trend and periodic features while suppressing noise. This binary choice directly implements the theoretical separation in SSA between a few dominant components capturing signal dynamics and numerous low-energy components representing noise, without introducing unnecessary complexity or requiring subjective threshold selection.

The training and prediction module are responsible for model training and iterative updating of network parameters. First, the loss between the theoretical output and the actual output is computed. Then, the Adam optimizer is used to update the model parameters, and convergence is accelerated through early stopping and learning rate scheduling mechanisms. As iterative training progresses, the model parameters are continuously adjusted to minimize the loss function, thereby gradually improving prediction performance.

3.2. Data preprocessing

The data preprocessing module prepares the input sequence for neural-network processing. The extracted trend sequence is divided into training and testing subsets according to a predefined ratio. Specifically, the first portion of the sequence is used as the training set:

while the remaining portion is used as the test set:

3.3. State prediction

The constructed CNN-BiLSTM neural network model is employed to process the time-series data and generate prediction results. First, a one-dimensional convolution (Conv1D) layer is applied to the input time-series signal to capture local temporal features. After convolution, a ReLU activation function is used to introduce nonlinearity and preserve the nonlinear characteristics of the signal. The resulting feature maps are then passed to a MaxPooling1D layer, which performs downsampling to reduce dimensionality and retain the most salient information.

The extracted feature maps are subsequently fed into Bidirectional Long Short-Term Memory (BiLSTM) layers. The BiLSTM architecture captures bidirectional temporal dependencies, enabling the network to learn relationships from both past and future contexts. To prevent overfitting, a Dropout layer is incorporated after the BiLSTM layer, randomly discarding a portion of neurons to enhance generalization capability. Another BiLSTM layer followed by a Dropout layer is then applied to further strengthen the model’s ability to learn complex temporal patterns. Finally, a fully connected layer maps the output of the BiLSTM network to the prediction space.

4.1. Experimental environment

Python was used as the programming language for the experiments, and the neural network algorithms, including LSTM and BiLSTM, were implemented using the TensorFlow framework. The experiments were conducted on a computation station with 11th Gen Intel(R) Core(TM) i7-11700 processor and 16 GB of RAM.

4.2. Dataset

The dataset used in this study consists of time-series monitoring signals collected from gas storage compressors. The experimental validation of the proposed BiLSTM model was conducted using monitoring signals collected from operational gas storage compressors. The collected dataset includes multiple types of time-series measurements, such as temperature-time signals, displacement-time signals, and pressure-time signals. Two representative displacement signals were selected for model validation and comparison: the motor non-drive-end terminal displacement (MNDTD) signal, as shown in Fig. 5(a), and the motor drive-end terminal displacement (MDTD) signal, as shown in Fig. 5(b). Both MNDTD and MDTD are time-series signals that reflect the displacement variations of the drive-end and non-drive-end terminals during compressor operation. The monitoring sequences contain approximately 130,000 data points per channel, collected over a continuous operational period at a sampling rate of 2 kHz, covering multiple operating phases including startup, steady-state operation, load variation, and shutdown. This diversity of signal types and operating phases demonstrates the richness of the dataset, and while this study focuses on displacement signals as the primary indicators of mechanical condition, the proposed method can be readily extended to other signal types.

4.3. Data preprocessing

The historical MNDTD and MDTD monitoring sequences contain approximately 130,000 data points. To improve computational efficiency while preserving representative system behavior, an 8,000-point segment was selected for analysis. Since the region around the 20,000th sampling point reflects a stable displacement state during compressor start-up, the segment beginning at the 19,000th point was extracted. The selection was to avoid the initial startup transients (approximately the first 15,000-18,000 points) and focus on the region where the compressor has reached relatively stable operating conditions, as evidenced by the consistent vibration amplitude range in Fig. 5(a). This ensures that the model learns patterns representative of normal operational behavior rather than transient startup dynamics. This interval avoids the noise-dominated regions at the beginning and end of the signal and retains the essential dynamic characteristics relevant for forecasting. The red circle in Fig. 5(a) marked the extracted segment MNDTD.

In K-means clustering, we set number of clusters $k=$ 2. Convergence is achieved when either: (a) the change in cluster centroids between consecutive iterations falls below 10^-4, or (b) the maximum number of iterations (set to 300) is reached. This ensures stable and reproducible clustering outcomes.

In preparing the experimental data, it is essential to determine an appropriate window length for the SSA decomposition, as the window size directly influences both the number and the characteristics of the decomposed components, thereby affecting the model’s ability to capture meaningful trend and periodic structures. $L=$ 100 is selected as it provides adequate resolution to capture cycle-related periodic behavior while maintaining statistical reliability $L<N/2$. This choice enables effective separation of the trend component from background noise, as validated by the subsequent prediction performance. Meanwhile, quantitative indicators such as Mean Squared Error (MSE) and Mean Absolute Percentage Error (MAPE) were employed to evaluate the discrepancy between the reconstructed and original signals.

The model handles different categories of input data in a consistent manner. When only a single sequence variable (e.g. MNDTD) is used, each prediction step yields a single output value. When multiple types of monitoring data are incorporated, such as displacement, temperature, and humidity, the model is capable of learning cross-dependencies among variables and producing multi-output predictions simultaneously. This enables the framework to adapt to both univariate and multivariate forecasting scenarios while leveraging the complementary information contained in heterogeneous monitoring signals.

Fig. 5Schematic diagram of original vibration signals for MNDTD and MDTD

a) MNDTD

b) MDTD

4.4. Experimental results analysis

The training process begins by feeding the data into the model. Before conducting full experiments, the model structure must be optimized. In this study, a controlled-variable approach is used to adjust and compare the effects of the time-window length and the number of BiLSTM layers on model performance.

First, with the network structure fixed at two BiLSTM layers, different window sizes (10, 20, and 50) were evaluated. The results are shown in Table 1. When the window size was set to 20, the model achieved the lowest MAPE and RMSE values on the test set, indicating the highest prediction accuracy. Therefore, a window size of 20 was selected for subsequent experiments.

The neural network input window of 20 time steps, selected through the empirical optimization in Table 1, also has a physical basis. With the same sampling rate of 10 kHz and rotational frequency of 50-100 Hz, a 20-point window corresponds to 10 % to 20 % of a full mechanical rotation. This window size captures local vibration patterns such as instantaneous acceleration, velocity changes, and short-term fluctuations within a fraction of a rotation. It provides sufficient temporal resolution to detect meaningful local dynamics while remaining short enough to enable rapid response to changing compressor conditions.

After determining a window size of 20, a data segment with a length of 1000 points was extracted to evaluate the effect of BiLSTM depth. Models with 2, 3, 4, and 5 BiLSTM layers were tested, with each BiLSTM layer followed by a Dropout layer with a dropout rate of 20 %. The results are presented in Table 2.

It can be observed that increasing the number of BiLSTM layers does not lead to significant improvements in MAPE, RMSE, or $R^2$. This indicates that deeper architectures bring limited benefits once sufficient temporal features have been learned and may even introduce a risk of overfitting. Therefore, two BiLSTM layers were selected as the final model configuration.

With a window size of 20 and a BiLSTM depth of two layers, the model training was conducted. The variation of the loss values during training is shown in Fig. 6. As the number of iterations increased, the loss values on both the training and validation sets steadily decreased and gradually converged to a stable level after approximately 20 iterations.

Fig. 6Iteration loss plot of the SSA-CNN-BiLSTM model

The prediction results after model training are shown in Fig. 7(a) illustrates the model’s performance on a sequence of 2,500 data points, while 7(b) provides a zoomed-in view of the first 100 data points for detailed comparison.

After 500 training iterations, the model achieved MAPE values below 2% on both the training and test sets, RMSE values below 0.5, and $R^2$ values above 0.8, demonstrating strong predictive performance. The results are shown in Table 3.

To further evaluate model performance, the proposed SSA-CNN-BiLSTM model was compared against three alternative approaches. These included a CNN-BiLSTM model without trend decomposition, a BiLSTM-only model, and a traditional LSTM model. In the comparative experiments, all CNNs we applied shared the same architecture: one convolutional layer, each followed by a max-pooling layer, with 64 filters per layer, a kernel size of 3, and a stride value of 1.

Since missing values or abnormal points may occur during the acquisition and transmission of experimental data, and noise components are also present in the signals, a moving average method is applied for noise reduction when SSA is not used to extract noise. The formula of the moving average method is as follows:

where $w$ denotes the window size and $i$ represents the index of the current data point.

Fig. 7Prediction results of the SSA-CNN-BiLSTM model

a)

b)

To apply the moving average method, an appropriate window size must first be selected. The window size determines how many neighboring points are used to compute each averaged value. For each point in the sequence, the average of the surrounding points (including itself) within the selected window is calculated and used as the smoothed value of that point. In this study, the window size is set to 3, meaning that every three points are averaged.

The training results of each model are presented in Table 4. The LSTM model achieves a lower MAPE on the training set (1.2826 %) than the CNN-BiLSTM and BiLSTM models. However, its $R^2$ value on the test set (0.7562) is relatively low, indicating insufficient generalization ability. Further verification shows that under complex operating conditions, the LSTM model struggles to consistently capture the nonlinear dynamic behavior of the signal. The CNN-BiLSTM and BiLSTM models exhibit relatively lower overall performance, but still demonstrate certain predictive capability. The CNN-BiLSTM model obtains a MAPE of 1.3473 % on the test set, while the BiLSTM model achieves 1.3307 %. Although their $R^2$ values are comparatively lower, both remain within an acceptable range, suggesting that these models retain adequate capability in feature extraction and trend representation. In contrast, the SSA-CNN-BiLSTM model achieves the lowest MAPE values on both the training and test sets (1.2391 % and 1.2414 %, respectively), indicating the smallest relative prediction error. Meanwhile, the model achieves the highest $R^2$ value on the test set (0.8159), showing the best agreement between predicted and actual values among all models. In addition, the RMSE of this model is also relatively low, further confirming its superior performance in complex time-series prediction tasks.

From the comparison results, it can be observed that the convolution operations in CNN enable the extraction of local patterns and non-stationary features, allowing the CNN-BiLSTM model to outperform the standalone BiLSTM model. This enhances the model’s feature representation capability and improves its understanding of local dependencies within the signal.

On the other hand, SSA plays a significant role in the overall model. By performing trend extraction and noise separation, SSA effectively removes most noise components while preserving essential nonlinear variations, resulting in a smoother and more informative signal. This process substantially improves the model’s sensitivity to the underlying data and enhances prediction accuracy.

The experimental results indicate that the proposed SSA-CNN-BiLSTM model can accurately and effectively predict the operating state of gas storage compressors.

The experimental design in Table 4 also serves as an ablation study to evaluate the individual contribution of each component. Removing SSA by comparing SSA-CNN-BiLSTM against CNN-BiLSTM increases test MAPE from 1.24 % to 1.35 %, a degradation of 0.11 percentage points. This confirms that SSA preprocessing plays a critical role in noise suppression and trend extraction. Adding CNN without SSA preprocessing, as seen by comparing CNN-BiLSTM with BiLSTM, yields minimal improvement with MAPE values of 1.35 % versus 1.33 %. This suggests that feature extraction alone is insufficient without proper denoising. Interestingly, BiLSTM achieves a MAPE of 1.33 %, which underperforms compared to standard LSTM at 1.27 %, indicating that bidirectional processing is not always beneficial for vibration time-series data. Nevertheless, the full SSA-CNN-BiLSTM model achieves the best overall performance with a MAPE of 1.24 %, demonstrating that the synergistic combination of all components rather than any single component is responsible for the observed improvements.

To rigorously validate the performance improvements of the proposed SSA-CNN-BiLSTM model, we conducted statistical significance testing using paired t-tests across 10 independent runs with different random initializations. As shown in Table 5, the proposed model achieves statistically significant improvements over the CNN-BiLSTM and BiLSTM baselines ($p<$ 0.01), confirming that the integration of SSA preprocessing and the CNN feature extractor provides genuine performance benefits beyond random variation. The improvement over the LSTM baseline, while positive, does not reach statistical significance at the 95 % confidence level ($p=$ 0.087), indicating that the performance difference is modest for this particular comparison. These results provide strong statistical evidence for the effectiveness of our proposed approach while also highlighting areas for future investigation.

To provide a more comprehensive performance context, we additionally compared our proposed model against Support Vector Regression (SVR), a representative non-deep learning method widely used for time-series prediction. The SVR model was implemented with an RBF kernel, and hyperparameters $(C, \gamma , \epsilon )$ were optimized via grid search with 5-fold cross-validation. The optimal configuration ($C=$100, $\gamma =$0.01, $\epsilon =$ 0.1) was then evaluated on the same test set used for the deep learning models.

The optimized SVR model achieved a test MAPE of 1.52 %, with corresponding RMSE of 0.61 and $R^2$ of 0.71. This performance is notably weaker than all deep learning approaches, which achieved MAPE values ranging from 1.24 % to 1.35 %. Specifically, the proposed SSA-CNN-BiLSTM model outperforms SVR by a substantial margin of 0.28 percentage points in MAPE, representing an 18.4 % relative improvement. These results confirm that deep learning methods are better suited for capturing the complex nonlinear patterns in compressor vibration signals, while also further validating the effectiveness of our proposed approach.

From a practical deployment perspective, the proposed SSA-CNN-BiLSTM model demonstrates strong potential for real-time industrial applications. With a total inference time of approximately 8.5 ms per prediction on standard hardware, the model can easily meet the real-time requirements of compressor condition monitoring, where predictions at 1-2 second intervals are sufficient. For long-term deployment, model accuracy can be maintained through periodic retraining (e.g., monthly) or performance-triggered updates when prediction errors exceed predefined thresholds. Regarding edge computing feasibility, the model's compact size (~1.2 MB) enables deployment on resource-constrained devices such as the Raspberry Pi 4 (inference time ~45 ms) or NVIDIA Jetson Nano (~12 ms), both well within practical latency limits. Future work will explore model quantization and pruning techniques to further optimize deployment on industrial edge devices, as well as incremental learning strategies for continuous model adaptation without full retraining.