# Least squares with SciML NonlinearSolve.jl - how to?

**URL:** <https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074>\
**Category:** General Usage\
**Tags:** nonlinearsolve, least-squares\
**Created:** [April 17, 2024, 2:25pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074 "2024-04-17T14:25:02Z")\
**Posts on this page:** 14\
**Page:** 1

<div class="post-metadata">

**Author:** ![Eben60](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/eben60/32/13475_2.png) [@Eben60](https://discourse.julialang.org/u/Eben60)\
**Post date:** [April 17, 2024, 2:25pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/1 "2024-04-17T14:25:02Z")

</div>

Trying to get how to use it - reading [Documentation](https://docs.sciml.ai/NonlinearSolve/stable/tutorials/getting_started/#Problem-Type-4:-Solving-Nonlinear-Least-Squares-Problems) unfortunately doesn’t help.

The use case - a model described by a formula, and a vector of (t, y) value pairs representing e.g. experimental data. Like the example problem in the [docs](https://docs.juliahub.com/General/LsqFit/stable/) of the LsqFit.jl package.

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [April 17, 2024, 3:14pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/2 "2024-04-17T15:14:31Z")

</div>

Can you show what you tried so far?

---

<div class="post-metadata">

**Author:** ![Eben60](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/eben60/32/13475_2.png) [@Eben60](https://discourse.julialang.org/u/Eben60)\
**Post date:** [April 17, 2024, 4:49pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/3 "2024-04-17T16:49:39Z")

</div>

I’ve tried to read the whole package documentation and analyze the example, spend all in all probably an hour, made different attempts to write an example of my own (which errored, and which I probably can’t reproduce now), and failed. The example problem in the documentation has a different structure from the (quite common) problem I’d like to be able to solve. I could find e.g. no explanation why ‘NonlinearFunction’ should be there.

There is an obvious mismatch between me and the docs, and you may consider the problem is not the documentation. On the other side, I had no such problems with the docs of some other SciML (sub-) packages, including `Problem Type 1: Solving Nonlinear Systems of Equations` from the same package.

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [April 17, 2024, 5:39pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/4 "2024-04-17T17:39:57Z")

</div>

```julia
using NonlinearSolve

function nlls!(du, u, p)
    du[1] = 2u[1] - 2
    du[2] = u[1] - 4u[2]
    du[3] = 0
end

u0 = [0.0, 0.0]
prob = NonlinearLeastSquaresProblem(
    NonlinearFunction(nlls!, resid_prototype = zeros(3)), u0)

solve(prob)

```

Here’s an example. What did you try? I’ll need a bit more feedback from you to figure out where the disconnect is here, like see what you misinterpreted from some example.

---

<div class="post-metadata">

**Author:** ![Eben60](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/eben60/32/13475_2.png) [@Eben60](https://discourse.julialang.org/u/Eben60)\
**Post date:** [April 17, 2024, 6:04pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/5 "2024-04-17T18:04:50Z")

</div>

E.g. this:

```julia
function nlls!(du, u, p)
    du[1] = 2u[1] - 2
    du[2] = (u[1] - 4u[2])^2 + 0.1
    du[3] = 0
end

```

`ERROR: LinearAlgebra.LAPACKException(2)`

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [April 17, 2024, 6:19pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/6 "2024-04-17T18:19:13Z")

</div>

@avikpal maybe we should default to pivoted QR?

---

<div class="post-metadata">

**Author:** ![Eben60](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/eben60/32/13475_2.png) [@Eben60](https://discourse.julialang.org/u/Eben60)\
**Post date:** [April 17, 2024, 6:25pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/7 "2024-04-17T18:25:42Z")

</div>

The function `nlls!` has a root `du.==0`. A typical least squares problem doesn’t. Why such an example?

---

<div class="post-metadata">

**Author:** ![avikpal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/avikpal/32/6550_2.png) [@avikpal](https://discourse.julialang.org/u/avikpal)\
**Post date:** [April 17, 2024, 6:31pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/8 "2024-04-17T18:31:05Z")

</div>

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

---

<div class="post-metadata">

**Author:** ![Eben60](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/eben60/32/13475_2.png) [@Eben60](https://discourse.julialang.org/u/Eben60)\
**Post date:** [April 17, 2024, 6:41pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/9 "2024-04-17T18:41:01Z")

</div>

to get an error, I’ve changed the function `nlls!`

```julia
du[2] = (u[1] - 4u[2])^2 + 0.1

```

---

<div class="post-metadata">

**Author:** ![avikpal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/avikpal/32/6550_2.png) [@avikpal](https://discourse.julialang.org/u/avikpal)\
**Post date:** [April 17, 2024, 6:52pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/10 "2024-04-17T18:52:05Z")

</div>

Probably, change the linear solve default right?

---

<div class="post-metadata">

**Author:** ![Paul\_Soderlind](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/paul_soderlind/32/1753_2.png) [@Paul\_Soderlind](https://discourse.julialang.org/u/Paul_Soderlind)\
**Post date:** [April 17, 2024, 7:01pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/11 "2024-04-17T19:01:53Z")

</div>

I use OSQP for things like that. See, for instance, this notebook for a restricted LS: [FinancialTheoryMSc/Ch08\_PerfEval.ipynb at master · PaulSoderlind/FinancialTheoryMSc · GitHub](https://github.com/PaulSoderlind/FinancialTheoryMSc/blob/master/Ch08_PerfEval.ipynb)

---

<div class="post-metadata">

**Author:** ![Eben60](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/eben60/32/13475_2.png) [@Eben60](https://discourse.julialang.org/u/Eben60)\
**Post date:** [April 18, 2024, 4:12pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/12 "2024-04-18T16:12:49Z")

</div>

Paul thank you for the suggestion. Actually there are quite a few packages out there, as discussed in [this thread](https://discourse.julialang.org/t/alternatives-to-lsqfit-for-nonlinear-curve-fitting/107685/10) as well as [here](https://discourse.julialang.org/t/comparing-non-linear-least-squares-solvers/104752).

For most cases, LsqFit.jl just does the job, however I got some issues with it, and @ChrisRackauckas recommended to avoid it.

What I need is just a package which I could use for simple _curve fitting_. If there is a SciML package doing what I need, that’d be my default choice.

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [April 18, 2024, 5:00pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/13 "2024-04-18T17:00:48Z")

</div>

> [@Eben60](#):
>
> What I need is just a package which I could use for simple _curve fitting_. If there is a SciML package doing what I need, that’d be my default choice.

There’s ongoing effort for that: [Migrate to SciML and extend? · Issue #35 · pjabardo/CurveFit.jl · GitHub](https://github.com/pjabardo/CurveFit.jl/issues/35). We’re close.

Your original problem here just had something in the defaults, which we can handle today. [Change the default QR to ColumnNorm by avik-pal · Pull Request #492 · SciML/LinearSolve.jl · GitHub](https://github.com/SciML/LinearSolve.jl/pull/492) May take a bit of discussion though.

---

<div class="post-metadata">

**Author:** ![Eben60](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/eben60/32/13475_2.png) [@Eben60](https://discourse.julialang.org/u/Eben60)\
**Post date:** [April 18, 2024, 9:18pm UTC](https://discourse.julialang.org/t/least-squares-with-sciml-nonlinearsolve-jl-how-to/113074/14 "2024-04-18T21:18:28Z")

</div>

OK, now I got it.

```julia
using NonlinearSolve, Plots

function curve_fit(model, u0, xdata, ydata, p)
    data = (xdata, ydata, p)

    function lossfn!(du, u, data)
        (xs, ys, p) = data   
        du .= model.(xs, Ref(u), Ref(p)) .- ys
        return nothing
    end

    prob = NonlinearLeastSquaresProblem(
        NonlinearFunction(lossfn!, resid_prototype = similar(ydata)), u0, data)
    sol = solve(prob)
    u = sol.u
    fit = model.(xdata, Ref(u), Ref(p))
    return (;sol, fit)
end

```

define a model:

```julia
function expmodel(x, u, t₀=0)
    y0 = u[1]
    a = u[2]
    τ = u[3]
    return y0 + a * exp(-(x-t₀)/τ)
end

```

synthesize data:

```julia
y0, a, τ, t₀ = 0.0, 1.0, 2.0, 0
xs = 0:0.1:3;
ys0 = expmodel.(xs, Ref([y0, a, τ]), Ref(t₀));
ys = @. ys0 + 0.02*randn()

```

Solve, plot:

```julia
(;sol, fit) = curve_fit(expmodel, [y0, a, τ], xs, ys, t₀)

plot(xs, [ys, fit])

```

 ![fit](https://global.discourse-cdn.com/julialang/original/3X/8/3/831ec37cdc7cc2059155bf129b9e9af64e93560c.png)  
and enjoy 🙂
