# How to speed up LinearModelFit, or use an alternative, for repeated use

I would like to implement a LOESS regression in Mathematica. My implementation depends on LinearModelFit, which is too slow to use on large datasets. (Can't use Fit because we need weights). Here is my current version:

loess[data_, order_, smoothing_] :=
Module[{d, weights, filtered, fData, fWeights, fit, x, weight, funcs},
weight[z_] := Piecewise[{{(1 - Abs[z]^3)^3, Abs[z] < 1}}];
funcs = Table[x^i, {i, 0, order}];
Table[
d = (data[[;; , 1]] - point[[1]])/smoothing;
weights = weight /@ d;
filtered =
Cases[Transpose[Join[Transpose[data], {weights}]],
a_ /; a[[3]] != 0];
fData = filtered[[;; , {1, 2}]];
fWeights = filtered[[;; , 3]];
fit = LinearModelFit[fData, funcs, x, Weights -> fWeights]
{point[[1]], fit[point[[1]]]}
, {point, data}]
]


With an example dataset:

data = Table[{x, 2*x^3 - 1.5 x^2 + RandomReal[{-0.05, 0.05}]}, {x, 0, 1, 0.001}];


running loess[data, 2, 0.1] takes about 7 seconds on my machine.

Is there a way to speed up this function? The accepted answer will also allow access to the slope of the LOESS regression; for example by returning {point[[1]], fit'[point[[1]]]} in the above implementation.

(Note that my definition of the smoothing parameter is non-standard.)

• I don't feel it's satisfying, but: Replacing LinearModelFit with NonlinearModelFit reduces the computation time on my machine from 13s to 9s. – Oscillon Sep 2 '16 at 20:16
• Reading your motivation for implementing LOESS regression, are you interested in any other algorithms to be applied that would be fast enough? – Anton Antonov Sep 2 '16 at 21:29
• I would certainly be interested in a completely different implementation of LOESS, if it were faster. Is that what you had in mind? If not, still curious what you are thinking of. – mdeceglie Sep 2 '16 at 21:44
• If Fit works faster, you can simply use data * weights instead of data for weighted fit. – Alexey Popkov Sep 3 '16 at 1:13
• You went to the trouble of creating filtered, fData and fWeights but then used data and weights in the fitting process. Typo? I get the same answer with both approaches with a factor of two speedup using the fData and fWeights. This would become more significant the larger the data set. – Jack LaVigne Sep 5 '16 at 2:51

As per discussion in the comments, here is an answer using Quantile Regression that produces a fitted curve for less than 0.5 seconds on my computer. (The function loess in the question took 8 seconds.)

Generate the data (as in the question):

data = Table[{x, 2*x^3 - 1.5 x^2 + RandomReal[{-0.05, 0.05}]}, {x, 0, 1, 0.001}];


Import["https://raw.githubusercontent.com/antononcube/MathematicaForPrediction/master/QuantileRegression.m"]


Do quantile regression with 12 B-spline basis nodes for the quantile 0.5:

AbsoluteTiming[
qfunc = QuantileRegression[data, 12, {0.5}][[1]];
]

(* {0.343834, Null} *)


Plot the data and the regression quantile curve:

ListPlot[{data, {#, qfunc[#]} & /@ data[[All, 1]]}]