Bidirectional gating recurrent neural network logging curve reconstruction method based on multistage wavelet decomposition

Abstract

The invention belongs to the technical field of oil and gas geophysical exploration and logging data processing, and particularly discloses a bidirectional gating recurrent neural network logging curve reconstruction method based on multistage wavelet decomposition. The method is characterized in that the reconstruction precision and the reservoir identification reliability are improved by mining logging data sequence dependence, multi-scale geologic features are captured in a layered mode, the interpretability of a model is enhanced by means of wavelet decomposition and natural relevance of stratum signals, and technical support is provided for oil-gas exploration and development. The method comprises the following steps: 1) preprocessing original logging data and constructing a data set; (2) an mWDN-BiGRU model framework is constructed; 3) model training and evaluation; the acoustic logging curve reconstruction method has the beneficial technical effects that in order to reconstruct the logging curve, the missing data of the acoustic logging curve is repaired by adopting a bidirectional gating recurrent neural network algorithm of multistage wavelet decomposition, so that the precision of the reconstructed data is improved, and the reservoir parameter interpretation precision and the sedimentary facies identification reliability are improved.

Description

Technical Field

[0001] This invention belongs to the field of oil and gas geophysical exploration and well logging data processing technology, specifically involving a well logging curve reconstruction method based on a two-way gated recurrent neural network with multi-level wavelet decomposition. Background Technology

[0002] Sonic logging, resistivity logging, and natural gamma logging are core data for identifying subsurface lithology, calculating reservoir properties, and performing detailed descriptions of oil and gas reservoirs. However, in actual operations, factors such as instrument malfunction, wellbore collapse, harsh drilling environments, or data loss often lead to distorted or severely missing logging curves in key well sections. Re-logging is costly and often not feasible for completed well sections. Therefore, logging curve reconstruction technology has become a key means to compensate for missing data and improve data quality.

[0003] Traditional curve reconstruction methods mainly include interpolation and empirical formula methods. Although these methods can fill in missing segments to some extent, they are mostly based on the assumption of linear or stationary signals, making it difficult to effectively capture the complex nonlinear relationships and local abrupt changes among underground strata parameters. In order to better explore the time-series dependence characteristics of well logging curves on depth sequences, recurrent neural networks and their variants, such as long short-term memory networks and gated recurrent units, have been gradually applied. However, traditional unidirectional networks cannot utilize "future" context information, resulting in insufficient perception of changes in strata interfaces and sedimentary cycles. Bidirectional recurrent neural networks solve this problem to some extent by combining historical and future information.

[0004] Nevertheless, existing models still have limitations when facing the core challenge of well logging curve reconstruction—decoupling of multi-scale features and noise suppression. Well logging curve data is a multi-scale signal, containing high-frequency noise such as local disturbances caused by wellbore collapse and low-frequency trends such as the overall changes in formation compaction. Conventional time-series modeling methods may be limited by high-frequency noise interference, low-frequency trend confusion, and difficulties in cross-scale modeling. Pure data-driven models are mainly based on time-domain features and are not sensitive enough to frequency-domain features, making them susceptible to high-frequency noise interference and thus difficult to accurately capture low-frequency geological trends. Summary of the Invention

[0005] The purpose of this invention is to provide a well logging curve reconstruction method based on a bidirectional gated recurrent unit network (mWDN-BiGRU) with multilevel wavelet decomposition, which solves the problems of decoupling multi-scale features and noise suppression in well logging curves. It incorporates frequency domain methods into the modeling task of well logging curve reconstruction. The multilevel wavelet decomposition module is used to separate high-frequency noise and low-frequency signals in the well logging curve, while the bidirectional gated recurrent network module is used for in-depth modeling of low-frequency principal components and optimized selection of high-frequency secondary components. The joint optimization of "frequency domain denoising + time-series modeling" achieves enhanced preservation of geological trends and extraction of favorable frequency band signals. The process of this invention is as follows: Figure 1 As shown.

[0006] The technical solution adopted in this invention is as follows: 1. A method for reconstructing well logging curves based on multi-level wavelet decomposition and a bidirectional gated recurrent neural network, the method comprising: Step 1: Raw logging data preprocessing and dataset construction: The raw logging data is subjected to quality checks and outlier cleaning. The raw logging data of different wells are divided into training set, validation set and test set according to the well number. The dataset is constructed by data transformation and normalization, dataset encapsulation, data iterator construction and data hyperparameter configuration. Step 2: Construct the mWDN-BiGRU model framework: First, extract the frequency information of the time series data by using the multi-level wavelet decomposition mWDN module. Connect the mWDN module with multiple GRU modules to form a parallel data flow structure, and set the basic network parameters and data flow dimensions, and configure the model parameters. Step 3, Model Training and Evaluation: Use the training set built in Step 1 to train the model built in Step 2, and build a training framework for deep temporal networks from three dimensions: loss function, optimizer, and evaluation metric. Step 4: Reconstructing logging curves for missing well sections: Verify the model's generalization performance using the test set, save the model parameter file, and use the trained model to reconstruct and predict logging curves for missing well sections.

[0007] Furthermore, the specific steps of step 1 are as follows: Step 1.1: Raw logging data inspection and cleaning. In order to transform the raw data of varying quality into a continuous, consistent, and complete standardized dataset, it is necessary to remove outliers caused by downhole environmental interference and instrument malfunctions in the raw logging data, ensure that the value range and distribution of the effective data of each logging curve conform to the actual geological and physical meaning, and correct systematic errors caused by instrument differences and changes in the measurement environment. Step 1.2: Dataset construction. Based on the number of effective wells in each region, the original data from the effective well sections is divided into training set, validation set and test set according to the well number. The purpose is to ensure that the data of the same well appears in only one dataset, so that the data between wells are independent of each other. The dataset is saved and read in three file formats: train.csv, valid.csv and test.csv. Step 1.3, Data Transformation and Normalization: To alleviate the problem of distribution skew, the resistivity RT curve is transformed into an RT_log curve using logarithmic transformation. The formula for logarithmic transformation is as follows: , In the formula, RT lg The resistivity curve after logarithmic transformation; RT The resistivity curve; This is the offset; To uniformly map logging curves of different dimensions to the [0,1] interval and eliminate the influence of feature scale differences on neural network training, the Min-Max normalization method was used to extract the transformed training set data, which was then stored as a JSON file in dictionary form. To prevent outlier maxima and minima, V max , V min Using the 1st and 99th quantiles, the formula for linear normalization is as follows: , In the formula, V norm The value is the normalized value; V log These are the original curve values; V max The maximum value of the curve; V min This is the minimum value of the curve; Step 1.4: Sequence data segmentation; the original dataset is represented as follows. D ∈R n×m Where n is the number of samples in a certain well dataset, and m is the number of features in the dataset, the dataset is divided based on the total amount of data and the target sequence length to obtain datasets corresponding to the sequence length. D i Where i = 1, 2...l, and l is the length of the target segmented sequence, then the array stacking function from the NumPy library in Python is used. np.stack Implement l datasets D 1, D 2, ..., D lStacking adds a sequence dimension, resulting in a three-dimensional tensor dataset with sequence dimensions. The dataset is expressed as follows: , In the formula, D new For a new three-dimensional tensor dataset; D 1, D 2, ..., D l The resulting sub-datasets; axis The axial direction parameter for the stacking operation; Step 1.5: Dataset Encapsulation and Data Iterator Construction. Based on the PyTorch framework, the Well_Dataset class for logging sequence datasets is defined by inheriting the abstract class torch.utils.data.Dataset to achieve dataset encapsulation. The dataset is initialized, receiving the segmented sequence sample data, and then returning the total number of samples in the dataset. The range of index values ​​is determined based on the returned total number of samples, and then a single training sample is returned based on the index and converted into a tensor. Again, based on the PyTorch framework, the torch.utils.data.DataLoader class is called, and Well_Dataset is passed in to initialize the data iterator. Step 1.6: Configure the data hyperparameters, including: standardization method, input sequence length, and sliding step size.

[0008] Furthermore, the specific steps of step 2 are as follows: Step 2.1: Construct the network framework. The network structure consists of five parts: an input layer, a multi-level wavelet decomposition network (mWDN) layer, a bidirectional recurrent gated neural network (BiGRU) layer, a fully connected neural network (FCN) layer, and an output layer. Layers 1 and 5 are set as the input and output layers, respectively. Layer 2 uses a multi-level wavelet decomposition network as the initial module to extract high-frequency information from the logging data and separate low-frequency trend information. Layer 3 uses a bidirectional gated recurrent neural network, which can extract the time-domain features of the logging curve along the depth sequence from two directions, allowing the model to consider the influence of adjacent values ​​during prediction. When used for high-frequency sequences, it processes the extracted frequency-domain information, thus better learning the time-domain and frequency-domain information in the logging curve. Layer 4 is a fully connected neural network composed of two linear layers. The mWDN-BiGRU model structure is as follows: Figure 2 As shown; Step 2.2: Implement the multi-level wavelet decomposition network model. First, define the weight matrix. W l ( i )and W h ( iThe formula is in the form of ) and is as follows: , In the formula, W l ( i )and W h ( i Let be the low-pass weight matrix and high-pass weight matrix of the i-th level wavelet decomposition; k The length of the filter coefficients; l 1. l 2、...、 l k These are the coefficients of the low-pass filter; h 1. h 2、...、 h k These are the coefficients of the high-pass filter; These are the filling elements in the matrix; To enable the model to perform multi-resolution analysis from the early stages of training, wavelet basis functions are used on the weight matrix. W l ( i )and W h ( i Initialization is performed by using the weight matrix as a filter to improve convergence efficiency and enhance interpretability. Then, multi-level decomposition is used to process the time series data, as shown in the following formula: , In the formula, a l ( i ), a h ( i ) represents the low-frequency and high-frequency intermediate variables of the i-th level wavelet decomposition; σ (·) is the activation function; W l (i) and W h (i) represents the low-pass weight matrix and high-pass weight matrix of the i-th level wavelet decomposition; x l (i-1) is the low-frequency approximation signal obtained from the (i-1)th level decomposition; b l (i) and b h (i) represents the trainable vector bias; The low-frequency approximate signal is obtained after the i-th level wavelet decomposition. x l ( i and high-frequency detail signals x h ( i ),right x l (i When performing the (i+1)th level decomposition, a sliding convolution is performed on it using a low-pass filter and a high-pass filter to obtain the intermediate variable sequence. The convolution formula is as follows: , In the formula, and Let n be the nth element of the low-frequency intermediate variable sequence and the nth element of the high-frequency intermediate variable sequence of the (i+1)th level wavelet decomposition; Let be the (n+k-1)th element of the i-th level low-frequency subsequence; K The length of the filter; k Internal index of the filter coefficients; l k and h k These are the low-pass filter coefficients and high-pass filter coefficients; To achieve effective compression of the sequence length, the intermediate variables are downsampled by half using average pooling. The average pooling formula is as follows: , In the formula, j This represents the length of the high- and low-frequency subsequences after the i-th level decomposition; and Let j be the j-th element of the low-frequency approximation signal and the j-th element of the high-frequency detail signal obtained from the i-th level wavelet decomposition. and Let be the 2j-th element and the (2j-1)-th element in the low-frequency intermediate variable sequence; and Let be the 2j-th element and the (2j-1)-th element in the high-frequency intermediate variable sequence; To prevent excessive deviation due to an increase in the number of training epochs, a regularization constraint needs to be added. The regularization constraint formula is as follows: , In the formula, W l ( i )and W h ( i Let be the low-pass weight matrix and high-pass weight matrix of the i-th level wavelet decomposition; , Weight matrix W l ( i ), W h ( i The initial state of ); L*、L ( θ ), Ll , L h The regularization loss is the final total loss function, the basic loss function, the low-frequency weight matrix, and the high-frequency weight matrix. α , β These are the regularization parameters for low and high frequencies; To iteratively update network parameters and optimize the objective function, the backpropagation (BP) algorithm is used. Gradient descent is employed to optimize the loss function, enabling the neural network to automatically learn the nonlinear mapping relationship from input features to the target curve from multi-source logging data. The update formula is as follows: , In the formula, W l ( i )and W h ( i Let be the low-pass weight matrix and high-pass weight matrix of the i-th level wavelet decomposition; , Weight matrix W l ( i ), W h ( i The initial state of ); α , β These are the regularization parameters for low and high frequencies; L The loss function; η The learning rate; During the forward propagation of the model, the final regularized loss value is compared with that of each wavelet decomposition. L l and L h The weighted calculations, stacking, and averaging are performed, and then the results are output together and used in the backpropagation optimization process of the loss function. Step 2.3: Configure the network skeleton, connect the mWDN module with multiple GRU modules to form a parallel data flow structure, and set the basic network parameters and data flow dimensions to reflect the dimensional changes of the input data during the forward propagation process of the network. Step 2.4: Configure the model hyperparameters, specifically including: wavelet decomposition level, wavelet basis function type, mWDN low-frequency matrix regularization parameter, mWDN high-frequency matrix regularization parameter, RNN hidden layer dimension, RNN hidden layer number, and FCN hidden layer dimension.

[0009] Furthermore, the specific steps of step 3 are as follows: Step 3.1: Select Mean Squared Error (MSE) as the loss function to constrain the global fitting accuracy. Its formula is as follows: , In the formula, y i The actual value; This is a predicted value; N The length of the logging sequence; Step 3.2: Adaptive Moment Estimation (Adam) is selected as the training optimizer. It combines momentum estimation with an adaptive learning rate strategy, dynamically adjusting the learning rate of each parameter by calculating the mean and variance of the gradient, resulting in more stable parameter updates and a more lenient initial learning rate setting. During training, a learning rate scheduling strategy combining warm-up and cosine annealing is employed. In the early stages of training, a linear increase in the learning rate is used to boost it, while in the later stages, a cosine annealing scheduling strategy is used to decrease the learning rate. The formulas for warm-up and cosine annealing are as follows: , In the formula, T warm To control the number of epochs for preheating and cosine annealing; T max This represents the maximum number of training rounds. t The number of iterations that have been completed for model training; η t Let be the learning rate at the t-th training round; η 0 represents the initial learning rate; η min Minimum learning rate; Using the validation set as the standard for early stopping, and by monitoring and evaluating metrics, the training is terminated early using the early stopping mechanism in order to maintain the optimal model state. Step 3.3: Select the Mean Absolute Error (MAE), Pearson correlation coefficient, and coefficient of determination R. 2 As evaluation indicators, to ensure a direct correlation between the assessment results and the logging curve values, all indicators are calculated under the original physical dimensions after inverse normalization. The formula for Mean Absolute Error (MAE) is as follows: , In the formula, y i The actual value; This is a predicted value; N The length of the logging sequence; The formula for the Pearson correlation coefficient is as follows: , In the formula, t This represents the total number of sampling points; m i and n i These are the sampled point values ​​of the true curve and the reconstructed curve; Coefficient of determination R 2 The value ranges from (-∞, 1]. The closer the value is to 1, the stronger the model's ability to explain data variations. A negative value indicates that the model's predictive performance is worse than directly using the mean as an estimate. The coefficient of determination R0 2 The formula is as follows: , In the formula, N The length of the logging sequence; The mean of the true curve; y i The actual value; This is a predicted value; Step 3.4: Configure the training hyperparameters, including: training batch size, initial learning rate, weight decay coefficient, gradient clipping ratio, and the number of epochs to wait for when the Early Stopping validation set loss does not improve.

[0010] Furthermore, the specific steps of step 4 are as follows: Step 4.1: Model Saving and Reloading. A file management system is established hierarchically according to the prediction target, training dataset, prediction model, and input sequence length. The training model is saved in a full-model saving mode. Data and training process parameters are saved using a state dictionary parameter, stored in the same directory as the model file, and saved as JSON files. Simultaneously, the command-line parameter management system in the argparse library is used to save all process parameters. During model reloading, the state dictionary containing data preprocessing parameters is loaded first, followed by the complete model file, reproducing the environment consistent with training. Step 4.2, Model Inference and Output: Input the preprocessed logging sequence data into the model to perform the forward propagation calculation process, directly output the predicted value of the target logging curve, evaluate the accuracy and analyze the effect of the prediction results on the verification well and test well, and construct the logging curve of the data missing well section.

[0011] In summary, due to the adoption of the above technical solution, the beneficial effects of the present invention are as follows: 1. In this invention, a bidirectional gated recurrent neural network algorithm based on artificial intelligence and deep learning is used to effectively mine the sequence dependency relationship of well logging data, accurately capture the curve change pattern, effectively repair missing data and correct distorted information, significantly improve the reliability and accuracy of reconstructed data, improve the accuracy of reservoir parameter interpretation and the reliability of sedimentary facies identification, and provide key technical support for reservoir evaluation and well location deployment in oil and gas exploration and development.

[0012] 2. In this invention, a multi-level wavelet decomposition mechanism is introduced and connected in parallel with multiple bidirectional gated recurrent neural networks to realize multi-scale signal decomposition and sequence representation of well logging data, forming a joint optimization model of "frequency domain noise reduction + time series modeling". This model can process geological features of different scales in layers, such as thick layer trends and thin layer abrupt changes. Dedicated network learning modules are used for different frequency sub-bands, which further improves the model's ability to capture multi-scale features.

[0013] 3. The physical meaning of wavelet decomposition is naturally related to the frequency domain characteristics of the formation signal. Therefore, the temporal importance of different input features can be analyzed, which enhances the interpretability of the model. Attached Figure Description

[0014] Figure 1 This is a flowchart illustrating the overall experimental process of the method of the present invention. Figure 2 A schematic diagram of a bidirectional recurrent neural network structure for three-level wavelet decomposition. Figure 3 A schematic diagram of the reconstruction results of the acoustic logging curve of well A. Figure 4 A schematic diagram of the reconstruction results of the acoustic logging curve of well B. Figure 5 Schematic diagram of the reconstructed C-wave logging curve of the test well. Detailed Implementation

[0015] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are only for explaining the invention and are not intended to limit the scope of protection of the invention.

[0016] A well logging curve reconstruction method based on multi-level wavelet decomposition and bidirectional gated recurrent neural network includes the following steps: Step 1: Raw logging data preprocessing and dataset construction: The raw logging data is subjected to quality checks and outlier cleaning. The raw logging data of different wells are divided into training set, validation set and test set according to the well number. The dataset is constructed through data transformation and normalization, dataset encapsulation, data iterator construction and data hyperparameter configuration.

[0017] Specifically, the details of step 1 are as follows: Step 1.1: Check and clean the raw logging data. Through logging software program identification and map observation, abnormal well sections in the raw logging data are removed to ensure that the effective data range and distribution of each logging curve conform to the actual geological and physical meaning. At the same time, the code of the detection range is compiled to correct the systematic errors caused by instrument differences and changes in the measurement environment. Step 1.2: Dataset construction. Based on the number of effective wells in each region, 100 wells are selected from the actual logging data. The dataset is divided into training set, validation set and test set according to the ratio of 80:10:10. The data of each well is ensured to be independent of each other. After the dataset of each region is divided, it is stored and read in three file formats: train.csv, valid.csv and test.csv. Step 1.3, Data Transformation and Normalization: The resistivity RT curve is transformed into an RT_log curve using logarithmic transformation to alleviate distribution skewness. The formula for logarithmic transformation is as follows: , In the formula, RT lg The resistivity curve after logarithmic transformation; RT The resistivity curve; This is the offset; The transformed training set data was extracted using the Min-Max normalization method and stored as a JSON file in dictionary form. This was done to prevent outliers from finding maxima and minima. V max , V min Using the 1st and 99th quantiles, the formula for linear normalization is as follows: , In the formula, V norm The value is the normalized value; V log These are the original curve values; V max The maximum value of the curve; V min This is the minimum value of the curve; Step 1.4: Sequence data segmentation. The original dataset can be represented as follows: D ∈R n×m Where n is the number of samples in a certain well dataset, and m is the number of features in the dataset, the dataset is divided based on the total amount of data and the target sequence length to obtain datasets corresponding to the sequence length. D i Where i = 1, 2...l, and l is the length of the target segmented sequence, the array stacking function in the NumPy library of Python is used. np.stack Implement l datasets D 1, D 2, ..., D l Stacking adds a sequence dimension, resulting in a three-dimensional tensor dataset with sequence dimensions. The dataset is expressed as follows: , In the formula, D new For a new three-dimensional tensor dataset; D 1, D 2, ..., D l The resulting sub-datasets; axis The axial direction parameter for the stacking operation; Step 1.5: Dataset encapsulation and data iterator construction. Based on the PyTorch framework, the well logging sequence dataset class Well_Dataset is defined by inheriting the abstract class torch.utils.data.Dataset to encapsulate the dataset, initialize the dataset, receive the segmented sequence sample data, and then return the total number of samples in the dataset. The range of index values ​​is determined based on the returned total number of samples, and then a single training sample is returned based on the index and converted into a tensor. Again, based on the PyTorch framework, the torch.utils.data.DataLoader class is called, and Well_Dataset is passed in to initialize the data iterator. Step 1.6: Configure data hyperparameters: Configure the standardization method scaler(str) to Min-Max normalization for data with well-defined boundaries; configure the input sequence length seq_len(int) to 32, which corresponds to a formation thickness of 4m when the sampling interval is 0.125m; configure the sliding step size stride(int) to 1 to maximize sample size or balance computational efficiency.

[0018] Step 2: Construct the mWDN-BiGRU model framework: First, extract the frequency information of the time series data by using the multi-level wavelet decomposition mWDN module. Connect the mWDN module with multiple GRU modules to form a parallel data flow structure, and set the basic network parameters and data flow dimensions, and configure the model parameters.

[0019] Specifically, the details of step 2 are as follows: Step 2.1: Construct the network framework. The network structure consists of five parts: an input layer, a multi-level wavelet decomposition network (mWDN) layer, a bidirectional recurrent gated neural network (BiGRU) layer, a fully connected neural network (FCN) layer, and an output layer. Layers 1 and 5 are the input and output layers, respectively. Layer 2 uses a multi-level wavelet decomposition network, layer 3 uses a bidirectional gated recurrent neural network, and layer 4 is a fully connected neural network consisting of two linear layers. Its model structure is as follows: Figure 2 As shown; Step 2.2: Implement the multi-level wavelet decomposition network model. First, define the weight matrix. W l ( i )and W h ( i The formula is in the form of ) and is as follows: , In the formula, W l ( i )and W h ( i Let be the low-pass weight matrix and high-pass weight matrix of the i-th level wavelet decomposition; k The length of the filter coefficients; l 1. l 2、...、 l k These are the coefficients of the low-pass filter; h 1. h 2、...、 h k These are the coefficients of the high-pass filter; Fill elements in the matrix; Using wavelet basis functions to weight matrix W l ( i )and W h ( i Initialize the wavelet basis functions by selecting db4 as the type. The wavelet basis functions are as follows: , In the formula, l and h These are the low-pass filter coefficients and high-pass filter coefficients; The weight matrix is ​​used as a filter, and then the time series data is processed through multi-level decomposition, as shown in the following formula: , In the formula, a l ( i ), ah ( i ) represents the low-frequency and high-frequency intermediate variables of the i-th level wavelet decomposition; σ ( · ) is the activation function; W l (i) and W h (i) represents the low-pass weight matrix and high-pass weight matrix of the i-th level wavelet decomposition; x l (i-1) is the low-frequency approximation signal obtained from the (i-1)th level decomposition; b l (i) and b h (i) is the trainable vector bias, with its initial value set to 0.0001; The low-frequency approximate signal is obtained after the i-th level wavelet decomposition. x l ( i and high-frequency detail signals x h ( i ),right x l ( i When performing the (i+1)th level decomposition, a sliding convolution is performed on it using a low-pass filter and a high-pass filter to obtain the intermediate variable sequence. The convolution formula is as follows: , In the formula, and Let n be the nth element of the low-frequency intermediate variable sequence and the nth element of the high-frequency intermediate variable sequence of the (i+1)th level wavelet decomposition; Let be the (n+k-1)th element of the i-th level low-frequency subsequence; K The length of the filter; k Internal index of the filter coefficients; l k and h k These are the low-pass filter coefficients and high-pass filter coefficients; intermediate variables a l ( i Half-sampling is performed using average pooling. The average pooling formula is as follows: , In the formula, j This represents the length of the high- and low-frequency subsequences after the i-th level decomposition; and Let j be the j-th element of the low-frequency approximation signal and the j-th element of the high-frequency detail signal obtained from the i-th level wavelet decomposition. and Let be the 2j-th element and the (2j-1)-th element in the low-frequency intermediate variable sequence; and Let be the 2j-th element and the (2j-1)-th element in the high-frequency intermediate variable sequence; Add regularization constraints; the regularization constraint formula is as follows: , In the formula, W l ( i )and W h ( i Let be the low-pass weight matrix and high-pass weight matrix of the i-th level wavelet decomposition; , Weight matrix W l ( i ), W h ( i The initial state of ); L *、 L ( θ ), L l , L h The regularization loss is the final total loss function, the basic loss function, the low-frequency weight matrix, and the high-frequency weight matrix. α , β These are the regularization parameters for low and high frequencies; Then, the network parameters are iteratively updated using the BP algorithm, with the following update formula: , In the formula, W l ( i )and W h ( i Let be the low-pass weight matrix and high-pass weight matrix of the i-th level wavelet decomposition; , Weight matrix W l ( i ), W h ( i The initial state of ); α , β These are the regularization parameters for low and high frequencies; L The loss function; η The learning rate; During the forward propagation of the model, the final regularized loss value is compared with that of each wavelet decomposition. L l and Lh The weighted calculations, stacking, and averaging are performed, and then the results are output together and used in the backpropagation optimization process of the loss function. Step 2.3: Configure basic network parameters and data flow dimensions. The initial input data dimension is (256, 128, 6), and its data dimension format can be understood as "number of samples × sequence length × number of features". In this embodiment, the multi-level wavelet decomposition has 3 levels. In the mWDN module, the feature dimension is downsampled by half through wavelet transform, and the data dimension is gradually compressed from (256, 128, 6) to (256, 64, 6), (256, 32, 6), and (256, 16, 6), and finally (256, 16, 6) is used as the current output dimension. In the GRU module, the sequence dependency is first learned through 4 layers of BiGRU, and its hidden dimension d The data dimension is expanded to (256, 16, 128) by setting the input dimension to 128. Then, the sequence length is restored to 128 through four interpolation upsampling operations in the nearest upsampling mode. At this point, the data dimension is (256, 128, 128). Finally, the results of the four interpolations are concatenated on a feature dimension dim of -1, expanding the data dimension to (256, 128, 128 × 4). In the FCN module, the feature dimension is mapped back to the initial input dimension d of 6 through a Linear layer, and the data dimension becomes (256, 128, 6). Then, the feature dimension is compressed to 1 through a Linear layer, resulting in the final output data dimension (256, 128, 1). Step 2.4: Configure model hyperparameters: Set the wavelet decomposition level to 3, the wavelet basis function type to db4, the empirical value of the low-frequency regularization parameter of the mWDN matrix to 0.01, the empirical value of the high-frequency regularization parameter of the mWDN matrix to 0.01, the RNN hidden layer dimension to a multiple of 2, with a value range of 32~512, and 128 is selected in this embodiment, the empirical value of the number of RNN hidden layers ranges from 1 to 6, and 2 is selected in this embodiment, the FCN hidden layer dimension is 1-4 times the hidden unit in the RNN layer, and 256 is selected in this embodiment.

[0020] Step 3, Model Training and Evaluation: Use the training set built in Step 1 to train the model built in Step 2, and construct the training framework of deep temporal network from three dimensions: loss function, optimizer, and evaluation metric.

[0021] Specifically, the details of step 3 are as follows: Step 3.1: Select MSE as the loss function, and its formula is as follows: , In the formula, y i For the true value, For predicted values, N The length of the logging sequence; Step 3.2: Select Adam as the training optimizer. During training, a learning rate scheduling strategy combining warm-up and cosine annealing is adopted. In the early stage of training, the learning rate is increased linearly. In the later stage of training, the learning rate is decreased using a cosine annealing scheduling strategy. The formulas for warm-up and cosine annealing are as follows: , In the formula, T warm To control the number of epochs for preheating and cosine annealing; T max This represents the maximum number of training rounds. t The number of iterations that have been completed for model training; η t Let be the learning rate at the t-th training round; η 0 represents the initial learning rate; η min Minimum learning rate; Monitoring and evaluation indicators, and adopting an early termination mechanism to end training in advance; Step 3.3: Select the Mean Absolute Error (MAE), Pearson correlation coefficient, and Coefficient of Determination (R²). 2 As evaluation metrics, all metrics are calculated under the original physical dimensions after inverse normalization. The formula for Mean Absolute Error (MAE) is as follows: , In the formula, y i The actual value; This is a predicted value; N The length of the logging sequence; The formula for the Pearson correlation coefficient is as follows: , In the formula, t This represents the total number of sampling points; m i and n i These are the sampled point values ​​of the true curve and the reconstructed curve; Coefficient of determination R 2 The range of values ​​is (-∞, 1]. The formula is as follows: , In the formula, N The length of the logging sequence; The mean of the true curve; y i The actual value; This is a predicted value; Step 3.4: Configure training hyperparameters: The batch size (batch_size) is 256, used to control the number of input samples for each training iteration. It needs to match the GPU memory and be a power of 2 to ensure gradient stability. The initial learning rate (lr) is 0.001, used for the step size of model parameter updates. Based on the optimizer type and experience, the weight decay coefficient (wdecay) is set to 1e-4, which prevents overfitting by adding a penalty term to the weights. The gradient clipping ratio (clip) is set to 0.1, used to limit the maximum range of gradients and avoid gradient explosion during training. The Early Stopping patience value (patience) is set to 5, used for the number of epochs to continue training when the validation set evaluation metrics show no improvement, preventing model overfitting.

[0022] Step 4: Reconstructing logging curves for missing well sections: Verify the model's generalization performance using the test set, save the model parameter file, and use the trained model to reconstruct and predict logging curves for missing well sections.

[0023] Specifically, the details of step 4 are as follows: Step 4.1, Model Saving and Reloading: Establish a file management system hierarchically according to the prediction target, training dataset, prediction model, and input sequence length. Save the training model in full-model saving mode, and store it in the file naming format of "{model name}-val_mae_{prediction accuracy}.pth". Use the state dictionary parameters to save data and training process parameters, store them in the same directory as the model file, and use the same naming method, storing them in JSON file format. At the same time, use the command line parameter management system in the argparse library to save all process parameters. During model reloading, first load the state dictionary containing data preprocessing parameters, and then load the complete model file to reproduce the same environment as during training. Step 4.2, Model Inference and Output: Input the preprocessed logging sequence data into the model to perform the forward propagation calculation process, and directly output the predicted value of the target acoustic transit time AC curve. The prediction results are used to evaluate the accuracy and analyze the effect on the verification well and test well to construct the logging curve of the data missing well section.

[0024] Example 1: This example discloses a method for reconstructing well logging curves based on a bidirectional gated recurrent neural network using multi-level wavelet decomposition. The method is performed according to the following steps: 1. This embodiment uses the original sandstone and mudstone logging dataset from Shengli Oilfield as the data source, consisting of 100 wells. The dataset is divided into training wells, validation wells, and test wells in a split ratio of 80:10:10. The input features only include four curves: DEN, CNL, GR, and CAL, with AC used as the prediction curve. In the network model, the number of network layers is set to 5: an input layer, a multi-level wavelet decomposition layer, a bidirectional gated recurrent neural network, a fully connected neural network layer, and an output layer. The wavelet decomposition level is 3, the bidirectional gated recurrent neural network layer is 2, and the hidden layer dimension is 32. For data prediction, a seq2seq method is used, the loss function is MSE, and the optimizer is Adam with a learning rate of 0.0001. Training ends when the loss on the validation set does not decrease for 10 consecutive rounds.

[0025] 2. The method proposed in this invention is used for training and predictive reconstruction of model parameters. The final reconstruction result is as follows: Figure 3 As shown in the figure; the first line represents depth, the second line represents logging lithology, the third line represents the wellbore caliber (CAL) and natural gamma ray (GR) curves, the fourth line represents the neutron density (CNL) and density (DEN) curves, the fifth line represents the interpretation conclusions, and the sixth line represents the comparison between the measured acoustic curves and the reconstructed acoustic curves.

[0026] 3. Validation set well A, the mean absolute error (MAE) of this well under the mWDN-BiGRU model is 29.727, R0 2 The correlation coefficient was 0.798, and the Pearson correlation coefficient was 0.897. Figure 3 Using the 2910-2930m well section as the demonstration section, the prediction accuracy results of the mWDN-BiGRU model showed good performance, and the overall acoustic reconstruction curve and the measured curve had a good matching effect.

[0027] 4. Validation set well B, the mean absolute error (MAE) of this well under the mWDN-BiGRU model is 23.054, R0 2 The correlation coefficient is 0.756, and the Pearson correlation coefficient is 0.895. Figure 4 The mWDN-BiGRU model demonstrated good performance in the prediction results of the upper water layer (2786-2797m) and the lower oil layer (2839-2843m). In the middle section of the large mudstone interlayer (2797-2839m), the mWDN-BiGRU model significantly improved the problem of error accumulation with sequence in this type of layer, and had better adaptability. Overall, the acoustic reconstruction curve and the measured curve had a good matching effect.

[0028] 5. Test set well C, the mean absolute error (MAE) of this well under the mWDN-BiGRU model is 13.819, R 2The correlation coefficient was 0.789, and the Pearson correlation coefficient was 0.921. Figure 5 This model can not only accurately capture the spatial distribution pattern of reservoir physical properties, but also achieve high-precision reproduction of the target curve morphology, demonstrating the superiority of the proposed model in well logging curve reconstruction tasks.

[0029] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions and improvements made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

Claims

1. A method for reconstructing well logging curves based on a bidirectional gated recurrent neural network using multi-level wavelet decomposition, characterized in that, The implementation steps are as follows: Step 1: Raw logging data preprocessing and dataset construction: The raw logging data is subjected to quality checks and outlier cleaning. The raw logging data of different wells are divided into training set, validation set and test set according to the well number. The dataset is constructed by data transformation and normalization, dataset encapsulation, data iterator construction and data hyperparameter configuration. Step 2: Construct the mWDN-BiGRU model framework: First, extract the frequency information of the time series data by using the multi-level wavelet decomposition mWDN module. Connect the mWDN module with multiple GRU modules to form a parallel data flow structure, and set the basic network parameters and data flow dimensions, and configure the model parameters. Step 3, Model Training and Evaluation: Use the training set built in Step 1 to train the model built in Step 2, and build a training framework for deep temporal networks from three dimensions: loss function, optimizer, and evaluation metric. Step 4: Reconstructing logging curves for missing well sections: Verify the model's generalization performance using the test set, save the model parameter file, and use the trained model to reconstruct and predict logging curves for missing well sections.

2. The method for reconstructing well logging curves based on a bidirectional gated recurrent neural network using multi-level wavelet decomposition as described in claim 1, characterized in that, Step 1 includes the following steps: Step 1.1: Check and clean the raw logging data. Remove outliers caused by downhole environmental interference and instrument malfunctions from the raw logging data. Ensure that the value range and distribution of the effective data of each logging curve conform to the actual geological and physical meaning. At the same time, correct systematic errors caused by instrument differences and changes in the measurement environment. Step 1.2: Dataset construction. Based on the number of effective wells in each region, the original logging data from the effective well sections is divided into training set, validation set and test set according to the well number. This ensures that the data from the same well appears in only one of the datasets. The datasets are saved and read in three file formats: train.csv, valid.csv and test.csv. Step 1.3, Data Conversion and Normalization: The resistivity RT curve is converted into an RT_log curve using logarithmic transformation. The formula for logarithmic transformation is as follows: ， In the formula, RT lg The resistivity curve after logarithmic transformation; RT The resistivity curve; This is the offset; The transformed training set data was extracted using the Min-Max normalization method and stored as a JSON file in dictionary form. The formula for linear normalization is as follows: ， In the formula, V norm The value is the normalized value; V log These are the original curve values; V max The maximum value of the curve; V min This is the minimum value of the curve; Step 1.4: Sequence data segmentation; the original dataset is represented as follows. D ∈R n×m Where n is the number of samples in a certain well dataset, and m is the number of features in the dataset, the dataset is divided based on the total amount of data and the target sequence length to obtain datasets corresponding to the sequence length. D i Where i = 1, 2...l, and l is the length of the target segmented sequence, then the array stacking function from the NumPy library in Python is used. np.stack Implement l datasets D 1, D 2, ..., D l Stacking adds sequence dimensions, resulting in a 3D tensor dataset with sequence dimensions. The dataset is expressed as follows: ， In the formula, D new For a new three-dimensional tensor dataset; D 1, D 2, ..., D l The resulting sub-datasets; axis The axial direction parameter for the stacking operation; Step 1.5: Dataset Encapsulation and Data Iterator Construction. Based on the PyTorch framework, the Well_Dataset class for logging sequence datasets is defined by inheriting the abstract class torch.utils.data.Dataset to encapsulate the dataset. The dataset is initialized, receiving the segmented sequence sample data and then returning the total number of samples in the dataset. The range of index values ​​is determined based on the returned total number of samples, and then a single training sample is returned based on the index and converted into a tensor. Again, based on the PyTorch framework, the torch.utils.data.DataLoader class is called, using Well_Dataset to initialize the data iterator. Step 1.6: Configure the data hyperparameters, including: standardization method, input sequence length, and sliding step size.

3. The method for reconstructing well logging curves based on a bidirectional gated recurrent neural network using multi-level wavelet decomposition as described in claim 1, characterized in that, Step 2 includes the following steps: Step 2.1: Construct the network framework. The network structure consists of five parts: an input layer, a multi-level wavelet decomposition network layer, a bidirectional recurrent gated neural network layer, a fully connected neural network layer, and an output layer. The first and fifth layers are set as the input layer and the output layer, respectively. The second layer uses a multi-level wavelet decomposition network, the third layer uses a bidirectional gated recurrent neural network, and the fourth layer is a fully connected neural network composed of two linear layers. Step 2.2: Implement the multi-level wavelet decomposition network model. First, define the weight matrix. W l ( i )and W h ( i The formula is in the form of ) and is as follows: ， In the formula, W l ( i )and W h ( i Let be the low-pass weight matrix and high-pass weight matrix of the i-th level wavelet decomposition; k The length of the filter coefficients; l 1. l 2、...、 l k These are the coefficients of the low-pass filter; h 1. h 2、...、 h k These are the coefficients of the high-pass filter; Fill elements in the matrix; Using wavelet basis functions to weight matrix W l ( i )and W h ( i The process involves initializing the data, using the weight matrix as a filter, and then processing the time series data through multi-level decomposition, as shown in the following formula: ， In the formula, a l ( i ), a h ( i ) represents the low-frequency and high-frequency intermediate variables of the i-th level wavelet decomposition; σ(·) is the activation function; W l (i) and W h (i) represents the low-pass weight matrix and high-pass weight matrix of the i-th level wavelet decomposition; x l (i-1) is the low-frequency approximation signal obtained from the (i-1)th level decomposition; b l (i) and b h (i) represents the trainable vector bias; The low-frequency approximate signal is obtained after the i-th level wavelet decomposition. x l ( i and high-frequency detail signals x h ( i ),right x l ( i When performing the (i+1)th level decomposition, a sliding convolution is performed on it using a low-pass filter and a high-pass filter to obtain the intermediate variable sequence. The convolution formula is as follows: ， In the formula, and Let n be the nth element of the low-frequency intermediate variable sequence and the nth element of the high-frequency intermediate variable sequence of the (i+1)th level wavelet decomposition; Let be the (n+k-1)th element of the i-th level low-frequency subsequence; K The length of the filter; k Internal index of the filter coefficients; l k and h k These are the low-pass filter coefficients and high-pass filter coefficients; The intermediate variable is downsampled using average pooling, as shown in the following formula: ， In the formula, j This represents the length of the high- and low-frequency subsequences after the i-th level decomposition; and Let j be the j-th element of the low-frequency approximation signal and the j-th element of the high-frequency detail signal obtained from the i-th level wavelet decomposition. and Let be the 2j-th element and the (2j-1)-th element in the low-frequency intermediate variable sequence; and Let be the 2j-th element and the (2j-1)-th element in the high-frequency intermediate variable sequence; Then, a regularization constraint is added. The regularization constraint formula is as follows: ， In the formula, W l ( i )and W h ( i Let be the low-pass weight matrix and high-pass weight matrix of the i-th level wavelet decomposition; , Weight matrix W l ( i ), W h ( i The initial state of ); L* , L ( θ ), L l , L h The regularization loss is the final total loss function, the basic loss function, the low-frequency weight matrix, and the high-frequency weight matrix. α , β These are the regularization parameters for low and high frequencies; Then, the network parameters are iteratively updated using the backpropagation (BP) algorithm, with the following update formula: ， In the formula, W l ( i )and W h ( i Let be the low-pass weight matrix and high-pass weight matrix of the i-th level wavelet decomposition; , Weight matrix W l ( i ), W h ( i The initial state of ); α , β These are the regularization parameters for low and high frequencies; L The loss function; η The learning rate; During the forward propagation of the model, the final regularized loss value is compared with that of each wavelet decomposition. L l and L h The weighted calculations, stacking, and averaging are performed, and then the results are output together and used in the backpropagation optimization process of the loss function. Step 2.3: Configure the network skeleton, connect the mWDN module with multiple GRU modules to form a parallel data flow structure, and set the basic network parameters and data flow dimensions; Step 2.4: Configure the model hyperparameters, specifically including: wavelet decomposition level, wavelet basis function type, mWDN low-frequency matrix regularization parameter, mWDN high-frequency matrix regularization parameter, RNN hidden layer dimension, RNN hidden layer number, and FCN hidden layer dimension.

4. The method for reconstructing well logging curves based on a bidirectional gated recurrent neural network using multi-level wavelet decomposition as described in claim 1, characterized in that, Step 3 includes the following steps: Step 3.1: Select MSE as the loss function, and its formula is as follows: ， In the formula, y i The actual value; This is a predicted value; N The length of the logging sequence; Step 3.2: Select Adam as the training optimizer. During training, a learning rate scheduling strategy combining preheating and cosine annealing is adopted. In the early stage of training, the learning rate is increased linearly. In the later stage of training, the learning rate is decreased using a cosine annealing scheduling strategy. The preheating and cosine annealing formulas are as follows: ， In the formula, T warm To control the number of epochs for preheating and cosine annealing; T max This represents the maximum number of training rounds. t The number of iterations that have been completed for model training; η t Let be the learning rate at the t-th training round; η 0 represents the initial learning rate; η min Minimum learning rate; Simultaneously monitor and evaluate indicators, and adopt an early stop mechanism to terminate training in advance; Step 3.3: Select the Mean Absolute Error (MAE), Pearson correlation coefficient, and Coefficient of Determination (R²). 2 As evaluation metrics, all metrics are calculated under the original physical dimensions after inverse normalization. The formula for Mean Absolute Error (MAE) is as follows: ， In the formula, y i The actual value; This is a predicted value; N The length of the logging sequence; The formula for the Pearson correlation coefficient is as follows: ， In the formula, t This represents the total number of sampling points; m i and n i These are the sampled point values ​​of the true curve and the reconstructed curve; Coefficient of determination R 2 The formula is as follows: ， In the formula, N The length of the logging sequence; The mean of the true curve; y i The actual value; This is a predicted value; Step 3.4: Configure the training hyperparameters, including: training batch size, initial learning rate, weight decay coefficient, gradient pruning ratio, and the number of epochs to wait for when the Early Stopping validation set evaluation metric shows no improvement.

5. The method for reconstructing well logging curves based on a bidirectional gated recurrent neural network using multi-level wavelet decomposition as described in claim 1, characterized in that, Step 4 includes the following steps: Step 4.1: Model Saving and Reloading. A file management system is established hierarchically according to the prediction target, training dataset, prediction model, and input sequence length. The training model is saved in a full-model saving mode. Data and training process parameters are saved using a state dictionary parameter, stored in the same directory as the model file, and saved as JSON files. Simultaneously, the command-line parameter management system in the argparse library is used to save all process parameters. During model reloading, the state dictionary containing data preprocessing parameters is loaded first, followed by the complete model file, reproducing the environment consistent with training. Step 4.2, Model Inference and Output: Input the preprocessed logging sequence data into the model to perform the forward propagation calculation process, directly output the predicted value of the target logging curve, evaluate the accuracy and analyze the effect of the prediction results on the verification well and test well, and construct the logging curve of the data missing well section.

Images