By Jonathan Leal, Rafael Jerónimo, Fabian Rada, Reinaldo Viloria and Rocky Roden Published with permission: First Break Volume 37, September 2019
A new approach has been applied to discriminate Net Reservoir using multi-attribute seismic analysis at single sample resolution, complemented by bivariate statistical analysis from petrophysical well logs. The combination of these techniques was used to calibrate the multi-attribute analysis to ground truth, thereby ensuring an accurate representation of the reservoir static properties and reducing the uncertainty related to reservoir distribution and storage capacity. Geographically, the study area is located in the south of Mexico. The reservoir rock consists of sandstones from the Upper Miocene age in a slope fan environment.
The first method in the process was the application of Principal Component Analysis (PCA), which was employed to identify the most prominent attributes for detecting lithological changes that might be associated with the Net Reservoir. The second method was the application of the Kohonen Self-Organizing Map (SOM) Neural Network Classification at voxel scale (i.e., sample rate and bin size dimensions from seismic data), instead of using waveform shape classification. The sample-level analysis revealed significant new information from different seismic attributes, providing greater insights into the characteristics of the reservoir distribution in a shaly sandstone. The third method was a data analysis technique based on contingency tables and Chi-Square test, which revealed relationships between two categorical variables (SOM volume neurons and Net Reservoir). Finally, a comparison between a SOM of simultaneous seismic inversion attributes and traditional attributes classification was made corroborating the delineated prospective areas. The authors found the SOM classification results are beneficial to the refinement of the sedimentary model in a way that more accurately identified the lateral and vertical distribution of the facies of economic interest, enabling decisions for new well locations and reducing the uncertainty associated with field exploitation. However, the Lithological Contrast SOM results from traditional attributes showed a better level of detail compared with seismic inversion SOM.
Self-Organizing Maps (SOM) is an unsupervised neural network – a form of machine learning – that has been used in multi-attribute seismic analysis to extract more information from the seismic response than would be practical using only single attributes. The most common use is in automated facies mapping. It is expected that every neuron or group of neurons can be associated with a single depositional environment, the reservoir´s lateral and vertical extension, porosity changes or fluid content (Marroquín et al., 2009). Of course, the SOM results must be calibrated with available well logs. In this paper, the authors generated petrophysical labels to apply statistical validation techniques between well logs and SOM results. Based on the application of PCA to a larger set of attributes, a smaller, distilled set of attributes were classified using the SOM process to identify lithological changes in the reservoir (Roden et al., 2015).
A bivariate statistical approach was then conducted to reveal the relationship between two categorical variables: the individual neurons comprising the SOM classification volume and Net Reservoir determined from petrophysical properties (percentage of occurrence of each neuron versus Net Reservoir).
The Chi-Square test compares the behavior of the observed frequencies (Agresti, 2002) for each SOM neuron lithological contrast against the Net Reservoir variable (grouped in “Net Reservoir” and “no reservoir” categories). Additional data analysis was conducted to determine which neurons responded to the presence of hydrocarbons using box plots showing Water Saturation, Clay Volume, and Effective Porosity as Net Pay indicators. The combination of these methods demonstrated an effective means of identifying the approximate region of the reservoir.
About the Study Area
The reservoir rock consists of sandstones from the Upper Miocene age in a slope fan environment. These sandstones correspond to channel facies, and slope lobes constituted mainly of quartz and potassium feldspars cemented in calcareous material of medium maturity. The submarine slope fans were deposited at the beginning of the deceleration of the relative sea-level fall, and consist of complex deposits associated with gravitational mass movements.
Stratigraphy and Sedimentology
The stratigraphic chart comprises tertiary terrigenous rocks from Upper Miocene to Holocene. The litho-stratigraphic units are described in Table 1.
Figure 1 (left) shows the facies distribution map of the sequence, corresponding to the first platform-basin system established in the region. The two dashed lines – one red and one dark brown – represent the platform edge at different times according to several regional integrated studies in the area. The predominant direction of contribution for studied Field W is south-north, which is consistent with the current regional sedimentary model. The field covers an area of approximately 46 km2 and is located in facies of distributary channels northeast of the main channel. The reservoir is also well-classified and consolidated in clay matrix, and it is thought that this texture corresponds to the middle portion of the turbidite system. The observed electrofacies logs of the reservoir are box-shaped in wells W-2, W-4, W-5, and W-6 derived from gamma ray logs and associated with facies of distributary channels that exhibit the highest average porosity. In contrast, wells W-3 and W-1 are different – associated with lobular facies – according to gamma ray logs. In Figure 1 (right), a sedimentary scheme of submarine fans proposed by Muti (1978).
The Stieber model was used to classify Clay Volume (VCL). The Effective Porosity (PIGN) was obtained using the Neutron-Density model and non-clay water intergranular Water Saturation (SUWI) was determined to have a salinity of 45,000 ppm using the Simandoux model. Petrophysical cut-off values used to distinguish Net Reservoir and Net Pay estimations were 0.45, 0.10 and 0.65, respectively.
The reservoir rock corresponds to sands with Net Pay thickness ranging from 9-12 m, porosity between 18-25%, average permeability of 8-15 mD, and Water Saturation of approximately 25%. The initial pressure was 790 kg / cm2 with the current pressure is 516 kg/cm2. The main problems affecting productivity in this volumetric reservoir are pressure drop, being the mechanism of displacement the rock-fluid expansion, and gas in solution. Additionally, there are sanding problems and asphaltene precipitation.
Multidisciplinary information was collected and validated to carry out seismic multi-attribute analysis. Static and dynamic characterization studies were conducted in the study area, revealing the most relevant reservoir characteristics and yielding a better sense of the proposed drilling locations. At present, six wells have been drilled.
The original available seismic volume and associated gathers employed in the generation of multiple attributes and for simultaneous inversion were determined to be of adequate quality. At target depth, the dominant frequency approaches 14 Hz, and the interval velocity is close to 3,300 m/s. Therefore, the vertical seismic resolution is 58 m. The production sand has an average thickness of 13 m, so it cannot be resolved with conventional seismic amplitude data.
Principal Component Analysis (PCA)
Principal Component Analysis (PCA) is one of the most common descriptive statistics procedures used to synthesize the information contained in a set of variables (volumes of seismic attributes) and to reduce the dimensionality of a problem. Applied to a collection of seismic attributes, PCA can be used to identify the seismic attributes that have the greatest “contribution,” based on the extent of their relative variance to a region of interest. Attributes identified through the use of PCA are responsive to specific geological features, e.g., lithological contrast, fracture zones, among others. The output of PCA is an Eigen spectrum that quantifies the relative contribution or energy of each seismic attribute to the studied characteristic.
PCA Applied for Lithological Contrast Detection
The PCA process was applied to the following attributes to identify the most significant attributes to the region to detect lithological contrasts at the depth of interest: Thin Bed Indicator, Envelope, Instantaneous Frequency, Imaginary Part, Relative Acoustic Impedance, Sweetness, Amplitude, and Real Part. Of the entire seismic volume, only the voxels in a time window (seismic samples) delimited by the horizon of interest were analyzed, specifically 56 milliseconds above and 32 milliseconds below the horizon. The results are shown for each principal component. In this case, the criterion used for the selected attributes were those whose maximum percentage contribution to the principle component was greater than or equal to 80%. Using this selection technique, the first five principal components were reviewed in the Eigen spectrum. In the end, six (6) attributes of the first two principal components were selected (Figure 2).
Figure 2. PCA results for lithological contrast detection.
Simultaneous Classification of Seismic Attributes Using a Self-Organizing Maps (SOM) Neural Network (Voxel Scale)
The SOM method is an unsupervised classification process in that the network is trained from the input data alone. A SOM consists of components (vectors) called neurons or classes and input vectors that have a position on the map. The values are compared employing neurons that are capable of detecting groupings through training (machine learning) and mapping. The SOM process non-linearly maps the neurons to a two dimensional, hexagonal or rectangular grid. SOM describes a mapping of a larger space to a smaller one. The procedure for locating a vector from the data space on the map is to find the neuron with the vector of weights (smaller metric distance) closer to the vector of the data space. (The subject of this analysis accounted for seismic samples located within the time window covering several samples above and below the target horizon throughout the study area). It is important to classify attributes that have the same common interpretive use, such as lithological indicators, fault delineation, among others. The SOM revealed patterns and identified natural organizational structures present in the data that are difficult to detect in any other way (Roden et al., 2015), since the SOM classification used in this study is applied on individual samples (using sample rate and bin size from seismic data, Figure 2, lower right box), detecting features below conventional seismic resolution, in contrast with traditional wavelet-based classification methods.
SOM Classification for Lithological Contrast Detection
The following six attributes were input to the SOM process with 25 classes (5 X 5) stipulated as the desired output: Envelope, Hilbert, Relative Acoustic Impedance, Sweetness, Amplitude, and Real Part.
As in the PCA analysis, the SOM was delimited to seismic samples (voxels) in a time window following the horizon of interest, specifically 56 milliseconds above to 32 milliseconds below. The resulting SOM classification volume was examined with several visualization and statistical analysis techniques to associate SOM classification patterns with reservoir rock.
3D and Plan Views
One way of identifying patterns or trends coherent with the sedimentary model of the area is visualizing all samples grouped by each neuron in 3D and plan views using stratal-slicing technique throughout the reservoir. The Kohonen SOM and the 2D Colormap in Figure 3 (lower right) ensures that the characteristics of neighboring neurons are similar. The upper part of Figure 3 shows groupings classified by all 5x5 (25) neurons comprising the neural network, while in the lower part there are groupings interpreted to be associated with the reservoir classified by a few neurons that are consistent with the regional sedimentary model, i.e., neurons N12, N13, N16, N17, N22, and N23.
Figure 3. Plan view with geological significance preliminary geobodies from Lithological Contrast SOM. Below: only neurons associated with reservoir are shown.
Vertical Seismic Section Showing Lithological Contrast SOM
The observed lithology in the reservoir sand is predominantly made up of clay sandstone. A discrete log for Net Reservoir was generated to calibrate the results of the Lithological Contrast SOM, using cut-off values according to Clay Volume and Effective Porosity. Figure 4 shows the SOM classification of Lithological Contrast with available well data and plan view. The samples grouped by neurons N17, N21, and N22 match with Net Reservoir discrete logs. It is notable that only the well W-3 (minor producer) intersected the samples grouped by the neuron N17 (light blue). The rest of the wells only intersected neurons N21 and N22. It is important to note that these features are not observed on the conventional seismic amplitude data (wiggle traces).
Figure 4. Vertical section composed by the SOM of Lithological Contrast, Amplitude attribute (wiggle), and Net Reservoir discrete property along wells.
Stratigraphic Well Section
A cross-section containing the wells (Figure 5) shows logs of Gamma Ray, Clay Volume, perforations, resistivity, Effective Porosity, Net Reservoir with lithological contrast SOM classification, and Net Pay.
The results of SOM were compared by observation with discrete well log data, relating specific neurons to the reservoir. At target zone depth, only the neurons N16, N17, N21, and N22 are present. It is noteworthy that only W-3 well (minor producer) intersect clusters formed by neuron N17 (light blue). The rest of the wells intersect neurons N16, N21, N22, and N23.
Traditionally, Vertical Proportion Curves (VPC) are qualitative and quantitative tools used by some sedimentologists to define succession, division, and variability of sedimentary sequences from well data, since logs describe vertical and lateral evolution of facies (Viloria et al., 2002). A VPC can be modeled as an accumulative histogram where the bars represent the facies proportion present at a given level in a stratigraphic unit. As part of the quality control and revision of the SOM classification volume for Lithological Contrasts, this statistical technique was used to identify whether in the stratigraphic unit or in the window of interest, a certain degree of succession and vertical distribution of specific neurons observed could be related to the reservoir.
The main objective of this statistical method is to identify how specific neurons are vertically concentrated along one or more logs. As an illustration of the technique, a diagram of the stratigraphic grid is shown in Figure 6. The VPC was extracted from the whole 3D grid of SOM classification volume for Lithological Contrast, and detection was generated by counting the occurrence among the 25 neurons or classes in each stratigraphic layer in the VPC extracted from the grid. The VPC of SOM neurons exhibits remarkable slowly-varying characteristics indicative of geologic depositional patterns. The reservoir top corresponds to stratigraphic layer No. 16. In the VPC on the right, only neurons N16, N17, N21, and N22 are present. These neurons have a higher percentage occurrence relative to all 25 classes from the top of the target sand downwards. Corroborating the statistics, these same neural classes appear in the map view in Figure 3 and the vertical section shown in Figure 4. The stratigraphic well section in Figure 5 also supports the statistical results. It is important to note that these neurons also detected seismic samples above the top of the sand top, although in a lesser proportion. This effect is consistent with the existence of layers with similar lithological characteristics, which can be seen from the well logs.
Figure 6. Vertical proportion Curve to identify neurons related to reservoir rock.
Bivariate Statistical Analysis Cross Tabs
The first step in this methodology is a bivariate analysis through cross-tabs (contingency table) to determine if two categorical variables are related based on observing the extent to which the occurrence of one variable is repeated in the categories of the second. Given that one variable is analyzed in terms of another, a distinction must be made between dependent and independent variables. With cross tabs analysis, the possibilities are extended to (in addition to frequency analyzes for each variable, separately) the analyses of the joint frequencies or those in which the analysis unit nature is defined by the combination of two variables.
The result was obtained by extracting the SOM classification volume along wells paths and constructing a discrete well log with two categories: “Net Reservoir” and “not reservoir.” The distinction between “Net Reservoir” and “not reservoir” simply means that the dependent variable might have a hydrocarbon storage capacity or not. In this case, the dependent variable corresponds to neurons of SOM classification for Lithological Contrast volume. It is of ordinal type, since it has an established internal order, and the change from one category to another is not the same. The neurons go from N1 to N25, organized in rows. The independent variable is Net Reservoir, which is also an ordinal type variable. In this tab, the values organized in rows correspond to neurons from the SOM classification volume for Lithological Contrast, and in the columns are discrete states of the “Net Reservoir” and “not reservoir” count for each neuron. Table 2 shows that the highest Net Reservoir counts are associated with neurons N21 and N22 at 47.0% and 28.2% respectively. Conversely, lower counts of Net Reservoir are associated with neurons N17 (8.9%), N16 (7.8%) and N23 (8.0%).
Table 2. Cross Tab for Lithological Contrast SOM versus Net reservoir.
Neuron N21 was detected at reservoir depth in wells W-2 (producer), W-4 (abandoned for technical reasons during drilling), W-5 (producer) and W-6 (producer). N21 showed higher percentages of occurrence in Net Reservoir, so this neuron could be identified as indicating the highest storage capacity. N22 was present in wells W-1 and W-6 at target sand depth but also detected in wells W-2, W-4 and W-5 in clay-sandy bodies overlying the highest quality zone in the reservoir. N22 was also detected in the upper section of target sand horizontally navigated by the W-6 well, which has no petrophysical evaluation. N17 was only detected in well W-3, a minor producer of oil, which was sedimentologically cataloged as lobular facies and had the lowest reservoir rock quality. N16 was detected in a very small proportion in wells W-4 (abandoned for technical reasons during drilling) and W-5 (producer). Finally, N23 was only detected towards the top of the sand in well W-6, and in clayey layers overlying it in the other wells. This is consistent with the observed percentage of 8% Net Reservoir, as shown in Table 2.
Chi-Square Independence Hypothesis Testing
After applying the cross-tab evaluation, this classified information was the basis of a Chi-Square goodness-of-fit test to assess the independence or determine the association between two categorical variables: Net Reservoir and SOM neurons. That is, it aims to highlight the absence of a relationship between the variables. The Chi-Square test compared the behavior of the observed frequencies for each Lithological Contrast neuron with respect to the Net Reservoir variable (grouped in “Net Reservoir” and “no reservoir”), and with the theoretically expected frequency distribution when the hypothesis is null.
As a starting point, the null hypothesis formulation was that the Lithological Contrast SOM neuron occurrences are independent of the presence of Net Reservoir. If the calculated Chi-Square value is equal to or greater than a certain critical theoretical value, the null hypothesis must be rejected. Consequently, the alternative hypothesis must be accepted. Observe the results in Table 3 where the calculated Chi-Square is greater than the theoretical critical value (296 ≥ 9.4, with four degrees of freedom and 5% confidence level), so the null hypothesis of the independence of Net Pay with SOM neurons is rejected, leaving a relationship between Net Reservoir and Lithological Contrast SOM variables.
The test does not report a goodness of fit magnitude (substantial, moderate or poor), however. To measure the degree of correlation between both variables, Pearson’s Phi (φ) and Cramer’s V (ν) measures were computed. Pearson’s φ coefficient was estimated from Eq. 1.1.
where X2: Chi-Square and n : No. of cases
Additionally, Cramer’s V was estimated using Eq. 1.2.
In both cases, values near zero indicate a poor or weak relationship while values close to one indicate a strong relation. The authors obtained values for φ, and Cramer´s ν equals to 0.559 (Table 3). Based on this result, we can interpret a moderate relation between both variables.
Table 3. Calculated and theoretical Chi-Square values and its correlation measures.
Box-and-whisker plots were constructed to compare and understand the behavior of petrophysical properties for the range that each neuron intersects the well paths in the SOM volume. Also, these quantify which neurons of interest respond to Net Reservoir and Net Pay properties (Figure 7). Five descriptive measures are shown for a box-and-whisker plot of each property:
• Median (thick black horizontal line)
• First quartile (lower limit of the box)
• Third quartile (upper limit of the box)
• Maximum value (upper end of the whisker)
• Minimum value (lower end of the whisker)
The graphs provide information about data dispersion, i.e., the longer the box and whiskers, the greater the dispersion and also data symmetry. If the median is relatively centered of the box, the distribution is symmetrical. If, on the contrary, it approaches the first or third quartile, the distribution could be skewed to these quartiles, respectively. Finally, these graphs identify outlier observations that depart from the rest of the data in an unusual way (these are represented by dots and asterisks as less or more distant from the data center). Horizontal dashed green line is the cut-off value for Effective Porosity (PIGN >0.10) while the dashed blue line represents the cut-off value for Clay Volume (VCL>0.45) and, dashed beige line is cut-off value for Water Saturation (SUWI<0.65).
Based on these data and the resulting analysis, it can be inferred that neurons N16, N17, N21, N22, and N23 respond positively to Net Reservoir. Of these neurons, the most valuable predictors are N21 and N22 since they present lower clay content in comparison with neurons N16 and N23 and associated higher Effective Porosity shown by neurons N16, N17, and N23 (Figure 7a). Neurons N21 and N22 are ascertained to represent the best reservoir rock quality. Finally, neuron N23 (Figure 7b), can be associated with rock lending itself with storage capacity, but clayey and with high Water Saturation, which allows discarding it as a significant neuron. It is important to note that this analysis was conducted by accounting for the simultaneous occurrence of the petrophysical values (VCL, PIGN, and SUWI) on the neurons initially intersected (Figure 7a), and then on the portion of the neurons that pass Net Reservoir cut-off values (Figure 7b), and finally on the portion of the neurons that pass net-pay cut-off values (Figure 7c). For all these petrophysical reasons, the neurons to be considered as a reference to estimate the lateral and vertical distribution of Net Reservoir associated with the target sand are in order of importance, N21, N22, N16, and N17.
Figure 7. Comparison between neurons according to petrophysical properties: VCL (Clay Volume), PIGN (Effective Porosity) and SUWI (Water Saturation). a) SOM neurons for lithological contrast detection, b) Those that pass Net Reservoir cut-off and c) Those that pass Net Pay cut-off.
Simultaneous Seismic Inversion
During this study, a simultaneous prestack inversion was performed using 3D seismic data and sonic logs, in order to estimate seismic petrophysical attributes as Acoustic Impedance (Zp), Shear Impedance (Zs), Density (Rho), as well as P&S-wave velocities, among others. They are commonly used as lithology indicators, possible fluids, and geomechanical properties. Figure 8a shows a scatter plot from well data of seismic attributes Lambda Rho and Mu Rho ratio versus Clay Volume (VCL) and as discriminator Vp/Vs ratio (Vp/Vs). The target sand corresponds to low Vp/Vs and Lambda/Mu values (circled in the figure). Another discriminator in the reservoir was S-wave impedance (Zs) (Figure 8b). From this, seismic inversion attributes were selected for classification by SOM neural network analysis. These attributes were Vp/Vs ratio, Lambda Rho/Mu Rho ratio, and Zs.
Figure 8. Scatter plots: a) Lambda Rho and Mu Rho ratio versus VCL and Vp/Vs y b) Zs versus VCL and Vp/Vs.
Self-Organizing Map (SOM) Comparison
Figure 9 is a plan view of neuron-extracted geobodies associated with the sand reservoir. In the upper part, a SOM classification for Lithological Contrast detection obtained from six traditional seismic attributes is shown; and in the lower part, a different SOM classification for Lithological Contrast detection was obtained from three attributes of simultaneous inversion. Both results are very similar. The selection of SOM classification neurons from inversion attributes was done through spatial pattern recognition, i.e., identifying geometry/shape of the clusters related to each of 25 neurons congruent with the sedimentary model, and by using a stratigraphic section for wells that includes both SOM classifications tracks.
Figure 9. Plan view of neurons with geological meaning. Up: SOM Classification from traditional attributes. Down: SOM Classification from simultaneous inversion attributes.
Figure 10 shows a well section that includes a track for Net Reservoir and Net Pay classification along with SOM classifications from traditional attributes and a second SOM from simultaneous inversion attributes defined from SOM volumes and well paths intersection. In fact, only the neurons numbers with geological meaning are shown.
Figure 10. Well section showing the target zone with tracks for discrete logs from Net Reservoir, Net Pay and both SOM classifications.
Discussion and Conclusions
Principal Component Analysis (PCA) identified the most significant seismic attributes to be classified by Self-Organizing Maps (SOM) neural network at single-sample basis to detect features associated with lithological contrast and recognize lateral and vertical extension in the reservoir. The interpretation of SOM classification volumes was supported by multidisciplinary sources (geological, petrophysical, and dynamic data). In this way, the clusters detected by certain neurons became the inputs for geobody interpretation. The statistical analysis and visualization techniques enabled the estimation of Net Reservoir for each neuron. Finally, the extension of reservoir rock geobodies derived from SOM classification of traditional attributes was corroborated by the SOM acting on simultaneous inversion attributes. Both multi-attribute machine learning analysis of traditional attributes and attributes of seismic inversion enable refinement of the sedimentary model to reveal more precisely the lateral and vertical distribution of facies. However, the Lithological Contrast SOM results from traditional attributes showed a better level of detail compared with seismic inversion SOM.
Collectively, the workflow may reduce uncertainty in proposing new drilling locations. Additionally, this methodology might be applied using specific attributes to identify faults and fracture zones, identify absorption phenomena, porosity changes, and direct hydrocarbon indicator features, and determine reservoir characteristics.
The authors thank Pemex and Oil and Gas Optimization for providing software and technical resources. Thanks also are extended to Geophysical Insights for the research and development of the Paradise® AI workbench and the machine learning applications used in this paper. Finally, thank Reinaldo Michelena, María Jerónimo, Tom Smith, and Hal Green for review of the manuscript.
Agresti, A., 2002, Categorical Data Analysis: John Wiley & Sons.
Marroquín I., J.J. Brault and B. Hart, 2009, A visual data mining methodology to conduct seismic facies analysis: Part 2 – Application to 3D seismic data: Geophysics, 1, 13-23.
Roden R., T. Smith and D. Sacrey, 2015, Geologic pattern recognition from seismic attributes: Principal component analysis and self-organizing maps: Interpretation, 4, 59-83.
Viloria R. and M. Taheri, 2002, Metodología para la Integración de la Interpretación Sedimentológica en el Modelaje Estocástico de Facies Sedimentarias, (INT-ID-9973, 2002). Technical Report INTEVEP-PDVSA.
Seismic attributes can be both powerful and challenging to incorporate into interpretation and analysis. Recent developments with machine learning have added new capabilities to multi-attribute seismic analysis. In 2018, Geophysical Insights conducted a proof of concept on 100 square miles of multi-client 3D data jointly owned by Geophysical Pursuit, Inc. (GPI) and Fairfield Geotechnologies (FFG) in the Denver-Julesburg Basin (DJ). The purpose of the study was to evaluate the effectiveness of a machine learning workflow to improve resolution within the reservoir intervals of the Niobrara and Codell formations, the primary targets for development in this portion of the basin.
The seismic data are from Phase 5 of the GPI/Fairfield Niobrara program in northern Colorado. A preliminary workflow which included synthetics, horizon picking and correlation of 28 wells was completed. The seismic volume was re-sampled from 2 ms to 1 ms. Detailed well time-depth charts were created for the Top Niobrara, Niobrara A, B and C benches, Fort Hays and Codell intervals. The interpretations, along with the seismic volume, were loaded into the Paradise® machine learning application, and two suites of attributes were generated, instantaneous and geometric. The first step in the machine learning workflow is Principal Component Analysis (PCA). PCA is a method of identifying attributes that have the greatest contribution to the data and that quantifies the relative contribution of each. PCA aids in the selection of which attributes are appropriate to use in a Self-Organizing Map (SOM). In this case, 15 instantaneous attribute volumes, plus the parent amplitude volume, were used in the PCA and eight were selected to use in SOMs. The SOM is a neural network-based machine learning process that is applied to multiple attribute volumes simultaneously. The SOM produces a non-linear classification of the data in a designated time or depth window.
For this study, a 60-ms interval that encompasses the Niobrara and Codell formations was evaluated using several SOM topologies. One of the main drilling targets, the B chalk, is approximately 30 feet thick; making horizontal well planning and execution a challenge for operators. An 8 X 8 SOM applied to 1 ms seismic data improves the stratigraphic resolution of the B bench. The neuron classification also images small but significant structural variations within the chalk bench. These variations correlate visually with the geometric curvature attributes. This improved resolution allows for precise well planning for horizontals within the bench. The 25 foot thick C bench and the 17 to 25 foot thick Codell are also seismically resolved via SOM analysis. Petrophysical analyses from wireline logs run in seven wells within the survey by Digital Formation; together with additional results from SOMs show the capability to differentiate a high TOC upper unit within the A marl which presents an additional exploration target. Utilizing 2D color maps and geobodies extracted from the SOMs combined with petrophysical results allows calculation of reserves for the individual reservoir units as well as the recently identified high TOC target within the A marl.
The results show that a multi-attribute machine learning workflow improves the seismic resolution within the Niobrara reservoirs of the DJ Basin and results can be utilized in both exploration and development.
Introduction and preliminary work
The Denver-Julesburg Basin is an asymmetrical foreland basin that covers approximately 70,000 square miles over parts of Colorado, Wyoming, Kansas and Nebraska. The basin has over 47,000 oil and gas wells with a production history that dates back to 1881 (Higley, 2015). In 2009, operators in the Wattenberg field began to drill and complete horizontal wells in the chalk benches of the Niobrara formation and within the Codell sandstone. As of October 2018, approximately 9500 horizontal wells have been drilled and completed within Colorado and Wyoming in the Niobrara and Codell formations (shaleprofile.com/2019/01/29/niobrara-co-wy-update-through-october-2018).
The transition to horizontal drilling necessitated the acquisition of modern, 3D seismic data (long offset, wide azimuth) to properly image the complex faulting and fracturing within the basin. In 2011, Geophysical Pursuit, Inc., in partnership with the former Geokinetics Inc., embarked on a multi-year, multi-client seismic program that ultimately resulted in the acquisition of 1580 square miles of contiguous 3D seismic data. In 2018, Geophysical Pursuit, Inc. (GPI) and joint-venture partner Fairfield Geotechnologies (FFG) provided Geophysical Insights with seismic data in the Denver-Julesburg Basin to conduct a proof of concept evaluation of the effectiveness of a machine learning workflow to improve resolution within the reservoir intervals of the Niobrara and Codell formations, currently the primary targets for development in this portion of the basin. The GPI/FFG seismic data analyzed are 100 square miles from the Niobrara Phase 5 multi-client 3D program in northern Colorado (Figure 1). Prior to the machine learning workflow, a preliminary interpretation workflow was carried out, that included synthetics, horizon picking and well correlation on 28 public wells with digital data. The seismic volume was resampled from 2 ms to 1 ms. Time depth charts were made with detailed well ties for the Top Niobrara, Niobrara A, B, and C benches, Fort Hays and Codell. The interpretations, along with the re-sampled seismic amplitude volume, were loaded into the Paradise® machine learning application. The machine learning software has several options for computing seismic attributes, and two suites were selected for the study: standard instantaneous attributes and geometric attributes from the AASPI (Attribute Assisted Seismic Processing and Interpretation) consortium (http://mcee.ou.edu/aaspi/).
Figure 1: Map of GPI FFG multi-client program and study area outline
Geologic Setting of the Niobrara and Surrounding Formations
The Niobrara formation is late Cretaceous in age and was deposited in the Western Interior Seaway (Kaufmann, 1977). The Niobrara is subdivided into the basal Fort Hays limestone and the Smoky Hill member. The Smoky Hill member is further subdivided into three subunits informally termed Niobrara A, B, and C. These units consist of fractured chalk benches which are primary reservoirs with marls and shales between the benches which comprise source rocks and secondary reservoir targets. (Figure 2). The Niobrara unconformably overlies the Codell sandstone and is overlain by the Sharon Springs member of the Pierre shale.
The Codell is also late Cretaceous in age, and unconformably underlies the Fort Hays member of the Niobrara formation. In general, the Codell thins from north to south due to erosional truncation (Sterling, Bottjer and Smith, 2016). In the study area, the thickness of the Codell ranges from 18 to 25 feet. Lewis (2013) inferred an eastern provenance for the Codell with a limited area of deposition or subsequent erosion through much of the DJ Basin. Based upon geochemical analyses, Sterling and others (2016) state that hydrocarbons produced from the Codell are sourced from the Niobrara, primarily the C marl, and the thermal maturity provides evidence of migration into the Codell. The same study found that oil produced from the Niobrara C chalk was generated in-situ.
Figure 2 (Sonnenberg, 2015) shows a generalized stratigraphic column and a structure map for the Niobrara in the DJ Basin along with an outline of the DJ basin and the location of the Wattenberg Field within which the study area is contained.
Figure 2: Outline of the DJ Basin with Niobrara structure contours and generalized stratigraphic column that shows the source rock and reservoir intervals for late Cretaceous units in the basin (from Sonnenberg, 2015).
Figure 3 shows the structural setting of the Niobrara in the study area, as well as types of fractures which can be expected to provide storage capacity and permeability for reservoirs within the chalk benches (Friedman and others, 1992). The study area covers approximately 100 square miles and shows large antiforms on the western edge. The area is normally faulted with most faults trending northeast to southwest. The Top Niobrara time structure also shows extensive small-scale structural relief which is visualized in a curvature attribute volume as shown in Figure 4. This implies that a significant amount of fracturing is present within the Niobrara.
Figure 3: Gross structure of the Niobrara in the study area in seismic two-way travel time. Insets from Friedman and others, 1992, showing predicted fracture types from structural elements. Area shown is approximately 100 square miles.
Figure 4: Most positive curvature, K1 on top Niobrara. The faulting and fractures are complex with both NE-SW and NW-SE trends apparent. Area shown is approximately 100 square miles. Seismic data provided courtesy of GPI and FFG.
Meissner and others (1984) and Landon and others (2001) have stated that the Niobrara formation kerogen is Type-II and oil-prone. Landon and others, and Finn and Johnson (2005) have also stated that the DJ basin contains the richest Niobrara source rocks with TOC contents reaching eight weight percent. Niobrara petroleum production is dependent on fractures in the hard, brittle, carbonate-rich zones. These zones are overlain and/or interbedded with soft, ductile marine shales that inhibit migration and seal the hydrocarbons in the fractured zones.
Why Utilize Machine Learning?
In the study area, the Niobrara to Greenhorn section is represented in approximately 60 milliseconds of two-way travel time in the seismic data. Figure 5 shows an amplitude section through a well within the study area. Figure 6 is an index map of wells used in the study with the Anderson 11-2 well highlighted in red. It is apparent that the top Niobrara is a well resolved positive amplitude or peak which can be picked on either a normal amplitude section or an instantaneous phase display. The individual units within the Niobrara A bench, A marl, B bench, B marl, C bench, C marl, Fort Hays and Codell present a significant challenge for an interpreter to resolve using only one or two attributes. The use of simultaneous multiple seismic attributes holds promise to resolve thin beds and a machine learning approach is one methodology which has been documented to successfully resolve stratigraphy below tuning (Roden and others, 2015, Santogrossi, 2017).
Figure 5: Amplitude section shows the approximately 60 milliseconds between marked horizons which contain the Niobrara and Codell reservoirs. Trace spacing is 110 feet, vertical scale is two-way time in seconds. Seismic data are shown courtesy of GPI and FFG.
Figure 6: Index map of vertical wells used in study. The dashed lines connect well names to well locations. Wells were obtained from the Colorado Oil and Gas Conservation Commission public database.
Machine Learning Data Preparation
The Niobrara Phase 5 3D data used for this study consisted of a 32-bit seismic amplitude volume that covers approximately 100 square miles. The survey contained 5.118 seconds of data with a bin spacing of 110 feet. Machine learning classifications benefit from sharper natural clusters of information through one level of finer trace sampling. Machine learned seismic resolution also benefits from sample-by-sample classification when compared to conventional wavelet analysis. Therefore, the data were upsampled to 1 ms from its original 2 ms interval by Geophysical Insights. The 1 ms amplitude data were used for seismic attribute generation.
Focus should be placed on the time interval that encompasses the geologic units of interest. The time interval selected for this study was 0.5 seconds to 2.2 seconds.
A total of 44 digital wells were obtained, 40 of which were within the seismic survey.
Classification by Principal Component Analysis (PCA)
Multi-dimensional analysis and multi-attribute analysis go hand in hand. Because individuals are grounded in three-dimensional space, it is difficult to visualize what data in a higher number dimensional space looks like. Fortunately, mathematics doesn’t have this limitation and the results can be easily understood with conventional 2D and 3D viewers.
Working with multiple instantaneous or geometric seismic attributes generates tremendous volumes of data. These volumes contain huge numbers of data points which may be highly continuous, greatly redundant, and/or noisy. (Coleou et al., 2003). Principal Component Analysis (PCA) is a linear technique for data reduction which maintains the variation associated with the larger data sets (Guo and others, 2009; Haykin, 2009; Roden and others, 2015). PCA has the ability to separate attribute types by frequency, distribution, and even character. PCA technology is used to determine which attributes to use and which may be ignored due to their very low impact on neural network solutions.
Figure 7 illustrates the analysis of a data cluster in two directions offset by 90 degrees. The first principal component (eigenvector 1) analyses the data cluster along the longest axis. The second principal component (eigenvector 2) analyses the data cluster variations perpendicular to the first principal component. As stated in the diagram, each eigenvector is associated with an eigenvalue which shows how much variance is in the data.
Figure 7: 2 attribute data set demonstrating the concept of PCA
Eigenvectors and eigenvalues from inline 1683 were consistently used for Principal Component Analysis because line 1683 bisected the deepest well in the study area. The entire pre-Niobrara, Niobrara, Codell, and post-Niobrara depositional events were present in the borehole.
PCA results for the first two eigenvectors for the interval Top Niobrara to Top Greenhorn are shown in Figure 8. Results show the most significant attributes in the first eigenvector are Sweetness, Envelope, and Relative Acoustic Impedance; each contributes approximately 60% of the maximum value for the eigenvector. PCA results for the second eigenvector show Thin Bed and Instantaneous Frequency are the most significant attributes. Figure 9 shows instantaneous attributes from the first eigenvector (sweetness) and second eigenvector (thin bed indicator) extracted near the B chalk of the Niobrara. The table shown in Figure 9 lists the instantaneous attributes that PCA indicated contain the most significance in the survey and the eigenvector associated with the attribute. This selection of attributes comprises a ‘recipe’ for input to the Self-Organizing Maps for the interval Niobrara to Greenhorn.
Figure 8: Eigenvalue charts for Eigenvectors 1 and 2 from PCA for Top Niobrara to Top Greenhorn. Attributes that contribute more than 50% of the maximum were selected for input to SOM
Figure 9: Instantaneous attributes near the Niobrara B chalk. These are prominent attributes in Eigenvectors 1 and 2. On the right of the figure is a list of eight selected attributes for SOM analysis. Seismic data is shown courtesy of GPI and FFG.
Teuvo Kohonen, a Finnish mathematician, invented the concepts of Self-organizing Maps (SOM) in 1982 (Kohonen, T., 2001). Self-Organizing Maps employ the use of unsupervised neural networks to reduce very high dimensions of data to a scale that can be easily visualized (Roden and others, 2015). Another important aspect of SOMs is that every seismic sample is used as input to classification as opposed to wavelet-based classification.
Figures 10 and 11 illustrate classification by SOM. Within the 3D seismic survey, samples are first organized into attribute points with similar properties called natural clusters in attribute space. Within each cluster new, empty, multi-attribute samples, named neurons, are introduced. The SOM neurons will seek out natural clusters of like characteristics in the seismic data and produce a 2D mesh that can be illustrated with a two- dimensional color map.
Figure 10: Example SOM classification of two attributes into 4 clusters (neurons)
In other words, the neurons “learn” the characteristics of a data cluster through an iterative process (epochs) of cooperative then competitive training. When the learning is completed each unique cluster is assigned to a neuron number and each seismic sample is now classified (Smith, 2016).
Figure 11: Illustration of how SOM works with 3D seismic volumes
Note that the two-dimensional color map in Figure 11 shows an 8X8 topology. Topology is important. The finer the topology of the two-dimensional color map the finer the data clusters associated with each neuron become. For example: an 8X8 topology distributes 64 neurons throughout an attribute set, while a 12X12 topology distributes 144 neurons. Finer topologies help to refine variations in lithologies, porosity, and other reservoir characteristics. Although there is no theoretical limit to a two-dimensional map topology, experience has shown that there is a practical limit to the number of neuron topologies for geological resolution. Conversely, a coarser neuron topology is associated with much larger data clusters and helps to define structural features. For the Niobrara project an 8X8 topology appeared to give the best stratigraphic resolution for instantaneous attributes and a 5X5 topology resolved the geometric attributes most effectively.
SOM Results for the Survey and their Interpretation
The SOM topology selected to best resolve the sub-Niobrara stratigraphy from the eight instantaneous attributes is an 8X8 hexagonal which yields 64 individual neurons. The SOM interval selected was Top Niobrara to Top Greenhorn. The next sequence of figures highlights the improved resolution provided by the SOM when compared to the original amplitude data. Figure 12 shows a north-south inline through the survey and through the Rotharmel 11-33 well which was one of the wells selected for petrophysical analysis. The original amplitude data is shown along with the SOM result for the interval.
Figure 12: North-South inline showing the original amplitude data (upper) and the 8X8 SOM result (lower) from Top Niobrara through Greenhorn horizons. Seismic data is shown courtesy of GPI and FFG.
The next image, Figure 13, zooms into the SOM and highlights the correlation with lithology from petrophysical analysis. The B chalk is noted by a stacked pattern of yellow-red-yellow neurons, with the red representing the maximum carbonate content within the middle of the chalk bench.
Figure 13: 8X8 Instantaneous SOM through Rotharmel 11-33 with well log composite. The B bench, highlighted in green on the wellbore, ties the yellow-red-yellow sequence of neurons. Seismic data is shown courtesy of GPI and FFG.
One can see on the SOM the sweet spot within the B chalk and that there is a fair amount of small-scale structural relief present. These results aid in the resolution of structural offset within the reservoir away from well control which is critical for staying in a 20 to 30 foot zone when drilling horizontally. Each classified sample is 1 ms in thickness which converted to depth equates to roughly 7 feet.
Figure 14 shows the K2 curvature attribute co-rendered with the SOM results in vertical sections. The Rotharmel 11-33 is at the intersection of the vertical sections. The curvature is extracted at the middle of the B chalk and shows good agreement with the SOM. The entire B bench is represented by only 5-6 ms of seismic data.
Figure 14: Most negative curvature, K2 rendered at the middle of the B chalk. Vertical sections are an 8X8 instantaneous SOM Top Niobrara to Top Greenhorn. Seismic data is shown courtesy of GPI and FFG.
A Marl Results
Seven wells within the survey were sent to a third party for petrophysical analysis (Figure 15). The analysis identified zones of interest within the Niobrara marls which are typically considered source rocks. The calculations show a high TOC zone in the upper A marl which the analysis identifies as shale pay (Figure 16). A seismic cross-section of the 8X8 instantaneous SOM (Figure 16) through the three wells depicted shows that this zone is well imaged. The neurons can be isolated and volumetric calculations derived from the representative neurons.
Figure 15: Index map for wells used in petrophysical analysis (in red)
Figure 16: Petrophysical results and SOM for three wells in the study area. The TOC curve (Track 12) and Shale pay curve (Track 10), highlighted in yellow, indicate the Upper A marl is both a rich source rock and a potential shale reservoir. Seismic data is shown courtesy of GPI and FFG.
The Codell sandstone in general and within the study area shows more heterogeneity in reservoir properties than the Niobrara chalk benches. The petrophysical analysis on 7 wells shows net pay ranging from zero feet to three feet. The gross thickness ranges from 17 feet to 25 feet. The SOM results reflect this heterogeneity, resolve the Codell gross interval throughout most of the study area, and thus, can be useful for horizontal well planning.
Figures 17 and 18 shows inline 60 through a well with the Top Niobrara to Greenhorn 8X8 SOM results. The 2D color map has been manipulated to emphasize the lower interval from approximately base Niobrara through the Codell. Figure 18 zooms into the well and shows the specific neurons associated with the Codell interval. Figures 19 shows a N-S traverse through four wells again with the Codell interval highlighted through use of a 2D color map. The western and southwest areas of the survey show a much more continuous character to the classification with only two neurons representing the Codell interval (6 and 48). Figure 20 shows both the N-S traverse and a crossline through the anomaly.
Figure 17: Instantaneous 8X8 SOM, Top Niobrara to Greenhorn. Seismic data is shown courtesy of GPI and FFG.
Figure 18: Detailed look at the Codell portion of the SOM at the Haythorn 4-12 with GR in background. The 2D color map shows how neurons can be isolated to show a specific stratigraphic interval. Seismic data is shown courtesy of GPI and FFG.
Figure 19: Traverse through 4 wells in the western part of the study area showing the isolation of the Codell sandstone within the SOM. The south west part of the line shows the Codell being represented by only 2 neurons (6 and 48). The colormap can be interrogated to determine which attributes contribute to any given neuron. Seismic data is shown courtesy of GPI and FFG.
Figure 20: View of the SW Codell anomaly where the neuron stacking pattern changes to two neurons only (6 and 47). Seismic data is shown courtesy of GPI and FFG.
Figure 21: 3D view of neurons isolated from the SOM in the Codell interval. The areas where red is prominent and continuous show the extent of Codell represented by neurons 6 and 47 only. Also, an area in the eastern part of the study is outlined. The Codell is not represented in this area by the six neurons highlighted in the 2D color map. Seismic data is shown courtesy of GPI and FFG.
Unfortunately, vertical well control was not available through this southwestern anomaly. To examine the extent of individual neurons within the SOM at Codell level, the next image, Figure 21, shows a 3D view of the isolated Codell neurons. The southwest anomaly is apparent as well as similar anomalies in the northern portion of the survey. What is also immediately apparent is that in the east-central portion of the survey, the Codell is not represented by the six neurons (6,7,47, 48, 55, 56) previously used to isolate it within the volume. Figure 22 takes a closer look at the SOM results through this area and also utilizes the original amplitude data. Both the SOM and the amplitude data show a change in character throughout the entire section, but the SOM results only change significantly in the lower Niobrara to Greenhorn portion of the interval.
The machine learning application has a feature in which individual neurons can be queried for statistics on how individual seismic attributes contribute to the cluster which makes up the neuron. Queries were done on all of the neurons within the Codell and shown are the results for neuron 6 which is one of 2 neurons characteristic of the southwestern Codell anomaly and on neuron 61in the area where the SOM changes significantly in Figure 23. Neuron 6 has equal contributions from Instantaneous Frequency, Hilbert, Thin Bed, and Relative Acoustic Impedance. Neuron 61 shows Instantaneous Q as the top attribute which is consistent with the interpretation of the section being structurally disturbed or highly fractured.
Figure 22: West-East crossline through two wells showing the SOM and amplitude data through the blank area from Figure 23. The seismic character and classification results differ significantly in this portion of the survey for the lower Niobrara, Fort Hays and Codell. This area is interpreted to be highly fractured. Seismic data is shown courtesy of GPI and FFG.
Figure 23: Example of attribute details for individual neurons (6 and 61). This shows the contribution of individual attributes to the neuron.
The machine learning workflow can be applied to geometric attributes. PCA and SOM need to be run separately from the instantaneous attributes since PCA assumes a Gaussian distribution of the attributes. This assumption doesn’t hold for geometric attributes but the SOM process does not assume any distribution and thus still finds patterns in the data. To produce a structural SOM, four attributes were selected from PCA: Curvature_K1, Similarity, Energy Ratio, Texture Entropy, and Texture Homogeneity. These were combined with the original amplitude data to generate SOMs from the Top Niobrara to Top Greenhorn interval. Several SOM topologies were generated with geometric attributes and a 5X5 yielded good results. Figure 24 shows the geometrical SOM results at the Top Niobrara, B bench, and Codell. The Top Niobrara level shows major faults, but not nearly as much structural disturbance as the mid-Niobrara B bench or the Codell level. The eastern part of the survey where the instantaneous classification changed also shows significant differences between the B bench and Codell and agrees with the interpretation that this is a highly fractured area for the lower Niobrara and Codell. The B bench appears more structurally disrupted than the Top Niobrara but shows fewer areal changes compared to Codell. Pressure and production data could help confirm how these features relate to reservoir quality.
Figure 24: 5X5 Structural SOM at 3 levels. There are significant changes both vertically and areally
Seismic multi-attribute analysis has always held the promise of improving interpretations via the integration of attributes which respond to subsurface conditions such as stratigraphy, lithology, faulting, fracturing, fluids, pressure, etc. Machine learning augments traditional interpretation and attribute analysis by utilizing attribute space to simultaneously classify suites of attributes into sample based, high dimension clusters that are subsequently visualized and further interpreted in the 3D seismic survey. 2D colormaps aid in their interpretation and visualization.
In the DJ Basin, we have resolved the primary reservoir targets, the Niobrara chalk benches and the Codell formation, represented within approximately 60 ms of data in two-way time, to the level of one to five neurons which is approximately 7 to 35 feet in thickness. Structural SOM classifications with a suite of geometric attributes better image the complex faulting and fracturing and its variations throughout the reservoir interval. The classification volumes are designed to aid in drilling target identification, reserves calculations and horizontal well planning.
The authors would like to thank their colleagues at Geophysical Insights for their valuable insight and suggestions and Digital Formation for the petrophysical analysis. We also thank Geophysical Pursuit, Inc. and Fairfield Geotechnologies for use of their data and permission to publish this paper.
Coleou, T., M. Poupon, and A. Kostia, 2003, Unsupervised seismic facies classification: A review and comparison of techniques and implementation: The Leading Edge, 22, 942–953, doi: 10.1190/1.1623635.
Finn, T. M. and Johnson, R. C., 2005, Niobrara Total Petroleum System in the Southwestern Wyoming Province, Chapter 6 of Petroleum Systems and Geologic Assessment of Oil and Gas in the Southwestern Wyoming Province, Wyoming, Colorado, and Utah, USGS Southwestern Wyoming Province Assessment Team, U.S. Geological Survey Digital Data Series DDS–69–D.
Guo, H., K. J. Marfurt, and J. Liu, 2009, Principal component spectral analysis: Geophysics, 74, no. 4, p. 35–43.
Kauffman, E.G., 1977, Geological and biological overview— Western Interior Cretaceous Basin, in Kauffman, E.G., ed., Cretaceous facies, faunas, and paleoenvironments across the Western Interior Basin: The Mountain Geologist, v. 14, nos. 3 and 4, p. 75–99.
Kohonen, T., 2001, Self-organizing maps: Third extended addition, Springer, Series in Information Services, Vol. 30.
Landon, S.M., Longman, M.W., and Luneau, B.A., 2001, Hydrocarbon source rock potential of the Upper Cretaceous Niobrara Formation, Western Interior Seaway of the Rocky Mountain region: The Mountain Geologist, v. 38, no. 1, p. 1–18.
Lewis, R.K., 2013, Stratigraphy and Depositional Environments of the Late Cretaceous (Late Turonian) Codell Sandstone and Juana Lopez Member of the Carlile Shale, Southeast Colorado: Colorado School of Mines MS Thesis, 190 p.
Longman, M.W., Luneau, B.A., and Landon, S.M., 1998, Nature and distribution of Niobrara lithologies in the Cretaceous Western Interior Seaway of the Rocky Mountain Region: The Mountain Geologist, v. 35, no. 4, p. 137–170.
Luneau, B., Longman, M., Kaufman, P., and Landon, S., 2011, Stratigraphy and Petrophysical Characteristics of the Niobrara Formation in the Denver Basin, Colorado and Wyoming, AAPG Search and Discovery Article #50469.
Meissner, F.F., Woodward, J., and Clayton, J.L., 1984, Stratigraphic relationships and distribution of source rocks in the greater Rocky Mountain region, in Woodward, J., Meissner, F.F., and Clayton, J.L., eds., Hydrocarbon source rocks of the greater Rocky Mountain region: Rocky Mountain Association of Geologists Guidebook, p. 1–34.
Molenaar, C.M., and Rice, D.D., 1988, Cretaceous rocks of the Western Interior Basin, in Sloss, L.L., ed., Sedimentary cover-North American craton, U.S.: Geological Society of America, The Geology of North America, v. D–2, p. 77–82.
Roden, R., Smith, T., and Sacrey, D., 2015, Geologic pattern recognition from seismic attributes: Principal component analysis and self-organizing maps, Interpretation, Vol. 3, No. 4, p. SAE59-SAE83.
Santogrossi, P., 2017, Classification/Corroboration of Facies Architecture in the Eagle Ford Group: A Case Study in Thin Bed Resolution, URTeC 2696775, doi 10.15530-urtec-2017-<2696775>.
Smith, T., 2016, Why SOM is an Appealing Learning Machine, Internal Geophysical Insights Paper.
Sonnenberg, S.A., 2015. Geologic Factors Controlling Production in the Codell Sandstone, Wattenberg Field, Colorado. URTeC Paper 2145312 presented at the Unconventional Resources Technology Conference, San Antonio, TX, July 20-22.
Sonnenberg, S.A., 2015. New reserves in an old field, the Niobrara/Codell resource plays in the Wattenberg Field, Denver Basin, Colorado. EAGE First Break, v. 33, p. 55-62.
Sterling, R., Bottjer, R. and Smith, K., 2016, Codell SS, A review of the Northern DJ oil resource play Laramie County, WY and Weld, County, CO, AAPG Search and Discovery Article #10754.
If you’re new to Machine Learning, let’s start at the top. The whole field of artificial intelligence is broken up into two categories – Strong AI and Narrow AI.
Strong AI is coming up with a robot that looks and behaves like a person. Narrow AI, or “neural networks” attempt to duplicate the brain’s neurological processes that have been perfected over millions of years of biological development.
Machine Learning is a subset of Narrow AI that does pattern classification. It’s an engine – an algorithm that learns without explicit programming. It learns from the data. What does that mean? Given one set of data, it’s going to come up with an answer. But given a different set of data, it will come up with something different.
A Self-Organizing Map is a type of neural network that adjusts to training data. However, it makes no assumptions about the characteristics of the data. So, if you look at the whole field of artificial intelligence, and then we look at machine learning as a subset of that, there’s two parts: supervised neural networks and unsupervised neural networks. Unsupervised is where you feed it the data and say “you go figure it out.” In supervised neural networks, you give it both the data and the right answer. Some examples of supervised neural networks would be convolutional neural networks and deep learning algorithms. Convolutional is a more classical type of a supervised neural network, where for every data sample, we know the answer.
Here’s a classical example of a supervised neural network: Your uncle just passed away and gave you the canning operations in Cordova, Alaska. You go there and observe the employees taking the fish off the conveyor and manually sorting them by type – buckets for eels and buckets for flounder and so forth. Can you use AI (machine learning) to do something more efficient? Perhaps have those employees do something more productive? Absolutely! As the eels come along, you weigh them, you take a picture of them, you see what the scales are, general texture, you get some idea about the general shape of them. There’s three properties already. You continue running eels through and maybe get up to four or five properties, including measurements, etc. The neural network is then trained on eels. Then, you do the same thing with all the flounder. There are going to be variations, of course, but in attribute space, of those four or five properties that we made for each one, they’re going to wind up in a different cluster in attribute space. And that’s how we tell the difference between eels and flounder. Everything else that you can’t classify very well, you don’t know. All of that goes into the algorithm. That’s the difference between supervised neural networks and unsupervised neural networks.
At Geophysical Insights, we believe we should be able to query our seismic data for information with learning machines just as effortlessly and with as much reliability as we query the web for the nearest gas station.
The application of machine learning to classify seismic attributes at single sample resolution is producing results that reveal more reservoir characterization information than is available from traditional interpretation methods. Two consequences of applying machine learning with several appropriately chosen seismic attributes include the discrimination of thin beds that are below conventional seismic tuning and the identification of Direct Hydrocarbon Indicators (DHIs). These capabilities enable a higher resolution interpretation of reservoirs and stratigraphy. An explanation of the machine learning methodology and its application to thin beds and DHIs is described briefly in this paper.
Machine Learning Methodology
Taking advantage of today’s computing technology, visualization techniques, and an understanding of machine learning on seismic data, Self-Organizing Maps (SOMs) (Kohonen, 2001), efficiently distills multiple seismic attributes into classification and probability volumes (Smith and Taner, 2010). When applied on a multi-attribute seismic sample basis, SOM is a powerful nonlinear cluster analysis and pattern recognition machine learning approach that helps interpreters identify patterns in their data that can relate to inherent geologic characteristics and different aspects of their data. SOM analysis, which is an unsupervised neural network application, when properly applied has been able to reveal both thin beds and DHIs in appropriate geologic settings. Figure 1 illustrates a single seismic amplitude trace and seven different seismic attributes computed from the amplitude data. All of these traces are displayed in a wiggle-trace variable area format. This display represents 100 ms vertically and each horizontal scale line represents a sample (4 ms). Each of these attributes are at different scales and in some cases vastly different scales. It is evident from this Figure that each of the attributes measures a different component of the total acoustic energy at every sample. SOM identifies clusters where different combinations of attributes congregate to reveal significant information about the natural groupings that are difficult to view any other way. The self-organizing property of SOM identifies and classifies natural clusters.
The SOM machine learning process is graphically presented in Figure 2. How large an area to select is dependent on the size of the geologic feature to be interpreted. For thin beds and DHIs, usually, a relatively thin zone of 50-250 ms around the anomalies is selected over a reasonable areal extent to provide sufficient data points for the SOM analysis. The selection of the seismic attributes is usually based on principal component analysis (PCA) and an interpreter’s knowledge of appropriate attributes for the area. Experience with SOM analysis has indicated that six to ten instantaneous seismic attributes are usually selected for thin beds and DHIs, depending on the geologic setting and data quality. In Figure 2 ten attributes are employed and all the data points from every sample from these attributes in the zone to be analyzed are placed in attribute space where they are normalized to put on the same scale. The SOM process employs cooperative and competitive learning techniques to identify the natural patterns or clusters in the data. Each pattern is identified by a neuron that sorts through the data in attribute space during the SOM training process of self-organization. In Figure 2 after training is completed, 64 winning neurons have identified 64 patterns or clusters in attribute space with an 8X8 neuron network. The SOM results are nonlinearly mapped back to a neuron topology map (2D colormap) where interpreters can select the winning neurons from the 2D colormap and identify in the 3D volume where the patterns and clusters occur for thin beds and DHIs.
Figure 1. Wiggle-trace variable area display format of a 100 ms window of seismic data with the amplitude trace and seven associated traces of attributes. Each attribute trace is at a different scale and each horizontal scale line is separated by the sample interval of 4 ms. If all these traces were employed in a SOM analysis, each red circle along a timing line indicates samples that would be input as a multi-attribute sample.
Figure 2. Display of SOM workflow where selected volume and data points from ten associated seismic attributes are input into Attribute Space. These data points are scaled and analyzed by the SOM process to identify 64 patterns by associated winning neurons. These neurons are nonlinearly mapped back to a 2D colormap where interpreters identify neurons and visually view the location of the patterns in the 3D survey.
In addition to the resultant classification volume, a probability volume is also generated which is a measure of the Euclidean distance from a data point to its associated winning neuron in attribute space (Roden et al., 2015). The winning neuron identifies a specific cluster or pattern. It has been discovered that a low classification probability corresponds to areas that are quite anomalous as opposed to high probability zones that relate to regional and common events in the data. Low probability anomalies identified by the SOM process are quite often associated with DHI characteristics.
Discriminating Thin Beds
The conventionally accepted definition of the tuning thickness (vertical resolution) is a bed that is ¼ wavelength in thickness, for which reflections from its upper and lower surfaces interfere and interference is constructive where the interface contrasts are of opposite polarity, often resulting in an exceptionally strong reflection (Sheriff, 2002). Several authors have described approaches to measure below tuning or thin beds usually employing various scaling techniques of amplitude or inversion data (MeckelandNath, 1977; Neidell and Poggiagliolmi, 1977; Schramm et al., 1977; Brown et al., 1986; and Connolly, 2007). However, these various techniques to determine thin beds have limitations and require assumptions that may not be met consistently (Simm, 2009). The application of SOM machine learning utilizing a multi-attribute classification has enabled the identification of thin beds and stratigraphy below tuning in a systematic and consistent manner as represented in the following case study.
The Eagle Ford Shale is a well-known unconventional resource play in Texas. Operators in this play must account for changing stratigraphy and facies to properly locate horizontal wells for optimum perforation intervals. The Eagle Ford stratigraphy is often associated with thin beds and facies well below conventional seismic resolution that change both vertically and laterally. This Eagle Ford case study contains 216 mi² (560 km²) of enhanced 3D PSTM data processed at a 2 ms sample interval. Conventional vertical resolution (tuning thickness) is 100-150 feet (30-45 meters) depending on the location within the Eagle Ford unit. In this study, over 300 wells were available for correlation including 23 type logs, 249 horizontal borehole montages, 9 vertical calibration wells with tops, logs, and time-depth corrections. Also available were five cores for which X-ray diffraction and saturation information was available. Well information was incorporated to corroborate the SOM results.
Figure 3. Resolution comparison between conventional seismic display and a Paradise® multi-attribute Self Organizing Map (SOM) classification: (a) Seismic amplitude profile through the 6 Well; and (b) the SOM results of the same profile identifying the Eagle Ford Group that comprises 26 sample based neuron clusters which are calibrated to facies and systems tracts. The 2D colormap displays the associated winning neuron cluster colors.
Ten instantaneous seismic attributes prominent in a Principal Component Analysis were selected for SOM. SOM training was conducted on a set of trial harvest lines and the successful line was then used to classify the entire survey. Figure 3a is a seismic amplitude line in color raster and wiggle-trace variable area formats across the location of Well 6 (V and 1H). The Figure shows the Austin Chalk-Eagle Ford Group-Buda stratigraphic interval represented by roughly 2.5 peak/trough cycles of a seismic trace. In the Eagle Ford, the amplitude appears continuous, yet any details are obscured because of the resolution limitations of the amplitude data where conventional tuning is 100-150 feet (30-45 meters). Figure 3b displays the equivalent line showing results of the SOM analysis. Sixty-four neurons were employed to identify 64 patterns in the data as seen on the associated 2D colormap. A seismic interval from 10 ms below the Buda to 100 ms above the Buda or near the top of the Austin Chalk was chosen for the SOM analysis. Shown clearly comparing Figures 3a and 3b is the resolution improvement provided by the SOM analysis over the seismic amplitude. The results reveal non-layer cake facies bands that include details in the Eagle Ford’s basal clay-rich shale, high resistivity and low resistivity Eagle Ford Shale objectives, the Eagle Ford ash, and the Upper Eagle Ford marl, which are overlain disconformably by the Austin Chalk (disconformity is a break in a sedimentary sequence that does not involve a difference in bedding angles). This interval of roughly 28 ms (or 14 samples) amounts to some 26 of the 64 SOM winning neurons to illuminate the various systems tracts within the Eagle Ford Group for this survey.
Adjacent to the SOM results at Well 6 are similar results at a nearby well. Figure 4a displays a zoomed-in vertical line display of the SOM results through Well 8(V) with the winning neurons identified. Figure 4b denotes the associated well log curves from well 8 and the correlative neuron associations. Winning neurons 63 and 64 are associated with the low resistivity Eagle Ford shale unit and neurons 53, 54, and 60 denote the high resistivity and more desirable Eagle Ford unit. The expanded display of the SOM results in Figure 4a denotes a low resistivity gold thin bed that is identified by a single neuron (#55) and is only one sample thick (2 ms). Shown here is clear evidence of consistent results between Wells 6 and 8 that lends itself to stratigraphic facies interpretation.
Figure 4. View of SOM results correlated with Well 8 logs: (a) Expanded view of SOM profile through Well 8; and (b) Well 8 logs with Eagle Ford units correlated with the SOM results of (a). Note that the resolution is down to the sample level of 10-12 feet and is illustrated by winning neuron 55 in (a). The dashed white line in (a) represents the top of the Eagle Ford and base of the Austin Chalk.
Over this survey area, 16 different winning neurons represent the various facies present in the Eagle Ford Shale over a 14 ms window (70-84 feet/21-26 meters). The facies of the entire Eagle Ford Group which includes the Basal Clay shale, Eagle Ford Shale, and Eagle Ford Marl, are defined by 26 different winning neurons over 28 ms (210-252 feet/64-77 meters). Individual facies units are as thin as one sample interval of 2 ms (10-12 feet/3-4 meters). These results of a single SOM classification are corroborated at multiple wells across the survey area.
Revealing Direct Hydrocarbon Indicators (DHIs)
The accurate interpretation of seismic DHI characteristics has proven to significantly improve the drilling success rates in the appropriate geologic setting where there is also adequate seismic data quality (Roden et al., 2012; Rudolph and Goulding, 2017). Specifically, DHIs are anomalies due to the presence of hydrocarbons induced by changes in rock physics properties (P and S wave velocities, and density). Typically, anomalies stand out as the difference in a hydrocarbon-filled reservoir in relation to the encasing rock or the brine portion of the reservoir. DHI characteristics are usually associated with anomalous seismic responses in a trapping configuration, such as structural traps, stratigraphic traps, or a combination of both. These include bright spots, flat spots, amplitude conformance to structure, etc. DHI anomalies are often compared to models, similar events, background trends, proven productive anomalies, and geologic features. DHI indicators can also be located below presumed trapped hydrocarbons where shadow zones or velocity pull-down effects may be present. These DHI effects can even be present and dispersed in the sediment column in the form of gas chimneys or clouds.
As described above there are numerous DHI characteristics and all of which should be evaluated in a consistent and systematic approach for any prospect or project. Forrest et al. (2010) and Roden et al. (2012) have identified the top four DHI characteristics, as related to commercially successful wells in a Class 3 AVO setting, based on an industry-wide database of almost 300 wells. These DHI characteristics include:
Anomaly conformance to structure
Phase or character change at the downdip edge of the anomaly
Anomaly consistency in the mapped target area
Figure 5. From the top of the producing reservoir: a) time structure map in contours with an amplitude overlay in color and b) SOM classification with low probability less than 1% denoted by white areas. The yellow line in b) represents the downdip edge of the high amplitude zone picked from a).
These DHI characteristics will be identified in the following case study by a multi-attribute SOM analysis. The case study is an offshore oil/gas field in 470 feet (145 meters) of water on the Louisiana continental shelf of the Gulf of Mexico. The field has two producing wells that were drilled on the upthrown side of a normal fault and into an amplitude anomaly. The normally pressured reservoir is approximately 100 feet (30 meters) thick and contains oil and gas. The hydrocarbon filled sandstone reservoir has low impedance compared to the encasing shales, indicative of a Class 3 AVO environment. The SOM analyzed a 170 ms window surrounding the reservoir. Applying the (SOM) multi-attribute analysis, a group of eight seismic attributes was selected based on Principal Component Analysis that would best expose Direct Hydrocarbon Indicators (DHIs). A neuron network of 5X5 (25 neurons) was employed. Figure 5a displays a time structure map as denoted by the contours with an amplitude overlay (color) from the mapped top of the reservoir in this field. The horizon at the top of the reservoir was picked on a trough (low impedance) on zero phase seismic data (SEG normal polarity). Figure 5a indicates there is a relatively good amplitude conformance to structure based on the amplitude as noted by the general agreement of the time contour and the red/green amplitude break (see color bar insert). Figure 5b is a display of classification probability from the SOM analysis at the top of the reservoir at the same scale as Figure 5a. This indicates that the top of the reservoir exhibits an anomalous response from the SOM analysis. Ordinary classifications such as those of green and yellow-greens are shown around the reservoir (see colormap insert). However, within the reservoir and in several other areas, anomalous classification of low probability below 1% are colored white. In comparing Figure 5a and 5b it is apparent that the low probability area corresponds closely to the amplitude conformance to structure as denoted by the yellow outline in Figure 5b. This confirms the identification of the productive area with low probability and proves the efficacy of this SOM approach. The consistency of the low probability SOM response in the field is another positive DHI indicator. In fact, the probabilities as low as .01% still produce a consistent response over the field indicating how significant evaluating low probability anomalies is critical in the interpretation of DHI characteristics.
Figure 6. North-south vertical profile 9411 through the middle of the field: a) stacked seismic amplitude display with the field location marked, b) SOM classification with 25 neurons (5X5) indicated by the 2D colormap over a 170 ms window, and c) three neurons highlighting the reservoir above the oil/water and gas/oil contacts and the hydrocarbon contacts (flat spots) as marked by winning neurons 15, 20, and 25. The expanded insets in the right column denote the details from the SOM results at the downdip edge of the field (see dashed black boxes on left).
This field contains an oil phase with a gas cap and before drilling, there were hints of possible flat spots suggesting hydrocarbon contacts on the seismic data, but the evidence was inconsistent and not definitive. Figure 6 displays a north-south vertical inline seismic profile through the middle of the field with its location denoted in Figure 5. Figure 6a exhibits the initial stacked amplitude data with the location of the field marked. Figure 6b denotes the SOM analysis results of this same vertical inline which incorporates the eight instantaneous attributes listed with the 5X5 neuron matrix in Figure 5. The associated 2D colormap in Figure 6b denotes the 25 natural patterns or clusters identified from the SOM process (see colormap insert). It is apparent in Figure 6b that the reservoir and portions of both the gas/oil contact and the oil/water contact are easily identified as opposed to the conventional seismic data shown in Figure 6a. This is more easily seen in Figure 6c where the 2D colormap indicates that the neurons highlighted in gray (winning neurons 20 and 25) are defining the hydrocarbon-bearing portions of the reservoir above the hydrocarbon contacts and the flat spots interpreted as hydrocarbon contacts are designated by the rust-colored winning neuron 15. The location of the reservoir and hydrocarbon contacts are corroborated by well control (not shown). Note that the southern edge of the reservoir is revealed in the enlargements in the column on the right (from the dashed black box on left). A character change at the downdip edge of the anomaly where the oil contact thins out is easily identified compared to the line of amplitude change. Downdip of the field is another undrilled anomaly defined by the SOM analysis that exhibits similar DHI characteristics identified by the same neurons.
A seismic multi-attribute analysis utilizing Self-Organizing-Maps is a machine learning approach that distills information from numerous attributes on a sample-by-sample basis to provide a more accurate assessment of thin beds and DHI characteristics than conventional methods. We have shown in these cases that the results readily correlate to conventional well log data. The SOM classification process takes advantage of natural patterns in multiple seismic attributes space that is not restricted to the resolution limits of conventional amplitude data. This process enables interpreters to produce higher resolution interpretations of reservoirs and stratigraphy. Another advantage of a SOM analysis is the generation of classification probability where low probability anomalies are often associated with DHIs. SOM analyses in appropriate geologic settings can improve confidence in interpreting DHI characteristics and more clearly define reservoir edges and thin bed components.
Smith, T. and M. T. Taner, 2010, Natural clusters in multi-attribute seismics found with self-organizing maps: Extended Abstracts, Robinson-Treitel Spring Symposium by GSH/SEG, March 10-11, 2010, Houston, Tx.
Rocky R. Roden has extensive knowledge of modern geoscience technical approaches (past Chairman-The Leading Edge Editorial Board). As former Chief Geophysicist and Director of Applied Technology for Repsol-YPF, his role comprised advising corporate officers, geoscientists, and managers on interpretation, strategy and technical analysis for exploration and development in offices in the U.S., Argentina, Spain, Egypt, Bolivia, Ecuador, Peru, Brazil, Venezuela, Malaysia, and Indonesia. He has been involved in the technical and economic evaluation of Gulf of Mexico lease sales, farmouts worldwide, and bid rounds in South America, Europe, and the Far East. Previous work experience includes exploration and development at Maxus Energy, Pogo Producing, Decca Survey, and Texaco. He holds a B.S. in Oceanographic Technology-Geology from Lamar University and a M.S. in Geological and Geophysical Oceanography from Texas A&M University.
Patricia Santogrossi is a Consultant to Geophysical Insights, producer of the Paradise multi-attribute analysis software platform. Formerly, she was a Leading Reservoir Geoscientist and Non-operated Projects Manager with Statoil USA E & P. In this role Ms. Santogrossi was engaged for nearly nine years in Gulf of Mexico business development, corporate integration, prospect maturation, and multiple appraisal projects in the deep and ultra-deepwater Gulf of Mexico. Ms. Santogrossi has previously worked with domestic and international Shell Companies, Marathon Oil Company, and Arco/Vastar Resources in research, exploration, leasehold and field appraisal as well as staff development. She has also been Chief Geologist for Chroma Energy, who possessed proprietary 3D voxel multi-attribute visualization technology, and for Knowledge Reservoir, a reservoir characterization and simulation firm that specialized in Deepwater project evaluations. She is a longtime member of SEPM, AAPG, GCSSEPM, HGS and SEG and has held various elected and appointed positions in many industry organizations. Ms. Santogrossi holds an M.S. in Geology from the University of Illinois, Champaign-Urbana.
In conventional geological settings, oil companies routinely evaluate prospects for their drilling portfolio where the process of interpreting seismic amplitude anomalies as Direct Hydrocarbon Indicators (DHIs) plays an important role. DHIs are an acoustic response owing to the presence of hydrocarbons and can have a significant impact on prospect risking and determining well locations (Roden et al., 2005; Fahmy 2006; Forrest et al., 2010; Roden et al., 2012; Rudolph and Goulding, 2017). DHI anomalies are caused by changes in rock physics properties (P and S wave velocities, and density) typically of the hydrocarbon-filled reservoir in relation to the encasing rock or the brine portion of the reservoir. Examples of DHIs include bright spots, flat spots, character/phase change at a projected oil or gas/water contact, amplitude conformance to structure, and an appropriate amplitude variation with offset on gathers. Many uncertainties should be considered and analyzed in the process of assigning a probability of success and resource estimate range before including a seismic amplitude anomaly prospect in an oil company’s prospect portfolio (Roden et al., 2012).
Seismic amplitude anomalies that are DHIs have played a major role in oil and gas exploration since the early 1970s (Hilterman, 2001). The technology and methods to identify and risk seismic amplitude anomalies have advanced considerably over the years, especially with the use of AVO (Amplitude vs. Offset) and improved acquisition and processing seismic technology (Roden et al., 2014). The proper evaluation of seismic direct hydrocarbon indicators for appropriate geologic settings has proven to have a significant impact on risking prospects. Rudolph and Goulding (2017) indicate from an ExxonMobil database of prospects that DHI-based prospects had over twice the success rate of non-DHI prospects on both a geologic and economic basis. In an industry-wide database of DHI prospects from around the world, Roden et al. (2012) indicate that when a prospect has a >20% DHI Index, a measure of the risk associated with DHI characteristics, almost all the wells were successful. Even with the use of advanced seismic technology and well-equipped interpretation workstations, the interpretation of DHI characteristics is not always easy or straightforward.
A key technology employed in evaluating potential DHI features is seismic attributes. Seismic attributes are any measurable property of seismic data including stacked or prestack data. Seismic attributes can be computed on a trace, multiple traces, on an entire volume, over isolated windows, on a horizon, and in either time or depth. There are hundreds of seismic attributes generated in our industry (Brown, 2004; Chen and Sidney, 1997; Chopra and Marfurt, 2007; Taner, 2003) and can be generally categorized as instantaneous, geometric, AVO, seismic inversion, and spectral decomposition attributes. The instantaneous, AVO, and inversion attributes are typically utilized to highlight and identify DHI features. For example, amplitude envelope, average energy, and sweetness are good instantaneous attributes to display how amplitudes stand out above the background, potentially identifying a bright spot and a potential hydrocarbon accumulation. AVO attributes such as intercept times gradient, fluid factor, Lambda/Mu/Rho and far offset-minus near offset-times the far offset can help to identify hydrocarbon-bearing reservoirs (Roden et al., 2014). However, not all amplitude anomalies are DHIs and interpreting numerous seismic attributes can be complicated and at times confusing. In addition, it is almost impossible for geoscientists to understand how numerous seismic attributes (>3) interrelate.
Over the last few years, machine learning has evolved to help interpreters handle numerous and large volumes of data (e.g. seismic attributes) and help to understand how these different types of data relate to each other. Machine learning uses computer algorithms that iteratively learn from the data and independently adapt to produce reliable, repeatable results. We incorporate a machine learning workflow where principal component analysis (PCA) and self-organizing maps (SOM) are employed to analyze combinations of seismic attributes for meaningful patterns that correspond to direct hydrocarbon indicators. A machine learning multi-attribute approach with the proper input parameters can help interpreters to more efficiently and accurately evaluate DHIs and help reduce risk in prospects and projects.
Table 1 Most important DHI characteristics for AVO classes 2 and 3 as determined by Forrest et al. (2010) and Roden et al. (2012)
DHI characteristics are usually associated with anomalous seismic responses in a trapping configuration: structural traps, stratigraphic traps, or a combination of both. These include bright spots, flat spots, amplitude conformance to structure, etc. DHI anomalies are also compared to other features such as models, similar events, background trends, proven productive anomalies, and geologic features. DHI indicators can also be located below presumed trapped hydrocarbons where shadow zones or velocity pull-down effects may be present. DHI effects can even be present dispersed in the sediment column in the form of gas chimneys or clouds. Forrest et al. (2010) and Roden et al. (2012) have documented the most important DHI characteristics based on well success rates in an industry-wide database of DHI prospects. Other than the amplitude strength above background (bright spots), Table 1 lists these DHI characteristics by AVO classes 2 and 3. These two AVO classes (Rutherford and Williams, 1989) relate to the amplitude with offset response from the top of gas sands which represent the specific geologic settings where most DHI characteristics are found. Therefore, the application of machine learning employing seismic multi-attribute analysis may help to clearly define DHI characteristics and assist the interpreter in making a more accurate assessment of prospect risk.
Table 2 Most important Class 3 DHI characteristics as denoted by Forrest et al. (2010) and Roden et al. (2012) and a designation of typical instantaneous attributes that identify these characteristics. Not all instantaneous attributes by themselves are conducive to identifying the top DHI characteristics.
Multi-attribute machine learning workflow
With the goal of identifying DHI characteristics, an interpreter must determine the specific attributes to employ in a machine learning workflow. A geoscientist can select appropriate attributes based on their previous knowledge and experience to define DHIs in a specific play or trend. Table 2 lists several common instantaneous attributes and the associated stacked seismic data DHI characteristics they tend to identify. These relationships are of course subjective and depend on the geologic setting and data quality. Class 3 DHIs are usually interpreted on full stack volumes and/or offset/angle volumes and their associated derivative products. Class 2 DHIs are typically interpreted on offset/angle volumes (especially far offset/angle volumes), gathers, and their associated derivative products including various types of crossplots. The relationships between attributes and DHI characteristics can be variable depending on the geologic setting and the seismic data quality. If it is unclear which attributes to employ, principal component analysis (PCA) can assist interpreters. PCA is a linear mathematical technique to reduce a large set of variables (seismic attributes) to a smaller set that still contains most of the variation of independent information in the larger set. In other words, PCA helps to determine the most meaningful seismic attributes. The first principal component accounts for as much of the variability in the data as possible and each succeeding component (orthogonal to each preceding) accounts for as much of the remaining variability. Given a set of seismic attributes generated from the same original volume, PCA can identify combinations of attributes producing the largest variability in the data suggesting these combinations of attributes that will better identify specific geologic features of interest and in this case specific DHI characteristics. Even though the first principal component represents the largest linear attribute combinations that best represents the variability of the bulk of the data, it may not identify specific features of interest to the interpreter. The interpreter should also evaluate succeeding principal components because they may be associated with DHI characteristics not identified with the first principal component. In fact, the top contributing seismic attributes from the first few principal components, when combined, often produce the best results for DHI delineation. In other words, PCA is a tool that, employed in an interpretation workflow with a geoscientist’s knowledge of DHI related attributes, can give direction to meaningful seismic attributes and improve interpretation results. It is logical, therefore, that a PCA evaluation may provide important information on appropriate seismic attributes to take into a self-organizing map generation.
After appropriate seismic attributes have been selected, the next level of interpretation requires pattern recognition and classification of often subtle information embedded in the seismic attributes. Taking advantage of today’s computing technology, visualization techniques, and understanding of appropriate parameters, self-organizing maps (SOMs) efficiently distill multiple seismic attributes into classification and probability volumes (Smith and Taner, 2010; Roden et al., 2015). Developed by Kohonen in 1982 (Kohonen, 2001), SOM is a powerful non-linear cluster analysis and pattern recognition approach that helps interpreters to identify patterns in their data that can relate to geologic features and DHI characteristics. The samples for each of the selected seismic attributes from the desired window in a 3D survey are placed in attribute space where they are normalized or standardized to the same scale. Also in attribute space are neurons, which are points in space that start at random locations and train from the attribute data and mathematically hunt for natural clusters of information in the seismic data. After the SOM analysis, each neuron will have identified a natural cluster as a pattern in the data. These clusters reveal significant information about the classification structure of natural groups that are difficult to view any other way. In addition to the resultant classification volume, a probability volume is also generated which is a measure of the Euclidean distance from a data point to its associated winning neuron (Roden et al., 2015). The winning neuron is the one that is nearest to the data point in attribute space. It has been discovered that a low classification probability corresponds to areas that are quite anomalous as opposed to high probability zones that relate to regional and common events in the data.
To interpret the SOM classification results, each neuron is displayed in a 2D color map. Highlighting a neuron or combination of neurons in a 2D color map identifies their associated natural clusters or patterns in the survey because each seismic attribute data point retains its physical location in the 3D survey. The identification of these patterns in the data enables interpreters to define geology not easily interpreted from conventional seismic amplitude displays alone. These visual cues are facilitated by an interactive workstation environment.
Low probability anomalies
After the SOM process and natural clusters have been identified, Roden et al. (2015) describe the calculation of a classification probability. This probability estimates the probable certainty that a winning neuron classification is successful. The classification probability ranges from zero to 100% and is based on goodness of fit of the Euclidean distances between the multi-attribute data points and their associated winning neuron. Those areas in the survey where the classification probability is low correspond to areas where no winning neurons fit the data very well. In other words, anomalous regions in the survey are noted by low probability. DHI characteristics are often associated with low classification probabilities because they are anomalous features that are usually not widespread throughout the survey.
SOM analysis for Class 3 DHI characteristics
A class 3 geologic setting is associated with low acoustic impedance reservoirs that are relatively unconsolidated. These reservoirs typically have porosities greater than 25% and velocities less than 2700 m/sec. The following DHI characteristics are identified by multi-attribute SOM analyses in an offshore Gulf of Mexico class 3 setting. This location is associated with a shallow oil and gas field (approximately 1200 m) in a water depth of 140 m that displayed a high seismic amplitude response. Two producing wells with approximately 30 m of pay each were drilled in this field on the upthrown side of an east-west trending normal fault. Before these wells were drilled, operators had drilled seven unsuccessful wells in the area based on prominent seismic amplitudes that were either wet or low saturation gas. Therefore, the goal was to identify as many DHI characteristics as possible to reduce risk and accurately define the field and to develop SOM analysis approaches that can help to identify other possible prospective targets in the area.
Initially, 20 instantaneous seismic attributes were run through a PCA analysis in a zone 20ms above and 150 ms below the top of the mapped producing reservoir. Based on these PCA results, various combinations of attributes were employed in different SOM analyses with neuron counts from 3X3, 5X5, 8X8, 10X10, and 12X12 employed for each set of seismic attributes. It is from this machine learning multi-attribute interpretation workflow that the results defining different DHI characteristics were interpreted and described below. All of the figures associated with this example are from a SOM analysis with a 5X5 neuron count and employed the instantaneous attributes listed below.
Relative Acoustic Impedance
Cosine of Instantaneous Phase
Final Raw Migration
Figure 1 From the top of the producing reservoir: a) time structure map in contours with an amplitude overlay in colour and b) SOM classification with low probability less than 1% denoted by white areas. The yellow line in b) represents with downdip edge of the high amplitude zone designated in a).
Figure 1a displays a time structure map as denoted by the contours with an amplitude overlay (color) from the mapped top of the reservoir in this field. The horizon at the top of the reservoir was picked on a trough (low impedance) on zero phase seismic data (SEG normal polarity). Figure 1a indicates that there is a relatively good amplitude conformance to structure based on the amplitude. Figure 1b is a display of classification probability from the SOM analysis at the top of the reservoir at the same scale as Figure 1a. This indicates that the top of this reservoir exhibits an anomalous response from the SOM analysis where any data points with a probability of less than 1% are displayed in the white areas. In comparing Figure 1a and 1b it is apparent that the low probability area corresponds closely to the amplitude conformance to structure as denoted by the yellow outline in Figure 1b. This confirms the identification of the productive area with low probability and proves the efficacy of this SOM approach. The consistency of the low probability SOM response in the field is another positive DHI indicator. In fact, the probabilities as low as .01% still produce a consistent response over the field indicating how the evaluation of low probability anomalies is critical in the interpretation of DHI characteristics.
This field contains oil with a gas cap and before drilling, there were hints of possible flat spots suggesting hydrocarbon contacts on the seismic data, but the evidence was inconsistent and not definitive. Figure 2 displays a north-south vertical inline profile through the middle of the field and its location is denoted in Figure 1. Figure 2a exhibits the initial stacked amplitude data with the location of the field annotated. Figure 2b denotes the SOM analysis results of this same vertical inline 9411 which incorporated the eight instantaneous attributes listed above in a 5X5 neuron matrix. The associated 2D color map in Figure 2b denotes the 25 natural patterns or clusters identified from the SOM process. It is apparent in this figure that the reservoir and portions of the gas/oil contact and the oil/water contact are easily identified. This is more easily seen in Figure 2c where the 2D color map indicates that the neurons highlighted in grey (20 and 25) are defining the hydrocarbon-bearing portions of the reservoir above the hydrocarbon contacts and the flat spots interpreted as hydrocarbon contacts are designated by the rust-colored neuron (15). The location of the reservoir and hydrocarbon contacts are corroborated by well control. The southern edge of the reservoir is revealed in the enlargements of the column displays on the right. Downdip of the field is another undrilled anomaly defined by the SOM analysis that exhibits similar DHI characteristics identified by the same neurons.
Figure 2 North-south vertical profile 9411 through the middle of the field: a) stacked seismic amplitude display with the field location designated, b) SOM classification with 25 neurons indicated by the 2D colour map over a 170 ms window, and c) three neurons highlighting the reservoir above the oil/water and gas/oil contacts and the hydrocarbon contacts (flat spots). The expanded insets denote the details from the SOM results at the downdip edge of the field.
Figure 3 West-east vertical profile 3183 through the field: a) stacked seismic amplitude display denoting tie with line 9411, b) SOM classification with 25 neurons indicated by the 2D colour map, and c) three neurons highlighting the gas/oil and oil/water contacts and the hydrocarbon contacts (flat spots). The expanded insets clearly display the edge of the field in the SOM classifications.
West to east crossline 3179 over the field is displayed in Figure 3 and with it the location designated in Figure 1. The stacked seismic amplitude display of Figure 3a indicates that its tie with inline 9411 is located in the updip portion of the reservoir where there is an apparent gas/oil contact. Figure 3b exhibits the SOM results of this west-east line utilizing 25 neurons as designated by the 2D color map. Similar to Figure 2b, Figure 3b indicates that the SOM analysis has clearly defined the reservoir by the grey neurons (20 and 25) and the hydrocarbon contacts in the rust-colored neuron (15). Towards the west, the rust-colored neuron (15) denotes the oil/water contact as defined by the flat spot on this crossline. Figure 3c displays only neurons 15, 20, and 25 to clearly define the reservoir, its relationship above the hydrocarbon contacts, and the contacts themselves. The three enlargements on the left are added for detail.
What is very evident from the SOM results in both Figures 2 and 3 is a clear character change and definition of the downdip edges of the reservoir. The downdip edge definition of an interpreted trap is an important DHI characteristic that is clearly defined by the SOM analysis in this field. The expanded insets in Figures 2 and 3 indicate that the SOM results are producing higher resolution results than the amplitude data alone and the edge terminations of the field are easily interpreted. These results substantiate that the SOM process with the appropriate set of seismic attributes can exhibit thin beds better than conventional amplitude data.
SOM analysis for Class 2 DHI characteristics
A class 2 geologic setting contains reservoirs more consolidated than class 3 and the acoustic impedance of the reservoirs are about equal to the encasing sediments. Typical porosities range from 15 to 25% and velocities 2700-3600 m/sec for these reservoirs. In class 2 settings, AVO attributes play a larger role in the evaluation of DHI characteristics than in class 3 (Roden et al., 2014). This example is located onshore Texas and targets Eocene sands at approximately 1830 m deep. The initial well B was drilled just downthrown on a small southwest-northeast regional fault, with a subsequent well drilled on the upthrown side (Well A). The reservoirs in the wells are approximately 5-m thick and composed of thinly laminated sands. The tops of these sands produce a class 2 AVO response with near zero amplitude on the near offsets and an increase in negative amplitude with offset (SEG normal polarity).
Figure 4 Time structure map at the top of the producing Eocene reservoir.
The goal of the multi-attribute analysis was to determine the full extent of the reservoirs revealed by any DHIs because both wells were performing much better than the size of the amplitude anomaly indicated from the stack and far offset seismic data. Figure 4 is a time-structure map from the top of the Eocene reservoir. This map indicates that both wells are located in stratigraphic traps with Well A situated on southeast dip and Well B located on the northwest dip that terminates into the regional fault. The defined anomaly conformance to downdip closure cannot be seen in the Well A reservoir because the areal extent of the reservoir is in a north-south channel and the downdip conformance location is very narrow. In the Well B reservoir, the downdip edge of the reservoir actually terminates into the fault so an interpretation of the downip conformance cannot be determined. The updip portion of the reservoir at Well B actually thins out towards the south-east forming an updip seal for the stratigraphic trap. The Well B reservoir was interpreted to have a stacked data amplitude anomaly of approximately 70 acres and the Well A reservoir was determined to only have an amplitude anomaly of only about 34 acres (Figure 5a).
Figure 5 At the top of the Eocene reservoir: a) stacked seismic amplitude, b) SOM classification with 64 neurons, and c) same classification as the middle display with low probability of less than 30% designated by the white areas.
Figure 6 North-south arbitrary line through Wells A and B with the location designated in Figure 4: a) stacked seismic amplitude display, b) SOM classification with 64 neurons indicated by the 2D colour map, and c) SOM classification with only four neurons in grey highlighting both the reservoirs associated with the wells.
The gathers associated with the 3D PSTM survey over this area were conditioned and employed in the generation of very specific AVO attributes conducive to the identification of class 2 AVO anomalies in this geologic setting. The ten AVO attributes used for the SOM analysis were selected from a PCA analysis of 18 AVO attributes. The AVO attributes that were selected for the SOM analysis are listed below:
Far – Near
Shuey 2 term approximation – Intercept
Shuey 2 term approximation – Gradient
Shuey 2 term approximation – 1/2 (Intercept + Gradient)
Shuey 2 term approximation – 1/2 (Intercept – Gradient)
Shuey 3 term approximation – Intercept
Shuey 3 term approximation – Gradient
Shuey 3 term approximation – 1/2 (Intercept + Gradient)
Several different neuron counts were generated with these ten AVO attributes and the results in the associated figures are from the 8X8 (64 neurons) count. Figure 5b displays the SOM results from the top of the Eocene reservoirs. The associated 2D color map indicates that neurons 47, 58, 62, and 63 are defining the reservoirs drilled by the two wells. Comparing the areal distribution of the amplitude defined reservoirs in 5a to the SOM defined reservoirs in Figure 5b indicates that the later is larger. In fact, the Well A amplitude defined area of 34 acres is compared to approximately 95 acres as denoted by four neurons in Figure 5b. The Well B amplitude defined reservoir area was determined to be 70 acres, whereas, the SOM defined area was determined to be approximately 200 acres. The SOM defined areal distributions were determined to be consistent with engineering and pressure data in the two wells. The anomaly consistency in the mapped target area is evident in Figure 5b and is better in defining the extent of the producing reservoirs than amplitudes.
Figure 5c displays the SOM results of 5b. However, less than 30% of the low classification probability results are displayed in white. It denotes that the core of the reservoirs at each of the well locations reveal low probability. Low probability is defining anomalous zones based on the ten AVO attributes run in the SOM classification process.
Figure 7 Northeast-southwest inline 2109 through Well B with location designated in Figure 4: a) stacked seismic amplitude display, b) SOM classification with 64 neurons as denoted by the 2D colour map, and c) SOM classification with only four grey neurons highlighting the reservoir at Well B. The expanded insets display the updip edges of the reservoir with the SOM results clearly defining the updip seal edge of the field.
Figure 6 is a north-south arbitrary line running through both Wells A and B with its location denoted in Figure 4. Figure 6a is the conventional stacked seismic amplitude display of this line. Figure 6b displays the SOM results and the reservoirs at both wells defined by neurons 47, 58, 62, and 63. In Figure 6c only these four neurons are turned on defining the location of the reservoirs on this line. The four neurons are clearly defining the field and the southern downdip limits of the reservoir associated with Well A and the updip limits of the reservoir at Well B where the sands are thinning out to the south. Figure 7 is northwest-southeast inline 2109 with its location denoted in Figure 4. Figure 7a is the stacked amplitude display and Figure 7b displays the SOM results defining the limits of the Well B reservoir as it terminates at the fault to the northwest. Figure 7c with only the four neurons defining the reservoir displayed indicates the thinning out of the reservoir updip much more clearly than with amplitudes alone. The insets of Figures 7b and 7c illustrate the details in the updip portion of the reservoir defined by the SOM process. The SOM analysis incorporates ten AVO attributes and is not limited by conventional amplitude/frequency limitations of thickness and areal distribution. The AVO attributes selected for this SOM analysis are specifically designed to bring out the appropriate AVO observations for a class 2 setting. It is clear from these results that the AVO attributes in this SOM analysis are clearly distinguishing the anomalous areas associated with the producing reservoirs from the equivalent events and zones outside these stratigraphic traps.
For more than 40 years seismic amplitudes have been employed to interpret DHIs in an attempt to better define prospects and fields. There are dozens of DHI characteristics associated primarily with class 2 and 3 geologic settings. Hundreds of seismic attributes have been developed in an effort to derive more information from the original seismic amplitude data and further improve DHI interpretations. A machine learning workflow incorporating seismic attributes, PCA, and SOM, has been proven to produce excellent results in the interpretation of DHIs. This machine learning workflow was applied to data in class 2 and 3 reservoirs in an effort to interpret the most important DHI characteristics as defined by a worldwide industry database. The SOM analysis employing instantaneous attributes in a class 3 setting successfully identified the top DHI characteristics and especially those defining edge effects and hydrocarbon contacts (flat spots). AVO attributes conducive to providing information in class 2 settings incorporated in a SOM analysis allowed the interpretation of DHI characteristics that better defined the areal extent of the producing reservoirs than amplitudes by clearly denoting the stratigraphic trap edges.
Low SOM classification probabilities have been proven to help identify DHI characteristics. These low probabilities relate to data regions where the attributes are very different from the data points of all of the attributes in the SOM analysis and their associated winning neurons, which has defined a natural cluster or pattern in the data. Anomalous regions in the data, for example, DHI characteristics, are noted by low probability. This analytical approach of defining low probabilities proved to be helpful in identifying DHI characteristics in both class 2 and 3 settings.
An important observation in these two case studies is that the use of appropriate seismic attributes in a SOM analysis can not only identify DHI characteristics not initially interpreted but can also increase or decrease confidence in already identified characteristics. This multi-attribute machine learning workflow provides a methodology to produce more accurate identification of DHI characteristics and a better risk assessment of a geoscientist’s interpretation.
The authors would like to thank the staff of Geophysical Insights for the research and development of the machine learning applications used in this paper. We would also like to thank the Rose & Associates DHI consortium, which has provided extremely valuable information on DHI characteristics. The seismic data in the offshore case study is courtesy of Petroleum Geo-Services. Thanks also go to Deborah Sacrey and Mike Dunn for reviewing the paper. Finally, we would like to thank Tom Smith for reviewing this paper and for the inspiration to push the boundaries of interpretation technology.
Hilterman, F.J. . Seismic amplitude interpretation. Distinguished instructor short course, SEG/EAGE. Kohonen, T. . Self Organizing Maps. Third extended addition, Springer Series in Information Services, Vol. 30.
Smith, T. and Taner, M.T. . Natural clusters in multi-attribute seismics found with self-organizing maps. Extended Abstracts, Robinson-Treitel Spring Symposium by GSH/SEG, March 10-11, 2010, Houston, Tx.
Taner, M.T. . Attributes revisited. http://www.rocksolidimages.com/pdf/attrib_revisited.htm, accessed 13 August 2013.
ROCKY RODEN owns his own consulting company, Rocky Ridge Resources Inc., and works with several oil companies on technical and prospect evaluation issues. He also is a principal in the Rose and Associates DHI Risk Analysis Consortium and was Chief Consulting Geophysicist with Seismic Micro-technology. He is a proven oil finder (36 years in the industry) with extensive knowledge of modern geoscience technical approaches (past Chairman – The Leading Edge Editorial Board). As Chief Geophysicist and Director of Applied Technology for Repsol-YPF, his role comprised advising corporate officers, geoscientists, and managers on interpretation, strategy and technical analysis for exploration and development in offices in the U.S., Argentina, Spain, Egypt, Bolivia, Ecuador, Peru, Brazil, Venezuela, Malaysia, and Indonesia. He has been involved in the technical and economic evaluation of Gulf of Mexico lease sales, farmouts worldwide, and bid rounds in South America, Europe, and the Far East. Previous work experience includes exploration and development at Maxus Energy, Pogo Producing, Decca Survey, and Texaco. He holds a BS in Oceanographic Technology-Geology from Lamar University and a MS in Geological and Geophysical Oceanography from Texas A&M University. Rocky is a member of SEG, AAPG, HGS, GSH, EAGE, and SIPES.
CHINGWEN CHEN, PH.D. received an M.S. (2007) and a Ph.D. (2011) in Geophysics from the University of Houston, studying global seismology. After graduation, she joined the industry as a geophysicist with Noble Energy where she supported both exploration and development projects. Dr. Chen has a great passion for quantitative seismic interpretation, and more specifically rock physics, seismic imaging and multi-seismic attribute analysis. She later joined Geophysical Insights as a Senior Geophysicist, where the application of machine learning techniques became a focus of her work. Since 2015, her primary interest has been in increasing the efficiency of seismic interpretation.