## 3 Methods
### Progressive Gaussian Decomposition  
The proposed PGD adopts an iterative approach to estimate the most appropriate number of Gaussian models for decomposing the individual waveform. It initiates an iteration according to the number of detectable peaks within the selected signal range, and the number of Gaussian models is gradually increased. The original peaks (OPs) were detected using findpeaks, a MATLAB peak detection function [41]. The amplitude, center, and width (\(A_{OPi}\), \(H_{OPi}\) and \(\sigma_{OPi}\), respectively) of each peak are used as the initial values for the first iteration. In this study, multiple Gaussian curve fitting was implemented using fit, a MATLAB curve fitting toolbox function [42]. Consequently, the number of initially approximated Gaussian components is equal to the number of OPs, and the estimated peaks (EPs) of each component are derived. The fitness of an approximated Gaussian mixture model is typically evaluated using residual-based measures, such as the error ratio [31]. In this study, \(\mathrm{R}^{2}\) is used to measure the residual between the original waveform (\(y\)) and approximated Gaussian mixture model (\(\hat{y}\)):  
\[R^{2}=1-\frac{\sum_{i=1}^{w}(y_{i}-\hat{y}_{i})^{2}}{\sum_{i=1}^{w}\bigl{(}y_ {i}-\mu_{y}\bigr{)}^{2}}. \tag{6}\]  
However, the evaluation using residual measures may ignore and miss components corresponding to relatively small amplitude peaks (e.g., deep bottom return). In view of this, another criterion is suggested in addition to the residual measure. If the Gaussian curve fitting is successfully performed, the OPs and their corresponding EPs should be located in similar positions. Therefore, we devised a new criterion, i.e., the time difference (\(\Delta t_{OP}\)) between each OP and corresponding (closest) EP, as follows:  
\[\Delta t_{OPj}=\underset{1\leq k\leq n}{min}\Big{|}t_{OPj}-t_{EPk}\Big{|}, \tag{7}\]  
where \(n=m+r-1\) represents the number of estimated Gaussian components, \(m\) is the number of OPs, and \(r\) represents the number of iterations. If all \(\Delta t_{OP}\) values are less than the temporal threshold (\(\uptau\)) and \(\mathrm{R}^{2}>0.95\), the decomposition is deemed successful,  
Figure 5: Example of signal range selection.  
and the process is terminated. Otherwise, if any \(\Delta t_{OP}\) is greater than \(\tau\) or R\({}^{2}\)\(\leqq\) 0.95, the decomposition is considered incomplete, and further iteration is performed. Because Gaussian curve fitting is dependent on the initial value, deciding the initial value to apply in the next iteration in addition to OPs is critical. If the EPs located far from the OPs are estimated, it means that the more dominant Gaussian components than the components corresponding to the OPs have not been undetected in the peak detection. To refer the undetected component to the initial values, the most distant EP from the OPs, i.e., the EP least related to the OPs, is selected as the PP and added as an initial value to the next iteration. In the rth iteration, r EPs are selected as PPs according to the order of those temporally farthest from the OPs. The time difference (\(\Delta t_{GP}\)) between an EP and the closest OP serves as the criterion for selecting the PP and is as follows:  
\[\Delta t_{EP_{k}}=\underset{1\leq j\leq m}{\min}\Big{|}t_{EP_{k}}-t_{OP_{j}} \Big{|} \tag{8}\]  
The PGD process can be expressed as the following steps:  
Step i. \(r=1\); \(S_{initial}=S_{OP}=\{OP_{1},\ OP_{2},\ \cdots,\ OP_{m}