How to decide on Theil-Sen outliers? After using the Theil-Sen estimator to obtain a robust linear regression (in my case 2D), I would like to decide on a set of outlier points.
Is there a natural way to do this given the Theil-Sen solution? Or is it simply taking all points whose residual is beyond some arbitrary cutoff?
 A: I assume that you are familiar with the notion of breakdown point of an estimator. A similar concept exists for outlier identification rules (see [3]). 


*

*(a) The breakdown point of the Theil-Sen estimator at 2D data is $1-\frac{1}{\sqrt{2}}$. 

*(b) Furthermore, the 2D Theil-Sen estimator is residual admissible (meaning that it only depends on the data through the vector of residuals of the fit). 


Denote $\pmb e=\{e_i\}_{i=1}^n$ the $n$-vector of fitted residuals. Because of (a) and (b) flagging as outliers the observations with residual larger than the $q_h=1-(1-\frac{1}{\sqrt{2}})$ quantile of  $|\pmb e|=\{|e_i|\}_{i=1}^n$ will enssure that your outlier identification rule has the same breakdown point that the estimator it is based on. Using a higher value of $q_h$ will decrease the breakdown point (of your outlier identification rule) while using a lower value of $q_h$ will not increase it (the breakdown of your outlier identification rule is bounded from below by the breakdown point of the fit you use to identify the outliers).
Increasing $q_h$ will reduce the risk of misclassifying non-outlying observations as outliers (setting $q_h=1$ will make sure that no non-outlying observations are miss-classified as outliers). But you could obtain a much better result by using one step reweighing (1). With one step re weighting you can set the asymptotic risk of misclassifying non-outlying observations as outliers to any small value $\epsilon$ without affecting the breakdown point of your outlier identification rule (though you will increase the minimum distance outliers have to be from the bulk of the data to be identifiable, this distance has no bearing on the notion of breakdown point of an outlier identification rule). As a cost, you will need to add assumptions about the distribution of the vector of residuals. 
In any case there is a range of admissible trade-off between the two type of risks (misclassifying non-outlying data points as outliers and misclassifying outliers as non-outlying data points) for any value of $q_h$ above a threshold $q_h^*$ corresponding to the breakdown point of the initial estimator. For outlier identification rule based on the Theil-Sen estimator, $q_h^*=1-(1-\frac{1}{\sqrt{2}})$.   
Using the Theil-Sen estimator to find outlier is sub-optimal from a statistical point of view. You will get better trade-off terms as well as also the choice of having more robustness to outliers by using a more modern method such as FastLTS. FastLTS also includes a re-weighting step, but is based on more robust initial estimates than the Theil-Sen fit (so that the $q_h^*$ of FastLTS can be as high as $\approx 0.5$). See 2 for a recent review. 


*

*[1] P. Cizek (2010). Reweighted Least Trimmed Squares: An Alternative to One-Step Estimators. 

*[2] M. Hubert, P. J. Rousseeuw, and S. Van Aelst (2008). High-Breakdown Robust Multivariate Methods. 

*[3] C. Becker and U. Gather (1999). 
The Masking Breakdown Point of Multivariate Outlier Identification Rules.

