Three weeks ago I replied to a question in the Wolfram Community site about finding a function fit to points of non-linear data.

Here is the data:

I proposed a different way of doing the non-linear fitting using Quantile regression with B-splines. The advantages of the approach are that it does not need guessing of the fit functions and it requires less experimentation.

Load the package QuantileRegression.m:

`Needs["QuantileRegression`"]`

After several experiments with the number of knots we can find a regression quantile using 3d order B-splines of 18 knots that approximates the data well.

qfunc = Simplify[QuantileRegression[data, 18, {0.5}, InterpolationOrder -> 2][[1]]];

```
```

`qfGr = Plot[qfunc[x], {x, Min[data[[All, 1]]], Max[data[[All, 1]]]}, PlotPoints -> 700];`

Show[{dGr, qfGr}, Frame -> True, ImageSize -> 700]

Here is the relative error estimate:

ListPlot[Map[{#[[1]], (qfunc[#1[[1]]] - #1[[2]])/#1[[2]]} &, data],

Filling -> Axis, Frame -> True, ImageSize -> 700]

I would like to point out that if we put large number of knots we can get aliasing errors effect:

qfunc = Simplify[

QuantileRegression[data, 40, {0.5}, InterpolationOrder -> 2][[1]]];

```
```

`qfGr = Plot[qfunc[x], {x, Min[data[[All, 1]]], Max[data[[All, 1]]]}, PlotPoints -> 700];`

Show[{dGr, qfGr}, Frame -> True, ImageSize -> 700]

If we select a non-uniform grid of knots according to gradients of the data we get a good approximation with one attempt.

In[107]:= knots =

Join[

Rescale[Range[0, 1, 0.2], {0, 1}, {Min[data[[All, 1]]],

Max[data[[All, 1]]] 3/4}],

Rescale[Range[0, 1, 0.1], {0, 1}, {Max[data[[All, 1]]] 3/4,

Max[data[[All, 1]]]}]

]

```
```Out[107]= {0., 3.75, 7.5, 11.25, 15., 18.75, 18.75, 19.375, 20., 20.625, 21.25, 21.875, 22.5, 23.125, 23.75, 24.375, 25.}

`In[108]:= qfunc = QuantileRegression[data, knots, {0.5}][[1]];`

The grid lines in the plot below are made with the knots.

qfGr = Plot[qfunc[x], {x, Min[data[[All, 1]]], Max[data[[All, 1]]]}];

Show[{dGr, qfGr}, GridLines -> {knots, None},

GridLinesStyle -> Directive[GrayLevel[0.8], Dashed], Frame -> True,

ImageSize -> 700]

Here are the relative errors:

ListPlot[Map[{#[[1]], (qfunc[#1[[1]]] - #1[[2]])/#1[[2]]} &, data],

Filling -> Axis, Frame -> True, ImageSize -> 700]

Here is the approximation function after using PiecewiseExpand and Simplify:

qfunc = Evaluate[Simplify[PiecewiseExpand[qfunc[[1]]]]] &;

qfunc

### Follow-up

If it is desired to make the fitting with a certain, known model function then the function `QuantileRegressionFit`

can be used. `QuantileRegressionFit`

is provided by the same package, QuantileRegression.m, used in the blog post.

This older blog post of mine discusses the use of `QuantileRegressionFit`

and points to a PDF with further details (theory and *Mathematica* commands).