Next Article in Journal
Introducing Uncertainty in Risk Calculation along Roads Using a Simple Stochastic Approach
Previous Article in Journal
Lifecycle of an Intermontane Plio-Pleistocene Fluvial Valley of the Northern Apennines: From Marine-Driven Incision to Tectonic Segmentation and Infill
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Enhancing Paleoreef Reservoir Characterization through Machine Learning and Multi-Attribute Seismic Analysis: Silurian Reef Examples from the Michigan Basin

1
Mewbourne College of Earth and Energy, School of Geosciences, University of Oklahoma Norman Campus, Norman, OK 73019, USA
2
Consumers Energy, Traverse City, Jackson, MI 49686, USA
*
Author to whom correspondence should be addressed.
Submission received: 31 December 2020 / Revised: 25 February 2021 / Accepted: 16 March 2021 / Published: 19 March 2021
(This article belongs to the Special Issue Integrated Stratigraphy of Carbonate Platforms)

Abstract

:
Historically, Silurian pinnacle reef complexes in the Michigan Basin have been largely identified using 2D seismic with very little research on the reservoir characterization of these reefs using 3D seismic data. By incorporating a high-resolution 3D dataset constrained by a well-studied and data-rich paleoreef reservoir, the Puttygut reef, seismic attributes were correlated to petrophysical properties through machine learning and self-organizing maps (SOMs). A suite of structural and frequency-based attributes was calculated from pre-stack time migrated (PSTM) seismic data, with only a subset of them selected as SOM inputs. Structural attributes enhanced details in the reef but frequency attributes were overall more useful for correlating with reservoir quality. A strong relationship between certain combination percentages of attributes and certain sections of the reef with porosity and permeability was found after the SOM results were compared to wireline log and core analysis data. Areas with high permeability and porosity correlated with the average frequency and spectral decomposition at 29 and 81 Hz. Areas with high porosity and varying permeability correlated with the average frequency and spectral decomposition at 29, 57, and 81 Hz. Areas with intermediate porosity correlated with the average frequency and spectral decomposition at 29 and 57 Hz. The efficacy of the procedure was then demonstrated on two nearby reefs with very similar results.

1. Introduction

The focus of this study was to characterize features in ancient reef reservoirs that are at and below conventional vertical seismic resolution by creating a machine learning and multi-attribute analysis workflow. A seismic attribute, first defined by Taner and Sheriff 1977 [1], is a measurement derived from seismic data, most commonly based on time, amplitude, frequency, and/or attenuation. The primary use for attributes is that they can help an interpreter visualize features, relationships or patterns that might go unnoticed [2,3]. An ideal area to develop and test such a workflow was in the Michigan Basin on Silurian Niagara-Lower Salina pinnacle reefs. The reefs form a sublinear trend on the northern and southern rims of the basin at depths between 900 and 2200 m [4]. The basin encompasses the entire lower peninsula and parts of the upper peninsula of Michigan, as well as adjacent areas of eastern Wisconsin, northeastern Illinois, northern Indiana, northwestern Ohio, and western Ontario [5]. As of 2019, these reefs have produced over 500 million barrels of oil and over 3 trillion cubic feet of natural gas [6]. In addition to what has already been produced, there is still great potential for secondary oil recovery through water, CO2 or gas injection. Michigan reef reservoirs have 26% average primary recovery and 12.5% average secondary recovery with only around 5% of the discovered fields having undergone secondary recovery efforts [7]. In addition to the production related motivation to better characterize these reef reservoirs, they have also been used for gas storage and CO2 sequestration over the past 50 years [8]. The majority of 3D seismic acquired and processed over these reefs in the past 20 years have been for gas storage and CO2 sequestration projects, with some of the gas storage reefs having core-rich datasets which are utilized in this study for correlation to seismic attributes. These are the first data sets of their kind where both 3D seismic and multiple cores are available within the same reef complex.
This study centers on paleoreefs located along the Southern reef trend in St Clair and Macomb Counties of Michigan. In addition to the apriori geologic knowledge (detailed core descriptions and conventional core analysis) the geophysical data included 3D pre-stack and post-stack time migrated seismic data over multiple reefs with well logs, and in some cases very detailed core descriptions and conventional core analysis data (personal communication, Matt Rine 2019). The Puttygut reef was chosen as the initial reef for workflow development, as it had the most complete and well-rounded dataset, with the most detailed core analysis of the reefs available. The Puttygut gas storage reef complex is approximately 150 acres in areal extent and 325 feet in subsurface height with a gas storage capacity of 14.6 billion cubic feet of gas (BCFG) ((personal communication, Matt Rine 2019). Figure 1A illustrates the area of interest located in Michigan. Within the inset (Figure 1B) are the three reefs used in this study, Puttygut, Ira, and Lenox. Zooming in further to the initial reef of interest, Puttygut, is Figure 1C which is a two-way time depth map of the Guelph Formation structure (referred to herein by driller’s terminology as the Brown Niagaran) which is comprised of the reef proper facies (refer to [9] for facies model). In general, well logs in these reefs generally reveal areas of 10% porosity or greater and regions that also have permeabilities of 50 mD or greater on average. These porous and permeable zones are what make these reefs world class reservoirs for both exploration and gas storage. These areas are relatively thin, however, sometimes only reaching thicknesses of 5 feet. In general, well log data generally exhibit good vertical resolution (in feet), but the lateral coverage provided is extremely limited (in inches). On the other hand, seismic data provide lateral coverage which is magnitudes greater but smaller lateral resolution owing to the Fresnel zones and bin size limitations [10]. In addition, seismic data are more economic to acquire and process when compared to drilling and logging, as well. The vertical resolution of the Puttygut reef seismic data at the reef level depth is approximately 60 to 80 feet. This conventional Rayleigh vertical seismic resolution was determined using the equation: VR = λ/4 where λ = (v)/f where λ is wavelength, v is velocity, and f is the dominant frequency [10]. In addition, the peak frequency was between 40–50 Hz and the interval velocity range at the depth of interest was 11,970–12,000 ft/s. Even with conventional seismic attribute analysis workflows it can be difficult to resolve the thin high porosity and permeability features within the reefs. Therefore, machine learning and the multi-attribute analysis was utilized in this study to drastically improve the level of detail seen in the seismic data. Seismic data can be used to create a suite of attributes that can, in turn, be used as inputs for a machine learning algorithm that creates self-organizing maps (SOMs). Attributes generate results at the wavelet scale, which can cause difficulty when comparing them side by side. SOMs, on the other hand, have the potential to analyze these attributes together all the way down to the sample interval of the data [11]. Once at that scale, SOMs then characterize data points into clusters based on their similarity. The SOM process greatly increases the level of detail seen in the data and consequently allows for more apt comparisons to the well log data.
It is important to increase the level of detail seen in seismic data to advance our understanding of the subsurface. With machine learning and multi-attribute analysis, it is possible to re-explore pre-existing seismic data through a new lens and see features on a much finer scale than previously possible. Machine learning can open new doors to exploring our geologic past and to increase our understanding of smaller scale events and processes. It can aid us in identifying and characterizing features of interest that were previously hidden in the data, and it is possible to use this method on any geologic regime. More practically, the machine learning method can be used to increase the level of detail seen in reservoirs which makes the task of characterization much easier. With the increased level of detail within the reservoir from the machine learning workflow, utilizing the reef for gas storage or carbon sequestration will be more precise, accurate, and easier overall. The workflow from this research was developed in an attempt to increase confidence when interpreting seismic data and the lower risk associated with planning and drilling future gas storage wells. Overall, if you have seismic data, a way to calculate attributes, and a method of machine learning, you will be able to explore that seismic data with increased detail and clarity which was not possible before with just the data and attributes. From this study, a relationship was found between the SOM results and the porosity and permeability within the Puttygut reef. Furthermore, applying this workflow on the Ira and Lenox reefs confirmed the repeatability of this method by finding very similar results.

2. Regional Geology

The Michigan Basin is a circular, intracratonic basin that covers an area of ~316,000 km2 shown in Figure 1 [5]. The reefs in this basin are Wenlock (Homerian) in age [8] and form a sublinear trend on the northern and southern rims at depth between 900 and 2200 m [4]. The Niagaran Guelph Formation reefs are unconformably overlain by, and laterally encased in, the cyclic limestones and evaporites of the Silurian Salina Group, and underlain by the Lockport Dolomite Formation [12]. See Figure 1D for a generalized stratigraphic scheme from Rine 2019 [8].
A summary of the time-step basin wide depositional evolution of the Michigan Basin developed by Rine 2019 [8] can be generalized as follows: Through a time of sea level oscillations the first formation to be deposited in the basin is the Lockport formation. Eventually, a sea level fall exposes the Lockport and ends its deposition. Following is a sea level rise which coincides with the initiation of the bioherm mounds and marks the beginning of the lower Guelph Formation. Next, the sea level falls again which exposes the bioherm mounds and causes significant wasting. That is what creates the flat tops seen today. A rapid sea level rise follows initiating the pinnacle reef growth in the upper Guelph Formation. Tidal flat facies begin to be deposited as the basin experiences restriction. The relative sea level fall in the basin results in the A-1 Evaporite deposition. The sea level rises, falls, and rises again resulting in the deposition of the lower A-1 Carbonate, the rabbit ear anhydrite (REA), and the upper A-1 Carbonate, respectively. Moreover, at this point in time, the reefs are completely covered in sediment. The A-2 Evaporite section is then deposited with another fall in sea level. Finally, through another rise and fall in sea level, the A-2 Carbonate and the B-Salt are deposited, respectively, which marks the end of the relevant sequence for this study. While the primary depositional model for these reefs is relatively well known, the post-depositional diagenetic history of these reefs is still relatively unknown. Since the reefs in SE Michigan have been 100% dolomitized and are partially salt plugged, the primary depositional facies are not the only factors that contribute to the petrophysical properties, hence the need for utilizing seismic attribute analysis for better understanding reservoir quality distribution spatially.
The reservoir interval for these reef reservoirs extends from the top of the A-1 Carbonate (Ruff Formation) and continues through the Brown Niagaran (Guelph Formation) down to the base of the reservoir or Gray Niagaran (Lockport Dolomite Formation). While these formation’s thicknesses vary widely throughout the reef, at well P-201 the thickness of the A-1 Carbonate is ~44 feet, the tidal flat cap is ~43 feet thick, the reef core is 115 feet thick (with a reservoir interval which is 105 feet thick), and the bioherm is 151 feet thick (with a reservoir interval which is 38 feet thick). See Figure 1 for locations of well P-201 and the formations. Since these formations are all comprised of dolomite, they do not exhibit good contrasts in acoustic impedance and are near or below seismic resolution. This general homogeneous lithology affected acoustic isotropy and contributes to the difficulty in determining where to pick these horizons by inspection of the seismic data alone. Figure 4 shows an inline (east-west direction) through the Puttygut reef. The seismic horizons were selected by an expert geologist (M. Rine) in the area with substantial apriori geologic and structural knowledge of the Michigan Basin in order to achieve the most accurate horizons picks by correlating the seismic with the available well data and cross sections. The A-2 Carbonate and Clinton Formation are easily visible owing to their acoustic impedance contrast to their surroundings. However, the A-2 Anhydrite, A-1 Carbonate, Brown, and Gray Formations are not so easy to resolve. These horizons are detectable at points along the line but are not confidently chosen without the prerequisite geologic and structural knowledge in the area. Unfortunately, the horizons which are more difficult to pick are the horizons that are of most interest in reservoir quality. The reef core and bioherm between the A-1 Carbonate and the Gray Formation are commonly the areas where the best reservoirs are found, as these areas of high permeability and porosity tend to exist in the skeletal wackestone and reef boundstone of the reef core, as well as the foreslope reef talus conglomerates that flank the windward (eastern) side of the pinnacles. The thrombolitic bindstone facies within the A-1 Carbonate that caps the reef complex in the reef crest position is also an area of high permeability and porosity [8].

3. Materials and Methods

The initial workflow consisted of running a large suite of attributes on the pre-stack 3D volume around the Puttygut reef. Geologically, we are primarily interested in improving the mapping of continuous areas of high porosity and permeability found in the sub-seismic formations that comprise the reef core and bioherm. As these areas are seismically difficult to detect, areas of interest with high permeability and porosity are initially chosen using core and well log data. There is one limitation with this approach, however, as the core and well log data are in depth and the seismic data are in time. There are no reliable sonic or density logs within the reef proper, only sonic logs in nearby off-reef wells. Therefore, producing a reliable well tie within the reef is particularly challenging without making assumptions that drastically affect the reliability of the results. The workflow for picking sub-seismic horizons was as follows: (1) Pick formation tops from geophysical logs and cores; (2) tie the wells to the seismic using a VSP from the P-201 well (see Figure 1C for well location) and synthetic seismograms generated from wells with DT and RHOB logs outside of the survey; (3) shift the wells to time using key horizons (i.e., A-2 Carbonate, A-2 Salt, A-1 Carbonate, Clinton. See Figure 2 and Figure 3 for the VSP results. See Supplementary Figures S1–S3 for uncropped versions); (4) pick sub-seismic horizons that cross-cut reflectors based on well control (i.e., Brown Niagaran, Bioherm, Gray Niagaran). From there, the smaller, sub-seismic horizons were mapped through apriori structural geologic knowledge and using the pseudo-tied wells to aid in picking locations. While this workflow is not ideal, there is a very high level of confidence with the accuracy of these provided horizons.
Attributes were generated and evaluated by the interpreter to determine how effective they were at displaying minute details within the target reef. Attributes that exceled at revealing the internal, less detectable aspects, of the reef were then selected for multi-attribute analysis in the form of self-organizing maps (SOMs) to achieve sub-seismic resolution.
The theory of using SOMs and multi-attribute analysis to image and interpret sub-seismic resolution has been successfully tested and proven in the past [11]. Most commonly SOMs were calculated over unconventional plays, such as in the Denver-Julesburg Basin and in the Eagle Ford shale [13,14,15]. They have also been tested in other areas, such as for interpreting direct hydrocarbon indicator (DHI) characteristics, classifying carbonate facies, and recognizing geologic patterns with varying degrees of success [16,17,18]. A SOM was explored in detail by Teuvo Kohonen [19] in his book where a very mathematical description of SOMs can be found. Kohonen describes SOMs as “nonparametric regressions” where a number of ordered discrete reference vectors are fit to a distribution of vectorial input samples. In order to approximate continuous functions though, the reference vectors are used to define nodes or neurons, in a hypothetical elastic network. Local interactions between these neurons along the neural network create that essence of elasticity. By setting up the map in this way, the neurons develop into specific detectors of signals in the input space. These detectors are formed in the map in a meaningful order as if falling into place on a feature coordinate system. These resultant feature maps are used for the preprocessing of patterns for recognition [19]. Stated simply, SOMs require some number “m” of input characteristics and a desired number “n” of output “detectors”. These output detectors are then fit to the input characteristics such that they both attempt to classify the data and are influenced by the other detectors. What is made in the end is a map which is organized in a meaningful manner that can be used for potential pattern recognition. These inputs can be anything from variations of textures in an image, to the variations in the quality of life between nations (an example created by Kohonen 1995 [19]), to the different attributes that can be created using seismic data. As mentioned previously, these SOMs analyze the input data at the sample interval level.
In machine learning, for a desired anthropologically reasoned outcome, rather than choose characteristics randomly, it is optimal to take advantage of the interpreter’s knowledge when choosing SOM attributes and parameters. In this case, two categories were decided upon when initially analyzing single attributes: Structural based attributes and frequency-based attributes. The structural based attributes were subjectively graded on how well they illuminated the internal structure of the reefs. The frequency-based attributes were similarly graded on how well they distinguished unique sections of the reef. Highlighting the sub-seismic internal structure of the reef and identifying different sections of the reef via frequency content potentially create the key to isolating areas of higher porosity and permeability. Narrowing down the attributes to a smaller set is also useful in the SOM stages later in the workflow. It is not necessarily beneficial to have more attributes when conducting multi-attribute analysis as many attributes are calculated in similar ways and can show similar features that can be unnecessary [3]. If all the attributes that were calculated are used in one SOM then there would be a great amount of unnecessary redundancy. In addition, using more data in the calculations will increase the amount of time required to complete a SOM. With these criteria in mind, the chosen set of attributes are defined below:
Structural attributes: The following structural defining attributes were chosen since they all excelled at defining some form of structure (internal or boundary) and since they are all calculated in reasonably different ways, which is important later during the multi-attribute analysis step.
  • Energy ratio similarity (ERS): An attribute which is designed for edge-detection. ERS is similar to other coherence attributes except that it is calculated along the structure and that it uses analytic traces rather than just the real traces. ERS is the ratio of the energy of the Karhunen-Loève filtered data over the total unfiltered energy of the input data within the analysis window [20]. For the purposes of this research, ERS discriminated best the reef margin as well as some internal features that other attributes were unable to detect. The idea for reef edge and internal structure detection by the use of coherence or similarity attributes has been tested extensively in the past. Chopra and Marfurt 2007a [20] and Skirius et al. 1999 [21] successfully used ERS to aid in carbonate reef interpretations.
  • Aberrancy: Defined as the third derivative of structure, it is the measure of the asymmetry of a curve about its normal [22]. Aberrancy measures the lateral changes in the curvature along a surface and is also able to detect faults that are below seismic resolution [23]. Qi and Marfurt 2018 [23] calculated aberrancy on a data volume over the Barnett shale gas reservoir in Fort Worth Basin, Texas and were able to identify small karst features unseen with coherency. Aberrancy was able to define smaller and more internal reef structures than other structural attributes for the purpose of this study.
  • Curvature: Defines how positively or negatively curving a feature is. Curvature is a second derivative calculation of the phase of the seismic waveform. If a feature has a positive curvature then it has an antiformal expression, but if a feature has a negative curvature it has a more synformal expression. Chopra and Marfurt 2007b [24], Duan et al., 2010 [25], and Wallet 2016 [26] successfully used a variety of curvature attributes to aid in their interpretations.
  • Gray level co-occurrence matrix (GLCM): GLCM is a statistical textural attribute that essentially compares a sample point within an analysis window to the other sample points within the same window. Textural attributes aid an interpreter comparable to the use of a similarity attribute. Eichkitz et al. (February 2015) [27] used a combination of a GLCM-based energy attribute calculated in four different directions in a comprehensive workflow to identify areas with high and low degrees of anisotropy. Eichkitz et al. (March 2015) [28] later used this attribute to help extract more information on channel interiors while using a coherence attribute to distinguish the edges of the channels. GLCM was useful for identifying additional features within the reef core that were unnoticed by the other structural attributes.
Frequency attributes: The following frequency attributes were chosen based on how well they differentiated the internal reef structure.
  • Cosine of Instantaneous phase: Instantaneous phase attribute provides the interpreter with the phase of the wavelet at a certain location in the data. Yet, instantaneous phase is cyclical in nature, wrapping around from −180 to +180°, and is also mathematically discontinuous. The cosine of the phase is implemented in order to remove those discontinuities [29]. Taking the cosine of the instantaneous phase also removes the cyclical nature of the attribute, creating a scale from −1 to +1, which is advantageous later on when being used in a SOM. Sukmono et al., 2006 [30] and Huang et al., 2011 [31] incorporated the generation of cosine of instantaneous phase (or just instantaneous phase) in order to aid in carbonate reef interpretations.
  • Average frequency: Instantaneous frequency is a measure of the frequency of the wavelet of the seismic trace at a given location in the data [29]. Toelle and Ganshin 2018 [32] found a relationship between the instantaneous frequency of seismic data and the amount of porosity present in reefs located in the northern Michigan Basin. They found that the instantaneous frequency decreased from 73 to 45 Hz as the porosity increased from 5% to 20%. They interpreted this decrease in frequency to be associated with frequency attenuation from the increase in porosity. Huang et al., 2011 [31] and Sarhan 2017 [33] also used instantaneous frequency to aid in carbonate interpretations. Unfortunately, instantaneous frequency tends to have spikes and noise that can mask important information. One way to counter this interference is to take the average of that instantaneous frequency. Consequently, the reduction in noise and spikes in this attribute from averaging outweigh the potential loss of resolution [34]. The average frequency from the Puttygut volume produced pockets of low frequency near wells with higher porosities and sections of high frequency in other areas.
  • Spectral decomposition (continuous wavelet transform (CWT)): This is an attribute that separates the amplitude seismic data cube into different frequency cubes. Although there are more than one way to calculate spectral decomposition, the CWT is a linear form which is based on the short time Fourier transform (STFT). Nejad et al., 2009 [35] and Saadatinejad et al., 2012 [36] used spectral decomposition successfully to aid in reef interpretation and exploration.
One of the ultimate goals of this research is to attempt to map/identify areas of high permeability and porosity within the reef reservoirs using machine learning and multi-attribute analysis. Yet, these two characteristics of reservoirs are challenging to unequivocally map with absolute accuracy [37,38,39,40,41]. These past studies have one thing in common when attempting to map or predict permeability or porosity: The attributes chosen are mostly derived from or related to the frequency and/or phase of the seismic data. In order to determine if the frequency-based attribute approach at identifying permeability and porosity is effective for the study area of this research, all attributes will be looked at initially on the provided high confidence horizons that represent the reef lithofacies. These detailed horizons make it possible to accurately estimate the location of high permeability and porosity zones surrounding wells that penetrate the reef. From there, it should be possible to work out from each well into the reef to determine how much variability in permeability and porosity there is. Four value ranges were set to differentiate “great” from “poor” permeability and porosity. The ranges for porosity are: 0–3% is poor, 3–6% is intermediate, 6–10% is good, and greater than 10% is great. The ranges for permeability are: <1 mD is poor, 1–10 mD is intermediate, 10–50 mD is good, and greater than 50 mD is great.

4. Puttygut Data

The Puttygut reef reservoir has a recent vintage 3D seismic acquisition and processing and seven cored wells with petrophysical properties and facies descriptions but has a relatively older vintage wireline log dataset. These log suites are lacking sonic and density logs which are needed for tying the wells to seismic. From the Puttygut seismic survey, two processed volumes were generated in 2019 (one is pre-stack time migrated and the other is post-stack time migrated). The inline and crossline spacing is 440 feet, with shot and receiver spacing of 110 feet, a final bin size of 55 by 55 feet, and a sample rate of 1 ms. The data is zero-phase with American polarity, where the positive impedance contrast is recorded as a peak. There are 30 wells in the area drilled within and around the reef (Figure 1C). All 30 wells had at least a gamma ray log, 14 wells had neutron porosity logs, two had resistivity logs, one well had an SNP log, and one well had a DT log. There were also 11 wells that had core data that were compiled and transferred into the LAS format. Eleven wells had permeability and porosity data, nine wells had oil and water saturation data, and four wells had horizontal permeability measured at 90° to the maximum permeability direction. In addition to all this well and 3D seismic data, there is a vertical seismic profile (VSP) at well P-201.
Generating well ties in the Puttygut seismic area was an iterative process. Owing to the lack of geophysical logs (sonic and density) within the Puttygut reef complex proper, fortunately a vertical seismic profile (VSP) was conducted in conjunction with the Puttygut 3D seismic acquisition in the P-201 wellbore. Since the reservoir interval (A-1 Carbonate to Gray Niagaran) was perforated throughout, a packer was set just above the A-1 Carbonate resulting in the A-2 Anhydrite being the deepest formation recorded during the VSP survey. That VSP survey (see Figure 2 and Figure 3 for the VSP results. See Supplementary Figures S1–S3 for uncropped versions) was used to generate time-depth tables for the P-201 and those time-depths tables were applied to the surrounding reef crest wells. In the off-reef position, the Buszek well (PN 23843) had a sonic log that was used to generate a synthetic seismogram. The synthetic seismogram tied well to the 3D seismic (Figure 4), and thus was used to generate a time-depth table for that well and was applied to the surrounding off-reef wells. Since the VSP was unable to be conducted within the reef proper, the top (A-1 Carbonate/Brown Niagaran) and base (Gray Niagaran) horizons of the reef reservoir had to be interpreted by comparing the seismic to structural cross sections. Furthermore, the lack of acoustic impedance on the reef crest between the A-2 Carbonate (dolomite), A-2 Anhydrite (anhydrite), A-1 Carbonate (dolomite), Brown Niagaran (dolomite), and Gray Niagaran (dolomite), results in a lack of strong reflectors to correlate horizons from off-reef to on-reef. Therefore, while those horizons were picked by closely correlating with structural cross-sections from the available well control, there is some uncertainty associated with them.

5. Results

A wide range of single attributes were calculated and compared by the interpreter to decide which attributes best showed features of interest in the reef. Therefore, after grading according to the criteria previously described, the following attributes were deemed most worthy. Four structural attributes (energy ratio similarity, negative curvature, aberrancy, and GLCM homogeneity) and four frequency attributes (cosine of instantaneous phase, average frequency, continuous wavelet transform (CWT), spectral decomposition at frequencies of 29, 57, and 81 Hz) were extracted from the seismic volume. The approximate thickness ranges of these frequencies that correspond to velocities at the reef interval are 107, 55, and 38 feet for 29, 57, and 81 Hz, respectively.
The following section details observations of the single attributes results extracted onto horizons picked by an expert geologist that works in the Michigan Basin. The single attributes that were calculated are extracted onto the Brown horizon slice and initial observations were recorded. The Brown horizon is the top of the Guelph reef complex (the top of the Reef Core, Reef Apron, and Reef Talus geobodies), which made it the most important horizon for displaying the attribute maps which is the most representative of the morphology of the Guelph reef complex. Then, six total SOMs were calculated. All six SOM calculations were output with a 5 × 5 neuron grid. Three of the six SOM runs used the aforementioned structural attributes in addition to the PSTM amplitude, while the other three used the aforementioned frequency attributes without the PSTM amplitude. The amplitude was removed from the frequency SOM runs since it tended to dominate the results. Following removal, there was a greater variability in which the attribute contributed the most for each cluster of data. The first two (one structure and one frequency) SOM runs were restricted between the A-2 Carbonate and Clinton layers and utilized all the inlines and crosslines. The next two (one structure and one frequency) were restricted between the same two horizons but the inlines and crosslines were restricted to isolate the reef from the surrounding geology. The final two SOM runs (one structure and one frequency) used the same inline and crossline restrictions as the second set of SOMs but were further confined vertically to be between the A-1 Carbonate and the Gray layers. The goal of this arbitrary section is to identify areas of interest along well bores where there are permeability and porosity data. Four type wells that had the most complete suite of porosity and permeability data were chosen for this task: P-106, P-201, P-102, and P-103.

Single Attributes on Brown

Figure 5 reveals the structural and frequency attributes extracted onto the Brown horizon in Petrel. The light pink outline in each image represents the reef extent. Figure 5A is a two-way time (TWT) map of the Brown horizon. The reef complex is a clear structural high compared to the off-reef areas, with a slightly higher structure observed in the northern section of the reef complex compared to the southern. The slope is gentler on the western edge of the reef and steeper on the eastern edge, which corresponds to the asymmetric reef model proposed by Rine et al., 2017 [9]. Figure 5B shows pre-stack time migrated (PSTM) amplitude data. There is a consistent negative amplitude that composes most of the reef top, followed by a ring of positive amplitude that contours the reef outline. Figure 5C shows the energy ratio similarity (ERS). The yellow arrows filled with black are showing areas of lower similarity along the eastern side of the reef, as well as one spot in the southern portion. The yellow lines outline some sinuous features, potentially talus deposits, that surround the main reef body. Figure 5D shows the negative curvature attribute. Outlined in yellow is a connected feature that runs throughout the reef body. There is also a broad trend of negative curvature that wraps around the reef outline itself everywhere except for the southern tip. The area that surrounds the reef complex is comprised of alternating carbonate (A-1 Carbonate) and halite (A-2 Salt) layers, which have a much greater acoustic impedance than the reef complex interior which is comprised of only carbonate. Figure 5E shows aberrancy. An interconnected feature is present here as well, but it takes up more of the reef. Within that interconnected feature there are small pockets of little to no aberrancy. At the bioherm level, these features might be associated with coral growth patterns within the reef. Figure 5F shows GLCM homogeneity. The yellow circle shows an area of high homogeneity which is concentrated in the northern section of the reef. There is also a low homogeneous zone surrounded by a tight ring of high homogeneity on the reef outline. The tight ring of homogeneity may be areas in which reef debris or talus collected as the reef was exposed and eroded during the sea-level fall. The yellow arrows filled with black are pointing to areas surrounding the reef which are high homogeneity concentrations. These concentrations might indicate where smaller bioherm mounds initially formed but were abandoned/drowned as the sea level rose during the time of reef growth. Figure 5G shows the average frequency. The black shapes outline areas of particularly low average frequency in the main reef body, approximately 20 Hz on average. The rest of the reef reveals a higher average frequency, approximately 70 Hz on average. The black arrows are pointing to areas of anomalously high or low average frequency which lie outside of the reef body. The varying frequencies found within the reef core potentially indicate a transition from a higher porosity and permeability to a lower porosity and permeability lithofacies. Figure 5H shows that the cosine of the instantaneous phase accentuates what the PSTM amplitude shows. There is a negative cosine value in the reef body, surrounded by a positive ring that hugs the reef body. The attribute also shows that there is another ring with a negative cosine value that tightly curves the inner features to the east and more broadly curves to the west. Outside of that ring lies mostly positive values speckled with negatives. Figure 5I–K displays the CWT spectral decomposition at frequencies of 29, 57, and 81 Hz, respectively. Here, 29 Hz corresponds to approximately 107-foot thicknesses, 57 Hz to 55-foot thicknesses, and 81 Hz to 38-foot thicknesses. Figure 5I shows that the CWT spectral decomposition at 29 Hz shows the main reef body has a lower presence of 29 Hz relative to its surroundings everywhere except where the yellow circle is in the southern portion of the reef. The presence of 29 Hz in the southern portion of the reef might indicate that it is slightly thicker than the rest of the reef. Figure 5J is the CWT spectral decomposition at 57 Hz. Around the reef there is a ring of higher presence of 57 Hz. The southern half of the ring exhibits a higher presence than the northern half. Most of the reef body is near the void of 57 Hz, except where the yellow lines are outlining. There appears to be small shapes or features that have some presence of 57 Hz within the reef itself. Going from 29 to 57 Hz, we begin to pick up on thinner features. In this case, we go from frequencies which are associated with 107- to 55-foot thicknesses. Seeing more detail within the reef itself as the frequency increases could imply that the internal growth patterns of the reef happened on a smaller scale. The ring surrounding the reef seen with the 57 Hz spectral decomposition may be associated with an area of reef debris that can be found surrounding the reef core. Figure 5K is the CWT spectral decomposition at 81 Hz. There is a seemingly random distribution at first glance within the reef outline but upon closer inspection, areas where the attribute displayed similar values are usually connected. The yellow lines outline the more drastic regions where the higher presence occurred within the reef. Looking out from the reef the same pattern can be seen. The highly interconnected areas seen with the 81 Hz spectral decomposition could be associated with the more complex internal coral growth patterns found within the reef core. Figure 5L is a traditional RGB co-blending of the three spectral decomposition volumes displayed at 430 ms. Through this co-blending the same internal features, with perhaps a little more detail, can be better visualized in order to interpret each spectral decomposition image individually. Yet, using the spectral decomposition volumes in a SOM, there is potentially more detail that can be extracted. Each of these attributes brings something different to the table that will become more important when combining them into a SOM. In addition, all SOM results output the data from the wavelet scale to the voxel scale.
After comparing the set of SOM runs, the third SOM run with the frequency attributes was next used. Figure 6 shows a comparison of the second and third SOM runs. The second SOM run was created with more input data than the third. Some notable differences are noted. The white circles show an area in the second SOM characterized by 3 neurons but in the third SOM that same area is characterized by at least 5 neurons. The white rectangle shows an area characterized by a few neurons in the second SOM and a significant amount more in the third SOM. Finally, the arrows point to areas where the general trends and patterns are the same. The gold and purple neurons in the second SOM are similar to the brown and green neurons in the third SOM. However, there are still some slight changes in these features, as well. These differences are caused by the varying amount of input data into each SOM. The third SOM was used due to the higher level of detail seen in areas of interest within the reef.

6. Discussion

The third SOM run with the frequency attributes was chosen to explore the reef further. The input data and parameters for the third SOM were best suited to achieve the goal of identifying areas of poor to great permeability and porosity. Similar graded ranges as used previously were used for both porosity and permeability. Horizons with high levels of confidence were used to bridge the gap between the well logs and the time volumes. A qualitative interpretation on the efficacy of the frequency SOM runs to identify sub-seismic areas of varying permeability and porosity was conducted using the available data. Figure 7 shows the A-A’ cross section from Figure 1. The top image in Figure 7 reveals the third SOM run from the frequency attributes with the wiggle traces of the seismic amplitude overlain. Clearly the SOM results are able to define and characterize the reef on a much finer scale. The size of the voxels for this particular dataset were 55 by 55 feet laterally and 1 ms vertically.
Areas of varying permeability and porosity were located on the well logs and then the neurons associated with those locations were noted. At well P-106, the best-case scenario was found with great permeability and porosity on average. Neurons associated with this case are numbers 19 and 23. Two areas of greater than 10% (good to great) porosity were found with widely varying permeability at wells P-201 and P-102 that correspond to neuron 19. Finally, well P-103 had locations where the porosity is roughly the same, but with vastly different permeability were found. The two neurons associated with this scenario are 10 and 25. The bottom image of Figure 7 shows these aforementioned neurons which are associated with the desirable porosity and permeability in blue, while the rest of the neurons are shown in grey. The well logs are also displayed over the SOM. The permeability log (left) has a logarithmic scale ranging from 1 to 1000 mD and the porosity log (right) has a linear scale ranging from 0% to 20%. Blue areas of each respective log are associated with greater than 50 mD permeability or greater than 10% porosity. The green areas are associated with 10–50 mD permeability or 6–10% porosity. Yellow areas are associated with 1–10 mD permeability or 3–6% porosity. Finally, the red areas are associated with <1 mD permeability or 0–3% porosity. Neuron 23 was only present at the one location where there is the highest values of permeability and porosity. Neuron 23 is also the only neuron with a high percentage of CWT spectral decomposition at 81 Hz at 40.6% contribution. The other top contributors for this neuron are the average frequency at 35.3% and CWT spectral decomposition at 29 Hz with a 22.4% contribution. Neuron 19 can be found where there was a wide range of permeability but intermediate to great porosity. Neuron 19 is also found next to neuron 23 at well P-106 where the best properties were found. The top three attribute contribution percentages for neuron 19 are CWT spectral decomposition at 29 Hz with 39.3%, average frequency with 29.9%, and CWT spectral decomposition at 57 Hz with 24.3%. The two neurons found at well P-103 may seem confusing at first, but upon closer inspection the attribute contribution percentages for neurons 10 and 25 are almost identical. The porosity in this well has slight variations from good to great through the whole section, while the permeability is more highly varied and does not seem to follow any specific trend. That information implies that these two neurons are mostly driven by the porosity and not as affected by the permeability. The blends of neurons 10 and 25 have attribute contribution percentages of ~38% for CWT spectral decomposition at 29 Hz, ~31–35% for CWT spectral decomposition at 57 Hz, ~14–15% for average frequency, and ~11% for CWT spectral decomposition at 81 Hz. The banding seems to come from the 3.2% difference in the cosine of instantaneous phase attribute contribution. There were also areas of undesirable rock properties that were found to be associated with two other neurons, 18 and 20. These two neurons only appear at well locations associated with poor porosity and permeability conditions. These two neurons are dominated by two attributes. Neuron 18 has a 51.2% contribution from CWT spectral decomposition at 29 Hz and a 30.5% contribution from average frequency. Neuron 20 has a 44.5% contribution from average frequency and a 39.3% contribution from CWT spectral decomposition at 29 Hz. Essentially, areas of high permeability and porosity are found where CWT spectral decomposition at 81 Hz, average frequency, and CWT spectral decomposition at 29 Hz are the highest contributing attributes. Areas that have high porosity, but widely varying permeability, were found to be characterized by two sets of attributes: (1) CWT spectral decomposition at 29 Hz, average frequency, and CWT spectral decomposition at 57 Hz. (2) CWT spectral decomposition at 29 Hz, CWT spectral decomposition at 57 Hz, average frequency, CWT spectral decomposition at 81 Hz. Finally, areas of poor porosity and permeability were characterized by attribute combinations dominated by average frequency CWT spectral decomposition at 29 Hz.
The next phase of this research applies the workflow developed over the Puttygut reef and is applied to other nearby reef complexes to test the repeatability and robustness of the method. These reefs have less porosity and permeability control, so this is a qualitative study in the end. Figure 1 shows the locations of the reef complexes in relation to each other. The following discussion will consist of a description of each reef, the parameters used to create the attributes and SOMs, and the process by which the neurons were filtered and selected. Then, an interpretation of the results and supplementary evidence of well information to confirm the claims follows.
The input data for the Ira and Lenox SOMs were the average frequency, cosine of instantaneous phase, CWT spectral decomposition at 29, 57, and 81 Hz. All these frequency attributes were created from the 3D PSTM amplitude volumes. Horizons for each reef were also provided from the A-2 Carbonate down to the Clinton. Only the A-1 Carbonate, the Brown, and the Gray horizons were used for this exercise. These attributes and horizons were imported, and the SOMs were created. The boundaries for each SOM were slightly different but followed the same principles as with the third Puttygut SOM run. They were vertically confined between the A-1 Carbonate and Gray horizons and confined laterally by a set of inlines and crosslines.
Once each SOM was completed, the resultant neurons were filtered based on their attribute contribution percentages. If there were any blends similar to the successful blends developed from the Puttygut trials then they were kept, if not they were discarded. There were no exact matches with the attribute percentages in neurons calculated at Ira, but there were three results that were very close. Lenox was primarily composed of three neurons within the reef body, so the filtering and selection process was very short. Again, there were no exact matches at this reef either, but two of the three neurons that were found within the reef were very close. The three blends of attributes that were developed at Puttygut and looked for in the Ira and Lenox reefs were: (1) Areas of >50 mD permeability and >10% porosity had a blend of CWT spectral decomposition at 81 Hz, average frequency, and CWT spectral decomposition at 29 Hz in order from the highest contribution to the lowest. (2) Areas with 6% to >10% porosity had a blend of CWT spectral decomposition at 29 and 57 Hz in higher percentages and average frequency and CWT spectral decomposition at 81 Hz in lower, but still significant, percentages. (3) Areas with porosity ranging from 3% to >10% had a blend of CWT spectral decomposition 29, average frequency, and CWT spectral decomposition at 57 Hz in order from the highest contribution to the lowest. A close proxy to each of these cases was found at the Ira reef, while proxies to cases one and two were found at the Lenox reef. The attribute blends were surprisingly similar to what was found at the Puttygut reef.
Ira 3D is located to the south of Puttygut. The reef has an area of 100 acres and is approximately 190 feet in height with a capacity of 6.3 BCFG. The seismic data that were used as the base for all the attribute and SOM analysis was the PSTM amplitude data shot in 2018. The inline and crossline spacing was 440 feet, the shot and receiver spacing was 110 feet, and bin size was 55 by 55 feet, and the sample rate was 1 ms. All the frequency-based attributes were calculated from the PSTM amplitude seismic volume. These attributes are the cosine of the instantaneous phase, the average frequency, and the CWT spectral decomposition at 29, 57, and 81 Hz. The input parameters for the SOM mimicked the Puttygut SOM as close as possible, with the lateral extent confined to the reef as much as possible, and the vertical extent between the same two horizons (the A-1 Carbonate and the Gray horizon). After filtering through the results of the SOM, three neurons were found that had attribute contribution percentages very similar to those found at Puttygut when looking for areas of varying permeability and porosity. The right two images of Figure 8 show crossline 58 at the Ira reef SOM which goes through five wells in the area. The top image is with all the neurons displayed, while the bottom image is with the neurons that are similar to the ones found at Puttygut displayed. Neuron 6 had the highest contributions from CWT spectral decomposition at 57 Hz of 48.6%, average frequency of 33.9%, and CWT spectral decomposition at 29 Hz of 17%. Neuron 16 has the highest contributions from CWT spectral decomposition at 57 and 81 Hz, average frequency, and CWT spectral decomposition at 29 Hz of 36.2%, 26.2%, 23.4%, and 13.4%, respectively. Finally, neuron 23 has the highest contributions from average frequency of 35.3% and CWT spectral decomposition at 29 and 81 Hz of 34.9% and 28.4%, respectively. When comparing the neuron attribute blends of Ira to Puttygut, it is observed that neuron 6 from Ira is most similar to the third case at Puttygut where there is a wider range of porosities from intermediate to great but changes in permeability are not captured. Neuron 16 is most similar to case two from Puttygut where there are areas of 6% to >10% porosity, but again changes in permeability are not well captured. Finally, neuron 23 is most similar to case one at Puttygut where there are areas of >50 mD permeability and >10% porosity present. From the images on the right side of Figure 8 we see that neuron 6, which had a similar attribute contribution percentage as case three from Puttygut, was not very abundant overall. There is a higher concentration of neurons 16 and 23 within the Ira reef. Neuron 16 is localized near the center of the reef, while neuron 23 is seen more in larger clumps on the outskirts of the reef, as well as in a single location in the center of the reef, in between a section of neuron 16. The vertical white lines on the crossline image represent the locations of wells that penetrate the reef. Three of the five wells land within areas of higher porosity and permeability (ranges mentioned earlier) indicated by the active neurons, with the other two wells to the north only coming close. Those three wells, I-201, I-111, and I-108 also happened to be the three highest ranked wells in the field. The other two wells that did not penetrate the neurons but only came close, I-104 and I-202, were ranked ninth and twelfth out of twelve in the field.
Lenox 3D is located to the west and slightly south of Puttygut. The reef has an area of 35 acres and is approximately 260 feet tall with a capacity of 3.2 BCFG. The seismic data that were used as the base for all the attribute and SOM analysis was the PSTM amplitude data shot in 2018. The inline and crossline spacing was 660 feet, the shot and receiver spacing was 110 feet, and bin size was 55 by 55 feet, and the sample rate was 1 ms. All the frequency-based attributes were calculated from the PSTM amplitude seismic volume. These attributes are the cosine of the instantaneous phase, the average frequency, and the CWT spectral decomposition at 29, 57, and 81 Hz. The input parameters for the SOM mimicked the Puttygut SOM as close as possible, with the lateral extent confined to the reef as much as possible, and the vertical extent between the same two horizons (the A-1 Carbonate and the Gray horizon). The filtering process for the Lenox reef was much easier, as there were only three neurons that were within the reef body, one of which did not match any of the attribute trends or relationships from Puttygut. That left two that were relatively similar to the cases discovered at Puttygut. Neuron 23 had the highest contributions from CWT spectral decomposition at 81 and 29 Hz, and average frequency of 34.1%, 30%, and 25.3%, respectively. In addition, neuron 24 which had the highest contributions from CWT spectral decomposition at 57 and 29 Hz, average frequency, and CWT spectral decomposition at 81 Hz of 29.1%, 23.5%, 23.4%, and 15.7%, respectively. It is very clear from the left image in Figure 8 that the dominant neurons within this reef are numbers 23 and 24. When comparing the neuron attribute blends of Lenox to Puttygut it can be seen that neuron 23 from Lenox is most similar to case one at Puttygut where there are areas of great permeability and porosity present, and neuron 24 at Lenox is most similar to case two from Puttygut where there are areas of good to great porosity, but changes in permeability are not well captured. The left image in Figure 8 displays the frequency SOM through the center of the reef at inline 56. It is very clear that neuron 24 is most of the reef interior, while neuron 23 makes up more of the reef boundary. There is some slight intrusion of another neuron, number 9, in the reef core, but its attribute contribution percentages did not closely resemble anything from the three Puttygut cases, so it was excluded from further analysis.
These results hold promise for future applications of machine learning to carbonates. The relationship found between the neurons and the porosity and permeability information demonstrates that this method can successfully locate areas of desirable rock properties. There are, however, some limitations with the current workflow. The main limitations being that the workflow was created from one reef, the Puttygut reef. Subsequently, that workflow was only tested on two other nearby reefs, the Ira and Lenox reefs. Reefs reservoirs are created by many different processes throughout the world, so it would be ideal to test this developed workflow in a new location and then to tailor it slightly to suit the needs of the geology in the area. Those changes to the workflow could include, but are not limited to, increasing or decreasing the amount of input data into the SOM, altering or changing the input attributes or potentially looking for different attribute contribution percentages in order to represent areas of higher porosity and permeability.

7. Conclusions and Future Work

From this work, there does seem to be a trend or relationship between porosity and permeability in the seismic/well data and the SOM calculations (see Figure 7). When comparing the well logs at the Puttygut reef to the SOM results, areas with >50 mD permeability and >10% porosity present coincide with attribute blends with high percentages of CWT spectral decomposition at 81 Hz, average frequency, and CWT spectral decomposition at 29 Hz. The core profile for well P-106 associates this depth with a heavily karsted reservoir zone. Other attribute combinations were also found that were mainly affected by porosity and not so much by permeability. Areas of 6% to >10% porosity were associated with attribute blends of CWT spectral decomposition at 29 and 57 Hz in higher percentages and average frequency and CWT spectral decomposition at 81 Hz in lower, but still significant, percentages. The relationship between those four attributes and the porosity was discovered at well P-103 and stretched from the A-1 Carbonate down almost to the base of the bioherm layer. These two neurons exhibit a similar vertical extent in multiple locations throughout the reef, creating pillars of desirable porosity. Finally, areas with a wider range of porosity, from 3% to >10%, had attribute blends in decreasing percentages of CWT spectral decomposition 29, average frequency, and CWT spectral decomposition at 57 Hz. These locations were near the upper section of the bioherm. Areas with attribute blends dominated with only CWT spectral decomposition at 29 Hz and average frequency should be avoided as they were associated with areas of lower porosity and permeability.
A Puttygut reef model created by Garrett 2016 [42] was compared to the four best porosity and permeability case neurons from the Puttygut SOM from this study. It is possible to make a qualitative comparison between the chosen neurons and the geology of the reef. The attributes associated with neuron 23 seem to correspond to areas where the reef talus is present. Neuron 19 corresponds to areas where the reef core is located. Neurons 10 and 25 correspond to areas where the reef core and bioherm are located. Garrett 2016 [42] also created the same Puttygut model populated with calculated porosity and permeability values. It can be seen that areas of higher porosity and permeability also correlate well with the chosen neurons from the SOM created from this workflow.
After repeating this study at the Ira and Lenox reefs, it became apparent that there was a relationship between certain attribute contribution percentages and rock properties. While there is no case between the attribute contribution percentages that are exactly the same, there is most definitely a consistent trend/relationship when looking solely at the highest contributing attribute combinations. For the case of >10% porosity and >50 mD permeability, CWT spectral decomposition at 81 Hz, average frequency, and CWT spectral decomposition at 29 Hz were always the highest contributors by a significant margin. In the second case where there is a range of 6% to >10% porosity but no noticeable effect on permeability, every attribute significantly contributes except for the cosine of the instantaneous phase. Finally, in the case where a wider range of 3% to >10% porosity are found, CWT spectral decomposition at 29 Hz, average frequency, and CWT spectral decomposition at 57 Hz are found to be the most prevalent contributors. From the results shown, it seems that the workflow created at Puttygut was able to translate over the Ira and Lenox wells. The similar attribute blends being found at similar geological locations within each reef increase the confidence that this workflow was able to accurately distinguish areas of higher porosity and potentially permeability.
In conclusion, this described SOM set up and search for neurons where the attribute contributions are heavily skewed towards CWT spectral decomposition at 81 Hz, average frequency, and CWT spectral decomposition at 29 Hz, provide the best quantification of permeability and porosity. While this workflow was only tested in a localized area, it is an initial demonstration on how the multi-dimensional, multi-attribute analysis can be applied to reefs. Reefs are formed by different processes depending on where they are geographically located. As such, different attributes or workflow modifications may be needed as this method is applied to other reefs around the world.

Supplementary Materials

The following are available online at https://0-www-mdpi-com.brum.beds.ac.uk/2076-3263/11/3/142/s1. Figure S1: L-plot at the P-201 well. Figure S2: Inline 91 versus the 8–110 Hz corridor stack and crossline 66 versus the 8–110 Hz corridor stack. Figure S3: OVSP overlain on crossline 66.

Author Contributions

Conceptualization, C.B. and H.B.; methodology, C.B. and H.B.; software, H.B.; geologic investigation, M.R.; formal analysis, C.B.; investigation, C.B.; resources, H.B.; data curation, M.R.; writing—original draft preparation, C.B.; writing—review and editing, H.B., M.R., and J.P.; supervision, H.B.; funding acquisition, H.B. All authors have read and agreed to the published version of the manuscript.

Funding

Bedle’s startup funds helped support this project and were provided by the University of Oklahoma.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Restrictions apply to the availability of these data. Data were obtained from Consumer’s Energy and are available from Matt Rine with the permission of Consumer’s Energy.

Acknowledgments

I would like to give a big thank you to Julian Chenin for helping me get started in the machine learning world. I would like to thank the SDA at the University of Oklahoma for being a great and friendly community that helped move this project forward. I would also like to acknowledge the AASPI consortium at the University of Oklahoma, Geophysical Insights, and Schlumberger who all provided the software that was used for this work.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Taner, M.T.; Sheriff, R.E. Application of amplitude, frequency, and other attributes to stratigraphic and hydrocarbon determination. In Applications to Hydrocarbon Exploration; Payton, C.E., Ed.; American Association of Petroleum Geologists: Tulsa, OK, USA, 1977; pp. 301–327. [Google Scholar]
  2. Brown, A.R. Interpretation of Three-Dimensional Seismic Data; Society of Exploration Geophysicists and American Association of Petroleum Geologists: Tulsa, OK, USA, 2011; p. 247. [Google Scholar] [CrossRef] [Green Version]
  3. Pigott, J.D.; Kang, M.-H.; Han, H.-C. First Order Seismic Attributes for Clastic Seismic Facies Inter-pretation: Examples from the East China Sea. J. Asian Earth Sci. 2013, 10, 34–54. [Google Scholar] [CrossRef]
  4. Sears, S.O.; Lucia, F.J. Reef-growth model for Silurian pinnacle reefs, northern Michigan reef trend. Geology 1979, 7, 299–302. [Google Scholar] [CrossRef]
  5. Catacosinos, P.; Daniels, P.; Harrison, W. Structure, stratigraphy, and petroleum geology of the Michigan Basin; Interior cratonic basins. AAPG Mem. 1991, 51, 561–601. [Google Scholar]
  6. Michigan Department of Environmental Quality. Michigan Office of Oil, Gas, and Minerals. Oil and Gas Miner Database. 2020. Available online: http://www.michigan.gov/deq (accessed on 1 February 2019).
  7. Brock, T.; Elenbaas, D.; Schmidt, G. Secondary recovery in Michigan reefs summary results, recovery achieved in engineering applications for Michigan: Society of Petroleum Engineers. In Proceedings of the Society of Petroleum Engineers—Michigan Technical Conference, Mt. Pleasant, MI, USA, 11–12 April 1995; p. 30. [Google Scholar]
  8. Rine, M.J.; McLaughlin, P.I.; Bancroft, A.M.; Harrison, W.B.; Kuglitsch, J.; Caruthers, A.H.; Ramezani, J.; Kaczmarek, S.E.; Emsbo, P. Linked Silurian carbon cycle perturbations, bursts of pinnacle reef growth, extreme sea-level oscillations, and evaporite deposition (Michigan Basin, USA). Palaeogeogr. Palaeoclim. Palaeoecol. 2020, 554, 109806. [Google Scholar] [CrossRef]
  9. Rine, M.J.; Garrett, J.D.; Kaczmarek, S.E. A New Facies Architecture Model for the Silurian Niagaran Pinnacle Reef Complexes of The Michigan Basin. SEPM Spec. Publ. 2017, 109, 70–86. [Google Scholar]
  10. Sheriff, R.E. Seismic resolution: A key element, Geophysical Corner. AAPG Explor. 1997, 18, 44–51. [Google Scholar]
  11. Roden, R.; Smith, T.A.; Santogrossi, P.; Sacrey, D.; Jones, G. Seismic interpretation below tuning with multi-attribute analysis. Lead. Edge 2017, 36, 330–339. [Google Scholar] [CrossRef]
  12. Mesolella, K.; Robinson, J.; McCormick, L.; Ormiston, A. Cyclic deposition of Silurian carbonates and evaporites in the Michigan Basin. Am. Assoc. Pet. Geol. Bull. 1974, 58, 34–62. [Google Scholar]
  13. Laudon, C.; Stanley, S.; Santogrossi, P. Machine Learning Applied to 3-D Seismic Data from the Denver-Julesburg Basin Improves Stratigraphic Resolution in the Niobrara. In Proceedings of the 7th Unconventional Resources Technology Conference, Brisbane, Australia, 18–19 November 2019; p. 337. [Google Scholar]
  14. Roden, R.; Sacrey, D. Seismic Pattern Recognition in Shale Resource Plays; Hart Energy: Houston, TX, USA, 2015; Available online: https://www.geoinsights.com/wp-content/uploads/2020/04/Eand-P-Seismic-Pattern-Recognition-in-Shale-Resource-Plays-1.pdf (accessed on 31 December 2020).
  15. Santogrossi, P. Classification/Corroboration of Facies Architecture in the Eagle Ford Group: A Case Study in Thin Bed Resolution. In Proceedings of the 4th Unconventional Resources Technology Conference, Austin, TX, USA, 24–26 July 2017; p. 2696775. [Google Scholar]
  16. Roden, R.; Chen, C.W. Interpretation of DHI characteristics with machine learning. First Break 2017, 35, 55–63. [Google Scholar] [CrossRef]
  17. Roden, R.; Santogrossi, P. Significant advancements in seismic reservoir characterization with machine learning. First 2017, 3, 14–19. [Google Scholar]
  18. Roden, R.; Smith, T.; Sacrey, D. Geologic pattern recognition from seismic attributes: Principal component analysis and self-organizing maps. Interpretation 2015, 3, 59–83. [Google Scholar] [CrossRef]
  19. Kohonen, T. Self-Organizing Maps; Springer: Berlin/Heidelberg, Germany, 1995. [Google Scholar]
  20. Chopra, S.; Marfurt, K.J. Seismic Attributes for Prospect Identification and Reservoir Characterization; Society of Exploration Geophysicists/European Association of Geoscientists and Engineers: Tulsa, OK, USA; Houten, The Netherlands, 2007. [Google Scholar]
  21. Skirius, C.; Nissen, S.; Haskell, N.; Marfurt, K.; Hadley, S.; Ternes, D.; Michel, K.; Reglar, I.; D’Amico, D.; Deliencourt, F.; et al. 3-D seismic attributes applied to carbonates. Lead. Edge 1999, 18, 384–393. [Google Scholar] [CrossRef]
  22. Schot, S. Aberrancy: Geometry of the third derivative. Math. Mag. 1978, 51, 259–275. [Google Scholar] [CrossRef]
  23. Qi, X.; Marfurt, K.; Popovici, A.M.; Fomel, S. Volumetric aberrancy to map subtle faults and flexures. Interpretation 2018, 6, 349–365. [Google Scholar] [CrossRef]
  24. Chopra, S.; Marfurt, K.J. Volumetric curvature attributes add value to 3D seismic data interpretation. Lead. Edge 2007, 26, 856–867. [Google Scholar] [CrossRef]
  25. Duan, Y.; Tao, J.; Zeng, L.; Bi, M.; Feng, B. Using Curvature Attributes to Identify a Carbonate Paleokarst System in the Ordos Basin. In Proceedings of the 72nd EAGE Conference and Exhibition Incorporating SPE EUROPEC, Barcelon, Spain, 14–17 June 2010. [Google Scholar]
  26. Wallet, B.C. Attribute expression of channel forms in a hybrid carbonate turbidite formation. Interpretation 2016, 4, 75–86. [Google Scholar] [CrossRef] [Green Version]
  27. Eichkitz, C.G.; Davies, J.; Amtmann, J.; Schreilechner, M.G.; De Groot, P. Grey level co-occurrence matrix and its application to seismic data. First Break 2015, 33, 71–77. [Google Scholar] [CrossRef]
  28. Eichkitz, C.G.; Schreilechner, M.G.; de Groot, P.; Amtmann, J. Mapping directional variations in seismic character using GLCM based attributes. Interpretation 2015, 3, 13–23. [Google Scholar] [CrossRef]
  29. Barnes, A. Handbook of Poststack Seismic Attributes; Society of Exploration Geophysicists: Tulsa, OK, USA, 2016; pp. 47–48. [Google Scholar]
  30. Sukmono, S.; Santoso, D.; Samodra, A.; Waluyo, W.; Tjiptoharsono, S. Integrating seismic attributes for reservoir characterization in Melandong Field, Indonesia. Lead. Edge 2006, 25, 532–538. [Google Scholar] [CrossRef]
  31. Huang, Y.-P.; Geng, J.-H.; Zhong, G.-F.; Guo, T.-L.; Pu, Y.; Ding, K.-Y.; Ma, J.-Q. Seismic attribute extraction based on HHT and its application in a marine carbonate area. Appl. Geophys. 2011, 8, 125–133. [Google Scholar] [CrossRef]
  32. Toelle, B.; Ganshin, Y.V. Porosity Characterization in a Silurian Reef, Northern Michigan Basin, Using Azimuthal Seismic Data and Potential Impacts for Enhanced Oil Recovery; Special Paper 531; The Geological Society of America: Boulder, CO, USA, 2018; pp. 157–173. [Google Scholar]
  33. Sarhan, M.A. The efficiency of seismic attributes to differentiate between massive and non-massive carbonate successions for hydrocarbon exploration activity. NRIAG J. Astron. Geophys. 2017, 6, 311–325. [Google Scholar] [CrossRef] [Green Version]
  34. Barnes, A.E. Weighted average seismic attributes. Geophysiscs 2000, 65, 275–285. [Google Scholar] [CrossRef]
  35. Nejad, M.R.S.; Hassani, H.; Khorasani, M.M.; Javaherian, A.; Sokooti, M.R. The investigation of the spectral decomposition application in detecting reef reservoir on Abadan Plain, Iran. Aust. J. Basin Appl. Sci. 2009, 3, 866–874. [Google Scholar]
  36. Saadatinejad, M.R.; Javaherian, A.; Sarkarinejad, K. Investigation of the Various Spectral Decomposition Methods to Detect and Explore Hidden Complex Reef Reservoir Structures and Their Hydrocarbon Potentials in Northwestern Part of the Persian Gulf. Energy Explor. Exploit. 2012, 30, 867–887. [Google Scholar] [CrossRef]
  37. Goloshubin, G.; Dmitry, S.; Vjacheslav, V.; Gleb, T.; Monir, L. Reservoir Permeability from Seismic Attribute Analysis. Lead. Edge 2008, 27, 376–381. [Google Scholar] [CrossRef] [Green Version]
  38. Kozlov, E. Seismic signature of a permeable, dual-porosity layer. Geophysics 2007, 72, 281–291. [Google Scholar] [CrossRef]
  39. Pramanik, A.G.; Singh, V.; Vig, R.; Srivastava, A.K.; Tiwary, D.N. Estimation of effective porosity using geostatistics and multiattribute transforms: A case study. Geophysics 2004, 69, 352–372. [Google Scholar] [CrossRef]
  40. Iturrarán-Viveros, U. Smooth regression to estimate effective porosity using seismic attributes. J. Appl. Geophys. 2012, 76, 1–12. [Google Scholar] [CrossRef]
  41. Schuelke, J.S.; Quirein, J.A.; Sag, J.F.; Altany, D.A.; Hunt, P.E. Reservoir architecture and porosity distribution, Pegasus Field, West Texas—An integrated sequence stratigraphic-seismic attribute study using neural networks. Geo-Triad’98 1998, 142–145. [Google Scholar] [CrossRef]
  42. Garrett, J. A Sequence Stratigraphic Model for the Silurian A-1 Carbonate and Modeling of a Niagara-Lower Salina Reef Complex, Michigan Basin, USA. Master’s Thesis, Western Michigan University, Kalamazoo, MI, USA, 2016. [Google Scholar]
Figure 1. Image (A) is a Google Earth image of Michigan. The yellow star is the area of interest. Image (B) is a zoom in of the yellow star area showing the three reefs used in this study, Puttygut, Ira, and Lenox. Image (C) is a further zoom in of the Puttygut reef showing the two-way time (TWT) map of the Guelph (Brown) Formation. An A-A’ line is drawn through the four wells used in this study. There are also all wells in and near the reef marked by the white and black dots. The black dots represent the wells that were cored. Image (D) is a generalized stratigraphic scheme of the reefs in the southern trend of the Michigan Basin from [8]. In this region, the Clinton formation is refered to as the Cordell formation.
Figure 1. Image (A) is a Google Earth image of Michigan. The yellow star is the area of interest. Image (B) is a zoom in of the yellow star area showing the three reefs used in this study, Puttygut, Ira, and Lenox. Image (C) is a further zoom in of the Puttygut reef showing the two-way time (TWT) map of the Guelph (Brown) Formation. An A-A’ line is drawn through the four wells used in this study. There are also all wells in and near the reef marked by the white and black dots. The black dots represent the wells that were cored. Image (D) is a generalized stratigraphic scheme of the reefs in the southern trend of the Michigan Basin from [8]. In this region, the Clinton formation is refered to as the Cordell formation.
Geosciences 11 00142 g001
Figure 2. Correlation plot created from the vertical seismic profile (VSP) at well P-201 (see Figure 1C for location) through crossline 66 (north-south direction).
Figure 2. Correlation plot created from the vertical seismic profile (VSP) at well P-201 (see Figure 1C for location) through crossline 66 (north-south direction).
Geosciences 11 00142 g002
Figure 3. Acquired and processed VSP results overlain onto crossline 66 (north-south directions) of the Puttygut seismic volume.
Figure 3. Acquired and processed VSP results overlain onto crossline 66 (north-south directions) of the Puttygut seismic volume.
Geosciences 11 00142 g003
Figure 4. Top image displays inline 104 (east-west direction) without any horizon interpretations but with wells, well tops, and synthetic seismogram displayed. Bottom image has horizon picks added. Warmer colors indicate a negative impedance contrast, while cooler colors indicate a positive impedance contrast. The top image is the inline without interpretation, the bottom is with the horizon interpretation. The A-2 Anhydrite, A-1 Carbonate, Brown, and Gray Formations are very difficult to resolve and require a great deal of apriori geologic and structural knowledge of the area to pick accurately. The A-2 Carbonate and Clinton horizons are much easier to resolve. The A-2 Salt is shown on the flanks of the reef.
Figure 4. Top image displays inline 104 (east-west direction) without any horizon interpretations but with wells, well tops, and synthetic seismogram displayed. Bottom image has horizon picks added. Warmer colors indicate a negative impedance contrast, while cooler colors indicate a positive impedance contrast. The top image is the inline without interpretation, the bottom is with the horizon interpretation. The A-2 Anhydrite, A-1 Carbonate, Brown, and Gray Formations are very difficult to resolve and require a great deal of apriori geologic and structural knowledge of the area to pick accurately. The A-2 Carbonate and Clinton horizons are much easier to resolve. The A-2 Salt is shown on the flanks of the reef.
Geosciences 11 00142 g004
Figure 5. All images, excluding (L), were extracted onto the Brown horizon. (A) TWT map. (B) Pre-stack time migrated (PSTM) amplitude. (C1) Energy ratio similarity. (C2) Interpreted energy ratio similarity. (D1) Negative curvature. (D2) Interpreted negative curvature. (E1) Aberrancy. (E2) Interpreted aberrancy. (F1) gray level co-occurrence matrix (GLCM) homogeneity. (F2) Interpreted GLCM homogeneity. (G1) Average frequency. (G2) Interpreted average frequency. (H) Cosine of instantaneous phase. (I1) Continuous wavelet transform (CWT) spectral decomposition at 29 Hz. (I2) Interpreted CWT spectral decomposition at 29 Hz. (J1) CWT spectral decomposition at 57 Hz. (J2) Interpreted CWT spectral decomposition at 57 Hz. (K1) CWT spectral decomposition at 81 Hz. (K2) Interpreted CWT spectral decomposition at 81 Hz. (L) RGB (red-green-blue) blending of the three CWT spectral decomposition volumes displayed at 430 ms where red is 29 Hz, green is 57 Hz, and blue is 81 Hz.
Figure 5. All images, excluding (L), were extracted onto the Brown horizon. (A) TWT map. (B) Pre-stack time migrated (PSTM) amplitude. (C1) Energy ratio similarity. (C2) Interpreted energy ratio similarity. (D1) Negative curvature. (D2) Interpreted negative curvature. (E1) Aberrancy. (E2) Interpreted aberrancy. (F1) gray level co-occurrence matrix (GLCM) homogeneity. (F2) Interpreted GLCM homogeneity. (G1) Average frequency. (G2) Interpreted average frequency. (H) Cosine of instantaneous phase. (I1) Continuous wavelet transform (CWT) spectral decomposition at 29 Hz. (I2) Interpreted CWT spectral decomposition at 29 Hz. (J1) CWT spectral decomposition at 57 Hz. (J2) Interpreted CWT spectral decomposition at 57 Hz. (K1) CWT spectral decomposition at 81 Hz. (K2) Interpreted CWT spectral decomposition at 81 Hz. (L) RGB (red-green-blue) blending of the three CWT spectral decomposition volumes displayed at 430 ms where red is 29 Hz, green is 57 Hz, and blue is 81 Hz.
Geosciences 11 00142 g005
Figure 6. Self-organizing map (SOM) comparison between the second and third SOMs. The second SOM, on the left, was created with more input data, while the third SOM, on the right, was created with less. The white circle and rectangle are showing the same area of each SOM being characterized by different amounts of neurons. The SOM on the left created with more input data tended to use less neurons to characterize the areas. The arrows are pointing to places where the general trends were similar, with only slight differences in the amounts of neurons used for characterization. The third SOM was used for the remainder of the analysis.
Figure 6. Self-organizing map (SOM) comparison between the second and third SOMs. The second SOM, on the left, was created with more input data, while the third SOM, on the right, was created with less. The white circle and rectangle are showing the same area of each SOM being characterized by different amounts of neurons. The SOM on the left created with more input data tended to use less neurons to characterize the areas. The arrows are pointing to places where the general trends were similar, with only slight differences in the amounts of neurons used for characterization. The third SOM was used for the remainder of the analysis.
Geosciences 11 00142 g006
Figure 7. The top image shows the frequency SOM result with a random color bar with the A-1 Carbonate horizon, the Brown horizon, the Gray horizon, and the amplitude seismic wiggle traces overlain. The SOM results output the data into voxel sizes that were 55 by 55 feet laterally and 1 ms vertically. The bottom image shows the same SOM with a different color scale where the four neurons of interest are blue, while the rest are grey. The permeability (left) and porosity (right) logs are also now displayed. The permeability logs have a logarithmic scale from 1–1000 mD. The porosity logs have a linear scale from 1–20%.
Figure 7. The top image shows the frequency SOM result with a random color bar with the A-1 Carbonate horizon, the Brown horizon, the Gray horizon, and the amplitude seismic wiggle traces overlain. The SOM results output the data into voxel sizes that were 55 by 55 feet laterally and 1 ms vertically. The bottom image shows the same SOM with a different color scale where the four neurons of interest are blue, while the rest are grey. The permeability (left) and porosity (right) logs are also now displayed. The permeability logs have a logarithmic scale from 1–1000 mD. The porosity logs have a linear scale from 1–20%.
Geosciences 11 00142 g007
Figure 8. (A) shows the frequency SOM results for the Lenox reef at inline 56. The reef zone is dominated by neurons 23 and 24. (B,C) are images are of the Ira reef frequency SOM results at crossline 58. (B) is with all the neurons displayed, while (C) is only displaying the neurons of interest. The five vertical white lines are wells on crossline 58. Wells I-201, I-111, and I-108 are the top three ranked wells in the field, while wells I-104 and I202 were ranked ninth and twelfth out of twelve.
Figure 8. (A) shows the frequency SOM results for the Lenox reef at inline 56. The reef zone is dominated by neurons 23 and 24. (B,C) are images are of the Ira reef frequency SOM results at crossline 58. (B) is with all the neurons displayed, while (C) is only displaying the neurons of interest. The five vertical white lines are wells on crossline 58. Wells I-201, I-111, and I-108 are the top three ranked wells in the field, while wells I-104 and I202 were ranked ninth and twelfth out of twelve.
Geosciences 11 00142 g008
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Buist, C.; Bedle, H.; Rine, M.; Pigott, J. Enhancing Paleoreef Reservoir Characterization through Machine Learning and Multi-Attribute Seismic Analysis: Silurian Reef Examples from the Michigan Basin. Geosciences 2021, 11, 142. https://0-doi-org.brum.beds.ac.uk/10.3390/geosciences11030142

AMA Style

Buist C, Bedle H, Rine M, Pigott J. Enhancing Paleoreef Reservoir Characterization through Machine Learning and Multi-Attribute Seismic Analysis: Silurian Reef Examples from the Michigan Basin. Geosciences. 2021; 11(3):142. https://0-doi-org.brum.beds.ac.uk/10.3390/geosciences11030142

Chicago/Turabian Style

Buist, Carl, Heather Bedle, Matthew Rine, and John Pigott. 2021. "Enhancing Paleoreef Reservoir Characterization through Machine Learning and Multi-Attribute Seismic Analysis: Silurian Reef Examples from the Michigan Basin" Geosciences 11, no. 3: 142. https://0-doi-org.brum.beds.ac.uk/10.3390/geosciences11030142

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop