Created by David Briguglio

The ability to process, manipulate and draw insight from large datasets has always been necessary in geosciences, however, as data velocity, variety and volume in the field continues to increase, the skills of a Data Scientist are increasingly relevant. My objective here is to show that a combination of domain expertise and some basic coding ability can result in a powerful set of skills that will greatly enhance geo-scientific data analysis.

In this post, I’ll demonstrate a method for combining dimensionality reduction with unsupervised machine learning to attempt to out-perform lithology classifications from cuttings descriptions. Cuttings descriptions are important for geoscientists as they often provide the only continuous representation of real rock information from a well or drill-hole. However, these are notoriously poor data sources which can suffer from error such as caving. Furthermore, the interpretation of cuttings are highly subjective and can differ between well-site geologists. Conventional drill-core offers the most accurate representation of the rock record, however, it is often incomplete or not collected at all.

Wireline log data is a powerful, high dimension dataset that offers almost continuous representation of lithofacies and fluid properties in the sub-surface. If collected properly it offers a consistent and high resolution dataset, that is unbiased by subjective interpretation. In this post, I’m going to show how wireline data can be used to classify lithofacies, without prior interpretation.

Datasets 

A LAS file containing wireline log data, csv of cuttings descriptions and the well completion report for the Roc-1 exploration well were all downloaded from the Geoscience Australia’s open source petroleum well database; NOPIMS. Some minor editing of the csv cuttings descriptions and wireline data has been conducted so that they could be used in this project.

Import libraries 

Before getting started, we must import the required python libraries. All packages that I will use during this tutorial are off-the-shelf python implementations.

Clean up the wireline data

For this experiment we intend to use three logs found in the las file including; Gamma Ray (GR), TNPH and RHOZ. The Gamma Ray measures the natural radioactivity of a formation and offers a powerful method for discriminating lithofacies. Shales emit high gamma radiation because they often contain abundant potassium compared with compositionally mature sands or limestones that emit much lower gamma radiation. TNPH is a thermal neutron log that measures neutron porosity within a formation, and RHOZ provides a measure of formation density. Together all three logs provide ideal variables for distinguishing between different rock types in the sub-surface.

Before using the data, we have to make sure it is of high enough quality, otherwise any results are unlikely to make geological sense. This is where domain knowledge is important. We can plot a histogram of the Density values from the RHOZ log and look for outliers. As can be seen below, there are some anomalously high values that are unlikely to be real. Fortunately python makes it easy to clean up the data.

Filter the data

To de-spike the density log, we can easily replace values of over 2.9 g/cm3 with the average of all the density data. This will provide a clean version of the variable without the need to remove data.

The Roc-1 exploration well was drilled to a total depth of over 5000m. However, in order to target a zone that includes a potential reservoir and seal, we are going to filter the dataset down to all data deeper than 4200m.

Attach the cuttings interpretations to the wireline data 

In order to visualise the cuttings descriptions alongside the wireline logs, we need to do some data preparation. The cuttings interpretations have to be appended to the dataframe containing wireline data. However, the cuttings are discrete records with different sample depths. We can use the pandas method ‘merge_asof’ to merge the two datasets and assign the nearest depth value from the wireline data, to that of the cuttings. Then we can append these to the original wireline dataframe and either back fill or forward fill null cuttings entries so that we have a complete dataframe to work with.

Implement UMAP data dimensionality reduction 

Finally, we are in a position to implement a data reduction technique and visualise the interpretation of the cuttings descriptions in reduced dimensions. For this, we are going to apply a technique known as Uniform Manifold Approximation and Projection (UMAP) first described by McInnes et al. (2018). This is a powerful dimension reduction technique that constructs an initial high dimension graph layout of the data before re-projecting to a structurally consistent low dimension form (Coenen & Pearce n.d. Google PAIR). Unlike Principal Component Analysis (PCA), UMAP can handle non linear relationships, and it has also been shown to be computationally more efficient than t-SNE.

The two most important parameters for UMAP are ‘n_neighbors’ and ‘min_dist’ (UMAP Docs n.d.). The n_neighbors parameter prescribes the number of nearest neighbours of each data point that is used to construct the initial graph layout. This parameter controls the balance between global and local structure with low values allowing scrutiny of local neighbourhoods and higher values providing a representation of broad structure. The min_dist parameter prescribes the minimum distance between data points in low dimension space and therefore controls the degree of clustering observed in the final result.

A best practice approach would be to iterate on these parameters until an acceptable balance between broad and local structure was found, but for this example we will invoke the defaults.

Since we only have three logs that we intend to involve as variables, the dataset does not suffer from high dimensionality. However the relationship between these variables is complex, non-linear and it is always difficult to comprehend how multiple petrophysical responses can reflect changes to lithofacies or fluid content in the sub-surface.

Implement unsupervised learning

As can be seen in the figure above, the two dimensional UMAP space clearly distinguishes a series of groups that presumably have similar petrophysical properties. However, when coloured by Cuttings Interpretation, it seems that the three lithology classes are spread across the groups, implying that the well-site geologist has not captured these differences. This highlights a pit-fall of logging rock types based on colour and without the textural context that conventionally cut drill-core provides. It also shows that visual logging can fail to detect subtle differences in physical rock properties that might be useful during reservoir modelling, or rock physics and inversion studies.

For the next step in this post, we shall test whether an unsupervised machine learning method known as Spectral Clustering does a better job of representing each of the observed petrophysical domains that are implied by the UMAP representation.

Spectral Clustering works well for a range of different cluster shapes, such as those observed in the 2D UMAP visualisation shown above. This is because the algorithm uses a connectivity approach to classify data points, meaning that data are more likely to be grouped when they are either touching or directly adjacent to one another (Keerthana 2021, Analytics Vidhya). This differs to distance based clustering techniques such as k-means which rely on computing the distance of data points from randomly assigned centroids.

In the example below we use the embeddings calculated by the UMAP dimension reduction as input for the clustering algorithm.

Pairplot showing distribution of each modelled cluster.

As can be seen in the pairplot above, the clusters each capture unique petrophysical domains. Cluster 2 has gamma ray values of generally lower than 200, normally distributed RHOZ values of ~2.5-2.8 g/cm3 and neutron porosity characteristic of brine filled sediments buried to more than 3000m. Cluster 0 also has low GR, but comparatively higher TNPH indicating reasonable quality sands. Cluster 1 has gamma ray values of over 200, generally much higher RHOZ values and is probably representative of claystone. Perhaps the most interesting is Cluster 3 which is made up of the fewest data points and characterised by low gamma ray indicative of the most mature and clean sands in the dataset. However, it also hosts the lowest density and neutron porosity, suggesting that these sands may also host hydrocarbon.

Visualise the spectral clustering result in UMAP space

As can be seen below, by colouring the data points in UMAP space by the spectral clustering labels, it is clear that the main petrophysical domains have been captured. The method has out-performed the cuttings interpretation as far as identifying unique rock properties in the well.

Visualise a strip log to compare the cuttings interpretations with the predicted classes

For this final section, we shall visualise the cuttings interpretations and compare them with the predicted petrophysical classes generated by dimension reduction and spectral clustering. We will view these alongside a gamma ray log to provide some geological context.

We can write a function that can plot a log alongside a series of discrete classes and specify the start and end depths that we would like to inspect.

Down hole strip log showing gamma ray wireline log (left), conventional log of the cuttings (centre) and the log produced by the dimensionality reduction and clustering here (right). See below for legend.

Legend for the dimensionality reduction and clustering log.

Legend for the conventional cuttings log.

The predicted classes clearly provide more granularity and appear to be geologically reasonable. The cuttings interpretation imply that high gamma zones between ~4200-4400m are the same lithofacies as those between ~4500m -TD. However, the predicted classes discriminate between these, and this is likely to be more accurate since the gamma ray is generally much higher for the deeper interval. The predicted classes also seem to discriminate between sand end members. Cluster 3 is coloured in yellow and as mentioned earlier, these exhibit a potential hydrocarbon response from the logs.

Of course, the true test of this approach is to compare the results with real rock information. There is no conventional core for Roc-1, so we need to compare our results with available side wall core (SWC) information. Since the cuttings and predicted classes mostly agree for the reservoir interval (~4390-4500m), we aren’t going to investigate that section further. Instead, lets compare the results of the two main fine grained packages mentioned earlier with SWC photos and descriptions provided in the Roc-1 Well Completion Report.

SWC sample R11-C10, was recovered from 4299.12m. According to the photos and descriptions below it represents a Mudstone with rare mica and ‘red and black heavy minerals or organic matter.’ This is far from a Siltstone interpreted from the Cuttings descriptions. Sample R11-C27 was recovered from 4607.76m and apparently represents a black mudstone with some bioturbation and diagenetic pyrite. This might suggest an initially sub-oxic to oxic ocean floor depositional space that later deepened to become stratified. More importantly it is a different facies to the previous sample and validates the discriminations made by Dimension Reduction and Spectral clustering.

Comparison of the conventional cuttings log (left) with the new clustered log (right). Reviewing the available side wall cores from the ROC 1 borehole show a good agreement between the facies predicted from the clustering approach here, validating the method.