Results (four‑parameter logistic fit)

Model fitted (increasing inhibition, log‑IC₅₀ parametrisation for well‑behaved errors):

$$y = \text{Bottom} + \frac{\text{Top}-\text{Bottom}}{1+10^{\,h\,(\log_{10}\text{IC}{50}-\log{10}x)}}$$

Unweighted least squares on all 16 points (duplicates kept as separate observations), n = 16, df = 12.

Fit quality: residual SD = 1.8 % inhibition, SSE = 39.2, R² = 0.998. CIs are Wald/t intervals (t₀.₉₇₅,₁₂ = 2.179); the IC₅₀ interval is obtained by back‑transforming the log₁₀IC₅₀ interval, so it is asymmetric (4.23–5.32) as it should be.

Sanity checks - Bottom and Top are statistically indistinguishable from 0 % and 100 %. Fixing them (2‑parameter Hill fit) gives IC₅₀ = 4.6 nM (4.32–4.95), h = 1.03 ± 0.03; extra‑sum‑of‑squares F(2,12) = 0.32, p ≈ 0.73 → the 4PL is not needed statistically, but its IC₅₀ is unchanged, which is reassuring. - Hill slope ≈ 1 ⇒ consistent with simple single‑site (1:1) competitive inhibition; no evidence of cooperativity or steep off‑target behaviour. - Both asymptotes are defined by real data (top plateau reached at 100–300 nM, bottom at 0.1–0.3 nM), so the IC₅₀ is not extrapolated — the usual failure mode of 4PL fits.

Figure

I can't write files to your machine or attach a PNG here, so below is a text preview of exactly what the script plots; running fit_4pl.py writes dose_response_4pl.png (300 dpi) and .pdf.

y (% inhibition)
100 |                                ···o
 95 |                            ··o·
 90 |                          o··
 85 |                         ·o
 80 |                       ··
 75 |                      ·
 70 |                     o
 65 |                    ·o
 60 |                   ·
 55 |
 50 |- - - - - - - - -·        <- 50% (midpoint) at IC50 = 4.7 nM
 45 |                ·
 40 |               o
 35 |              ·
 30 |             ·
 25 |            ·
 20 |           ·
 15 |          o
 10 |       ···
  5 |  o··o
  0 +o--------------------------------------
     0.1       1        10       100    300   concentration (nM, log scale)
   o = replicate data,  · = 4PL fit,  dashed = IC50 marker

The real figure additionally shows the 95 % confidence band of the fitted curve (delta method) and an annotation box with IC₅₀, its CI, the Hill slope and R².

Fitting script — save as fit_4pl.py

#!/usr/bin/env python3
"""
fit_4pl.py -- Four-parameter logistic (4PL) fit of a dose-response inhibition curve.

Outputs
-------
  dose_response_4pl.png / .pdf : data + fit + 95% band + IC50 marker
  fit_4pl_results.csv          : parameter estimates, SEs, 95% CIs
  console                      : full fit report, constrained fit, F-test

Requires: numpy, scipy, matplotlib
Run:      python fit_4pl.py
"""
import numpy as np
from scipy.optimize import curve_fit
from scipy import stats
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt

# ----------------------------------------------------------------------
# 1. Data (duplicates as separate observations)
# ----------------------------------------------------------------------
conc = np.array([0.1, 0.1, 0.3, 0.3, 1, 1, 3, 3,
                 10, 10, 30, 30, 100, 100, 300, 300], dtype=float)   # nM
resp = np.array([2, 4, 5, 7, 14, 18, 38, 42,
                 66, 70, 87, 89, 96, 97, 98, 99], dtype=float)       # % inhibition

# ----------------------------------------------------------------------
# 2. Model: fit log10(IC50) rather than IC50 (better conditioned, symmetric errors)
# ----------------------------------------------------------------------
def fourPL(x, bottom, top, log_ic50, hill):
    """Increasing 4PL: y(0)=bottom, y(inf)=top, y(IC50)=(bottom+top)/2."""
    return bottom + (top - bottom) / (1.0 + 10.0 ** (hill * (log_ic50 - np.log10(x))))

def twoPL(x, log_ic50, hill):          # constrained: bottom=0, top=100
    return fourPL(x, 0.0, 100.0, log_ic50, hill)

# ----------------------------------------------------------------------
# 3. Fit
# ----------------------------------------------------------------------
p0 = [0.0, 100.0, np.log10(5.0), 1.0]
popt, pcov = curve_fit(fourPL, conc, resp, p0=p0, maxfev=20000)
bottom, top, log_ic50, hill = popt

n, p = len(resp), len(popt)
dof = n - p
resid = resp - fourPL(conc, *popt)
sse = float(resid @ resid)
s2 = sse / dof                                   # residual variance
perr = np.sqrt(np.diag(pcov))                    # asymptotic SEs
tcrit = stats.t.ppf(0.975, dof)
ci_lo, ci_hi = popt - tcrit * perr, popt + tcrit * perr

ic50 = 10 ** log_ic50
ic50_lo, ic50_hi = 10 ** ci_lo[2], 10 ** ci_hi[2]   # back-transformed -> asymmetric
r2 = 1 - sse / np.sum((resp - resp.mean()) ** 2)

names = ["Bottom (%)", "Top (%)", "log10 IC50", "Hill slope"]
print("\n=== 4PL fit (unweighted least squares) ===")
print(f"n = {n}, parameters = {p}, dof = {dof}")
for nm, v, e, lo, hi in zip(names, popt, perr, ci_lo, ci_hi):
    print(f"  {nm:>12s} = {v:9.4f}  +/- {e:6.4f}   95% CI [{lo:8.4f}, {hi:8.4f}]")
print(f"\n  IC50           = {ic50:.3g} nM   95% CI [{ic50_lo:.3g}, {ic50_hi:.3g}] nM")
print(f"  Hill slope     = {hill:.3f}        95% CI [{ci_lo[3]:.3f}, {ci_hi[3]:.3f}]")
print(f"  residual SD    = {np.sqrt(s2):.3f} % inhibition")
print(f"  SSE = {sse:.3f},  R^2 = {r2:.5f}")

# --- optional: profile-likelihood CI for log10 IC50 (more honest for small n) ---
def profile_ci_logic50(lo_span=0.4, npts=241):
    grid = np.linspace(log_ic50 - lo_span, log_ic50 + lo_span, npts)
    thr = sse * (1 + p / dof * stats.f.ppf(0.95, 1, dof))   # ~F-based threshold
    keep = []
    for g in grid:
        f = lambda x, b, t, h: fourPL(x, b, t, g, h)
        try:
            q, _ = curve_fit(f, conc, resp, p0=[bottom, top, hill], maxfev=20000)
            r = resp - f(conc, *q)
            if float(r @ r) <= thr:
                keep.append(g)
        except RuntimeError:
            pass
    return (10 ** min(keep), 10 ** max(keep)) if keep else (np.nan, np.nan)

pl_lo, pl_hi = profile_ci_logic50()
print(f"  profile-likelihood 95% CI for IC50: [{pl_lo:.3g}, {pl_hi:.3g}] nM")

# --- constrained fit (Bottom=0, Top=100) + extra-sum-of-squares F test ---
popt2, pcov2 = curve_fit(twoPL, conc, resp, p0=[log_ic50, hill], maxfev=20000)
resid2 = resp - twoPL(conc, *popt2)
sse2 = float(resid2 @ resid2)
F = ((sse2 - sse) / (p - 2)) / (sse / dof)
pval = 1 - stats.f.cdf(F, p - 2, dof)
print("\n=== Constrained fit (Bottom=0, Top=100) ===")
print(f"  IC50 = {10**popt2[0]:.3g} nM, Hill = {popt2[1]:.3f}, SSE = {sse2:.3f}")
print(f"  F({p-2},{dof}) = {F:.2f}, p = {pval:.3f}  "
      f"(p>0.05 -> plateaus not different from 0/100)")

# ----------------------------------------------------------------------
# 4. 95% confidence band for the curve (delta method)
# ----------------------------------------------------------------------
xfit = np.logspace(np.log10(conc.min() / 2), np.log10(conc.max() * 2), 400)
yfit = fourPL(xfit, *popt)

def jac(x, theta, h=1e-6):
    J = np.empty((len(x), len(theta)))
    for i in range(len(theta)):
        tp, tm = theta.copy(), theta.copy()
        tp[i] += h; tm[i] -= h
        J[:, i] = (fourPL(x, *tp) - fourPL(x, *tm)) / (2 * h)
    return J

J = jac(xfit, popt)
se_fit = np.sqrt(np.einsum("ij,jk,ik->i", J, pcov, J))
band = tcrit * se_fit
ymid = bottom + (top - bottom) / 2.0

# ----------------------------------------------------------------------
# 5. Figure
# ----------------------------------------------------------------------
fig, ax = plt.subplots(figsize=(6.2, 4.6))
ax.fill_between(xfit, yfit - band, yfit + band, color="tab:blue", alpha=0.15,
                lw=0, label="95% CI of fit")
ax.semilogx(xfit, yfit, "-", color="tab:blue", lw=2, label="4PL fit")
ax.semilogx(conc, resp, "o", mfc="none", mec="k", ms=7, mew=1.4,
            label="data (duplicates)")

# IC50 marker
ax.axvline(ic50, color="crimson", ls="--", lw=1.2)
ax.plot([ic50], [ymid], "D", color="crimson", ms=7, zorder=5)
ax.hlines(ymid, xfit[0], ic50, color="crimson", ls=":", lw=1.0)
ax.axvspan(ic50_lo, ic50_hi, color="crimson", alpha=0.10, lw=0)
ax.annotate(f"IC$_{{50}}$ = {ic50:.2f} nM\n95% CI {ic50_lo:.2f}–{ic50_hi:.2f} nM",
            xy=(ic50, ymid), xytext=(0.06, 0.72), textcoords="axes fraction",
            arrowprops=dict(arrowstyle="->", color="crimson", lw=1.2),
            color="crimson", fontsize=9)

ax.text(0.97, 0.06,
        f"Hill slope = {hill:.2f} (95% CI {ci_lo[3]:.2f}–{ci_hi[3]:.2f})\n"
        f"Top = {top:.1f}%, Bottom = {bottom:.1f}%\n"
        f"R$^2$ = {r2:.4f}, s = {np.sqrt(s2):.2f}%",
        transform=ax.transAxes, ha="right", va="bottom", fontsize=8.5,
        bbox=dict(fc="white", ec="0.7", alpha=0.9))

ax.set_xlabel("Concentration (nM)")
ax.set_ylabel("Inhibition (%)")
ax.set_title("Dose–response: four-parameter logistic fit")
ax.set_ylim(-5, 105)
ax.grid(True, which="both", ls=":", lw=0.5, alpha=0.6)
ax.legend(loc="upper left", frameon=False, fontsize=9)
fig.tight_layout()
fig.savefig("dose_response_4pl.png", dpi=300)
fig.savefig("dose_response_4pl.pdf")
print("\nSaved: dose_response_4pl.png / .pdf")

# ----------------------------------------------------------------------
# 6. Results table
# ----------------------------------------------------------------------
with open("fit_4pl_results.csv", "w") as fh:
    fh.write("parameter,estimate,std_error,ci95_low,ci95_high\n")
    for nm, v, e, lo, hi in zip(names, popt, perr, ci_lo, ci_hi):
        fh.write(f"{nm},{v:.6g},{e:.6g},{lo:.6g},{hi:.6g}\n")
    fh.write(f"IC50_nM,{ic50:.6g},,{ic50_lo:.6g},{ic50_hi:.6g}\n")
    fh.write(f"residual_SD,{np.sqrt(s2):.6g},,,\n")
    fh.write(f"R2,{r2:.6g},,,\n")
print("Saved: fit_4pl_results.csv")

Reporting sentence you can paste into a report: “Inhibition data (8 concentrations, duplicate) were fitted by unweighted nonlinear least squares to a four‑parameter logistic model. IC₅₀ = 4.7 nM (95 % CI 4.2–5.3 nM), Hill slope = 1.05 (95 % CI 0.93–1.18), Bottom = 0.9 % and Top = 100.0 % (R² = 0.998, residual SD 1.8 %). The slope of unity and well‑defined plateaus are consistent with simple reversible single‑site inhibition.”

One caveat worth keeping in mind: the CIs above assume constant, normally distributed error in % inhibition. If your assay noise grows in the mid‑range (common), refit with weights (sigma= in curve_fit) or bootstrap the residuals — with these data the IC₅₀ will move very little, but the interval may widen slightly.