# Fit Data Samples with a Robust Fit

I have a data from a sensor which the connection model of $$x$$ and $$y$$ is known:

For instance, in the case above, the model is linear.

The issue is how to handle outliers.
Specifically when there are many samples yet the model parameters is small.

Linear Least Squares fails in this case as the number of outliers can be up to 15-20%.

So the question is, what would be a robust fit which can handle large amount of data and a rate of up to 20% of outliers?

• Hi Royi! Good to see you. Commented Jan 23 at 22:59
• @robertbristow-johnson, Hi. Good to see you.
– Royi
Commented Jan 26 at 16:50

A simple and cheerful approach is to solve regular LMS for the sequence using different subsets of the total sample set, then pick the model that best satisfies some fitting criterion for the entire sample set.

Scales poorly to long sequences and large numbers of potential outliers.

• This is similar to RANSAC. Even less efficient, as in this case RANSAC will only need ~50-60 samples and the LS calculation will be negligible. I wanted something more bullet proof than iterative way to mark the outliers.
– Royi
Commented Jan 24 at 20:32

One way, which can handle ~30% of outliers is given by solving:

$$\arg \min_{\boldsymbol{p}} {\left\| \boldsymbol{X} \boldsymbol{p} - \boldsymbol{y} \right\|}_{1}$$

Using the $${L}^{1}$$ norm instead of the $${L}^{2}$$ norm means the fit is more robust to outliers.
The intuition is form the fact that the $${L}^{1}$$ norm minimizer is the median (See The Median Minimizes the Sum of Absolute Deviations (the $${\ell}_{1}$$ Norm)) while the $${L}^{2}$$ norm minimizer is the mean.

There are several ways to solve this:

1. Forming the problem as a Linear Programming problem.
2. Using the Sub Gradient method.
3. Using Iterative Reweighted Least Squares (IRLS).

The IRLS method is very quick to converge.
Given $$\boldsymbol{X} \in \mathbb{R}^{m \times n}$$ and $$n \ll m$$ then it can be solved pretty quickly as the most demanding step will be solving a linear system composed by an SPSD matrix of size $$n \times n$$.

For affine function, as above, it means solving a $$2 \times 2$$ system.

As can be seen, it works pretty well handling all the samples.
Solving this with efficient code is about ~100 [Mili Sec] for ~10,000 samples.
Run time could be even improved by a sub sampling.

The full Julia code is available on my StackExchange Signal Processing GitHub Repository (Look at the SignalProcessing\Q91788 folder).

• Wondering where the Convex.jl line is?
– Peter K.
Commented Jan 22 at 1:38
• @PeterK. , I load it at github.com/RoyiAvital/StackExchangeCodes/blob/…. Using it at github.com/RoyiAvital/StackExchangeCodes/blob/….
– Royi
Commented Jan 22 at 6:13
• @PeterK., At the graph you can't see it as the IRLS solution is on top of it. They are the same. I use Convex.jl to verify the implementation of the IRLS code.
– Royi
Commented Jan 22 at 6:15
• @PeterK., I created a tag named data-fitting. Could you make fit-data, data-fit and fitting-data its synonymous?
– Royi
Commented Jan 23 at 20:29
• No existing tag fit-data or fitting-data seem to be there, but I've merged data-fit onto data-fitting and made data-fit a synonym.
– Peter K.
Commented Jan 23 at 21:32