CamPetro

Curve Selection and Missing Data

On this page

Summary

A facies model is only as good as its inputs. Choosing the curves means picking logs that carry different information, dropping redundant ones, and deciding in advance what the model does when a curve is missing: return null (Null propagation), fill the gap, or use a smaller model. The curves must also be put on a common scale across wells before one model is applied to several of them.

Inputs and outputs

Item Units
Input Candidate curves for every well: Gamma ray, Bulk density, Neutron porosity, Compressional slowness, True formation resistivity, Photoelectric factor and derived curves such as Neutron-density separation mixed
Parameter the required curves, the missing-data strategy, the standardisation statistics none
Output The curve set, a standardised matrix for every well, and a Facies code that is null wherever a required curve is null none

Equations

Redundancy between two candidate curves is measured with the correlation coefficient of their values over the well:

\[ r_{ab} = \frac{\sum_i (a_i - \bar a)(b_i - \bar b)}{\sqrt{\sum_i (a_i - \bar a)^2}\,\sqrt{\sum_i (b_i - \bar b)^2}} \]

A pair with \(|r|\) near 1 carries one piece of information twice and doubles its weight in a distance-based method. Density porosity is a rewriting of bulk density, so \(|r| = 1\) and only one of them is kept.

Standardisation puts every curve on a common scale. With a mean \(\mu\) and standard deviation \(\sigma\) from the key well,

\[ z = \frac{x - \mu}{\sigma} \]

Resistivity enters as \(\log_{10}\). A well logged with a different tool is first matched to the key well, for example by shifting and scaling so that its 5th and 95th percentiles equal those of the key well, as in Percentile Picks:

\[ x' = \left(x - P_5^{B}\right)\frac{P_{95}^{A} - P_5^{A}}{P_{95}^{B} - P_5^{B}} + P_5^{A} \]

When a curve is missing at some depths, a distance-based method can use only the curves that are present. The squared distance from a sample to a class centre is averaged over the \(m\) available curves, so that samples with fewer curves stay comparable:

\[ d^{2} = \frac{1}{m}\sum_{j \in \text{present}} \left(z_j - \mu_j\right)^{2} \]

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

Five candidate curves in a synthetic well are gamma ray, bulk density, neutron porosity, log resistivity and density porosity. Density porosity and bulk density correlate at -1.00, so one is redundant. Gamma ray, neutron porosity and log resistivity are also strongly related, with correlations of 0.81 (gamma ray and neutron), -0.78 (gamma ray and resistivity) and -0.84 (neutron and resistivity), because here shale is high in gamma ray and neutron and low in resistivity. Standardisation across wells matters more than the choice among these: a nearest-centroid model trained on well A scores 0.740 on a well B whose density reads 0.10 g/cm³ high and whose gamma ray is scaled by 1.35, and 0.898 once well B has been matched to well A by percentiles. When the density is missing over the lower 40% of well B, the average over 30 simulated wells shows the accuracy inside the gap falling from 0.941 with all curves to 0.836 when only the other curves are used and 0.845 when the density is replaced by the well mean. The two fill-in strategies are about equal, because both lose the density information; the only safe alternative is a null result, which here covers 60% of the well.

Parameter guidance

Which logs. Use a log from each family that responds to a different property: gamma ray (clay and radioactive minerals), bulk density and photoelectric factor (matrix), neutron (hydrogen index), sonic (compaction and porosity) and resistivity (fluid and organic matter). Drop a curve that is an algebraic rewrite of another. Which transformation. Resistivity as log10, and photoelectric factor and gamma ray often as they are; test others on the key well. Required curves. List the curves without which the model must return null, and the curves that can be left out with a warning. A model that is trained on a four-curve set should not be run silently on three. Missing curves. In order of preference: acquire or repair the curve, for which see Washout Repair and Null Handling; run a reduced model that was trained and validated on the smaller curve set; return null. Filling with the well mean or a regression value hides the gap and gives an answer that looks complete. Standardisation across wells. Normalize the curves first, see Normalization, and then standardise with the key-well statistics, using the same values for every well. Check the result on the histograms of each curve, well by well, before the model is applied.

Worked example

Two synthetic wells, the second logged with a different tool. The page computes the correlation matrix of five candidate curves, scores a nearest-centroid model on the second well before and after percentile matching, and compares missing-curve strategies over 30 simulated wells:

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)


CURVES = ["GR", "RHOB", "NPHI", "logRt", "DPHI"]


def feature_matrix(logs):
    """GR, RHOB, NPHI, log10(Rt) and density porosity (limestone matrix), which is RHOB rewritten."""
    x = logs.copy()
    x[:, 3] = np.log10(x[:, 3])
    dphi = (2.71 - x[:, 1]) / 1.71
    return np.column_stack([x, dphi])


def nearest_centroid(z, centroids):
    """Assign each row to the nearest centroid using only the columns that are present in that row.
    The squared distance is averaged over the available columns, so rows with fewer curves stay comparable.
    A row with no curves at all gets -1 (no result)."""
    ok = ~np.isnan(z)
    diff = np.where(ok[:, None, :], z[:, None, :] - centroids[None, :, :], 0.0)
    d = (diff ** 2).sum(axis=2) / np.maximum(ok.sum(axis=1), 1)[:, None]
    lab = d.argmin(axis=1)
    lab[ok.sum(axis=1) == 0] = -1
    return lab


def tool_b(logs):
    """Well B was logged with a different tool: RHOB reads 0.10 high and GR is scaled by 1.35."""
    return logs * np.array([1.35, 1.0, 1.0, 1.0]) + np.array([0.0, 0.10, 0.0, 0.0])


logs_a, y_a = synthetic_well(n_beds=50, seed=21)
logs_b, y_b = synthetic_well(n_beds=50, seed=22)
logs_b = tool_b(logs_b)
xa, xb = feature_matrix(logs_a), feature_matrix(logs_b)

print("1. redundancy: correlation matrix of the five candidate curves in well A")
r = np.corrcoef(xa.T)
print("          " + " ".join(f"{c:>7s}" for c in CURVES))
for i, c in enumerate(CURVES):
    print(f"  {c:7s} " + " ".join(f"{v:7.2f}" for v in r[i]))

# model for steps 2 and 3: standardise with well A statistics, class centroids from the key well,
# using the four non-redundant curves (columns 0 to 3)
cols = [0, 1, 2, 3]
mu, sd = xa[:, cols].mean(axis=0), xa[:, cols].std(axis=0)
za = (xa[:, cols] - mu) / sd
cent = np.array([za[y_a == k].mean(axis=0) for k in range(4)])
p_a = np.percentile(xa[:, cols], [5, 95], axis=0)


def match_percentiles(x):
    """Shift and scale each curve so that its 5th and 95th percentiles equal those of the key well."""
    p = np.percentile(x, [5, 95], axis=0)
    return (x - p[0]) * (p_a[1] - p_a[0]) / (p[1] - p[0]) + p_a[0]


print("\n2. standardisation across wells: model trained on well A, applied to well B")
raw = nearest_centroid((xb[:, cols] - mu) / sd, cent)
fixed = nearest_centroid((match_percentiles(xb[:, cols]) - mu) / sd, cent)
print(f"  raw well B curves            accuracy {np.mean(raw == y_b):.3f}")
print(f"  percentile-matched curves    accuracy {np.mean(fixed == y_b):.3f}")

print("\n3. a missing curve: RHOB is absent over the lower 40% of the well (mean of 30 simulated wells B)")


def gap_test(seed, curve=1):
    lb, yb = synthetic_well(n_beds=50, seed=seed)
    z = (match_percentiles(feature_matrix(tool_b(lb))[:, cols]) - mu) / sd
    gap = np.arange(len(z)) >= int(0.6 * len(z))
    zg = z.copy()
    zg[gap, curve] = np.nan
    filled = zg.copy()
    filled[gap, curve] = 0.0                       # the well mean, which is zero once standardised
    out = {}
    for name, zz in (("all curves present", z), ("use the curves present", zg), ("fill with the mean", filled)):
        lab = nearest_centroid(zz, cent)
        out[name] = (np.mean(lab == yb), np.mean(lab[gap] == yb[gap]))
    return out


results = [gap_test(s) for s in range(100, 130)]
print(f"  {'strategy':24s} {'whole well':>10s} {'inside gap':>10s}")
for name in results[0]:
    a = np.mean([res[name][0] for res in results])
    g = np.mean([res[name][1] for res in results])
    print(f"  {name:24s} {a:10.3f} {g:10.3f}")
print("  null where any missing   no result (coverage 0.60) inside the gap")

Output

1. redundancy: correlation matrix of the five candidate curves in well A
               GR    RHOB    NPHI   logRt    DPHI
  GR         1.00    0.12    0.81   -0.78   -0.12
  RHOB       0.12    1.00   -0.11    0.10   -1.00
  NPHI       0.81   -0.11    1.00   -0.84    0.11
  logRt     -0.78    0.10   -0.84    1.00   -0.10
  DPHI      -0.12   -1.00    0.11   -0.10    1.00

2. standardisation across wells: model trained on well A, applied to well B
  raw well B curves            accuracy 0.740
  percentile-matched curves    accuracy 0.898

3. a missing curve: RHOB is absent over the lower 40% of the well (mean of 30 simulated wells B)
  strategy                 whole well inside gap
  all curves present            0.937      0.941
  use the curves present        0.895      0.836
  fill with the mean            0.899      0.845
  null where any missing   no result (coverage 0.60) inside the gap

Assumptions and limitations

  • The logs differ between wells because of the tools, the borehole and the calibration, and not because the rock is different. Normalization removes the first and keeps the second. If the rock differs, matching the percentiles removes the real signal.
  • The key well is representative. Its statistics define the scale for all the other wells, and an unusual key well biases every result.
  • The correlation matrix describes redundancy in this well. A pair that is redundant in one lithology may not be in another.
  • A missing curve is missing at random with respect to the facies. If the density is missing exactly in washouts, which are in shale, the missing data themselves carry information, and a null in the facies track has a geology to it.
  • A reduced model, or a fill value, is validated. The accuracy of a model trained with all curves and run on fewer curves is not known until it has been tested.

QC checks

  • The facies track is null where, and only where, a required curve is null. Nothing is computed from a placeholder value such as -999.25, and no gap has been filled in silence.
  • The curve statistics, mean and standard deviation, are the same in every well that uses the model, after normalization, and the histograms of the same curve in different wells overlap.
  • The proportion of null results by well is reported. A well that is 40% null has a different quality of result from one that is complete.
  • No two curves of the set have a correlation above about 0.9 in absolute value unless that is deliberate.
  • Removing one curve at a time from the set and rescoring against core shows which curves carry the result. A curve whose removal changes nothing is a candidate to drop.

Going Deeper

Curve selection looks like a technical detail and decides most of the outcome. Sets with a gamma ray, a density, a neutron and a resistivity are the working minimum of facies work. Redundancy is not harmful for tree-based classifiers in the way it is for distance-based ones, but it makes the importance scores hard to read. Null propagation is a matter of principle: a number computed from a filled value looks the same as one computed from data, and an interpreter further down the workflow cannot tell them apart, so the null is the honest result. Multi-well work adds a second problem, which is that tool differences between vintages are larger than the facies differences the model is supposed to see. Normalization is therefore a precondition, and a facies model that works in one well and fails in the next is more often a normalization problem than a model problem. Some workflows add the missing-data pattern itself as a feature, or train one model per curve subset.

References

References will be added once verified.

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.