Scaling of Parameters

[ghost]

Hi everyone,
I'm getting started with lmfit and stumbled across the issue of (auto)scaling parameters. According to an old issue the least_sq method does autoscaling, but for me, the fit itself and the estimation of parameter errors fails if the values differ by orders of magnitude, see minimal working example below.

Interestingly, it seem to work when using method="least_squares" for fitting.

import numpy as np
from lmfit.models import LinearModel

np.random.seed(1)

# real data
a = 1e-13
b = 0.1
y_err = 1e-14

x = np.linspace(0, 5, 8)
y = a*x+b + np.random.normal(0, y_err, x.shape)

model = LinearModel()
params = model.guess(y, x=x)
result = model.fit(y, params, x=x)

print(result.fit_report())

[newville]

@bleolta Actually, the fit does succeed. Evaluating the uncertainties does not.

That is because the value of one of the variables is extremely close to the machine precision estimate of 0, and the "noise" added to your data is extremely small, again close to machine precision. That means that the chi-square is also exceptionally small, and so that it is very hard to estimate uncertainties. If, for example, you changed a to 1.e-9, both the fit and estimation of the uncertainties would "work".
If you had set a to 1.e-12 and the noise scale to 1.e-12, both woould have worked.

This is a well-known feature of using floating point numbers. It is not something we are going to attempt to fix. It is always recommended to scale parameters to have magnitudes in the range of 1e-10 to 1e10.

[ghost]

Ah ok, that's indeed a good reason :)

So autoscaling is applied, but it's just too close to machine precision?

[newville]

@bleolta The leastsq method (MINPACK-1 lmdif.f or lmder.f depending on whether a Jacobian function is provided) does auto-scaling of parameter values so that the initial derivatives are on the same range of magnitudes.

Our code in eval_uncertainties does not do that.

But, that is probably not important here, where parameter values differ by 13 orders of magnitude. The advice is simple : Don't do that. scale values in your function or change units or whatever so that the variable parameters have magnitudes on the range [1-e6, 1e6].

You are also asking for the fit to match the data at the 1.e-13 level (because you added noise at that level), with then claimed the noise level to be 1.e-14. That's very close to machine precision and problem cases are likely to appear.

If the noise level in your data was truly better than 1.e-12, you probably would not be asking basic questions about fitting and floating point on GitHub ;)