TopolearnTurning measurements into model parameters (rates, half-lives, growth curves) is a daily scientific task.
optimize.curve_fit(model, xdata, ydata, p0=...) finds parameters that minimise the squared residuals between model(xdata, *params) and ydata.
The model function takes the independent variable first, then each parameter as a separate argument: def model(t, a, k): return a * np.exp(-k * t).
It returns popt (best parameters) and pcov (their covariance matrix). Standard errors are np.sqrt(np.diag(pcov)); large values or inf mean the parameter is poorly determined by the data.
Practical advice:
- Give a sensible starting guess
p0. The default is all ones, which fails for many non-linear models (exponentials, sigmoids). - Use
bounds=([lo...], [hi...])to keep parameters physical, for example a positive rate. - If you know measurement errors, pass
sigma=andabsolute_sigma=Truesopcovreflects them. - Plot residuals: a pattern means the model shape is wrong, whatever the fit statistics say.
curve_fit wraps optimize.least_squares, which you can call directly for more control: robust losses such as loss="soft_l1" reduce outlier influence, and it works with residual functions of any form. For straight lines, stats.linregress or np.polyfit are simpler.
Example
import numpy as np
from scipy import optimize
def decay(t, a, k): # x first, then parameters
return a * np.exp(-k * t)
rng = np.random.default_rng(0)
t = np.linspace(0, 5, 30)
y = decay(t, 10.0, 0.8) + rng.normal(0, 0.2, t.size)
popt, pcov = optimize.curve_fit(decay, t, y, p0=[5, 1], bounds=([0, 0], [np.inf, np.inf]))
perr = np.sqrt(np.diag(pcov))
print(popt.round(3)) # [10.094 0.818]
print(perr.round(3)) # [0.108 0.013] (standard errors)
print(round(np.log(2) / popt[1], 2)) # 0.85 (half-life)Remember
- The model signature is
f(x, p1, p2, ...) - Provide
p0andboundsfor non-linear models - Standard errors come from
sqrt(diag(pcov))
Watch out for
- Relying on the default p0 of all ones for an exponential or logistic model.