Supervised Classification
On this page
Summary
Supervised classification learns the relation between the logs and core-described facies from labelled depths and applies it to the rest of the well and to other wells. It is the right method when core description exists and the facies are geological, not just log-defined. The result has to be tested by Depth-block cross-validation, with whole depth blocks held out, because neighbouring samples are not independent and a random split overstates the skill.
Inputs and outputs
| Item | Units | |
|---|---|---|
| Input | Feature matrix of normalized logs at the depths of core-described facies, plus the core facies as labels | mixed |
| Parameter | classifier settings (here the number of neighbours k = 5), class weights | none |
| Output | Facies code for every sample in the wells to be classified, with a Confusion matrix from held-out depth blocks | none |
Equations
The training set is a set of pairs \((x_i, y_i)\), with \(x_i\) the standardised log vector at depth \(i\) and \(y_i\) the facies of the core description at that depth. The example uses a k-nearest-neighbour rule: the facies of a new sample is the one with the largest total weight among its \(k\) nearest training samples, with the weights \(w_c\) for class \(c\):
\(N_k(x)\) is the set of the \(k\) training samples nearest to \(x\). The standardisation is computed from the training rows only. With weights equal to 1 this is the plain rule; setting \(w_c\) in inverse proportion to the number of training samples of class \(c\) is the simplest remedy for Class imbalance. Decision trees, random forests, support vector machines and neural networks are the other usual choices; they differ in the model but the data handling is the same.
Skill is judged on data the model has not seen. In depth-block cross-validation the samples are split into \(F\) contiguous blocks, and each block in turn is predicted by a model trained on the other \(F-1\) blocks. The predictions are pooled into a Confusion matrix \(C\), with the measures
Single-value calculator
No calculator: this method is a procedure, not a single equation. The worked example below runs the whole procedure on synthetic logs.
Behavior
The synthetic core set has 961 samples from 60 beds: 325 shale, 383 sandstone, 214 limestone and only 39 dolomite. Shuffling the samples into five folds gives an accuracy of 0.971 and a balanced accuracy of 0.914, because the neighbours of every test sample, from the same bed, are in the training set. Holding out five contiguous depth blocks gives 0.912 and 0.780. The gap of six points in accuracy and 13 in balanced accuracy is leakage, and only the second figure estimates performance in a new well. The block confusion matrix shows where it fails: dolomite has a recall of 0.333 with 26 of 39 samples called limestone, while shale and sandstone are above 0.94. Weighting the votes in inverse proportion to the class frequency raises the dolomite recall to 0.538, but the overall accuracy falls from 0.912 to 0.865 and the balanced accuracy does not improve (0.775 against 0.780). Weighting moves errors between classes and does not create information that the logs do not hold.
Parameter guidance
Labels. Core description at a consistent scale, depth-shifted to the logs, and a facies scheme of a few classes that the logs can in principle separate. Merge classes that the core describers separate by texture alone. Features. The normalized logs that respond to the facies, and derived curves such as the neutron-density separation or the logarithm of resistivity; see Curve Selection and Missing Data. Add context with caution: depth, or the logs of the neighbouring samples, makes the classifier smoother and also easier to overfit. Validation. Hold out whole depth blocks, or whole wells. A leave-one-well-out test is the best estimate of performance in a new well. Choose settings (here k) on a split that is separate from the final test. Imbalance. Report recall by class and balanced accuracy, not only overall accuracy. Class weights, resampling and gathering more examples of the rare facies are the remedies, in order of preference the last. Applying the model. Apply it only to wells whose logs have been normalized to the training wells, and mark as null the samples that fall far from every training sample.
Worked example
A k-nearest-neighbour classifier in numpy, trained on a synthetic core-described well whose beds have large bed-to-bed offsets, as real beds do. The same data are scored with a random five-fold split and with five depth blocks, the confusion matrix is built by hand, and class weights are tried:
import numpy as np
NAMES = ["Shale", "Sandstone", "Limestone", "Dolomite"]
# facies means: GR (API), RHOB (g/cm3), NPHI (v/v, limestone scale), Rt (ohm.m). Illustrative values.
MEAN = np.array([[105.0, 2.55, 0.30, 3.0],
[45.0, 2.34, 0.17, 12.0],
[28.0, 2.60, 0.09, 35.0],
[28.0, 2.72, 0.07, 45.0]])
def synthetic_well(n_beds=40, seed=7, probs=(0.35, 0.30, 0.25, 0.10), bed_scale=1.0):
"""Beds of 8 to 25 samples. Each bed has its own offset, shared by all its samples, and each
sample adds independent noise. bed_scale multiplies the bed offsets. Resistivity is simulated in log10."""
rng = np.random.default_rng(seed)
truth, logs = [], []
for _ in range(n_beds):
f = rng.choice(4, p=probs)
n = rng.integers(8, 26)
bed = rng.normal(0, np.array([8.0, 0.04, 0.025, 0.20]) * bed_scale)
gr = MEAN[f, 0] + bed[0] + rng.normal(0, 12.0, n)
rhob = MEAN[f, 1] + bed[1] + rng.normal(0, 0.04, n)
nphi = MEAN[f, 2] + bed[2] + rng.normal(0, 0.025, n)
rt = 10 ** (np.log10(MEAN[f, 3]) + bed[3] + rng.normal(0, 0.15, n))
truth += [f] * n
logs.append(np.column_stack([gr, rhob, nphi, rt]))
return np.vstack(logs), np.array(truth)
def features(logs):
"""GR, RHOB, NPHI and log10(Rt)."""
x = logs.copy()
x[:, 3] = np.log10(x[:, 3])
return x
def knn_predict(x_train, y_train, x_test, k=5, class_weight=None):
"""k-nearest-neighbour vote. Scaling uses the training rows only, so nothing leaks from the test rows.
class_weight (one value per class) multiplies each neighbour's vote."""
mu, sd = x_train.mean(axis=0), x_train.std(axis=0)
a, b = (x_train - mu) / sd, (x_test - mu) / sd
d = ((b[:, None, :] - a[None, :, :]) ** 2).sum(axis=2)
nn = np.argsort(d, axis=1)[:, :k]
w = np.ones(4) if class_weight is None else np.asarray(class_weight, dtype=float)
votes = np.zeros((len(b), 4))
for c in range(4):
votes[:, c] = ((y_train[nn] == c) * w[c]).sum(axis=1)
return votes.argmax(axis=1)
def confusion(truth, pred, n=4):
cm = np.zeros((n, n), dtype=int)
np.add.at(cm, (truth, pred), 1)
return cm
def cross_validate(x, y, folds, **kw):
"""folds is an array giving the fold number of each sample. Returns the pooled confusion matrix."""
pred = np.empty_like(y)
for f in np.unique(folds):
test = folds == f
pred[test] = knn_predict(x[~test], y[~test], x[test], **kw)
return confusion(y, pred)
def report(cm, names=NAMES):
recall = np.diag(cm) / cm.sum(axis=1)
return np.trace(cm) / cm.sum(), recall.mean(), recall
logs, y = synthetic_well(n_beds=60, seed=11, probs=(0.35, 0.35, 0.22, 0.08), bed_scale=2.0)
x = features(logs)
n = len(y)
print(f"{n} core-described samples; class counts {np.bincount(y).tolist()} ({', '.join(NAMES)})")
rng = np.random.default_rng(0)
random_folds = rng.permutation(n) % 5 # samples shuffled into 5 folds: neighbours leak
block_folds = np.arange(n) * 5 // n # 5 contiguous depth blocks: whole beds held out
for label, folds in (("random 5-fold", random_folds), ("depth-block 5-fold", block_folds)):
acc, bal, rec = report(cross_validate(x, y, folds))
print(f"{label:20s} accuracy {acc:.3f} balanced accuracy {bal:.3f}")
cm = cross_validate(x, y, block_folds)
print("\nconfusion matrix, depth-block CV (rows true, columns predicted):")
print(" " + " ".join(f"{nm[:5]:>6s}" for nm in NAMES))
for k, nm in enumerate(NAMES):
print(f" {nm:10s}" + " ".join(f"{v:6d}" for v in cm[k]))
acc, bal, rec = report(cm)
prec = np.diag(cm) / np.maximum(cm.sum(axis=0), 1)
f1 = 2 * prec * rec / (prec + rec)
for k, nm in enumerate(NAMES):
print(f" {nm:10s} recall {rec[k]:.3f} precision {prec[k]:.3f} F1 {f1[k]:.3f}")
# class imbalance: weight each vote by the inverse class frequency
w = 1.0 / np.bincount(y)
w = w / w.mean()
cm_w = cross_validate(x, y, block_folds, class_weight=w)
acc_w, bal_w, rec_w = report(cm_w)
print(f"\nwith inverse-frequency vote weights: accuracy {acc_w:.3f}, balanced accuracy {bal_w:.3f}")
print(f"dolomite recall {rec[3]:.3f} -> {rec_w[3]:.3f}")
Output
961 core-described samples; class counts [325, 383, 214, 39] (Shale, Sandstone, Limestone, Dolomite)
random 5-fold accuracy 0.971 balanced accuracy 0.914
depth-block 5-fold accuracy 0.912 balanced accuracy 0.780
confusion matrix, depth-block CV (rows true, columns predicted):
Shale Sands Limes Dolom
Shale 316 9 0 0
Sandstone 0 360 23 0
Limestone 0 11 187 16
Dolomite 0 0 26 13
Shale recall 0.972 precision 1.000 F1 0.986
Sandstone recall 0.940 precision 0.947 F1 0.944
Limestone recall 0.874 precision 0.792 F1 0.831
Dolomite recall 0.333 precision 0.448 F1 0.382
with inverse-frequency vote weights: accuracy 0.865, balanced accuracy 0.775
dolomite recall 0.333 -> 0.538
Assumptions and limitations
- The core description is correct and at the scale of the logs. Describers disagree, and a thin lamination in the core is not seen by a log of one foot resolution.
- The training depths represent the facies, the logs and the wells to which the model is applied. A model trained on one well and applied to a well with another tool or another formation is an extrapolation.
- The log vector of a facies is stable from one bed to the next. If it is not, the model learns the beds it saw, which is what a random split rewards.
- The classes are separable on the logs chosen. If two classes have the same log response, no classifier separates them, and the confusion matrix shows it.
- The facies are a fixed scheme. A new facies in a new well is assigned to the nearest known one, with no warning, unless the distance to the training data is checked.
QC checks
- The reported accuracy comes from held-out depth blocks or wells, never from a random split of the samples or from the training data.
- The confusion matrix, with recall and precision for every facies, is shown. A rare facies with a recall near zero is flagged, however good the overall accuracy.
- Predicted facies proportions are close to the core proportions. Shifts of proportion in the application wells are examined.
- The standardisation and any other preprocessing are fitted on the training data only. Statistics from the whole well leak information from the test blocks.
- The result is geologically coherent down the well. A flickering prediction suggests overfitting or noise, and a depth filter is a remedy.
- Samples far from every training sample are flagged. A confident prediction there has no basis.
Going Deeper
Supervised facies classification became a public benchmark in 2016, when a contest on a Kansas gas field with nine facies compared teams on held-out wells. A common finding of that kind of work is that the gain of one machine-learning model over another is small compared with the gain from careful features, depth blocking and attention to the rare classes. Predictions from a model are probabilities once the vote is normalised, and carrying them forward gives the confidence of every sample. The class-imbalance problem is also geological: the facies that matter most, thin pay sands or a dolomite stringer, are often the rarest. The model learns what the core describers saw, so its limit is the quality of the description. When no core is available, labels made from a deterministic rule table or from named clusters can serve, with the caution that the classifier then reproduces those rules and does not add knowledge of the rock.
References
- Roberts, D.R. et al., 2017. Cross-validation strategies for data with temporal, spatial, hierarchical, or phylogenetic structure. Ecography, 40(8), 913–929.
- Dubois, M.K., Bohling, G.C. and Chakrabarti, S., 2007. Comparison of four approaches to a rock facies classification problem. Computers & Geosciences, 33(5), 599–617.
- Hall, B., 2016. Facies classification using machine learning. The Leading Edge, 35(10), 906–909.
Python reference implementation
Python reference implementation
The Python reference implementation is available to registered users with a verified email address. Register or sign in to view it.