Skip to content
Open in molab View on GitHub Download .py

How to estimate the uncertainty of a slope for static calibration or regression

A calibration experiment gives you pairs of readings: you command an input U and the instrument reports I. Plotting I against U should give a straight line, and the slope of that line is the sensitivity of the instrument — the number you actually want, with an uncertainty attached.

Here is one such series. We will do three things, in order:

  1. plot the measurements,
  2. attach error bars and fit a straight line,
  3. find the steepest and shallowest lines the data allow.
import matplotlib.pyplot as plt
import numpy as np

# Six calibration points: commanded voltage against reported current.
U = np.array([4.5, 6.0, 7.0, 8.0, 9.0, 10.0])  # voltage, V
I = np.array([109.5, 115.3, 129.2, 134.8, 142.1, 144.6])  # current, mA

du = 0.15  # standard uncertainty of the voltage reading, V

1. The measurements

Six points. There is a clear upward trend, but they do not lie exactly on a line — that scatter is the measurement uncertainty we need to quantify.

fig1, ax1 = plt.subplots(figsize=(7, 4.5))
ax1.plot(U, I, "o", markersize=7, color="tab:blue", label="measurements")
ax1.set_xlabel("$U$ [V]")
ax1.set_ylabel("$I$ [mA]")
ax1.set_title("Calibration measurements")
ax1.grid(alpha=0.3)
ax1.legend()
fig1

2. Error bars and a least-squares fit

Each reading carries an uncertainty, drawn as an error bar. Drag the slider to change the standard uncertainty u(I) of the current measurement and watch the picture change.

np.polyfit(U, I, 1) returns the least-squares slope and intercept, which is our best estimate \(m\).

m_best, intercept = np.polyfit(U, I, 1)
print(f"least-squares slope   m = {m_best:.2f} mA/V")
print(f"intercept                 = {intercept:.1f} mA")

least-squares slope   m = 6.96 mA/V
intercept                 = 77.6 mA

u_I = mo.ui.slider(
    start=1.0,
    stop=6.0,
    step=0.25,
    value=4.5,
    label="standard uncertainty of the current reading, u(I) [mA]",
)
u_I

3. The steepest and shallowest acceptable lines

Rather than trusting a formula, ask a geometric question: what is the range of slopes for which some straight line still passes through every error bar?

The lines that just touch the error bars at one end give the two limits. Each one ignores most of the points — that is the point. Only the pair of points that binds the answer matters; the rest only have to be consistent.

def slope_limits(U, I, du, di):
    """Steepest and shallowest lines that cross every error bar.

    The line b + m*x crosses box i when it is above the bottom-right corner
    and below the top-left corner, i.e. for some b,

        (I_i - di) - m (U_i + du)  <=  b  <=  (I_i + di) - m (U_i - du)

    A b exists for all boxes only when every box is compatible with every
    other, which reduces to a bound on m per pair i, j:

        (I_i - di) - (I_j + di)  <=  m * [ (U_i + du) - (U_j - du) ]

    Taking the largest lower bound and the smallest upper bound over all
    pairs gives the interval of acceptable slopes. A box never constrains
    itself, so i == j is excluded.
    """
    n = len(U)
    i, j = np.meshgrid(np.arange(n), np.arange(n), indexing="ij")
    other = i != j
    gap = (U[i] + du) - (U[j] - du)
    rise = (I[i] - di) - (I[j] + di)
    safe = np.where(gap == 0, 1.0, gap)

    # gap > 0 forces m to be at least rise/gap; take the strongest such bound.
    m_lo = float(np.where(other & (gap > 0), rise / safe, -np.inf).max())
    # gap < 0 forces m to be at most rise/gap; take the strongest such bound.
    m_hi = float(np.where(other & (gap < 0), rise / safe, np.inf).min())

    feasible = m_lo <= m_hi
    return (m_lo, m_hi) if feasible else (np.nan, np.nan)

def touching_intercept(U, I, du, di, m):
    """Largest intercept for which the line of slope m still crosses every bar."""
    return float(np.max((I - di) - m * (U + du)))

m_min, m_max = slope_limits(U, I, du, u_I.value)
feasible = np.isfinite(m_min) and np.isfinite(m_max)

if feasible:
    b_min = touching_intercept(U, I, du, u_I.value, m_min)
    b_max = touching_intercept(U, I, du, u_I.value, m_max)
else:
    b_min = b_max = np.nan

fig2, ax2 = plt.subplots(figsize=(7.5, 4.8))
ax2.errorbar(
    U, I, xerr=du, yerr=u_I.value, fmt="o", markersize=7,
    capsize=4, color="tab:blue", label="measurements $\\pm u$",
)

if feasible:
    grid = np.linspace(U.min() - du, U.max() + du, 200)
    ax2.plot(grid, m_best * grid + intercept, "k-", lw=2,
            label=f"least squares, $m={m_best:.2f}$")
    ax2.plot(grid, m_min * grid + b_min, color="tab:green", ls="--", lw=2,
            label=f"shallowest, $m_{{\\min}}={m_min:.2f}$")
    ax2.plot(grid, m_max * grid + b_max, color="tab:red", ls="--", lw=2,
            label=f"steepest, $m_{{\\max}}={m_max:.2f}$")
else:
    ax2.plot([], [], " ")
    ax2.text(
        0.5, 0.06,
        "no single line crosses every error bar\nat this uncertainty",
        transform=ax2.transAxes, ha="center", va="bottom",
        color="tab:red", fontsize=11, weight="bold",
    )

ax2.set_xlabel("$U$ [V]")
ax2.set_ylabel("$I$ [mA]")
ax2.set_title(f"Slope uncertainty, $u(I)={u_I.value:.2f}$ mA")
ax2.grid(alpha=0.3)
ax2.legend(fontsize=9, loc="upper left")
fig2

if feasible:
    u_m = (m_max - m_min) / 2
    report = f"""
With $u(I) = {u_I.value:.2f}$ mA the acceptable slopes run from
$m_{{\\min}} = {m_min:.2f}$ to $m_{{\\max}} = {m_max:.2f}$ mA/V, so

$$\\Delta m = \\frac{{m_{{\\max}} - m_{{\\min}}}}{{2}} = {u_m:.2f}\\ \\text{{mA/V}}$$

and the slope is quoted as $m = {m_best:.1f} \\pm {u_m:.1f}$ mA/V.

Shrink the slider: below about 2.2 mA **no** straight line passes through every
error bar, which is the method correctly telling you that the uncertainty you
claimed is too small to explain your own scatter.
"""
else:
    report = f"""
At $u(I) = {u_I.value:.2f}$ mA there is **no** straight line that passes through
every error bar — the points are more scattered than the quoted uncertainty
allows, so the slope is not bracketed at all. Widen the error bars until a line
fits.
"""
mo.md(report)

With \(u(I) = 4.50\) mA the acceptable slopes run from \(m_{\min} = 5.39\) to \(m_{\max} = 8.48\) mA/V, so

\[\Delta m = \frac{m_{\max} - m_{\min}}{2} = 1.54\ \text{mA/V}\]

and the slope is quoted as \(m = 7.0 \pm 1.5\) mA/V.

Shrink the slider: below about 2.2 mA no straight line passes through every error bar, which is the method correctly telling you that the uncertainty you claimed is too small to explain your own scatter.

Where the bound comes from

Write the candidate line as \(I = b + m\,U\). Reading \(I\) at \(U\) means interpolating within a rectangle of half-widths \(d_U\) in voltage and \(d_I\) in current. The line enters the rectangle at its left edge and leaves at its right edge, so it crosses that rectangle exactly when

\[ I_i - d_I \;\le\; b + m\,(U_i + d_U) \qquad\text{and}\qquad b + m\,(U_i - d_U) \;\le\; I_i + d_I , \]

that is, when some intercept \(b\) satisfies

\[ (I_i - d_I) - m\,(U_i + d_U) \;\le\; b \;\le\; (I_i + d_I) - m\,(U_i - d_U). \]

A single line must satisfy this for every point at once, so the best any intercept can do must sit above every lower bound and below every upper bound. Pairing a lower bound from point \(i\) with an upper bound from point \(j\) gives

\[ (I_i - d_I) - (I_j + d_I) \;\le\; m\,\big[(U_i + d_U) - (U_j - d_U)\big]. \]

When the bracket is positive this is a lower bound on the slope; when it is negative it is an upper bound. Collecting every pair and keeping the strongest of each kind,

\[ m_{\min} = \max_{i \ne j} \frac{(I_i - d_I) - (I_j + d_I)}{(U_i + d_U) - (U_j - d_U)}, \qquad m_{\max} = \min_{i \ne j} \frac{(I_i - d_I) - (I_j + d_I)}{(U_i + d_U) - (U_j - d_U)} \]

over the pairs where the denominator has the required sign. A pair \(i = j\) is excluded, since a box never constrains itself.

Note which points survive. Only the pair that produces the largest lower bound actually sets \(m_{\min}\), and only the pair producing the smallest upper bound sets \(m_{\max}\) — every other point merely has to stay consistent. That is the "ignore some of the points" step: the extreme slopes are pinned down by two readings, not by a fit to all six.

Reporting

The interval is not a statistical confidence interval, so it is quoted as a half-width rather than as a standard uncertainty, and with the correct number of significant figures:

\[ m = m_{\text{best}} \pm \frac{m_{\max} - m_{\min}}{2}, \]

both terms rounded to the same decimal place — \(7.3 \pm 1.9\) mA/V, not \(7.28 \pm 1.874\) mA/V.

Two cautions

  • This construction assumes a single line really is the right model. If the scatter is curved rather than random, no line will cross every error bar at any width, and you have found a hysteresis or a nonlinearity instead of an uncertainty.
  • The bound widens as the error bars grow, so quoting an uncertainty makes the method agree with it. Judge it against a second, independent method — a weighted least-squares fit with its own residual estimate, or a Monte Carlo propagation.