Noise Reduction Based on Improved Variational Mode Decomposition for Acoustic Emission Signal of Coal Failure

Featured Application

This study proposes an enhanced variational mode decomposition method using the beetle antennae search optimizer algorithm. Simulation signal analysis demonstrates the effectiveness and superiority of the proposed method. Additionally, the study highlights that the multifractal parameters of the acoustic emission frequency-domain spectrum can serve as a reference for the early warning of coal damage and instability.

Abstract

Acoustic emission (AE) signal processing and interpretation are essential in mining engineering to acquire source information about AE events. However, AE signals obtained from coal mine monitoring systems often contain nonlinear noise, limiting the effectiveness of conventional analysis methods. To address this issue, a novel denoising approach using enhanced variational mode decomposition (VMD) and fuzzy entropy is proposed in this study. The denoised AE signal’s spectral multifractal features are analyzed. The optimization algorithm based on VMD with a weighted frequency index is introduced to avoid mode mixing and outperform other decomposition methods. The characteristic parameter Δα of the AE spectral multifractal parameter serves as an early warning indicator of coal instability. These findings contribute to the accurate extraction of time–frequency features and provide insights for on-site AE signal processing.

1. Introduction

Acoustic emission (AE) is a crucial technique in structural health monitoring. Many studies have reported that the energy accumulated in coal and rock material under loading can be partly released in the form of AE [1,2,3]. Through the processing and analysis of the AE signal, the damage degree and damage mechanism inside the material can be revealed [4,5,6,7], which can serve as a reference for analyzing on-site AE signals, and further predict the structural life of coal mass [8].

However, the complex underground monitoring environment often introduces significant environmental noise during AE signal transmission. These noises have obvious nonlinear characteristics. These nonlinear noises hinder AE event identification and source localization and reduce the accuracy of coal and rock mass disaster prediction. Therefore, it is of great significance to develop a new noise reduction method for the structural life prediction of coal mass.

The AE signal can be viewed as a composition of multifrequency information [9]. By decomposing the AE signal, different frequency components can be obtained, allowing for the more precise extraction of time–frequency characteristics after noise frequency filtering. Numerous methods have been developed for AE noise reduction, including fast Fourier transform [10], wavelet decomposition [11], wavelet packet decomposition [12], empirical mode decomposition (EMD) [13], and local mean decomposition [14]. However, fast Fourier transform is limited to stationary signal analysis, providing only frequency components within a signal segment and demonstrating poor performance with nonlinear signals [15]. Wavelet decomposition and wavelet packet decomposition can better deal with nonstationary signals, but it is necessary to select a suitable wavelet basis before signal decomposition. Although the EMD and local mean decomposition methods are adaptive techniques without the aforementioned limitations, they still encounter mode mixing and end effects in practical applications.

Mode mixing refers to the blending and interference of different intrinsic mode functions (IMFs) during the decomposition process. It occurs when the characteristics of a signal cannot be accurately separated into distinct IMFs due to overlapping or mutual influence, complicating signal analysis and interpretation. End effects occur at the boundaries of the signal during decomposition, leading to distortions and instability in the generated IMFs. Insufficient data near the signal edges exacerbate these effects, especially at the beginning and end of the signal [16].

To address these challenges, variational mode decomposition (VMD) [17] has emerged as an adaptive, non-recursive signal processing approach. VMD overcomes mode mixing and end effects encountered in EMD while reducing nonstationarity in complex and highly nonlinear time series. In recent years, VMD has been gradually applied in the field of signal processing. For instance, An and Yang [18] employed VMD to denoise the vibration signal from a hydropower unit, utilizing approximate entropy [19] to filter out the useful components, demonstrating its superiority over wavelet analysis. Yang et al. [20] concluded that VMD outperforms EMD in studying the wind speed influence by analyzing wind turbine states. Zhang et al. [21] compared VMD and EMD in bearing defect feature extraction from rolling bearing simulation signals, affirming VMD’s accurate extraction of the principal bearing fault modes and superior performance over EMD. Zhang et al. [22] employed VMD to identify microseismic signals and blasting signals, successfully separating each component and accurately identifying the target signals. Xue et al. [23] utilized VMD for seismic signal analysis, demonstrating higher spectral and spatial resolutions in the instantaneous spectra compared with EMD.

However, in the aforementioned studies, the selection of mode number k in VMD was often determined empirically and conveniently. For example, An and Yang [18] applied VMD to analyze simulation and blasting signals, without providing a basis for choosing k. Yang et al. [20] determined k based on the number of decomposition layers in EMD, which is an empirical approach. Yang et al. [24] used correlation coefficients and energy ratios to determine k, but this method is not suitable to AE signals and lacks physical meaning. The accurate determination of k in VMD is crucial as it greatly impacts the decomposition results and can lead to poor outcomes. Therefore, selecting an appropriate k holds significant importance in VMD.

Several methods have been proposed for selecting mode number k, including correlation coefficient [25,26], kurtosis index [27,28], and entropy [29]. In this study, we combine the advantages of AE signal frequency characteristics and correlation coefficient to construct a weighted frequency index as the objective for optimizing VMD parameters. Additionally, we employ a novel optimization algorithm called beetle antennae search algorithm (BAS) [30] to optimize VMD parameters. BAS is a simple method that has the advantages of a simple-principle, parameter-efficient search algorithm that significantly reduces computational requirements compared with the particle swarm optimization algorithm [31].

To address the difficulty in filtering out noise from AE signals, this paper proposes an improved VMD method for AE signal analysis and noise reduction. The remaining sections of this paper are organized as follows: Section 2 introduces the fundamental theory of VMD and fuzzy entropy; Section 3 provides a brief overview of the improved VMD-related algorithms; Section 4 verifies the proposed method using simulation signals; Section 5 applies the proposed algorithm to reduce AE signal noise and analyzes the multifractal characteristics of AE signals.

2. Basic Theory

This section introduces the calculation process of VMD and fuzzy entropy. Although the decomposition effects of EMD and VMD are compared in this paper, due to space limitations, the calculation process of EMD is not given in this section. The relevant calculation process can be found in [13].

2.1. Variational Mode Decomposition

VMD is a novel non-recursive signal decomposition algorithm that is particularly well suited to the analysis and processing of nonlinear and nonstationary signals. VMD can be calculated as follows:

$$\left\{\begin{array}{c}mi{n}_{\left\{{\omega }_k\right\},\left\{u_k\right\}}\left\{{\displaystyle \sum }_{k=1}^K\Vert {\partial }_t\left[\left(\delta \left(t\right)+\frac{j}{\pi t}\right)u_k\left(t\right)\right]e^{-j{\omega }_{k}t}{\Vert }^2\right\}\\ {\displaystyle \sum }_{k=1}^K{u}_k=f\left(t\right)\end{array}.\right.$$

(1)

where $\left\{u_k\right\}=\left\{u_1,u_2,u_3,\cdots ,u_k\right\}$ represent mode components, with k being the mode number; $\left\{{\omega }_k\right\}=\left\{{\omega }_1,{\omega }_2,{\omega }_3,\cdots ,{\omega }_k\right\}$ corresponds to the center frequency of $\left\{u_k\right\}$; and $f\left(t\right)$ denotes the original signal to be decomposed. To compute the solution of Equation (1), penalty factor α and Lagrange multiplier λ(t) are introduced, transforming the constrained problem into an unconstrained variational problem. The augmented Lagrange, L, can be described as follows:

$$L\left(\left\{u_k\right\},\left\{{\omega }_k\right\},\lambda \right)=\alpha {{\displaystyle \sum }}_{k=1}^K\Vert {\partial }_t\left[\left(\delta \left(t\right)+\frac{j}{\pi t}\right)\ast u_k\left(t\right)\right]e^{-j{\omega }_{k}t}{\Vert }_2^2+\Vert f\left(t\right)-{{\displaystyle \sum }}_{k=1}^K{u}_k\left(t\right){\Vert }_2^2+\left(\lambda \left(t\right),f\left(t\right)-{{\displaystyle \sum }}_{k=1}^K{u}_k\left(t\right)\right).$$

(2)

where α is utilized to ensure signal reconstruction accuracy, while λ is employed to maintain the strictness of the constraint conditions. According to Equation (2), the Lagrange expressions are sought by iteratively updating $u_k^{n+1}$, ${\omega }_k^{n+1}$, and ${\lambda }^{n+1}$:

$${\hat{u}}_k^{n+1}\left(\omega \right)=\frac{\hat{f}\left(\omega \right)-{{\displaystyle \sum }}_{i\ne k}{\hat{u}}_i\left(\omega \right)+\frac{\hat{\lambda }\left(\omega \right)}{2}}{1+2\alpha {\left(\omega -{\omega }_k\right)}^2},$$

(3)

$${\omega }_k^{n+1}=\frac{{{\displaystyle \int }}_0^{\infty }\omega {\left|{\hat{u}}_k\left(\omega \right)\right|}^{2}d\omega }{{{\displaystyle \int }}_0^{\infty }{\left|{\hat{u}}_k\left(\omega \right)\right|}^{2}d\omega },$$

(4)

$${\hat{\lambda }}^{n+1}\left(\omega \right)={\hat{\lambda }}^n\left(\omega \right)+\tau \left[f\left(\omega \right)-{{\displaystyle \sum }}_{k=1}^K{\hat{u}}_k^{n+1}\left(w\right)\right].$$

(5)

We update with repeated iterations until the iteration stopping condition is satisfied:

$${{\displaystyle \sum }}_{k=1}^K\frac{\Vert {\hat{u}}_k^{n+1}-{\hat{u}}_k^n{\Vert }_2^2}{\Vert {\hat{u}}_k^n{\Vert }_2^2}<\epsilon .$$

(6)

where n represents the number of iterations and τ denotes the iteration step length.

2.2. Introduction of Fuzzy Entropy

In this paper, fuzzy entropy [32] is used to determine the signal-to-noise boundary. Fuzzy entropy is described as shown below.

For each time series [h(1), h(2), …, h(N)], we reconstruct U space by introducing non-negative integers m (m ≤ N − 2), and the reconstructed signal is as follows:

$$S^m\left[i\right]=\left\{h\left(i\right),h\left(i+1\right),\dots ,h\left(i+m-1\right)\right\}-h_0\left(i\right),$$

(7)

where i = 1, 2, …, N − m + 1, and $h_0\left(i\right)$ can be defined as follows:

$$h_0\left(i\right)=\frac{1}{m}{{\displaystyle \sum }}_{j=0}^{m-1}h\left(i+j\right).$$

(8)

The membership function is as follows:

$$A\left(x\right)=\left\{\begin{array}{c}1, x=0\\ e^{\left[-ln\left(2\right){\left(\frac{x}{r}\right)}^2\right]},x>0\end{array},\right.$$

(9)

where r is similarity tolerance. According to Equation (9), the membership function is defined as follows:

$$A_{ij}^m=\exp \left(-\ln \ln \left(2\right)·{\left(\frac{d_{ij}^m}{r}\right)}^2\right).$$

(10)

We define distance $d_{ij}^m$ between two vectors $S^m\left[i\right]$ and $S^m\left[j\right]$ as the maximum absolute difference between the two corresponding elements, and the distance is given by the following:

$$d_{ij}^m=\left(\left|h\left(i+p-1\right)-h_0\left(i\right)\right|-\left|h\left(j+p-1\right)-h_0\left(j\right)\right|\right),$$

(11)

For each i value, the average value is obtained as follows:

$$C_i^m\left(r\right)=\frac{1}{N-m}{{\displaystyle \sum }}_{j=1,j\ne i}^{N-m+1}{A}_{ij}^m,$$

(12)

$${\phi }^m\left(r\right)=\frac{1}{N-m+1}{{\displaystyle \sum }}_{i=1}^{N-m+1}{C}_i^m\left(r\right).$$

(13)

Therefore, the fuzzy entropy of the sequence is defined as follows:

$$FuzzyEn\left(m,r,N\right)=ln{\phi }^m\left(r\right)-ln{\phi }^{m+1}.$$

(14)

Based on the findings of previous research [33], we set m = 2, r = 0.2, and N = 1000.

3. FCI-VMD Method Based on BAS

From the above descriptions, there are four parameters that require determination before signal decomposition: mode number k, penalty parameter α, noise tolerance ξ, and tolerance of convergence criterion ε. Among these parameters, ξ and ε have minimal impact on decomposition. ξ is typically set to zero to ensure effective convergence of the algorithm, while ε is set to a small positive value to ensure reconstruction accuracy. Conversely, mode number k significantly influences the decomposition outcome. If k is set too high, the signal is over-decomposed, resulting in different modes for signals at the same frequency and causing mode mixing. On the other hand, setting k too low leads to under-decomposition, where the same mode exhibits varying frequency characteristics. Thus, selecting an appropriate value for k holds paramount importance in signal decomposition. To address this issue, this paper proposes the FCI-VMD algorithm to optimize the selection of k.

3.1. Weighted Frequency Index

The AE spectrum of coal failure carries valuable information about the progression of fracture damage, with the main frequency being a crucial feature in the frequency domain. Hence, it is essential to consider the impact of the main frequency on the selection of mode number k during AE signal decomposition.

After the original signal, $f\left(t\right)$, is decomposed using VMD, the decomposed component ${\chi }_k\left(t\right)$ undergoes the Hilbert transform:

$${\hat{\chi }}_k\left(t\right)=\frac{1}{\pi }{{\displaystyle \int }}_{-\infty }^{\infty }\frac{{\chi }_k\left(\tau \right)}{t-\tau }d\tau .$$

(15)

The signal is defined as follows:

$$z\left(t\right)={\chi }_k\left(t\right)+j{\hat{\chi }}_k\left(t\right)=\lambda \left(t\right)e^{j\theta \left(t\right)},$$

(16)

The instantaneous amplitude is as follows:

$$\lambda \left(t\right)=\sqrt{{\chi }_k^2\left(t\right)+{\hat{\chi }}_k^2\left(t\right)}.$$

(17)

The instantaneous frequency is as follows:

$$\rho \left(t\right)=\frac{1}{2\pi }\frac{d\theta \left(t\right)}{dt}.$$

(18)

The main frequency of the m mode is ${\rho }_m$, ($2<m\le K)$, and the weighted frequency index, FCI, is defined as follows:

$$FCI=\frac{FI}{C},$$

(19)

$$FI=\left|{\rho }_m-{\rho }_{m+1}\right|.$$

(20)

$$C=\frac{E\left[\left(x-\overline{x}\right)\left(y-\overline{y}\right)\right]}{E\left[{\left(x-\overline{x}\right)}^2\right]E\left[{\left(y-\overline{y}\right)}^2\right]}.$$

(21)

where ${\rho }_m$ is the main frequency of the modes and FI represents the main frequency distance between the two modes. x and y represent adjacent modal components.

From Equations (19)–(21), we can know that a smaller FI suggests a higher likelihood of mode mixing when the main frequencies of two signals are close. Conversely, if FI is too large, further decomposition of the signal is necessary. Hence, finding a suitable FI helps avoid these issues. In this study, the AE frequency range during coal sample failure was primarily distributed between 100 kHz and 300 kHz. The decomposition results were relatively stable when the minimum FI value exceeded 10 kHz, which was chosen as the threshold in this paper. The correlation coefficient, C, can assess the degree of overlap between modes and provide valuable insights. A higher correlation coefficient indicates a greater overlap between modes. Considering the advantages and limitations of both indices, this paper adopts the weighted frequency index as the objective function for optimizing the k value.

3.2. Beetle Antennae Search

The beetle antennae search algorithm (BAS) is an intelligent algorithm that imitates beetle foraging [30]. In nature, beetles rely on their antennae to detect chemical signals. If the right antennae detect a stronger smell compared with the left antennae, the beetle chooses to fly towards the right in the subsequent step; otherwise, it flies towards the left. By mimicking this behavior, BAS performs a random search in the vicinity based on an objective function to be optimized. Notably, BAS utilizes a single beetle, enabling it to achieve high search speed. The specific steps of the BAS algorithm are reported below.

Assuming that the position of the beetle in the D-dimensional space is $X_D=\left(x1,x2,x3,\dots ,xD\right)$, the antennae on the left and right sides of the beetle can be defined by the following formulas:

$$\left\{\begin{array}{c}{X}_r=X_i+l\overline{d}\\ X_l=X_i-l\overline{d}\end{array}\right.,$$

(22)

where $X_i$ is the current position of the beetle, $X_r$ and $X_l$ denote the positions of these two antennae, l is the sensing length of antennae, and $\overline{d}$ represents the orientation of the beetle in D space:

$$\overline{d}=\frac{rands\left(n,1\right)}{\Vert rands\left(n,1\right)\Vert },$$

(23)

where n is the domain dimension and $rands\left(·\right)$ represents a random function.

The update location of the beetle is described as follows:

$$X^{i+1}=X^i+{\delta }_i\overline{d}sign\left(\Phi \left(X_r\right)-\Phi \left(X_l\right)\right).$$

(24)

where sign is the sign function, $\Phi \left(·\right)$ is the objective function, ${\delta }_i$ is the step size, and i is the number of iterations. The beetle continuously updates the position with Equation (24) until the objective function finally satisfies the optimal solution.

3.3. Proposed Method

The FCI-VMD method proposed in this paper takes the weighted frequency index as the objective function. The signal–noise separation point is determined by the fuzzy entropy. Since BAS is an algorithm for finding the minimum value, according to Equation (24), the parameter optimization process can be described as follows:

$$\left\{\begin{array}{c}\Phi \left(X_r\right),\Phi \left(X_l\right)=\min \left\{FCI\right\}\\ k=3, 4,\dots , 10\end{array}\right..$$

(25)

where FCI is the objective function to be optimized and k is the parameter to be optimized for VMD. Generally speaking, k is an integer between 3 and 10. The specific steps are as follows:

(1)

Set the VMD initial parameters of the original signal, $f\left(t\right)$, and set the initial parameters of BAS, such as the range of parameter k, number of iterations, and iteration step size.

(2)

Decompose the AE signal using VMD, and calculate the $\Phi \left(X_r\right)$ and $\Phi \left(X_l\right)$ of all modes.

(3)

Determine if an iteration termination condition is reached: If $\Phi \left(X_r\right),\Phi \left(X_l\right)>FCI$, then k = k + 1, and continue the iteration; else, k = k − 1, and stop the iteration.

(4)

Decompose the signal again with the optimized parameters.

(5)

Determine the antennae for the decomposed signal by calculating the fuzzy entropy; then, reconstruct the effective signal.

4. Analysis of the Simulation Signal

To verify the effectiveness of the method proposed in this paper for denoising AE signals, simulation signals with different frequencies are used to verify the proposed method. The mathematical model of AE constructed by [34] is expressed as follows:

$$f\left(t\right)={{\displaystyle \sum }}_{i=1}^n{A}_i{e}^{\left[-P_i{\left(t-t_i\right)}^2\right]}\sin \left[2\pi {\rho }_m\left(t-t_i\right)\right].$$

(26)

where $A_i$, $P_i$, $t_i$, and ${\rho }_m$ are the amplitude, attenuation coefficient, delay time and dominant frequency of the superimposed signal of i, respectively. β is the number of superimposed signals in the model.

The AE signal is obtained by using the above model, and its parameters are shown in

Table 1. AE simulation parameters.

β	${\mathit{A}}_{\mathit{i}}$	${\mathit{P}}_{\mathit{i}}$	${\mathit{t}}_{\mathit{i}}$/ms	${\mathit{\rho }}_{\mathit{m}}$/kHz
1	2	6 × 108	0.4	70
2	2	8 × 108	0.6	40
3	2	7 × 108	0.8	110

.

To simulate the AE signal more realistically, white noise with signal-to-noise ratios (SNRs) of 5 dB, 8 dB, and 10 dB is added to the simulation signal. This section only specifically analyzes the AE simulation signal with an SNR of 10 dB (Figure 1).

4.1. Results of Simulation Signal Decomposition

After determining the value of k using the FCI-VMD algorithm, the simulation signal is decomposed using VMD with a selected number of decomposition layers, namely, k = 4 (see Figure 2). The spectrogram reveals that IMF2, IMF3, and IMF4 correspond to the main frequencies ${\rho }_3$ = 110 kHz, ${\rho }_2$ = 70 kHz, and ${\rho }_1$ = 40 kHz, respectively. In contrast, the IMF1 component exhibits a wider frequency band with high frequency and low energy, indicating it as the noise signal. To verify the effectiveness of the proposed optimization algorithm, we also decompose the signal with k = 5 using VMD, and the resulting components are shown in Figure 3. Notably, the main frequencies of the IMF2 and IMF3 components in Figure 3 are found to be relatively close, indicating the occurrence of mode mixing. Consequently, the optimal number of decomposition layers, k, for the signal is determined to be 4.

To evaluate the advantages of the FCI-VMD method, a comparison is made with two alternative approaches: EMD and EVMD. The decomposition results obtained using EMD are depicted in Figure 4. As illustrated in the figure, IMF2 to IMF5 all exhibit varying degrees of mode mixing, making it challenging to distinguish between effective modes and noise modes. On the other hand, EVMD (empirical variational mode decomposition) is a variational mode decomposition method in which mode number k is empirically determined. Figure 5 showcases the decomposition results obtained using EVMD, clearly indicating excessive decomposition of the signal.

4.2. Signal-to-Noise Separation with Fuzzy Entropy

Entropy serves as a metric for assessing the complexity of a signal. with higher entropy indicating greater complexity. Using VMD decomposition, the AE signal is separated into several components with distinct characteristic frequencies. A higher entropy value signifies a more chaotic nature of the signal. In comparison to the effective signal, the AE noise signal exhibits heightened complexity and uncertainty, resulting in higher entropy. Consequently, in this section, fuzzy entropy is employed to eliminate noise components from the AE signals.

The fuzzy entropy of different modes obtained using EMD, EVMD, and FCI-VMD is computed and illustrated in Figure 6. From Figure 6, it can be observed that EMD yields larger fuzzy entropy values for IMF1 and IMF2. VMD, on the other hand, produces a larger fuzzy entropy value for the IMF1 component. Notably, EVMD introduces more noise components due to signal over-decomposition. To remove noise, modes with fuzzy entropy exceeding 0.2 are eliminated, while the remaining IMF components are retained as optimal components. The denoised simulation AE signal is subsequently reconstructed based on these optimal components.

Figure 7 depicts the denoising results obtained with the EMD method applied to a noise-contaminated analog AE signal. It can be observed that the 110 kHz frequency signal is filtered out, leading to waveform distortion. This indicates that filtering the complex simulation AE signal using EMD decomposition compromises the integrity of the effective signal. Conversely, upon decomposing the simulation AE signal using FCI-VMD, the desirable frequencies of 40 kHz, 70 kHz, and 110 kHz are successfully extracted, allowing for the separation of noise and effective signals.

4.3. Evaluation Indicators

To compare the noise reduction effect of these three methods, the percentage of signal energy after noise reduction to the original signal energy (ESN) and the root mean square error (RMSE) are used as noise reduction evaluation indicators. ESN is defined as follows:

$$ESN=\frac{E^{\prime }}{E}.$$

(27)

where E is the total energy of the original signal and $E^{\prime }$ is the energy of the signal after noise reduction. The RMSE can be defined as follows:

$$RMSE=\sqrt{\frac{1}{N}{\left(f\left(t\right)-f{\left(t\right)}^{\prime }\right)}^2}.$$

(28)

The results of EMD, EVMD, and FCI-VMD denoising evaluation indicators for signals with different SNRs are shown in

Table 2. Evaluation index of three noise reduction methods.

Method	Metrics	Simulation AE Signal
		SNR = 5	SNR = 10	SNR = 8
EMD	ESN	0.7423	0.6490	0.7423
	RMSE	0.3254	0.2452	0.3154
EVMD	ESN	0.9047	0.9309	0.9122
	RMSE	0.0628	0.0942	0.0621
FCI-VMD	ESN	0.9273	0.9456	0.9742
	RMSE	0.0412	0.0878	0.0251

. It is generally believed that the smaller the RMSE and the larger the ESN are, the closer the denoised signal is to the real signal. Compared with EMD and EVMD, after denoising the noisy analog AE signal with the FCI-VMD method, the signal retains a higher proportion of energy, and the RMSE is also smaller. This shows that at different noise levels, FCI-VMD is better than EMD and EVMD in processing AE signals.

5. Experiments

Based on the above simulation signal analysis and verification of the effectiveness of the proposed optimization algorithm, we carried out a uniaxial compression test on coal samples to investigate the AE characteristics in different stress stages. Following the denoising of the AE signal, we further extracted the spectral characteristics using the Hilbert transform and explored the multifractal properties of the spectrogram.

5.1. Uniaxial Compression Test

The experimental system includes a loading system and an AE system. An MTS C45.104 testing machine is used in the loading system. The AE instrument utilizes a PCI-Express8 multi-channel AE system manufactured by American Physical Acoustics Company, along with Nano 30 miniature sensors. The sensor has a frequency range of 125~750 kHz, and it is equipped with a 2/4/6 voltage preamplifier that allows for a selectable range of 20, 40, and 60 dB. Before the test, the sampling rate is set to 1 Msample/s. The tests are displacement-controlled at a rate of 1 mm/min.

All the coal samples tested are taken from Wanglou Coal Mine, Shandong Province, China. All samples are prepared in accordance with the ISRM-suggested shape and size and are processed into cylinders with a diameter of 50 mm and a height of 100 mm. Five samples are prepared, numbered wl-2, wl-3, wl-5, wl-6, and wl-8. The stress–strain curve is shown in Figure 8, and the basic physical and mechanical parameters are shown in

Table 3. Physical and mechanical properties of coal samples.

Number	Density (g/cm3)	P-Wave Velocity (m/s)	Uniaxial Compressive Strength (MPa)	Young’s Modulus (GPa)
wl-2	1.45	2023.79	12.15	2.59
wl-3	1.37	2050.09	12.35	2.48
wl-5	1.37	1550.46	9.7	2.12
wl-6	1.32	1981.38	8.67	2.17
wl-8	1.37	1532.41	6.25	1.26

. It can be observed that the uniaxial compressive strength of specimens wl-5, wl-6, and wl-8 is relatively low. This is because coal is a material with significant heterogeneity, containing many inherent fractures and voids. Additionally, for coal samples with significant differences in porosity, the uniaxial compressive strength tends to vary significantly.

The stress–AE count–cumulative AE count curve of coal under uniaxial compression is presented in Figure 9. In the initial stages of low stress, the coal sample undergoes pore compaction, resulting in a relatively small number of AE counts. As the stress gradually increases, microcracks begin to form in the coal sample, leading to a continuous increase in AE counts. When the coal sample approaches its peak strength, a significant number of dense, high-amplitude signals are generated, indicating the occurrence of instability failure. The above observations demonstrate that AE signals can effectively reflect the internal damage and propagation of microcracks in coal. By processing and analyzing the collected AE signals, it becomes possible to assess the stress stage of coal damage, offering a valuable reference for early warning systems pertaining to coal damage and instability.

5.2. AE Signal Noise Reduction

To reduce the computational load associated with the high sampling frequency and large number of AE signal data, 100 sets of data are selected at equal intervals from the time sequence. Each data set consists of 1024 sampling points. The FCI-VMD algorithm is applied to decompose each set of signals, simultaneously calculating the fuzzy entropy of the resulting IMF components. Figure 10 illustrates the decomposition results and the corresponding fuzzy entropy values. In the figure, σ represents stress; σ_c represents peak stress; and σ/σ_c represents different stress stages. It can be observed that the optimal decomposition level for the AE signal is 5, and there is a significant variation in fuzzy entropy among the IMF components. The fuzzy entropy values are arranged in descending order, with IMF1 exhibiting the highest entropy, and IMF5, the lowest. This suggests that IMF1 contains a greater proportion of noise signal components, and each layer of IMF components contains noise modes. In this study, components with fuzzy entropy exceeding a certain threshold are filtered out, while the remaining components are reconstructed to obtain the effective signals. This approach avoids filtering out only IMF1 or selectively filtering out some effective signals, thereby preventing the presence of residual noise in the reconstructed signals. Based on statistical analysis of the fuzzy entropy values of each IMF component, it is observed that the IMF component displays relative discreteness when fuzzy entropy exceeds 0.015. Hence, a threshold of 0.015 is set in this study.

After filtering the noise signal and reconstructing the effective signal, the signal after noise reduction outcomes is shown in Figure 11.

Figure 11 demonstrates that the AE signal waveform becomes smoother, and the high-frequency noise is effectively eliminated after denoising using FCI-VMD. The improvement is particularly noticeable in the frequency domain, where the elimination of noise is evident. This indicates that the improved VMD approach not only filters out the noise signal but also preserves the integrity of the effective signal without distortion or loss.

5.3. Time–Frequency Feature Extraction of the AE Signal

Time–frequency analysis is a signal analysis method that transforms time-domain signals into time–frequency images. It describes the spectral components of the signal at different times, so it contains more information [35]. In this paper, the Hilbert transform is used to convert the denoised one-dimensional time-domain signal into a three-dimensional time–frequency image, and the spectrum of typical AE signals in the four stages is selected, as shown in Figure 12.

Figure 12 shows the distinct spectral characteristics of the AE signal in different stages of coal uniaxial compression. The main frequency range of AE during coal failure is predominantly between 100 and 300 kHz. During the initial 20% ${\sigma }_c$ of the test, the main frequency band is narrow, centered around approximately 296 kHz. As the sample continues to be loaded, the dominant frequency of AE spans between 100 and 300 kHz at 50% ${\sigma }_c$ and 80% ${\sigma }_c$, with a longer duration. When the sample reaches peak strength, the main frequency is approximately 92.5 kHz. At this time, the main frequency bandwidth decreases, and the amplitude increases. Therefore, the AE frequency bandwidth gradually widens during the coal loading process, but the main frequency bandwidth suddenly decreases, and the main frequency decreases when the coal is destroyed. There are differences in the spectral characteristics of AE in the four typical stress stages. Thus, we can try to extract the spectral characteristics of AE by using the multifractal spectrum.

5.4. AE Multifractal Spectrum

The fractal dimension can be used to describe nonlinear signals and images, and it is an effective tool to describe complexity and inhomogeneity. Complex systems generally have self-similar characteristics. In the critical instability stage of disordered media, there is implicit complexity in space and time, generally forming multiscale features and fractal structures. The spectral distribution of AE is relatively complex, and simple fractal dimensions alone cannot accurately describe the fractal characteristics of the acoustic emission spectrum. In this section, the box counting method [36] is used to calculate the multifractal spectrum of the AE spectrum. First, the spectrogram is converted into a binary image; then, the grid with side length r is divided into grids of the same size. Then, each grid with side length r continues to be divided into small squares with grids with side length α. If there are pixels in the small square, we define the square as black; otherwise, it is white. We define the number of small black squares in each square with side length r as $N_{ij}$; then, the probability measure is as follows:

$$P_{ij}\left(r\right)=\frac{N_{ij}}{{{\displaystyle \sum }}_{ij}{N}_{ij}}.$$

(29)

In the scale-free self-similar region, probability measure $P_{ij}$ has the following relationship:

$$P_{ij}\left(r\right)~r^{{\alpha }_{ij}}.$$

(30)

where ${\alpha }_{ij}$ is the singularity index. We denote the number of all lattices with the same singular index ${\alpha }_{ij}$ in the spectrogram by $N_{{\alpha }_{ij}}$; then, the following relationship also exists in the scale-free self-similar interval:

$$N_{{\alpha }_{ij}}\left(r\right)~r^{-f\left({\alpha }_{ij}\right)}.$$

(31)

where f(αij) represents the probability of a in all singularity index sets, also known as the fractal dimension, when the singularity index is ${\alpha }_{ij}$. The multifractal spectrum is formed by plotting all $\left[{\alpha }_{ij}, f\left({\alpha }_{ij}\right)\right]$ in the Cartesian coordinate system, which can reflect the inhomogeneity of the acoustic emission spectrum. The width of the multifractal spectrum of Δα is defined as follows:

$$\Delta \alpha ={\alpha }_{max}-{\alpha }_{min}.$$

(32)

Spectral width Δα reflects the complexity of the entire AE spectrogram. The larger Δα is, the more complex the AE spectrogram and the wider the frequency band are. The relationship between the fractal dimension of the small and large probability subsets can be obtained with the following calculation:

$$\Delta f=f\left({\alpha }_{max}\right)-f\left({\alpha }_{min}\right).$$

(33)

where $\Delta f$ reflects the morphological characteristics of the spectrogram. According to multifractal theory, the multifractal spectra of coal in the four stress loading stages are obtained, as shown in Figure 13.

Figure 13 presents the multifractal spectra of AE spectral images obtained from coal samples at different stress levels: 20%${\sigma }_c$, 50% ${\sigma }_c$, 80%${\sigma }_c$, and 100% ${\sigma }_c$. The results demonstrate noticeable variations in the spectral width, Δα, of the multifractal spectrum across the different stress stages while exhibiting similar patterns among different coal samples. In the middle stage of loading, multifractal spectral parameter Δα is the largest, and in the failure stage, Δα reaches the minimum. Specifically, during the intermediate loading stage, parameter Δα exhibits the highest values, whereas it reaches its minimum during the failure stage. This behavior primarily reflects the progressive evolution of cracks induced by the coal loading process.

Table 4. The multifractal spectral characteristic parameters of AE spectral images in different stages.

Sample	State	Characteristic Parameter of Multifractal Spectrum
		αmin	f (αmin)	αmax	f (αmax)	Δα	Δf
Wl-2	20% σc	2.0795	2.0795	2.2445	1.9967	0.1650	−0.0828
	50% σc	2.0778	2.0778	2.2473	1.9987	0.1695	−0.0791
	80% σc	2.0770	2.0770	2.2555	1.9960	0.1785	−0.0810
	100% σc	2.0786	2.0034	2.2188	1.9901	0.1402	−0.0751
Wl-3	20% σc	2.0823	2.0823	2.2179	2.0081	0.1356	−0.0742
	50% σc	2.0809	2.0809	2.2373	1.9987	0.1564	−0.0822
	80% σc	2.0790	2.0790	2.2561	1.9907	0.1771	−0.0883
	100% σc	2.0824	2.0824	2.2102	2.0126	0.1278	−0.0698

summarizes the multifractal spectral parameters of different coal samples. According to Table 4, in the initial loading stage (20% ${\sigma }_c$), Δα = 0.1650. During this stage, the coal sample experiences gradual closure of cracks, resulting in low acoustic emission signal amplitude, high main frequency, and narrow main frequency bandwidth. With the continuous loading of the coal sample, at 50%${\sigma }_c$ and 80% ${\sigma }_c$, the presence of microfractures and of a limited number of larger fractures becomes more prominent. Consequently, the main frequency bandwidth of AE widens, and Δα values of 0.1695 and 0.1785 are observed, respectively. These stages exhibit a gradual widening trend in Δα. When the stress reaches its peak strength, macroscopic cracks develop, leading to subsequent instability failure. At this time, the AE amplitude reaches its maximum; the main frequency decreases; the main frequency bandwidth narrows; and Δα measures 0.1402. The order of Δα across different loading stages is as follows: 80% ${\sigma }_c$ > 50% ${\sigma }_c$ > 20% ${\sigma }_c$ > 100% ${\sigma }_c$. This order aligns with the trend observed in the main frequency of AE, which first increases and then decreases throughout the entire loading stage. These findings indicate that Δα is related to the main frequency of AE and can serve as an early warning indicator of coal damage and instability.

Another multifractal parameter, $\Delta f$, reflects the overall difference degree of the AE spectrogram. $\Delta f$ < 0 indicates that the grayscale image is dominated by a smaller grayscale, which is due to the low proportion of the AE spectral distribution in the entire image. With the continuous loading of coal samples, the AE frequency band changes from widening to narrowing, and the spectral band is the narrowest in the destruction stage. $\Delta f$ reaches its maximum value.

6. Conclusions

To address the challenge of effectively filtering emitted noise, this paper proposes the FCI-VMD algorithm based on variational mode decomposition and applies it to the study of indoor sound emission signal separation and vibration characteristics. Additionally, multiple fractal spectra are utilized to extract AE spectral features. The main findings are summarized as follows:

• A weighted frequency index is proposed as the objective function to optimize VMD decomposition level k. The BAS algorithm is used to automatically search for the objective function and obtain the effective signal components. The results show that the proposed method can effectively avoid the mode mixing issues and overcome the need to determine the k value beforehand when employing VMD for AE signal decomposition.
• In comparison to EMD and EVMD, FCI-VMD proves more suitable for AE signal denoising. By utilizing fuzzy entropy as the basis for selecting the noise components in the four stress stages, the method avoids the poor noise reduction effect caused by filtering out only the components with more noise content.
• The multifractal parameter, Δα, exhibits a sharp decrease as coal damage occurs, aligning with the evolution pattern of the AE main frequency. Consequently, Δα can be employed as a quantitative analysis parameter for AE spectra. Moreover, Δα can serve as an early warning indicator for coal damage and instability.