Example 1 — Basic logistic fitting

Fits a logistic growth model to synthetic plate-reader data, plots observed vs fitted curves, and compares two models by AICc.

Full script: examples/scripts/02_nl_fitting.py

Setup

import numpy as np
import matplotlib.pyplot as plt
import pykinbiont
from pykinbiont import GrowthData, FitOptions, ModelSpec, MODEL_REGISTRY, fit

pykinbiont.configure("/path/to/KinBiont.jl")  # skip if using registry version

Generate synthetic data

def logistic(t, K=1.2, mu=0.5, N0=0.01):
    return K / (1 + ((K - N0) / N0) * np.exp(-mu * t))

rng   = np.random.default_rng(1)
times = np.linspace(0, 20, 80)

curves = np.stack([
    np.maximum(
        logistic(times, K=1.0 + 0.1 * i, mu=0.4 + 0.05 * i, N0=0.01)
        + rng.normal(0, 0.007, len(times)),
        1e-4,   # clip to avoid log(negative) errors
    )
    for i in range(4)
])

data = GrowthData(
    curves=curves,
    times=times,
    labels=["W1", "W2", "W3", "W4"],
)

Why clip to 1e-4?

KinBiont computes log(OD) internally to estimate growth rates. Any OD ≤ 0 raises a DomainError. Clipping Gaussian noise at a small positive floor prevents this.

Fit with logistic model

logistic_model = MODEL_REGISTRY["NL_logistic"]

spec = ModelSpec(
    models=[logistic_model],
    params=[[1.2, 0.01, 0.5]],    # [N_max, N0, mu]
    lower=[[0.3, 1e-3, 0.05]],
    upper=[[3.0, 0.05, 3.0]],
)
opts = FitOptions(smooth=True, smooth_method="rolling_avg", multistart=True, n_restart=20)

results = fit(data, spec, opts)

Inspect results

df = results.to_dataframe()
print(df[["label", "best_model", "aic", "param_1", "param_2", "param_3"]])

label	best_model	aic	param_1 (N_max)	param_2 (N0)	param_3 (mu)
W1	NL_logistic	−142	1.007	0.010	0.401
W2	NL_logistic	−138	1.103	0.010	0.449

Plot

fig, axes = plt.subplots(2, 2, figsize=(9, 6))

for ax, r in zip(axes.flat, results):
    idx = data.labels.index(r.label)
    ax.plot(data.times, data.curves[idx], ".", markersize=3, alpha=0.5, label="data")
    ax.plot(r.times, r.fitted_curve, "-", lw=2, label="logistic fit")
    ax.set_title(r.label)
    ax.set_xlabel("Time (h)")
    ax.set_ylabel("OD")
    ax.legend(fontsize=7)

plt.tight_layout()
plt.savefig("logistic_fits.png", dpi=150)

Two-model comparison

gompertz_model = MODEL_REGISTRY["NL_gompertz"]

spec_both = ModelSpec(
    models=[logistic_model, gompertz_model],
    params=[[1.2, 0.01, 0.5], [1.2, 0.5, 1.0]],
    lower=[[0.3, 1e-3, 0.05], [0.3, 0.0, 0.0]],
    upper=[[3.0, 0.05, 3.0],  [3.0, 3.0, 10.0]],
)

results_both = fit(data, spec_both, opts)
print(results_both.to_dataframe()[["label", "best_model", "aic"]])

The model with the lower AICc is selected automatically per curve.