Bio Applied Enzyme Kinetics

Enzyme Kinetics: Computational Patterns

When to Use

• Fitting Vmax and Km from substrate-velocity data (spectrophotometric or fluorescent enzyme assays)
• Determining inhibitor mechanism (competitive, noncompetitive, uncompetitive, mixed) and estimating Ki
• Detecting allosteric/cooperative binding (Hill coefficient n) e.g. hemoglobin-O2 saturation
• Comparing multi-substrate mechanisms (ping-pong vs sequential/ternary-complex) from initial-rate data
• Reporting kcat/Km catalytic efficiency with confidence intervals for enzyme characterization

Version Compatibility

Python ≥3.10, NumPy ≥1.24, SciPy ≥1.11 (scipy.optimize.curve_fit, scipy.stats.chi2). Matplotlib ≥3.7 for plots. No enzyme-specific package needed — this is pure curve-fitting on scipy/numpy.

Prerequisites

pip install numpy scipy matplotlib. You need initial-rate (v0) data at multiple substrate concentrations, ideally spanning ~0.1×Km to 10×Km, and (for inhibition studies) velocities at several inhibitor concentrations. Prior concept: nonlinear least squares fitting.

Goal: Extract Vmax, Km, kcat, and catalytic efficiency from substrate-velocity data with uncertainty estimates. Approach: Fit the Michaelis-Menten equation with curve_fit using data-derived initial guesses, then bootstrap residuals for a 95% confidence band.

import numpy as np
from scipy.optimize import curve_fit


def michaelis_menten(S, vmax, km):
    """Michaelis-Menten velocity: v = Vmax*[S] / (Km + [S])."""
    return (vmax * S) / (km + S)


def fit_michaelis_menten(substrate_uM, observed_v, enzyme_conc_uM=None):
    """Fit MM parameters with data-driven initial guesses.

    Returns dict with vmax, km, their 95% CIs, and (if enzyme_conc_uM given)
    kcat and catalytic efficiency (kcat/Km, in µM^-1 s^-1).
    """
    # Good initial guesses matter: bad p0 crashes curve_fit or converges to nonsense
    v_max_guess = max(observed_v) * 1.1
    half_max_idx = np.argmin(np.abs(observed_v - max(observed_v) / 2))
    km_guess = substrate_uM[half_max_idx]

    popt, pcov = curve_fit(
        michaelis_menten, substrate_uM, observed_v,
        p0=[v_max_guess, km_guess],
        bounds=([0, 0], [np.inf, np.inf]),
        maxfev=5000,
    )
    vmax_fit, km_fit = popt
    perr = np.sqrt(np.diag(pcov))  # standard errors
    ci95 = 1.96 * perr

    result = {"vmax": vmax_fit, "km": km_fit, "vmax_ci95": ci95[0], "km_ci95": ci95[1]}
    if enzyme_conc_uM is not None:
        kcat = vmax_fit / enzyme_conc_uM  # turnover number, s^-1
        result["kcat"] = kcat
        result["catalytic_efficiency"] = kcat / km_fit  # µM^-1 s^-1
    return result


# Substrate concentrations should span 0.1*Km to 10*Km for a well-constrained fit
substrate_uM = np.array([0.5, 1.0, 2.0, 3.0, 5.0, 7.5, 10.0, 15.0, 20.0, 30.0, 50.0])

Goal: Quantify fit uncertainty beyond the covariance-matrix standard errors. Approach: Parametric bootstrap — resample residuals around the fitted curve, refit, and take the 2.5th/97.5th percentile envelope.

def bootstrap_mm_ci(substrate_uM, observed_v, popt, n_boot=500, s_smooth=None):
    """Parametric bootstrap for MM fit confidence band.

    Assumes homoscedastic residuals (constant noise SD). If assay noise scales
    with velocity, pass sigma=observed_v to curve_fit instead (weighted LS).
    """
    if s_smooth is None:
        s_smooth = np.linspace(0, substrate_uM.max() * 1.1, 300)

    residuals = observed_v - michaelis_menten(substrate_uM, *popt)
    residual_sd = np.std(residuals)

    boot_curves = np.full((n_boot, len(s_smooth)), np.nan)
    boot_params = np.full((n_boot, 2), np.nan)

    for i in range(n_boot):
        v_boot = michaelis_menten(substrate_uM, *popt) + np.random.normal(0, residual_sd, len(substrate_uM))
        v_boot = np.clip(v_boot, 0, None)
        try:
            p_boot, _ = curve_fit(michaelis_menten, substrate_uM, v_boot,
                                   p0=popt, bounds=([0, 0], [np.inf, np.inf]), maxfev=2000)
            boot_curves[i] = michaelis_menten(s_smooth, *p_boot)
            boot_params[i] = p_boot
        except RuntimeError:
            continue  # fit failed to converge for this resample; leave as NaN

    ci_low = np.nanpercentile(boot_curves, 2.5, axis=0)
    ci_high = np.nanpercentile(boot_curves, 97.5, axis=0)
    return s_smooth, ci_low, ci_high, boot_params

Inhibition Models and Identification

Goal: Determine inhibitor mechanism and estimate Ki from dose-response data. Approach: Fit each candidate model (competitive/noncompetitive/uncompetitive/mixed), compare Vmax/Km shifts, and confirm with AIC or a Dixon plot.

def competitive(S, vmax, km, I, Ki):
    """Km increases with [I]; Vmax unchanged. Converging LB lines left of y-axis."""
    return (vmax * S) / (km * (1 + I / Ki) + S)


def uncompetitive(S, vmax, km, I, Ki):
    """Both Vmax and Km decrease by the same factor alpha. Parallel LB lines."""
    alpha = 1 + I / Ki
    return (vmax / alpha * S) / (km / alpha + S)


def noncompetitive(S, vmax, km, I, Ki):
    """Vmax decreases with [I]; Km unchanged. LB lines converge on x-axis."""
    return (vmax / (1 + I / Ki) * S) / (km + S)


def mixed_inhibition(S, vmax, km, I, Ki, alpha):
    """General case: alpha=1 -> noncompetitive; alpha -> inf -> competitive."""
    return (vmax * S) / (alpha * km * (1 + I / Ki) + S)


def estimate_ki_from_dixon(inhibitor_concs, apparent_kms):
    """Estimate Ki and true Km from a series of apparent-Km fits at different [I]
    (competitive inhibition: Km_app = Km * (1 + [I]/Ki), linear in [I]).
    """
    slope, intercept = np.polyfit(inhibitor_concs, apparent_kms, 1)
    km_true = intercept
    ki = intercept / slope
    return km_true, ki

Allosteric Cooperativity (Hill Equation)

Goal: Detect and quantify positive/negative cooperativity (e.g. hemoglobin-O2 binding). Approach: Fit the Hill equation; n > 1 means positive cooperativity, n < 1 negative, n = 1 reduces to Michaelis-Menten.

def hill_equation(S, vmax, k_half, n):
    """Hill equation. n > 1: positive cooperativity; n < 1: negative; n = 1: hyperbolic (MM)."""
    return (vmax * S**n) / (k_half**n + S**n)


def fit_hill(substrate_uM, observed_v):
    """Fit Hill equation; K_half is [S] at half-maximal velocity (= Km when n=1)."""
    popt_hill, pcov_hill = curve_fit(
        hill_equation, substrate_uM, observed_v,
        p0=[max(observed_v), np.median(substrate_uM), 1.5],
        bounds=([0, 0, 0.1], [np.inf, np.inf, 10]),
    )
    vmax_h, k_half, n_hill = popt_hill
    return vmax_h, k_half, n_hill

Multi-Substrate Kinetics (Ping-Pong vs Sequential)

def ping_pong(S_A, S_B, vmax, Km_A, Km_B):
    """Ping-pong bi-bi: parallel lines in double-reciprocal plot at varied [B]."""
    return vmax / (1 + Km_A / S_A + Km_B / S_B)


def sequential(S_A, S_B, vmax, Km_A, Km_B, Ki_A):
    """Sequential (ternary complex): intersecting lines in double-reciprocal plot."""
    return (vmax * S_A * S_B) / (Ki_A * Km_B + Km_B * S_A + Km_A * S_B + S_A * S_B)

Goodness of Fit

from scipy.stats import chi2


def evaluate_fit(substrate_uM, observed_v, popt, measurement_sd=None):
    """Report R^2 and, if per-point measurement SDs are known, chi-squared p-value."""
    residuals = observed_v - michaelis_menten(substrate_uM, *popt)
    ss_res = np.sum(residuals**2)
    ss_tot = np.sum((observed_v - observed_v.mean())**2)
    r_squared = 1 - ss_res / ss_tot

    out = {"r_squared": r_squared}
    if measurement_sd is not None:
        chi2_stat = np.sum((residuals / measurement_sd)**2)
        out["chi2_p"] = chi2.sf(chi2_stat, df=len(observed_v) - 2)
    return out

Pitfalls

• Bad initial guesses crash curve_fit: estimate p0 from the data — Vmax ≈ max(observed_v) × 1.1, Km ≈ [S] where v ≈ Vmax/2; never use default p0=[1,1] for kinetic data
• Bootstrap assumes homoscedastic residuals: enzyme kinetic noise is often proportional to velocity (heteroscedastic); in that case use weighted least squares (sigma=observed_v) in curve_fit
• Hill n is sensitive to data range: fitting the Hill equation requires data both below and above K_half; sparse data near the inflection point gives unstable n estimates
• Competitive vs mixed inhibition: mixed inhibition with large alpha looks like competitive — measure at multiple inhibitor concentrations and use a Dixon plot (apparent Km vs [I], or 1/v vs [I]) to distinguish
• kcat requires total enzyme concentration: convert Vmax to kcat only if [E]_total is from active-site titration, not total protein concentration (often contains inactive enzyme)
• pcov returns inf when the fit is underdetermined: too few data points relative to parameters, or parameters are not independently identifiable — add more substrate concentrations
• Never use Lineweaver-Burk for fitting: linearization amplifies error at low [S] (1/[S] blows up), badly biasing Vmax/Km; use it only for visual mechanism diagnosis, fit with nonlinear least squares

See Also

• bio-systems-biology-flux-balance-analysis — reaction-level kinetics in genome-scale metabolic models
• bio-chemoinformatics-admet-prediction — small-molecule inhibitor property prediction
• statistical-analysis — general nonlinear regression and bootstrap CI methodology
• scikit-learn — alternative optimization/regression backends if scipy fit fails to converge