# How to fit a function to measurements with error?

**URL:** <https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369>\
**Category:** Data\
**Tags:** question, curve-fitting, measurements\
**Created:** [July 27, 2021, 1:23pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369 "2021-07-27T13:23:29Z")\
**Posts on this page:** 13\
**Page:** 1

<div class="post-metadata">

**Author:** ![DaWa](https://avatars.discourse-cdn.com/v4/letter/d/46a35a/32.png) [@DaWa](https://discourse.julialang.org/u/DaWa)\
**Post date:** [July 27, 2021, 1:23pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/1 "2021-07-27T13:23:30Z")

</div>

Hi all,  
I am currently using Measurements.jl for error propagation and LsqFit.jl for fitting functions to data. Is there a simple way to fit a function to data with errors? It would be no problem to use an other package if that makes things easier.

Thanks in advance for your help.

---

<div class="post-metadata">

**Author:** ![nilshg](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nilshg/32/2283_2.png) [@nilshg](https://discourse.julialang.org/u/nilshg)\
**Post date:** [July 27, 2021, 1:29pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/2 "2021-07-27T13:29:40Z")

</div>

Thanks for heeding my advice and welcome to Discourse! Could you explain a bit more about what you mean by “measurements with error”? LsqFit solves Y=m(x,γ)+ϵ, where Y and x are your observations, so it naturally applies to data with errors.

---

<div class="post-metadata">

**Author:** ![DaWa](https://avatars.discourse-cdn.com/v4/letter/d/46a35a/32.png) [@DaWa](https://discourse.julialang.org/u/DaWa)\
**Post date:** [July 27, 2021, 1:36pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/3 "2021-07-27T13:36:04Z")

</div>

Hello nice to see you again!

currently I am storing my data in a dataframe like this

```julia
data= DataFrame("xdata" => 
        [
        9.7 ± 0.1,
        6.6 ± 0.1,
        ...

```

if I try to fit it using

`fit = curve_fit(lineareq, data.xdata,data.ydata, [1.,1.])`

I get an error message

MethodError: no method matching Float64(::Measurement{Float64})

---

<div class="post-metadata">

**Author:** ![nilshg](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nilshg/32/2283_2.png) [@nilshg](https://discourse.julialang.org/u/nilshg)\
**Post date:** [July 27, 2021, 2:44pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/4 "2021-07-27T14:44:53Z")

</div>

Ah I see - it looks like LsqFit is doing an explicit promotion to float somewhere. Not sure whether one should expect LsqFit to work with Measurements - maybe @pkofod can comment.

What do you expect to happen with the uncertainty in estimation? It looks like measurement error in the covariates (right hand side) of your equation to me, which traditionally in econometrics would just lead to attenuation bias. Or do you have different (known) degrees of uncertainties in different parts of your data set? In that case this discussion here might be of relevance, which shows the use of inverse variance weights in LsqFit:

> [@Weighted linear regression with confidence interval fitted to error bars](https://discourse.julialang.org/t/weighted-linear-regression-with-confidence-interval-fitted-to-error-bars/60743/6):
>
> As indicated in GLM.jl’s doc above, glm does not handle inverse-variance weighting. For this purpose, you may use LsqFit.jl. using DataFrames, LsqFit, Printf, Plots; gr() df = DataFrame(x = [0.0, 0.0669873, 0.25, 0.5, 0.75, 0.933013, 1.0], y = [0.223, 0.291, 0.393, 0.549, 0.73, 0.85, 0.896], u\_y = [0.023, 0.024, 0.027, 0.031, 0.037, 0.041, 0.043]) x, y = df.x, df.y wt = 1 ./ df.u\_y .^2 p0 = [0.5, 0.5] m(x, p) = p[1] .+ p[2] \* x # p: model parameters fit = curve\_fit(m…

---

<div class="post-metadata">

**Author:** ![DaWa](https://avatars.discourse-cdn.com/v4/letter/d/46a35a/32.png) [@DaWa](https://discourse.julialang.org/u/DaWa)\
**Post date:** [July 27, 2021, 4:45pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/5 "2021-07-27T16:45:21Z")

</div>

> [@nilshg](#):
>
> What do you expect to happen with the uncertainty in estimation? It looks like measurement error in the covariates (right hand side) of your equation to me, which traditionally in econometrics would just lead to attenuation bias. Or do you have different (known) degrees of uncertainties in different parts of your data set? In that case this discussion here might be of relevance, which shows the use of inverse variance weights in LsqFit:

Thanks again for the reply. The sample code you provided was very helpfull. I completely overlooked, that there is an option to weight your data. Sadly I dont know that much about statistics, so it is not easy to point out what I need.

I want something that is equal to using weighted mean instead of a mean.

Here in my linear regression example it seems I can use a weight wt = 1 ./ df.u\_y .^2 to weight my y values. It would be nice to also consider the x error and if possible automatically calculate the weights for more difficult functions. I know it is much to ask, but maybe there is something like this out there.

Ill read up on this so maybe I can clarify later on

---

<div class="post-metadata">

**Author:** ![jling](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jling/32/212909_2.png) [@jling](https://discourse.julialang.org/u/jling)\
**Post date:** [August 15, 2021, 12:58pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/6 "2021-08-15T12:58:52Z")

</div>

![image](https://global.discourse-cdn.com/julialang/original/3X/7/d/7dcd210bf6d3dc38c57f7b2382ff36abfca22432.png)

Polynomials work, but I agree having LsqFit just work would be great

---

<div class="post-metadata">

**Author:** ![giordano](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/giordano/32/2166_2.png) [@giordano](https://discourse.julialang.org/u/giordano)\
**Post date:** [August 15, 2021, 2:28pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/7 "2021-08-15T14:28:10Z")

</div>

> [@jling](#):
>
> but I agree having LsqFit just work would be great

Probably it needs to be taught how to weight the errors

---

<div class="post-metadata">

**Author:** ![jling](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jling/32/212909_2.png) [@jling](https://discourse.julialang.org/u/jling)\
**Post date:** [August 15, 2021, 3:11pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/8 "2021-08-15T15:11:33Z")

</div>

I think it will be tricky to teach it to do the right&sensible thing. Polynomials (i.e. linear regression) is non controversial right? Just propagate error through the matrix operations. But idk if `mse` minimization has some property saying “the error on fitted parameter is still 1 sigma”

---

<div class="post-metadata">

**Author:** ![sijo](https://avatars.discourse-cdn.com/v4/letter/s/da6949/32.png) [@sijo](https://discourse.julialang.org/u/sijo)\
**Post date:** [August 15, 2021, 7:00pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/9 "2021-08-15T19:00:15Z")

</div>

Not really what is being asked, but another approach would be a Bayesian inference: defining the MvNormal corresponding to each (x\_i,y\_i) data point, and drawing e.g. 100 samples from that distribution, then using MCMC to estimate the model parameters on 100n “data points”. Does anyone have an idea whether this should work?

(Even better would be to use the MVNormals directly in the model definition without sampling… I don’t know if that’s possible in one of Juli’a probabilistic DSLs.)

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [August 15, 2021, 8:00pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/10 "2021-08-15T20:00:13Z")

</div>

There is this package for dealing with outliers as well:

[https://github.com/fsobral/RAFF.jl](https://github.com/fsobral/RAFF.jl)

---

<div class="post-metadata">

**Author:** ![sijo](https://avatars.discourse-cdn.com/v4/letter/s/da6949/32.png) [@sijo](https://discourse.julialang.org/u/sijo)\
**Post date:** [August 15, 2021, 9:47pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/11 "2021-08-15T21:47:14Z")

</div>

Another (simpler) idea also based on samples from the distribution of each (x\_i,y\_i) data point: make an “ensemble fit” and calculate errors from the set of results. Basically what you would get with MonteCarloMeasurements.jl if it worked with LsqFit.jl (I assume it would fail with the same error as with Measurements.jl).

---

<div class="post-metadata">

**Author:** ![Fourier](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fourier/32/38176_2.png) [@Fourier](https://discourse.julialang.org/u/Fourier)\
**Post date:** [September 8, 2024, 7:51am UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/12 "2024-09-08T07:51:21Z")

</div>

Did someone find a better version yet? Something that also plots the uncertainty of the model given by the uncertainty in the model parameters? Maybe something from SciML, they seem to have things that go into that direction: [Optimization Under Uncertainty · Overview of Julia's SciML](https://docs.sciml.ai/Overview/stable/showcase/optimization_under_uncertainty/)

---

<div class="post-metadata">

**Author:** ![Fourier](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fourier/32/38176_2.png) [@Fourier](https://discourse.julialang.org/u/Fourier)\
**Post date:** [September 8, 2024, 8:28pm UTC](https://discourse.julialang.org/t/how-to-fit-a-function-to-measurements-with-error/65369/13 "2024-09-08T20:28:27Z")

</div>

I found something that works!

```julia
pkg> add https://github.com/chriswaudby/LsqFit.jl.git#MeasurementsExtension

```

For example

```julia
using LsqFit, Plots, Measurements
@. gauss(x,p) = p[1] + p[2] * exp(-(x-p[3])^2/p[4]^2)
x = range(-1,1,100)
yerr = gauss(x,[0,1,0,0.5]) + 0.1*randn(length(x)) .± 0.1
fit = curve_fit(gauss, x,yerr,Float64[0,1,0,1])
begin
	scatter(x,yerr)
	plot!(x,gauss(x,fit.param);ribbon = gauss(x,fit.param.+stderror(fit)) - gauss(x,fit.param.-stderror(fit)))
end
savefig("fit.png")

```

![fit](https://global.discourse-cdn.com/julialang/original/3X/e/b/eb44eb8c4d6bce39afe75c41594e30dac26e6ac0.png)
