https://orcid.org/0009-0003-4877-139X
Abstract
Seismic data is often contaminated with various low-frequency noises such as swell noise, instrumental noise, and tow noise. Typically, in the processing of seismic data, a high-pass filter is used to eliminate these unwanted low-frequency components. However, this filtering not only removes noise but also eliminates meaningful low-frequency signals essential for broadband seismic imaging. Numerous studies have been conducted to reconstruct low-frequency components by using machine learning techniques. However, most studies use convolutional neural networks and transform the data into image form rather than raw time-series data directly. In this study, we propose a low-frequency reconstruction method based on a recurrent neural network tailored to the sequence nature of seismic data. To improve the accuracy of the low-frequency component reconstruction network, we augment the feature information by incorporating envelope data, which provides additional aspects to seismic traces. Synthetic data tests confirm that our proposed method can reconstruct low-frequency components accurately. Furthermore, the validated network also successfully reconstructs low-frequency components in real dataset from the Arctic Sea.
1. Introduction
The low-frequency component of seismic data is an important factor to obtain high-fidelity seismic images. Furthermore, without low-frequency components of seismic data, a velocity model building using full waveform inversion (FWI) (Tarantola 1984) becomes challenging (Ovcharenko et al. 2022; Sun & Demanet 2020b; Zhang et al. 2017). Unfortunately, various noise are recorded during seismic exploration (Shin and Cha, 2009). Some are external noise sources caused by artificial factors such as ships and exploration facilities, and some are environmental noise caused by wind, waves, ocean currents, and wildlife. The wave-related noise that is the most problematic in marine seismic exploration has a frequency band of up to 10 Hz. Furthermore, various random noises are added to these low-frequency components (Cheng et al. 2015; Wang and Herrmann, 2016). Technical limitations of available sensor technology made research into the lowest frequency components of the seismic background field challenging for years, especially in the marine environment (Olofsson 2010). Afterwards, with the advent of high-sensitivity, broadband seismometers, research began to address the problem mechanically. However, challenges with resolution at low frequencies can arise due to the sensitivity and intricate design of the sensor, as well as the frequency attenuation characteristics of the hydrophone (Harris 1994). Therefore, in the processing of seismic data, unreliable low-frequency components are removed using a high-pass filter. However, this filtering also removes meaningful low-frequency components essential for the imaging process. To address the challenge posed by the absence of low frequencies in seismic imaging, Oldenburg et al. (1983) proposed extrapolating the low-frequency components of a single channel using a prediction filter, and Claerbout & Fomel (2008) reconstructed the low-frequency components through amplitude spectrum extrapolation using a multidimensional regression filter. Additionally, Li & Demanet (2015) attempted to estimate low frequencies based on the phase-tracking method. Chiu (2020) introduced a novel approach wherein a multi-channel prediction filter is derived from known higher frequencies, which is then utilized to extrapolate the low frequencies in both Radon domain and XT domain.
In studies focusing on reconstructing or utilizing existing low-frequency components, the calculations involved are often complicated, and there are limitations depending on the data acquisition environment (Zhan et al. 2015). Recently, various efforts have been made to address these challenges by integrating seismic data processing with machine learning techniques. Ovcharenko et al. (2019) proposed a deep CNN-based low-frequency reconstruction algorithm using the frequency map of the source-receiver as input. Sun & Demanet (2020a) proposed an algorithm for trace-by-trace low-frequency component prediction by training a CNN using time-domain data from numerical modeling. Fang et al. (2020) used a CNN architecture based on convolutional autoencoders to train the relationship between low-frequency data and high-frequency data to predict the low-frequency component. Kuijpers et al. (2021) designed a CNN-recurrent neural network (RNN) network to predict low-frequency components in shotgather. Recently, studies have been conducted to improve the applicability of field data of existing studies that relied on synthetic data. Cheng et al. (2024) proposed a low-frequency component prediction algorithm that applied the characteristics of field data in training based on a self-supervised framework. In addition, various methods of machine learning-based studies using seismic images are being conducted (Aharchaou and Baumstein, 2020; Hu et al. 2021; Jin et al. 2021). Like this, research is being conducted to successfully reconstruct low-frequency components by learning patterns in data using CNNs. However, CNNs face challenges in learning the sequential characteristics of the data, and as the size of the input data increases, the network becomes heavier, requiring more powerful hardware resources (Alzubaidi et al. 2021; Lee et al. 2024). Therefore, this study aims to utilize RNNs to extract sequential features from individual traces of the observed data to perform low-frequency reconstruction. By adopting this approach, the computational burden during the network training process was reduced. This choice enables us to capitalize on the inherent sequential characteristics of seismic data. Additionally, to enrich the learning process with diverse features, we incorporate transformed seismic data that unveil other pertinent characteristics of seismic data. Specifically, we introduce the envelope of seismic data as supplementary input data. The envelope, distinguished by its simplicity and prominent long-wavelength characteristics (Liu et al. 2020; Yang, 2017), serves to augment the training process effectively.
In this paper, we present an algorithm designed to reconstruct the missing low-frequency components of seismic data using a stacked-bidirectional long short-term memory (BiLSTM)-based trained network. We validate the accuracy of our proposed method through synthetic data test. Furthermore, we apply the trained network to field data to verify the efficacy of our low-frequency component reconstruction approach in practical application.
2. Methodology
2.1. Network architecture
This study uses BiLSTM to reconstruct low-frequency components. BiLSTM is a bidirectional learning module that extracts features from sequence data (Siami-Namini et al. 2019). It is suitable for handling long sequence data among RNNs because it simultaneously learns signals in both the forward and backward directions. It is also not constrained by the size of the input data, allowing it to be applied to dataset having different sizes.
The network architecture, illustrated in Fig. 1, is composed of an Input layer, BiLSTM layer, dropout layer, fully connected layer, and regression layer. The dropout layer plays a crucial role in preventing the training data from being overfitted.
The proposed stacked-BiLSTM architecture for estimating the low-frequency information. The number of hidden units is shown under each layer, while the numbers under the dropout layers is the dropout rate.
2.2. Generate the training data
This study aims to reconstruct the low-frequency components of seismic data acquired from diverse environments. To achieve this goal, it is essential to generate a substantial amount of varied data, encompassing different exploration environments, source characteristics, and amplitude features. Consequently, we developed and employed a convolution model to efficiently generate random training data. Figure 2 illustrates the process of generating the convolutional model for training purposes.
To create a convolutional model, a reflectivity model is required. Consequently, a random reflectivity model was generated with a maximum amplitude value ranging from −1 to 1.
Source wavelet for the convolution model, a Ricker wavelet with a dominant frequency ranging from 10 to 60 Hz was employed.
Convolve the reflectivity from (i) with the source wavelet from (ii), then normalize so that the maximum amplitude value lies between −1 and 1.
A random number within the range of 50% to 150% of the maximum value from the field data was selected, and this was applied as a weighting factor to the generated convolution model to incorporate the amplitude characteristics of the field data.
Proposed workflow for the recovery of low-frequency information with stacked-BiLSTM.
The amplitude-adjusted convolutional model is used as the label data, while the input data for training consists of low-frequency component-filtered data obtained through a high-pass filter. The high-pass filter used to generate the input data was designed as a ninth-order Butterworth high-pass filter with a cutoff frequency of 8 Hz. This configuration effectively removes low-frequency noise from the field data. We implement the filter using the “butter” command in MATLAB. Figure 3 displays the frequency spectrum of the Butterworth-designed high-pass filter whose dominant cutoff frequency is set to 8 Hz.
Frequency response of the designed Butterworth high-pass filter.
The input data for the proposed network comprises two time-series data: the high-pass filtered data and the envelope of high-pass filtered data. Extracting features for RNN training from seismic time-series data can be challenging due to the complexity of the features involved. Therefore, this study incorporates the envelope of the data, with its low-frequency components removed, as additional features. Huang et al. (2015) used envelope inversion to construct long-wavelength velocity models from data without low-frequency components. In this way, the envelope is frequently employed in geophysics to capture the low-frequency characteristics of time-series data (Wu et al. 2014). In this study, the envelope, known for its relative insensitivity to frequency variations, is utilized as an additional input. This approach aims to provide additional feature information to the network, mitigating the difficulty in learning complex characteristics from single time-series data, particularly high-pass filtered data. In Fig. 4, although the waveforms of the two Ricker signals are different, the envelopes of the two signals are similar. Therefore, by providing an envelope with characteristics similar to the long-wavelength features of the time-series data before high-pass filtering as additional input data, it is expected that the learning challenge will be alleviated.
Comparison the time-series data with the envelope data generated with a 30-Hz Ricker wavelet. (a) Ricker wavelet. (b) High-pass filtered Ricker wavelet. (c) Comparison of the envelope data.
2.3. Train strategy
The input dataset consists of high-pass filtered data and envelopes generated using a convolution model. The network is trained using 17 000 datasets for training and 2000 data sets for validation. In addition, the performance of the trained network is tested using 1000 test data samples generated through the same process as the training data. The hyperparameters and training environment are configured as specified in Tables 1 and 2. The convergence of the training and validation root mean-squared error (RMSE) curve depicted in Fig. 5 indicates that the proposed network is trained stably without overfitting. Figure 6 shows test results for data using sources with dominant frequencies of 20, 30, and 40 Hz. When comparing the original data, before the removal of low-frequency components, with the data where low-frequency components were reconstructed using the proposed network, some amplitude differences are observed. However, the phase information is successfully predicted. These amplitude discrepancies decreased as the frequency component of the data increased. Additionally, the frequency spectrum results confirm the reproduction of low-frequency components. The prediction performance is better for high-frequency data. This trend is also reflected in the analysis of prediction accuracy by dominant frequency of the sources used, as shown in Fig. 7. Lower dominant frequencies resulted in higher RMSE values, while an increase in dominant frequency corresponded to a convergence toward lower RMSE values.
RMSE curves along iterations for stacked-BiLSTM.
Comparison of the (a, d, g) time series, (b, e, h) frequency spectrum of trace, and (b, e, h) zoom in on the frequency spectrum of the trace 0–20 Hz at a (a)–(c) 20 Hz, (d)–(f) 30 Hz, and (g)–(i) 40 Hz.
RMSE for the test data according to the dominant frequency of the source wavelet.
Hyperparameter and training settings for the proposed stacked-BiLSTM neural network model.
| Hyperparameter | Value |
|---|---|
| Signal length of the input | 1s |
| Sampling interval | 0.1 × 10−2s |
| Number of epochs | 10 |
| Data shuffle | Per every epoch |
| Mini batch size | 128 |
| Optimizer | Adam |
| Learning rate | 0.1 × 10−3 |
| L2 regularization | 0.001 |
| Hyperparameter | Value |
|---|---|
| Signal length of the input | 1s |
| Sampling interval | 0.1 × 10−2s |
| Number of epochs | 10 |
| Data shuffle | Per every epoch |
| Mini batch size | 128 |
| Optimizer | Adam |
| Learning rate | 0.1 × 10−3 |
| L2 regularization | 0.001 |
Hyperparameter and training settings for the proposed stacked-BiLSTM neural network model.
| Hyperparameter | Value |
|---|---|
| Signal length of the input | 1s |
| Sampling interval | 0.1 × 10−2s |
| Number of epochs | 10 |
| Data shuffle | Per every epoch |
| Mini batch size | 128 |
| Optimizer | Adam |
| Learning rate | 0.1 × 10−3 |
| L2 regularization | 0.001 |
| Hyperparameter | Value |
|---|---|
| Signal length of the input | 1s |
| Sampling interval | 0.1 × 10−2s |
| Number of epochs | 10 |
| Data shuffle | Per every epoch |
| Mini batch size | 128 |
| Optimizer | Adam |
| Learning rate | 0.1 × 10−3 |
| L2 regularization | 0.001 |
The hardware specification used for testing.
| Experiment system parameter | |
|---|---|
| CPU | Intel(R) Core(TM) i9-10980XE CPU @ 3.00GHz 3.00 GHz |
| GPU | Dual NVIDIA GeForce RTX 3090 24GB |
| Operating system | Windows 11 64bit |
| MATLAB version | MATLAB R2022a |
| Lab environment | Deep Learning Toolbox |
| Experiment system parameter | |
|---|---|
| CPU | Intel(R) Core(TM) i9-10980XE CPU @ 3.00GHz 3.00 GHz |
| GPU | Dual NVIDIA GeForce RTX 3090 24GB |
| Operating system | Windows 11 64bit |
| MATLAB version | MATLAB R2022a |
| Lab environment | Deep Learning Toolbox |
The hardware specification used for testing.
| Experiment system parameter | |
|---|---|
| CPU | Intel(R) Core(TM) i9-10980XE CPU @ 3.00GHz 3.00 GHz |
| GPU | Dual NVIDIA GeForce RTX 3090 24GB |
| Operating system | Windows 11 64bit |
| MATLAB version | MATLAB R2022a |
| Lab environment | Deep Learning Toolbox |
| Experiment system parameter | |
|---|---|
| CPU | Intel(R) Core(TM) i9-10980XE CPU @ 3.00GHz 3.00 GHz |
| GPU | Dual NVIDIA GeForce RTX 3090 24GB |
| Operating system | Windows 11 64bit |
| MATLAB version | MATLAB R2022a |
| Lab environment | Deep Learning Toolbox |
3. Numerical example
To validate the accuracy of the trained network, two numerical experiments were conducted. The proposed network was applied to synthetic data simulating the characteristics of the Arctic data and to field data obtained from the Chukchi Sea.
3.1. Synthetic data test
To simulate the Arctic Sea environment with the presence of permafrost layers, a salt dome model with high-velocity anomalies is utilized (Aminzadeh et al. 1996). Permafrost layers are typically located beneath the seabed at depths of 50 to 600 m and exhibit P-wave velocities ranging from 3.4 s to 4.35 km/s depending on the proportion of ice (Taylor et al. 2013; Kang et al. 2021). To align with these characteristics, a modified salt dome model is used to construct a synthetic environment that mimics the acquisition conditions of the target field data (Fig. 8). In this environment, the wave propagation from the source is solved using the 2D acoustic wave equation in the finite-difference time domain. Additionally, perfectly matched layer (PML) boundary conditions (Hastings et al. 1996) are applied to eliminate unwanted artificial reflections at the modeling boundaries. The source used in the modeling considered the frequency range of training data, transmitting a Ricker wavelet with a center frequency of 10–60 Hz from the upper part of the salt dome model. The velocity model used for numerical examples has the following specifications: a resampled 4800 × 1200 size velocity model of the salt dome with a grid size of 0.5 m. The source is positioned at a depth of 6 m, and receivers, organized in a streamer with 120 channels spaced 12.5 m apart horizontally, are located at a depth of 6 m. The recording time is 1.5 seconds, and the sampling rate is 0.5×10−4s. The source applied Ricker wavelets with dominant frequencies ranging from 10 to 60 Hz with an interval of 1 Hz, obtaining a total of 51 synthetic datasets with different frequency ranges.
Comparison of the 2D salt dome synthetic model and the source and receiver position in this synthetic data set.
After removing the low-frequency components from the acquired synthetic data, the low-frequency components are restored by applying the proposed method, and the RMSEs for frequencies below 17 Hz of the restored low-frequency components and label data are compared. The analysis results reveal that, when applying the proposed network, the components lost due to the filter are reconstructed, resulting in lower RMSE values compared to the high-pass filtered data. Additionally, we see that as the dominant frequency of the source increased, the energy ratio of frequency components below 17 Hz decreased, leading to a relatively lower RMSE (Fig. 9). For a more detailed comparison, an analysis is performed on 20, 40, and 60 Hz data, encompassing the benchmarked data's frequency range, from the set of 60 generated examples.
RMSE trend of prediction accuracy by the dominant frequency of the source wavelet.
Figure 10 shows the comparison of normal low-frequency components in label, input with low-frequency components lost due to the high-pass filter, output reconstructed by low-frequency components through the proposed method, and the difference between label and reconstructed low-frequency components for frequencies of 20, 40, and 60 Hz. In Fig. 11, we zoom in on the initial seconds (0–0.5s) to examine the early signal characteristics. In the modeling results, the designed high-pass filter is applied to suppress low-frequency noise in Arctic data, and then the proposed method is used to reconstruct the low-frequency components. Figure 11 parts b, f, and j display a ringing effect around the early signal and major reflection signals due to the lack of frequency components below 8 Hz. By contrast, the label data (Fig. 11c, g, and k) shows no ringing effect and shows a clear signal as it contains frequency components below 8 Hz. This effect is more evident in the frequency domain (Figs. 12 and 13). The data with reconstructed low-frequency components using the proposed method (Fig. 12c, g, and k) show relatively less signal vibration and appear sharper since the frequency components below 8 Hz, present in Fig. 12b, f, and j, have been reconstructed. To assess the accuracy of the predicted low-frequency components, a comparison analysis with label containing the correct frequency components (Fig. 12a) confirms that the predicted low-frequency components are very similar to the low-frequency components before the high-pass filter. In particular, frequency components below 5 Hz are crucial in inversion process where seismic velocity analysis is used. Figure 13 shows the performance of the low-frequency reconstruction algorithm in reconstructing ultra-low-frequency components below 5 Hz. While the input data processed through the high-pass filter (Fig. 13b, f, and j) lacks low-frequency components below ∼8 Hz, the data reconstructed through the proposed network (Fig. 13c, g, and k) demonstrates the recovery of low-frequency components. The reconstructed low-frequency component is very similar to the correct answer (label data), as there is little difference between the label data and the predicted data.
Comparison of common-shot gathers: (a–d) results from data using a 20 Hz source wavelet, (e–h) results from data using a 40 Hz source wavelet, and (i–l) results from data using a 60 Hz source wavelet. (a, e, i) Raw data, (b, f, j) high-pass filtered data, (c, g, k) predicted data, and (d, h, l) the difference between raw and predicted data.
Zoomed-in comparison of common-shot gathers (0–0.5 s): (a–d) results from data using a 20 Hz source wavelet, (e–h) results from data using a 40 Hz source wavelet, and (i–l) results from data using a 60 Hz source wavelet. (a, e, i) Raw data, (b, f, j) high-pass filtered data, (c, g, k) predicted data, and (d, h, l) differences between raw and predicted data.
Comparison of 2D frequency analysis result: (a–d) results from data using a 20 Hz source wavelet, (e–h) results from data using a 40 Hz source wavelet, and (i–l) results from data using a 60 Hz source wavelet. (a, e, i) Raw data, (b, f, j) high-pass filtered data, (c, g, k) predicted data, and (d, h, l) the difference between raw and predicted data.
Zoomed-in comparison of 2D frequency analysis result (0–16 Hz): (a–d) results from data using a 20 Hz source wavelet, (e–h) results from data using a 40 Hz source wavelet, (i–l) results from data using a 60 Hz source wavelet. (a, e, i) Raw data, (b, f, j) high-pass filtered data, (c, g, k) predicted data, and (d, h, l) the difference between raw and predicted data.
For a more detailed analysis, we extract traces acquired at the 30th channel (offset = 375m) for the 20, 40, and 60 Hz data. We compare the data with low-frequency components, data with low-frequency components removed by high-pass filter, and the reconstructed data after applying the proposed algorithm. In Fig. 14, the waveform is distorted due to the absence of low-frequency components and the ringing effect caused by the low-cut filter. By applying the proposed algorithm, the previously removed low-frequency components were reproduced, and the distorted waveform was reconstructed to be similar to the raw data (Fig. 14a, c, and e). The reconstruction of low-frequency components is evident in the fast Fourier transform results shown in Figs. 14b, d, and f, where low-frequency components below 16 Hz are removed and then successfully reconstructed, resembling the original data. We further quantified the similarity using RMSE, comparing the data with removed and reconstructed low-frequency components to the raw data. The RMSE values between the trace with removed low-frequency components obtained through the high-pass filter and the raw data are 9.405. By contrast, the RMSE between the data with reconstructed low-frequency components through the proposed method and the raw data was 0.404. This indicates a higher similarity between the data with reconstructed low-frequency components and raw data, confirming the effectiveness of the proposed method in improving seismic data by reconstructing low-frequency components. The 1D frequency analysis results also showed consistent recovery of low-frequency components similar to the frequency components before removal. The consistent low-frequency reconstruction effect was observed for both 40 and 60 Hz, with relatively more accurate reconstruction results for the higher-frequency data compared to 20 Hz.
Comparison of trace and frequency spectrograms: (a) and (b) is result from 20 Hz data, (c) and (d) is result from 40 Hz data, (e) and (f) is result from 60 Hz data. (a, c, e) Time-domain data, and (b, d, f) Frequency-domain data.
3.2. Field data test
To evaluate the applicability of the proposed algorithm to field data, we first remove the contaminated low-frequency components from benchmarked data. Then, we apply an algorithm to reconstruct the low-frequency components. Figure 15 illustrates the procedure of the processing. To examine the low-frequency reconstruction attributes of the proposed algorithm, we extract common-shot gather data from a survey line. The raw field data displays minimal visibility of primary signals due to severe low-frequency component noises (Fig. 15a). Subsequently, we use high-pass filter, as employed during training, to eliminate the low-frequency components (Fig. 15b). While this filtering reduces overall noise, it results in abrupt oscillations in the initial signals causing distortion. To address this issue, we apply muting on the filtered data (Fig. 15c) to eliminate noisy signals before the direct arrival. On reconstructing the low-frequency components, we can find that both the phase and amplitude from the reconstructed data (shown in Fig. 15d) closely resemble those from the high-pass filtered data (shown in Fig. 15b). To examine the frequency characteristics of each dataset, we conduct an analysis in the F-K domain. As depicted in Fig. 16, noticeable differences are observed in the target frequency range of below 17 Hz, while showing minimal variations in upper frequency ranges. We further zoom in on the 0–17 Hz range, affected by the high-pass filter applied to the data (Fig. 17). In Fig. 17a, severe low-frequency contamination in the raw data hinders discernment of signals from the artificial source. Therefore, we remove such noise components (Fig. 17b and c). Throughout this process, to reconstruct meaningful low-frequency information that may have been lost, we employ the proposed network, resulting in the generation of coherent low-frequency components (Fig. 17d).
A procedure of field data processing: (a) a sample of raw data, (b) high-pass filtered data from (a), (c) data after muting noisy signals before the direct waves, and (d) data reconstructed low-frequency component through our proposed algorithm.
Comparison of 2D frequency analysis result: (a) raw data, (b) high-pass filtered data, (c) muted data from (b), and (d) low-frequency component reconstructed data.
Zoomed-in version of the results in Fig. 16 (0–16 Hz).
To conduct a more detailed analysis of the low-frequency reconstruction results obtained from the proposed method, we extract and compare the 29th trace (offset = 362.5 m) in Fig. 18. Figure 18 depicts the raw data contaminated with low-frequency noise, the high-pass filtered data after low-frequency noise removal, and the predicted data where the removed low-frequency component is reconstructed from the high-pass filtered data as input. In Fig. 18a, we see that the raw data is distorted by low-frequency noise, affecting the waveform. Although the distortion can be suppressed by applying a high-pass filter, a distortion such as a phase change can persists even after applying a high-pass filter due to the absence of low-frequency components (Fig. 18b). To address this problem, the proposed network compensates for the lost information in the removed low-frequency component, resulting in a reduction of waveform distortion. This effect is also discernible through the comparison of the frequency spectrograms (Fig. 18c and d). The high-pass filter successfully removes severe low-frequency contamination from the raw data, while the proposed network effectively reconstructs the attenuated low-frequency component.
Comparison of traces and frequency spectra: (a) and (b) results of time-domain data, and (c) and (d) results of frequency-domain data. (b) and (d) are enlarged versions from 0 to 20 Hz.
4. Conclusions
We introduce a method aimed at reconstructing low-frequency components from seismic data using RNNs. The filter designed to suppress low-frequency contamination can also eliminate meaningful low-frequency components. To address this, we propose a stacked-BiLSTM network as a reconstruction network. This network is trained using arbitrary signals with specific ranges of frequencies (10–60 Hz), enabling it to learn how to reconstruct the low-frequency components removed by the designed filter, as validated by numerical experiments. However, because this network is trained trace by trace, it is confirmed that some discontinuities between traces may be appeared when applied to shot gather. Nevertheless, despite this limitation, the network effectively reconstructs low-frequency components regardless of the source signature of the input data and can be applied without size constraints on the input data. The reconstructed seismic data, enriched with low-frequency components, can be utilized in seismic model building technique such as FWI, facilitating the construction of stable long-wavelength velocity models.
Acknowledgments
This research was supported by the Basic Research Project (grant no. 24–3313) of the Korea Institute of Geoscience and Mineral Resources (KIGAM) funded by the Ministry of Science, ICT and Future Planning of Korea and also was supported by Korea Institute of Marine Science & Technology Promotion (KIMST) funded by the Ministry of Oceans and Fisheries, Korea (RS-2023–00259633).
Conflict of interest statement
None declared.
Data availability
The field data used in this study cannot be shared. However, other data supporting the findings of this study may be available upon reasonable request to the corresponding author.
