Ensemble Deep Learning-Based Porosity Inversion From Seismic Attributes

Underground porosity is important in many earth sciences and engineering fields, including hydrocarbon reservoir characterization and geothermal energy production. Popular methods largely rely on the analysis of lithological core data, well logs, and seismic inversion methods. While these methods are reliable, they are also time-consuming, expensive, and difficult to implement. In addition, seismic inversion has nonlinearity, data dimensionality, and non-uniqueness issues. However, deep learning (DL) can provide a more flexible, efficient, and accurate capability by mapping directly from seismic attributes to underground porosity. Therefore, we trained several DL models with different optimization functions. In the training steps, we labelled every seismic attribute data point with its corresponding porosity derived from the well-logs. In contrast to popular ensemble techniques, we proposed a weighted prediction approach based on the strengths of each model. Testing results showed a coefficient of determination ( $\text{R}^{2}$ ) of 0.94345 and a Pearson’s correlation coefficient of 0.9725 between the actual model and the model of the proposed approach, versus 0.9681 and 0.9716 for the best single and popular ensemble models, respectively. Further, we tested the effectiveness of our method using real seismic data from the North Sea. With a Pearson’s correlation value of 0.9743, the inverted model ranges from 27 to 35%, compared to the reference model, which has an overall range of 20 to 33%. These results provide insights into the potential of the proposed method and its applicability to any other seismic volume to determine spatially varying underground porosity from seismic attributes directly.


I. INTRODUCTION
Porosity is a critical parameter in many applied geoscience and engineering fields for delineating and characterizing petroleum and geothermal energy reservoirs in order to further quantify reserves and drill production wells. Predicting and analyzing subsurface porosity enables understanding the properties of fluid flow [1], [2] and the mechanical behavior of rocks [3]. In addition, porosity is key in the assessment of subsurface contamination [4] and CO 2 sequestration [5]. Commonly, a variety of techniques are used to determine the seismic inversion suffer from strong nonlinearity, nonuniqueness, stability, and uncertainty. In addition, [15] and many other researchers noted that joint inversion increases the dimensionality of the data, making the problem very complex and requiring an exorbitant computational cost for evaluation.
However, innovative methods based on machine learning (ML) and DL, which have a strong ability to handle nonlinear features, are now proving to be very powerful in the field of earth sciences. Actually, ML and DL provide a collection of algorithms that imitate the ways of intelligent entities, learning from data to make predictions [16], to the extent that problems that an expert usually solves intuitively can also be automatically solved [17] using ML and DL. So far, by utilizing these algorithms, many researchers have found good results.
For instance, [18] proved the concept by using unsupervised ML and a colour feature blending technique to perform reservoir prediction from multiple seismic attributes. The results showed that the approach could help highlight interesting geological and hydrocarbon characteristics and improve traditional seismic interpretation techniques. Likewise, [19] used the Random Forest algorithm to deduce porosity from a number of seismic attributes, and the results proved the effectiveness of ML to characterize spatially varying porosity in a reservoir rock with quantification of the uncertainty. Reference [20] proposed applying DL to predict lithofacies from seismic data and found that the approach improved resolution in the presence of complex geological environments. In addition, [21] exploited the possibility of interpolating lithology between many wells by using geophysical inverse methods in a DL framework. Reference [22] discussed the possibility of using statistical and classification approaches of ML as a supporting workflow to the subsurface property characterization, integrating heterogeneous geophysical data that includes electromagnetic, seismic, gravity, and well data, which combines the approaches of both geophysical modelling and geophysical data inversion. Reference [23] evaluated the effectiveness of a multi-attribute ML method in the Denver-Julesburg basin, which resulted in enhanced seismic resolution. References [24] and [25] got much cleaner and more accurate interpreted geological structures than those obtained from post-stack structural seismic images, which were used to train a convolutional neural network to attenuate noise in marine seismic data. [26] used a multiresolution deep neural network to track seismic horizons. The results showed accurate predictions even in areas located far from known horizons. A state-of-the-art ML method for image processing successfully allowed [27] and [28] to delineate subtle salt domes and faults from 3D seismic images. Yet, [29] used synthetic seismic images to train a convolutional neural network to accurately estimate fault orientations in seismic images. Reference [30] discussed an automated approach to seismic interpretation by using DL.
As demonstrated by the works cited above, researchers have used single models to solve a specific problem. Thus, despite the success of applying these models individually, and since every algorithm may have its own ability to adapt to the training samples, meaning that, by using different algorithms, the user may obtain different models with different errors or accuracies, the reliability of the result may be affected. To address this issue, researchers developed the concept of ensemble learning, which allows for results that are more accurate with lower variance and errors. Thus, commonly, ensemble learning relies on bagging, stacking, or boosting approaches [31] to combine predictions from base models by voting and averaging them to improve the result. Although the latter has recently succeeded in many fields, no application of DL ensemble as a seismic inversion method has been reported, and in contrast to the common ensemble learningrelated techniques, no other approach has been performed.

II. DL POROSITY INVERSION A. DNN ALGORITHMS
In general, deep neural networks, referred to as ''DNNs,'' are designed to model a function that underpins the connections between predictor X and prediction Y in a flexible way that meets the user's strategy. Thus, many different learning techniques exist in practice. Some methods are used to handle classification problems alone, whereas others can only be used to solve regression scenarios.
The supervised regression learning approach was the focus of our study. Actually, for a specific supervised multiregression problem, the overall connection between X and Y is captured by This is a linear equation for predicting values of Y from those of x, which minimises the sum of the squared errors as follows: This expression gives the error for each predicted porosity. Considering i =1. . . n in (1), we can obtain n correspondences between multiple variables as follows: . .
Formulation (3) shows the link between several predictors and their associated predictions that our developed DNNs would try to replicate. Seismic inversion, on the other hand, can be viewed as a process of reconstructing a subsurface property model of interest, such as porosity, from the observed seismic data. Usually, the conventional seismic inversion methods are formulated as deterministic or stochastic processes, in which typical inversion workflows for reservoir properties such as porosity, permeability, lithology, and fluid saturations, the first step is the determination of elastic properties such as Pwave velocity, S-wave velocity, density, Poisson ratio, etc.
Next, one may apply established rock-physics models to convert those elastic properties into the desired reservoir property [32], [33], [34].
However, we designed our DNN-based inversion operator [35] to use multiple attributes from the observed seismic amplitudes as inputs to output a number of subsurface porosity models. As an optimization problem, the training process was defined as supervised learning by where i denotes the i th training sample, m i is the seismic amplitude model from which we generate seismic attribute data (d i ), and is the inversion operator used to predict the inverted porosity model ( ⌢ m i ) from the input attributes and the weight vector (W ).
We designed DNNs in a sequential architecture, using Python's Keras and TensorFlow libraries. We made up the baseline model with several hyper-parameters, including four dense layers in total with 500, 1000, 1000, and one node, respectively. Each hidden layer was activated by the rectified linear unit function (ReLU), after which a dropout probability of 0.1 was defined. A learning rate of 0.001 was defined to be used in each of the four selected optimizers, namely adaptive moment estimation (Adam) and SGD [36], as well as RMSProp and Nadam.
A batch size of 128 was used to fit the model. Both loss and evaluation metrics were used to track the performance of the models during the training and testing processes. Table 1 summarizes the trained networks' properties, which enable the networks to aggregate the inputs in a more sophisticated way to produce a rich prediction capability.

B. ENSEMBLE LEARNING
In this study, the ensemble learning approach was used effectively to improve the performance of the designed predictive models. Due to its ability to use insights from various models to aid in decision-making, ensemble learning has gained popularity recently [37]. Actually, ensemble learning is a remarkable technique that involves combining the performance capabilities of various models to provide better outcomes with a reduced variance while allowing for the avoidance of overfitting [38] and yielding a much superior model [39]. Because the ensemble model includes all the base models, even if one is incorrect, the other base models might be able to fix it.
There are three main methods of ensemble learning. One of them is the so-called ''boosting technique.'' The important feature of the boosting method is the idea of revising prediction errors. The algorithms are fitted and summed to the ensemble in order, with the second model attempting to correct the first model's predictions, the third model correcting the second model, etc. Simple averaging is used to combine the predictions of the poor models to obtain one strong model.
In the proposed approach, the contribution of each model is weighed proportionally to its performance, and Fig. 1 shows our workflow of a weighted ensemble approach to enhance porosity inversion from seismic attributes.
With our proposed ensemble approach, porosity inversion is a simple and low-cost operation that can be carried out after models have been trained. In the results and discussion section, we will compare the results of the DNN's porosity inversion for field data to those of one of the traditional seismic inversion methods.
The primary seismic data used to train the above-described DNNs is part of the F3 block seismic data in the offshore domain of the Netherlands in the North Sea. The size of the considered area ranges from 624 milliseconds to 924 milliseconds in depth and, laterally, from 400 to 600 inlines and 750 to 950 cross-lines. In total, we computed 26 seismic attributes to be used as porosity predictors, which were labeled by corresponding porosity and measured using the so-called horizon cube inversion method, developed by ''dGB Earth Sciences B.V.,'' which successfully blends well log attributes and seven horizons that are approved by professionals.
To boost model performance and accuracy by reducing data dimensionality, we performed feature selection based on the variance threshold method. This method removes all the features whose variance is zero, just like those with the same values in training samples. Further, we performed feature ranking using the random forest algorithm and finally got 14 attributes (Fig. 2), which are relevant for porosity prediction. Table 2 summarizes the final dataset used to train our DNN models.
The Inst., TB, HC, and SD suffixes in ''Attributes'' column of the table refer to instantaneous, thick beds, horizon cubes, and spectral decomposition, respectively. Seismic attributes are quantities derived from seismic data to highlight subtle information in the geological interpretation of the data. The documentation regarding the above-chosen attributes can be accessed online at the ''dGB Earth Sciences B.V.'' website.

2) OPTIMIZERS AND METRICS
Training a DNN is a challenging task. The model may look well-trained, but if it is unable to make accurate predictions, then it is useless. According to [40] and [41], this problem is commonly due to the challenges of fitting the algorithm during training steps. Overfitting occurs when the algorithm predicts well on training samples but not on new data. If there is a lack of data in the training steps, an under-fitted model cannot predict based on training samples. Therefore, a good generalization ability is critical for any predictive model. During the training steps, we must provide sufficient data so that the algorithm can efficiently learn from it. In addition, we have to set its parameters in such a way that its predictions and training labels are close to each other.
One of the powerful practices is using the cross-validation technique, which allows us to use every data point from the whole dataset in the training and testing steps [42], [43]. This popular technique usually results in less biased estimation of the DL-trained algorithm's skill. Since in supervised learning every training data point is labeled, meaning that the outcome is known, we are able to compare predictions with labels to the extent that we can change the parameters of the algorithm until both predictions and labels fit well [44]. The main goal here is to give the algorithm enough ability to successfully generalize and strengthen its interpretability.
As an optimization problem, the aim is to determine the DNN's efficiency at the lowest possible cost. The weights and biases of the model have to be updated through optimization functions to lower errors.
The Adaptive Momentum Estimation, also known as Adam [45], [46], and the Nesterov Accelerated Adaptive Moment Estimation, popularly known as Nadam [47], were two of the four optimizers that we used in this study. However, choosing an appropriate optimization function helps the algorithm produce optimal and rapid predictions.
An optimized model will also result in minimal errors and losses on the testing dataset, which represents the model's real-world performance. This testing data is a portion of the data that we have divided from the full training data. The algorithm won't observe or learn from it during the training process, in other words. Consequently, we predict using the testing data in order to assess the algorithm's performance.
Recently, it has been quite common to apply certain advanced means for dividing the complete dataset into training and testing samples. In this work, we applied a technique known as ''k-fold cross-validation. '' By allowing each data point from the primary data to be used separately in the training and validation processes, this data sampling approach divides the original data into several subsets. The parameter k represents how many subsets were extracted from the entire set of data. We set k to three so that the original data could be divided into three sets. The k-fold cross-validation approach is a trustworthy way of evaluating algorithm performance.
Comparative evaluation performances of all the trained DNNs using a 3-fold cross-validation approach are shown in Figs. 3, 4, 5, 6, and 7. In total, four models have been trained with different optimizers. The number of epochs was set to 30, and for all the models, the loss was decreasing with the number of epochs, while the R 2 that gives insights into the accuracy of the trained models was increasing with the number of epochs.
Since we have defined our porosity inversion problem to be a multiple regression task, according to [48], we can use the mean squared error (MSE) expressed below to measure the loss of the trained models. MSE is the squared average of the difference between the measured and inverted porosities.
Like in the estimate of classification model accuracy, the R 2 in the equation below gives a comprehensive measure of the regression model's precision by indicating the extent to which The NADAM optimizer-based model showed the lowest loss, followed by the ADAM, RMSP, and SGD optimizer-based models. After 11 epochs, the RMSP model tended to overfit. The early stopping of the hyperparameter helped to save its best performance on the 11 th . While a small gap can be seen between ADAM and Fig. 1. NADAM-based models, RMSP and SDG-based models show a considerable gap between them, which is why we did not want to display the SGD learning curve to keep the detail of the other curves visible. The NADAM optimizer-based model showed the highest R 2 , followed by the ADAM, RMSP, and SGD optimizer-based models, respectively. As for the loss, after 11 epochs, the RMSP model tended to overfit, but the early stopping of the hyperparameter helped to stop the training and save its best performance at the 11 th epoch. While a small gap can be seen between ADAM and NADAM-based models, RMSP and SDG-based models show a very large gap between them, so we did not display the SGD learning curve to keep the detail of the other curves visible.
the dependent variable can be predicted.  the NADAM optimizer-based model showed the lowest error, respectively, followed by ADAM, RMSP, and SGD optimizer-based models. After 11 epochs, the RMSP model was tending to overfit. The early stop of the hyper-parameter intervened by stopping the training and saving the best performance at the 11th epoch. A small gap can be seen between the ADAM and NADAM models, but the RMSP and SDG models show a considerable gap between them. The SGD learning curve is not displayed to keep the detail of the other curves visible.
In addition to that, we used the mean absolute error (MAE), expressed below, to measure the error of the trained DNNs. This function has the advantage of not considering variable directions. It calculates a linear score by equally weighting the average of all the individual differences.

MAE y,ŷ metric
The root mean squared error (RMSE), as shown by the following equation, is a quadratic scoring rule that is used to square the errors before they are averaged. For comparison with the MAE, the RMSE gives high weight to large errors. For this reason, it is useful where large errors are undesirable. In all of these equations, ⌢ y i and y i refer to the inverted and measured porosities, respectively.
As shown in Figs. 8, 9 and 10, during the testing step, the weight-averaging ensemble DNN resulted in the best findings, followed by the popular averaging-based ensemble. The NADAM DNN-based porosity inversion model resulted in superior results among all the DNNs taken individually. RMSE y,ŷ metric = 1 n samples

III. APPLICATION EXAMPLE A. FIELD DATA RESULTS
After training and testing the different DNNs, we performed porosity inversion using the pre-trained models. The data used herein is a subset cube cropped from F3 seismic data. In order to give insights into the distribution and quality of the inversion results, in this section we display and comparatively analyze several maps drawn from raw seismic data: the results from one of the conventional inversion methods, called horizon cube, and those of DNNs, namely, single DNN, common ensemble DNN, and the proposed weighted DNN. The maps highlight subtle differences between them in accordance with the performance power of each method. We consider the horizon cube result as our reference model.

B. PERFORMANCE DISCUSSION
In addition to the spatial distribution maps shown above, we computed Pearson correlation coefficients for all the inversion results in the testing steps, which we reported in Table 3. Further, using the trained models, we applied inversion to a new data set and compared each inverted model with the reference model in Table 4. Pearson correlation is a way of quantifying the relationship between inverted results and the measured, referred to as the reference model herein. Pearson's method outputs a range of coefficients that take values between minus one, indicating a negative linear correlation, which is perfect, and plus one, indicating a positive linear correlation, which is also perfect.
When the coefficient is zero, there is no linear correlation at all [34], and the closer the coefficients are to one, the stronger the relationship between the results and the reference porosity model. It is important to keep in mind that the geology of the chosen application interval is generally acknowledged to be composed of deposits from a significant fluvial-deltaic system, as stated in the TerraNubis-Data Info of the F3 Demo 2020 webpage. There are several seismic facies, such  However, it has been noted that the data in this area is noisy, thus before computing the necessary attributes for use with our trained models, we applied a dip-steered median filter to decrease noise. The region has a high porosity overall (20-33%), and the predicted values of porosity, which broadly range from 0.27 to 0.35, are pretty close to those already known in the area. This high porosity might be connected to the sand-and shale-containing deltaic package.
In conclusion, the testing phase demonstrated outstanding performance for both the weighted ensemble DNN and common ensemble methods. The results from the NADAM DNN-based porosity inversion were the best when all the DNNs were considered separately. The weighted ensemble model WEDNN-Inversion, followed by the common ensemble model EDNN-Inversion, displayed the strongest correlation with the reference model HC-Inversion during the real inversion application. While the SDNN-Inversion of the NADAM DNN-based porosity inversion was the best of all DNNs when viewed separately, it fell short of the performance of the common ensemble.
Although the proposed methodology is repeatable, we cannot guarantee that the trained models will perform on any new datasets that do not have the same properties as those used during their training. We recommend quickly retraining the models for the new dataset. We also suggest that, in order to improve the generalization ability of the new models, future research would focus on comprehending the topology and other model parameters.

IV. CONCLUSION
Conventional seismic inversion methods suffer from problems of strong nonlinearity, non-uniqueness, model stability, and data dimensionality, which increase the computational costs and affect the expected solutions. To address these problems, in this study, we have proposed the approach of using DL, which has proven to be powerful for working with non-linear futures. While many researchers previously used single DL methods as well as ensemble learning-based approaches for subsurface porosity inversion from multiple seismic attributes, despite these methods producing good results, we tried another approach that we called ''weighted'' ensemble DL-based method for subsurface porosity inversion from multiple seismic attributes.
The proposed method, when compared to the traditional ensemble method, provides a more accurate spatial distribution of the results with lower error and minizers the computational cost. Testing of the proposed method resulted in an R 2 of 0.94345 and a Pearson's correlation coefficient of 0.9725 between actual measurements and the model of the proposed approach, versus 0.9681 and 0.9716 for the best single and popular ensemble models, respectively. The application case on F3 seismic data proved the effectiveness of the suggested approach.
With a Pearson's correlation value of 0.9743, the inverted model ranges from 27 to 35%, compared to the reference model, which has an overall range of 20 to 33%. The trained model needs only fourteen seismic attributes from the original seismic data (reduced from 26 during the training stage) to make a porosity inversion. This study proved that a weighted ensemble learning approach could be used to more precisely infer any spatially changing reservoir property than porosity itself.