# Nonlinear curve fitting with weights with LsqFit

Hello,

I need to fit the exponential decay law to some data with weights. I tried to use LsqFit.jl, however the `standard_errors` function returns too low values. I discovered it by comparing my results with my friend and I also made the fit in R.

I am sorry, I did not think it would be useful as I just called the relevant functions, but anyway here it is:

``````using LsqFit

x = [0, 975, 1767, 2682, 5959]
y = [5490.0, 4546.0, 3679.0, 3134.0, 1451.0]
weights = [109.0, 102.0, 96.0, 89.0, 71.0]

@. model(t, p) = p[1]*(0.5)^(t/p[2])
fit  = curve_fit(model, x, y, weights, [0.1, 80])
err = round(standard_errors(fit)[2], sigdigits=1)
T = fit.param[2]

println("\$T +- \$err")
``````

Just to clarify, the `weights` variable contains the stderr of measurement of y.

Sorry, please provide a piece of code that is completely defined. I tried your code and get the message:

``````
julia> @. model(t, p) = p[1]*(0.5)^(t/p[2])
model (generic function with 1 method)

julia> fit  = curve_fit(model, x, y, weights, [0.1, 60])
ERROR: UndefVarError: x not defined
Stacktrace:
[1] top-level scope
@ REPL[3]:1
``````

What are x, y and the weights?

In particular:
"
Do your best to make your example self-contained (âminimal working exampleâ, MWE ), so that it runs (or gets to the error that you want help with) as is. This means including package loading (e.g. `using ThatPackage` ) and any data that the code operates on. If your data is large or proprietary, generate example data if possible and include that. This StackOverflow answer explains how a DataFrame can be turned into a string representation that can be copy/pasted.
"

This is probably the issue. In `LsqFit.jl`, the weights are assumed to be the inverse of the variances.

1 Like

I edited the code in my previous post.

I tried that also, but the error is still suspiciously low.

Using the inverse of the variances as weights, we obtain ` using LsqFit` :

``````p1 = 5541 +/- 86
p2 = 3135 +/- 110
``````

Which from the QC plot below does not seem unreasonable.
What values do the other software packages compute?

1 Like

Here is my take on it

``````using LsqFit

x = [0, 975, 1767, 2682, 5959]
y = [5490.0, 4546.0, 3679.0, 3134.0, 1451.0]
weights = [109.0, 102.0, 96.0, 89.0, 71.0]

@. model(t, p) = p[1]*(0.5)^(t/p[2])
fit  = curve_fit(model, x, y, [0.1, 80])
err = round(standard_errors(fit)[2], sigdigits=1)
T = fit.param[2]

println("\$T +- \$err")

@. model(t, p) = p[1]*(0.5)^(t/p[2])
fit  = curve_fit(model, x, y, weights, [0.1, 80])
err = round(standard_errors(fit)[2], sigdigits=1)
T = fit.param[2]

println("\$T +- \$err")

@. model(t, p) = p[1]*(0.5)^(t/p[2])
fit  = curve_fit(model, x, y, 1 ./ (weights .* weights), [0.1, 80])
err = round(standard_errors(fit)[2], sigdigits=1)
T = fit.param[2]

println("\$T +- \$err")
``````

yielding

``````3158.802386807851 +- 100.0
3170.773189838312 +- 0.1
3134.79231889205 +- 100.0
``````

It is much more clear to write:

``````Ďy = [109., 102., 96., 89., 71.]   # standard deviations of y
w = 1 ./ Ďy.^2                     # weights for LsqFit
``````
1 Like

I see the problem now, I used `1/Ď` instead of `1/Ď^2`. Thank you, for your help!