# Alternatives to LsqFit for nonlinear curve fitting?
**URL:** https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685
**Category:** General Usage
**Tags:** curve-fitting
**Created:** [December 15, 2023, 7:58pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685 "2023-12-15T19:58:11Z")
**Posts on this page:** 13
**Page:** 1
### Author:  [@PeX](https://discourse.julialang.org/u/PeX)
#### Post date: [December 15, 2023, 7:58pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/1 "2023-12-15T19:58:11Z")
Hi!
I’m trying to use Julia to make a curve fit with the functional form:
```julia
model = ( p[1] / x ) * exp( -p[2] / x )
```
to a data-set I have.
I tried LsqFit, however, if I’m not making really good initial guesses for the parameters `p[1]` and `p[2]`, the algorithm can’t find a good fit (or even remotely close one), so this kind of misses the point.
Are there alternatives to LsqFit that are more robust? Or is there a way to use LsqFit in a better way than simply calling ` curve_fit`?
Thanks !!
---
### Author:  [@lmiq](https://discourse.julialang.org/u/lmiq)
#### Post date: [December 15, 2023, 8:11pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/2 "2023-12-15T20:11:12Z")
Have you tried fitting the log?
---
### Author:  [@PeX](https://discourse.julialang.org/u/PeX)
#### Post date: [December 15, 2023, 8:24pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/3 "2023-12-15T20:24:07Z")
Do you mean to transfer x → log10(x) within the same model structure?
---
### Author:  [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)
#### Post date: [December 15, 2023, 8:55pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/4 "2023-12-15T20:55:41Z")
What are the ranges of `x` and `y` in the data to be fitted?
---
### Author:  [@PeX](https://discourse.julialang.org/u/PeX)
#### Post date: [December 15, 2023, 9:11pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/5 "2023-12-15T21:11:09Z")
` 1e-5 < x < 1e-3` and `0.6 < y < 0.8`
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [December 15, 2023, 9:13pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/6 "2023-12-15T21:13:26Z")
> [@PeX](#):
>
> I tried LsqFit, however, if I’m not making really good initial guesses for the parameters `p[1]` and `p[2]`, the algorithm can’t find a good fit (or even remotely close one), so this kind of misses the point.
I can’t reproduce your problem. For this particular model, I’m finding that LsqFit is pretty robustly finding the parameters. But maybe you are in a very different regime?
For example, I considered a ground-truth model `p = [1,1]`, with 30 noisy data points on x \in (0,7), given by the following code:
```julia
model(x, p) = ( p[1] / x ) * exp( -p[2] / x )
p = [1,1]
xdata = rand(30) * 7
ydata = model.(xdata, Ref(p)) .+ randn.() .* 0.02
```
which looks like:
Then, if I run `curve_fit` from various starting points, it always seems to recover p \approx [1,1]. Even a starting guess that is way off works:
```julia
fit = curve_fit((x, p) -> model.(x, Ref(p)), xdata, ydata, [12.0, 32.0])
fit.param
```
gives `[1.0053649510098448, 0.9865169899979824]`.
What does the data that you are trying to fit look like?
It might be worth plotting your least-square objective function as a function of `p` to see if it has multiple local minima for your data, or if it is badly conditioned (compute the condition number of the Hessian matrix at the optimum — ForwardDiff.jl or TaylorDiff.jl should have no trouble computing the Hessian of a simple two-parameter function like this).
---
### Author:  [@PeX](https://discourse.julialang.org/u/PeX)
#### Post date: [December 15, 2023, 9:17pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/7 "2023-12-15T21:17:42Z")
First of all, thank you for the help!
Second, you are all probably right, my data falls at the very initial section of the functional form (the very small values of x where the curve starts to incline in a semi-parabolic fashion). That’s probably why it is difficult to reproduce the fit with the algorithm
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [December 15, 2023, 9:44pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/9 "2023-12-15T21:44:15Z")
> [@PeX](#):
>
> Second, you are all probably right, my data falls at the very initial section of the functional form
Ah, I see. If I plot the sum-of-squares error (SSE) that you are minimizing, computed by
```julia
sse(model, xdata, ydata, p) = sum(((x,y),) -> abs2(y - model(x, p)), zip(xdata,ydata))
```
then for my example data above (x \in (0,7)), I get:
which looks like it has a single well-conditioned local minimum at ≈ the ground-truth p = [1,1]. This is why least-square fitting works so well for me.
In contrast, if I use x data just in (0,0.5), before the peak of the model, then the same plot looks like:
i.e. there is a whole curve in p space that has nearly the same SSE as the ground-truth p = [1,1]. I don’t know if these are all local minima or if it slopes _very_ shallowly down towards the ground-truth p, but at the very least the minimum is quite ill conditioned. In particular, we can compute the condition number of the Hessian using:
```julia
using ForwardDiff, LinearAlgebra
H = ForwardDiff.hessian(p -> sse(model, xdata, ydata, p), [1,1])
cond(H)
```
which gives `≈ 145` — the second derivative at the ground truth is 145x smaller in one direction than the other, which makes the location of the optimum quite sensitive to noise, and also makes optimization converge slowly at best.
That being said, I find that I can still get a decent approximation for the ground-truth optimum by lowering the tolerance on `curve_fit`, so that it runs for more iterations. For example,
```julia
fit = curve_fit((x, p) -> model.(x, Ref(p)), xdata, ydata, [12.0, 11.0]; maxIter=1000, show_trace=true, x_tol=0, g_tol=0)
@show fit.param
```
where I’ve forced it to run for 1000 iterations, converges to:
```julia
[1.2752929578813368, 1.091880904890356]
```
which is probably as close to the ground-truth p = [1,1] as the noise allows — again, because the minimum is badly conditioned, the noise in the problem is amplified to a large shift in in the minimum (unless you have a huge amount of data).
I’m guessing that as your x range shrinks, the condition number gets worse, so the problem is exacerbated. Still, you might still get acceptable results if you lower the tolerances / increase the iterations as I’ve done above.
(But get data from larger x if you can!)
---
### Author:  [@PeX](https://discourse.julialang.org/u/PeX)
#### Post date: [December 15, 2023, 9:58pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/10 "2023-12-15T21:58:43Z")
Fantastic response! Very enlightening and insightful! I definitely learned a lot.
Much appreciated!
---
### Author:  [@bjarthur](https://discourse.julialang.org/u/bjarthur)
#### Post date: [December 19, 2023, 3:59pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/11 "2023-12-19T15:59:41Z")
> [@PeX](#):
>
> Are there alternatives to LsqFit that are more robust?
so LsqFit.jl then is the go-to package for fitting equations like this in julia? no one who responded to this post so far has addressed this portion of the OP.
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [December 19, 2023, 4:28pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/12 "2023-12-19T16:28:03Z")
> [@bjarthur](#):
>
> so [LsqFit.jl](https://juliahub.com/ui/Packages/LsqFit) then is the go-to package for fitting equations like this in julia? no one who responded to this post so far has addressed this portion of the OP.
There are plenty of other optimization packages if you want to treat least-square fitting as a general optimization problem.
But if you want to take advantage of the special structure of least-square optimization, some options currently seem to be LsqFit.jl, [LeastSquaresOptim.jl](https://github.com/matthieugomez/LeastSquaresOptim.jl), [NLLSsolver.jl](https://github.com/ojwoodford/NLLSsolver.jl), and [ManOpt.jl](https://manoptjl.org/stable/solvers/LevenbergMarquardt/). (_Update:_ See replies for more packages.) See also [Comparing non-linear least squares solvers](https://discourse.julialang.org/t/comparing-non-linear-least-squares-solvers/104752)
In the particular problem raised in the present thread, however, the lack of robustness was intrinsic to the optimization problem being posed rather than a property of the algorithm, since the minimum was ill-conditioned.
---
### Author:  [@Tbl](https://discourse.julialang.org/u/Tbl)
#### Post date: [December 21, 2023, 4:29pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/13 "2023-12-21T16:29:05Z")
`LsqFit` is actually mainly an API for quick use of `Optim.jl` for fitting. I’d also recommend `GLM.jl` for anything that can be linearized, which is quite common.
---
### Author:  [@RomeoV](https://discourse.julialang.org/u/RomeoV)
#### Post date: [December 21, 2023, 4:48pm UTC](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/14 "2023-12-21T16:48:52Z")
I’ve found the pretty new NonlinearSolve.jl to be by far the best and easiest to use. [Nonlinear Least Squares Solvers · NonlinearSolve.jl](https://docs.sciml.ai/NonlinearSolve/stable/solvers/NonlinearLeastSquaresSolvers/)