Leveraged least trimmed absolute deviations

The design of regression models that are not affected by outliers is an important task which has been subject of numerous papers within the statistics community for the last decades. Prominent examples of robust regression models are least trimmed squares (LTS), where the k largest squared deviations are ignored, and least trimmed absolute deviations (LTA) which ignores the k largest absolute deviations. The numerical complexity of both models is driven by the number of binary variables and by the value k of ignored deviations. We introduce leveraged least trimmed absolute deviations (LLTA) which exploits that LTA is already immune against y-outliers. Therefore, LLTA has only to be guarded against outlying values in x, so-called leverage points, which can be computed beforehand, in contrast to y-outliers. Thus, while the mixed-integer formulations of LTS and LTA have as many binary variables as data points, LLTA only needs one binary variable per leverage point, resulting in a significant reduction of binary variables. Based on 11 data sets from the literature, we demonstrate that (1) LLTA’s prediction quality improves much faster than LTS and as fast as LTA for increasing values of k and (2) that LLTA solves the benchmark problems about 80 times faster than LTS and about five times faster than LTA, in median.


Machine learning, regression and optimization
Machine learning is a fast growing field of research that inherits and combines methods from statistics, computer science and optimization to tackle a vast variety of applications like fraud detection, recommender systems, predictive maintenance and autonomous driving (Marsland 2015). One subfield of machine learning is supervised learning, whose task is to train a function on labeled data. This stands in contrast to other applications, like anomaly detection or clustering, where no labels are available and which are therefore examples of unsupervised learning (Bishop 2006). The most popular examples of supervised learning are classification and regression tasks. While classification aims at assigning discrete values to data points (e.g., binary values for cancer detection), regression methods train functions that assign continuous numbers to data points (e.g., prediction of house prices). Commonly used candidate mappings are (piecewise) linear functions, splines, tree-based models and neural networks (Clark and Pregibon 2015;Goldberg et al. 2021;Krasko and Rebennack 2017;Micula and Micula 2012;Rebennack and Kallrath 2015;Rebennack and Krasko 2020;Specht 1991).
The training procedure involves the minimization of a so-called loss function that measures the distance of the observations to the corresponding predictions. Minimizing the loss function results typically in an unconstrained smooth optimization problem that is tackled by variants of the stochastic gradient descent method which is a lightweight modification of gradient descent where only parts of the gradient are evaluated in each iteration (Schmidt et al. 2017;Robbins and Monro 1951). Other optimization-related topics within machine learning are Bayesian optimization (Snoek et al. 2012) or the optimization of pretrained machine learning models in the feature space (Thebelt et al. 2020b).
Mixed-integer linear optimization (MILO) models involve linear terms in the decision variables as well as integrality restrictions (for some) of the decision variables (Jünger et al. 2009;Wolsey and Nemhauser 1999). A very rich class of optimization problems in practice can be modeled using MILO models. Current state-of-the-art solvers for general MILO models use so-called branch-and-cut algorithms. The idea of branch-and-cut algorithms is to repeatedly solve linear optimization problems (these are easy to solve) which are obtained by relaxing the integrality restrictions on the decision variables. The linear optimization problems are updated by additional restrictions on the relaxed variables in order to cut out fractional values. This is called branching. In the worst case, there are exponentially many such branches in the number of integer variables. The branching is accompanied by cutting planes whose goal is to cut away fractional solutions (without the need to execute the costly branching). Therefore, as a general ruleof-thumb, fewer integer variables lead to lower computational times (though this is not always true). We make use of this observation in this paper.
Dimitris Bertsimas was one of the first researchers to point out that recent advances in linear and quadratic mixed-integer optimization have been rarely noticed in the statistics and machine learning communities. This inspired him to publish a series of papers under the motto "Machine Learning under a Modern Optimization Lens" that are summarized in the eponymous book (Bertsimas and Dunn 2019). Bertsimas' assessment was confirmed by some of the most renowned researchers in the statistics and machine learning community, Trevor Hastie and Robert Tibshirani, who state (Hastie et al. 2017): In exciting new work, Bertsimas et al. (2016) showed that the classical best subset selection problem in regression modeling can be formulated as a mixed integer optimization (MIO) problem. Using recent advances in MIO algorithms, they demonstrated that best subset selection can now be solved at much larger problem sizes that what was thought possible in the statistics community.
This paper was heavily inspired by Bertsimas' observation that mixed-integer optimization is still relatively unknown but can be applied to many optimization problems in the context of machine learning and statistics. Therefore, we present Leveraged Least Trimmed Absolute Deviations (LLTA), a mixed-integer based robust regression model, whose main idea we explain now.

Motivation
Let f ∶ ℝ n → ℝ, be a linear candidate function whose parameter ∈ ℝ n+1 we want to determine optimally with respect to some labeled training data with x j ∈ ℝ n for j = 1, … , N . The most popular idea of obtaining such a function f is referred to as Ordinary Least Squares (OLS) and goes back to Legendre or Gauß (Stigler 1981) at the end of the 18th century. Let be the residual of f with respect to the jth data point, j = 1, … , N . Then, OLS computes by solving the unconstrained convex quadratic optimization problem OLS is computationally attractive as it possesses a closed-form solution. However, it is very sensitive with respect to outliers. To soften this sensitivity to outliers, an alternative approach is to minimize the 1 -norm of the residual vector r = (r 1, , … , r j, , … , r N, ) instead of the 2 -norm. This results in the Least Absolute Deviations (LAD), a problem that was stated in 1757 by Boscovich (Koenker and Bassett 1985), even before OLS. LAD results in the unconstrained convex piecewise linear optimization problem In contrast to (1), LAD does not have a closed-form solution but can be reformulated and solved as a linear (continuous) optimization problem (LP) or tackled by a subgradient-based method. When only estimating 0 and all i = 0 , then this leads to the so-called location model (Bassett 1991). For this case, an optimal estimator for 0 is simply the median of the sorted data points y (j) , j = 1, … , N.
A crucial property of LAD is its robustness against so-called y-outliers. This is a consequence of Theorem 1 which states that LAD is not affected by changes in y i for data points that do not lie directly on the regression line as long as the signs of the residuals are not reverted.
Theorem 1 (Dodge 1997) Suppose ⋆ is a minimizer of Then, ⋆ is also a minimizer of However, while being robust to y-outliers, LAD is still affected by leverage points, i.e., outliers in x. This is illustrated in Fig. 1.
Inspired by this observation, Rousseeuw (1984) proposed the Least Trimmed Squares (LTS) in 1984, whose formulation as an optimization problem is given by the mixed-integer nonlinear optimization problem (MINLP) with k ∈ ℕ and n 2 < k < n , where we use " ⋅ " whenever multiplying decision variables. By design, (3) minimizes the sum of squares while ignoring the k largest squared deviations.
Similar to the LTS, in 1999, the least trimmed sum of absolute deviations (LTA) is proposed by Hawkins and Olive (1999). LTA can be formulated as the MINLP To compute the LTA regression, Hawkins and Olive propose an enumeration algorithm over all possible subsets with k elements. This algorithm is of particular interest for the location model, as the LTA for a fixed subset is then obtained by evaluating the N − k + 1 subsets of ordered data points y (k) , y (k+1) , … , y (k+h−1) , for all k = 1, … , N − k + 1 (Bassett 1991;Tableman 1994).
Next to an enumeration algorithm, the LTA regression problem is solved by Flores (2011) via a tailored continuous global optimization algorithm for the cases that the intercept is zero, i.e., 0 = 0 . First, the MINLP is reformulated as an NLP by introducing the nonconvex constraint b 2 j − b j = 0 for continuous variables b j instead of the binary restriction on b j . The resulting continuous nonconvex global optimization problem is then solved by a tailored global optimization algorithm in the spirit of Lasserre (2001).
The book about "open problems in optimization and data analysis" contains a chapter which discusses the connection between optimization and statistical robust estimators in the context of LTA regression (Pardalos and Migdalas 2018). For the location model, Zioutas et al. present a MINLP model of type (4) and a MILP reformulation of the bilinear terms using standard techniques. The LTA for the location model is extended to take into account outliers violating the correlational structure of the data set via a two-level approach in Chatzinakos et al. (2016).
Nowadays, within the statistics and machine learning communities, the LTS and LTA are solved by applying heuristics since both the LTS and the LTA are considered to be intractable as being NP-hard (Bernholt 2006). However, during the period 1991-2015, due to algorithmic advances, mixed-integer linear optimization problem (MILP) solvers have experienced an average speedup factor of 780,000, cf. Bertsimas et al. 2016; Bixby 2012 and the references therein. These machine-independent advances have been accompanied by an impressive progress in hardware performance. Thus, many real-world applications that could not be solved in the 1980s or 1990s are now solvable to global optimality within seconds. Similarly, modern software packages like CPLEX and GUROBI can now also solve large-scale nonconvex mixed-integer quadratic optimization problems. Despite the latest solver developments, LTS and LTA can only be solved for medium-sized problem instances. Therefore, we introduce Leveraged Least Trimmed Absolute Deviations (LLTA), which is a two-step approach that trains a linear function on possibly infiltrated data. The two steps are: 1. Identify all leverage points. 2. Minimize the total absolute deviations and ignore the k ∈ ℕ largest deviations to data points that are leverage points, for some chosen k with n 2 < k < n.
Consider now Fig. 2. LTS needs 11 binary decision variables and k = 3 to achieve a reasonable fit (Fig. 2a). LTA yields the same result with 11 binary decision variables and k = 1 (Fig. 2b). However, LLTA produces the same high-quality fit for k = 1 using only one binary decision variable (Fig. 2c). These indicated advantages of LLTA compared to LTS and LTA are further examined in Sect. 3 after a formal introduction of LLTA in Sect. 2.

Statement of contributions
The unique contributions of this paper are: We 1. introduce Leveraged Least Trimmed Absolute Deviations (LLTA), 2. demonstrate that LLTA outperforms LTS with respect to regression-quality and computational speed, 3. show that the regression-quality of LLTA is comparable to LTA while being much faster in terms of run time, 4. first benchmark the LTS and LTA with current MIQP solvers (the LTS and LTA are only solved by heuristic methods in the literature, ignoring the recent progress in MIQP algorithms and software developments).
The remainder of this paper is organized as follows. In Sect. 2, we introduce LLTA. We provide the benchmarking of LLTA with LTA and LTS in Sect. 3 before we conclude with Sect. 4.

Leveraged least trimmed absolute deviations
We start noting that the complexity of LTS and LTA is governed by 1. the number of ignored data points k since there are N k subsets of length k among the N data points and N k grows exponentially in k for fixed N and k < N∕2, 2. and the number of binary variables since the search space also grows exponentially in the number of binary variables.
To mitigate the computational complexity resulting from the second point, we introduce Leveraged Least Trimmed Absolute Deviations (LLTA). LLTA is a twostep procedure. Let D = {1, … , N} .
1. Compute the index set O ⊊ D of leverage points; see Sect. 2.2 for details.

Solve the optimization problem
Since the classical 1 -regression LAD is immune to y-outliers, we only protect our regression function with respect to leverage points. This is achieved through the parameter k allowing the optimal fit to ignore k data points within the index set of leverage points O . Note that we utilize here that the set of leverage points can be computed beforehand, while this is not possible for the y-outliers because they are regression-function dependent. Therefore, we obtain 1. a significant reduction of binary decision variables of (5) compared to (3) and (4), because LLTA only introduces one binary decision variable for each leverage point instead of one binary decision variable for each data point. 2. The possibility to choose smaller values of k compared to LTS, since LAD is already immune with respect to y-outliers. (2001), there are two approaches to statistical modeling which are very different: In the Data Modeling Culture, a stochastic data model with an underlying distribution is assumed whose distribution is estimated from successive draws. Examples for that are given in Liu (1996) and Vanhatalo et al. (2009). In the Algorithmic Modeling Culture, a function y = f (x) is fitted to observed data where the data generating process remains a black box with no further distributional assumptions. Many successful methods from machine learning like deep neural networks or gradient boosted trees are treated in the spirit of the latter culture. In this introductory work, we made the conscious decision to perform the analysis of LLTA within the framework of Algorithmic Modeling Culture. This yields a clear uncluttered overview about the main ideas and may serve as starting point for extensions from both cultures. To introduce underlying stochastic assumptions and to apply distributional-free sensitivity analysis using methods like bootstrapping are then possible extensions, cf. Sect. 2.5.

Remark 2
We assume that the number of infiltrated data points, i.e., the number of possible outliers, is strictly smaller than N/2. This is a standard assumption that is also posed in LTS and LTA. Therefore, also k must not exceed N/2 and we focus on the better half of the residuals as elaborated in Sect. 2.4.

Epigraph reformulation
In order to implement LTS, LTA and LLTA in a modern mixed-integer optimization solver, we first have to apply some reformulations. An LTS model is trained by solving the mixed-integer quadratically-constraint quadratic optimization problem (MIQCQP)

3
Leveraged least trimmed absolute deviations where we have avoided the trilinear terms r 2 j, ⋅ b j in the objective function by using only bilinear and quadratic expressions in the objective functions and constraints, respectively. Specifically, in an optimal solution, An LTA-estimate is computed as an optimal solution of the mixed-integer quadratic optimization problem (MIQP) The absolute value term |r j, | in the objective function is modeled through two linear constraints, for every j ∈ D . This is possible because LTS is a minimization problem.
Finally, we compute a linear regression function for LLTA by minimizing the MIQP where we rewrite the absolute value terms like in the MIQP for the LTA above.
The prefactors 1 N 2 and 1 N in the three formulations above do not affect the optimal solutions, but are added to enhance numerical stability. Because GUROBI version 9.0 is capable of dealing with MIQCQPs and MIQPs, there is no need to reformulate the optimization problems as MILPs. In this way, we avoid the introduction of Big-M constraints which are known to yield notoriously weak relaxations.

Computation of leverage points
Let x 1 , … , x j , … , x N ∈ ℝ n be the data points and q 0.
This leads us to the definition of the index set of all leverage points The definition of a leverage point given in this paper is more specific than commonly defined in existing literature where a leverage point is usually described as "data point that has an extreme value for one of the explanatory variables" (Dodge 1997).

Remark 3
The outlier tolerance 1.5 might be treated as an hyper parameter t of LLTA which influences the number of binary variables as well as the prediction quality. However, for the remainder of this work, we set t = 1.5 , which coincides with the definition of an outlier for boxplots (Tukey 1977).

Remark 4
Note that based on this definition of leverage points, there might be a tendency to identify more leverage points for high-dimensional data sets since the probability mass within a multidimensional probability distribution tends "to move away from its center" (van Handel 2014). The combination of domain knowledge and the selection of a tailored problem-specific outlier detection method (Hodge and Austin 2004) probably yields the best definition of leverage points for the problem at hand.

Choosing the number of outliers k
The choice of k may have a significant influence on the computed regression functions for LTS, LTA as well as LLTA. Quite generally, methods developed for LTS and LTA to choose k can also be applied toward LLTA.
By inspecting the optimization models (3), (4) and (5), we observe that their objective functions are monotone decreasing in the number k, i.e., allowing more outliers leads to a better fitting regression function. At the same time, increasing the number of outliers beyond the actual number of outliers in the data set implies loss of information. Consequently, k should not be chosen "too large." For practical problems, one would compute regression functions for different numbers of k and choose the k heuristically which seems to yield a good compromise between outlier detection and regression quality-one might choose k such that there is a significant improvement in the fit compared to k − 1 and where k + 1 yields only some minor improvement.
To choose k optimally, one would need an (objective) function quantifying both the regression fit and the "loss" from excluding potentially useful data.

Performance evaluation
Classical performance measures, like the root-mean-square error (RMSE) or meanabsolute error (MAE), are not suitable to measure the quality of statistical models in the presence of outliers because they evaluate the residuals for all data points. In contrast, a good robust model ignores some data points for being outliers. We evaluate the performance of the models by sorting the absolute residuals in ascending order, i.e., and computing the trimmed MAE on the better half of all residuals where ⌊N∕2⌋ denotes the floor function of N/2.
To motivate tMAE as performance metric, consider the following synthetic example where we have 30 "good" data points, ten "x-outliers," i.e., leverage points and ten "y-outliers." The scatter plot and regression lines for LLTA and LTS with k = 10 are depicted in Fig. 3.
It is not surprising to see that LTS is affected by the outliers, while LLTA yields a very good approximation of the "good" data points. However, how can we measure that in a multidimensional setting where visual inspections are not that easy to perform? We notice first that the classical performance measures RMSE and MAE are not suited since LTS outperforms LLTA in both, RMSE (3258 vs. 3660) and MAE (1471 vs. 1560), despite the fact that LTS' fit is obviously worse for this illustrative and synthetic example with 20 outliers. If we knew which of the data points are the "good" ones, we could just evaluate RMSE and MAE on these points. Unfortunately, we do not know that in real-world applications; otherwise, we would not have to immunize the regression function against outliers. Even worse, up to half of the data could be infiltrated and, in fact, we have an outlier rate of 40% in this example.
Let us take a look at the empirical distribution of the residuals of both methods which are illustrated in Fig. 4.
We recognize that the distributions are right-skewed-while the majority of all residuals are rather small, there exist some outliers, i.e., some data points show large residuals. This is not a problem per se since our goal is to design robust methods that ignore outliers on purpose. However, we should then also ignore these residuals when calculating our performance metric. Since this might affect up to 50% of our residuals, we decided to use trimmed MAE as performance metric for LLTA. If an upper bound b on the relative share of outliers is available, an obvious adjustment of tMAE would be to evaluate the residuals not only on the better half of the residuals but on (1 − b) ⋅ 100 % of the data.

Uncertainty measurement
The quantification of prediction uncertainty is important as it tries to determine the trustworthiness of a particular prediction or even of the statistical model in general. In particular, if decision making is based on good predictions, an uncertainty measure for the prediction is crucial.
A prominent example for model predictions, embedded in a "predict-tell" cycle, is Bayesian optimization. In Bayesian optimization, in each iteration an acquisition function is optimized, taking into account the prediction of a surrogate model s d Fig. 4 Empirical distribution of residuals for LLTA and LTS as well as the model uncertainty (Pelikan et al. 1999). The predominant surrogate models are Gaussian Processes which rely on a normal distribution assumption also yielding an uncertainty estimate.
However, more recent distributional-free approaches to Bayesian optimization also work with gradient tree ensemble methods as surrogate models and distancebased uncertainty measures (Thebelt et al. 2020a). A distance-based uncertainty measure does not assume any probability distribution in the data and is therefore also applicable to LLTA. The main idea is to measure the distance of an x-value, whose y-value is to be predicted, to the set of existing training data since the statistical model might have bad extrapolation properties. Distance-based uncertainty measures have a nice intuitive interpretation, but might provide misleading information for large datasets where especially the 2 -norm shows counterintuitive behavior (Aggarwal et al. 2001). Distance-based uncertainty measures that use the 1 -norm might be a useful uncertainty measure for LLTA.
Next to distance-based uncertainty measures, the second (by far more popular) approach to uncertainty estimation without distributional assumptions is bootstrapping (Diaconis and Efron 1983). Koenker and Hallock use bootstrapping in their famous work (Koenker and Hallock 2001) to quantify uncertainty in quantile regression which is a very popular approach to robust regression. The main idea of bootstrapping is to train statistical models on random samples of the training data and to compare its statistical properties. If models trained on different samples tend to vary much, then this might be an indication of a high model uncertainty. Therefore, using bootstrapping, or one of its many variants, is recommended as uncertainty measure for LLTA.

Computational results
We perform computational tests, comparing the two approaches LTS and LTA from the literature to the new model LLTA. All models are implemented in GUROBI 9.0 via its Python-API, and they are solved on a standard desktop computer possessing four cores each with 2.71 GHz and 16 GB RAM.
We start with an example on the body-brain data set in Sect. 3.1 with the goal to illustrate the usage of LLTA and to give an indication regarding the strengths and weaknesses of the examined models. We then use 11 instances from the literature in Sect. 3.2 to demonstrate the computational differences of the models resulting from LTS, LTA and LLTA.

Comparison with LTS and LTA on the body-brain data set
To illustrate the new model LLTA, we compare it to LTS and LTA based on the "Brain and Body Weights" dataset (Rousseeuw and Leroy 1987;Weisberg 1985). This dataset contains the following average body and brain weights in kg in the format (body, brain) for 65 animals which are depicted in Fig. 5a, b using logarithmic scales.
For the body-brain data set, we have n = 1 . Its quartiles are given by q 0.25 = 0.75 and q 0.75 = 60 which results in an interquartile range of iqr = 59.25 . We compute i.e., we have 13 leverage points with the x-values 465, 11700, 2547, 187.1, 521, 529, 207, 6654, 9400, 87000, 192, 250, 160. Due to the presence of these leverage points, OLS and LAD are heavily affected such that there is need for a robust statistical model. Therefore, we train LTS, LTA and LLTA on the data set for different values of k. In Fig. 6, we observe that the asymptotic trimmed mean absolute errors tMAE of the residuals of the three models is comparable, whereas LTA and LLTA obtain a better score for k < 25.
In addition to the good statistical performance of LLTA, we observe in Fig. 7 that LLTA outperforms LTS and LTA with respect to the number of visited nodes in the branch-and-bound tree and with respect to the run time (at a time limit of 600 s). LLTA visits at most 151 nodes and solves most optimization problems within the root node. In contrast, LTS and LTA visit up to 155,572 and 1,071,225 nodes, respectively. Regarding the run time, LLTA needs at most 0.04 s to solve , 6, 7, 8, 9, 12, 13, 15, 16, 26, 28, 52, 61}, any of the instances in contrast to LTA whose run time increases up to 17 seconds. In turn, LTA is much better than LTS, where an optimality certificate cannot be computed within the time limit. As such, we obtain a maximum speedup of 425 (or 99.8%) of LLTA compared to LTA and of 15,000 (or 99.99%) of LLTA to LTS. For some values k ≥ 16 , LTS does not manage to close the optimality gap within the time limit, as depicted in Fig. 8.

Datasets from the literature
We extracted 11 datasets from the existing literature on robust regression. Table 1 summarizes some relevant information about the datasets we use for benchmarking.
We run LTS, LTA and LLTA for all data sets at a time limit of 600 seconds for all k ∈ {0, … , ⌊N∕2 − 1⌋} for LTS and LTA as well as k ∈ {0, … , |O|} for LLTA. The results are depicted in Figs. 9,10,11,12,13,14,15,16,17,18 and 19, where we see the trimmed MAE, the number of visited branch-and-bound nodes and the run time for each method. Among the 11 datasets, we compare 118 different regression functions. 1  Siegel's exact fit example data 9 2 Rousseeuw and Leroy (1987) 10 Steam usage data (excerpt) 25 9 Norman and Draper (1981) 11 Modified data on wood specific gravity 20 6 Rousseeuw and Leroy (1987) Consider now Figs. 9a, 10a, 11a, 12a, 13a, 14a, 15a, 16a, 17a, 18a and 19a, where the tMAE is shown for different values of k. Among the 118 regression functions, the performance of LTS is never better than the one of LTA. This is because LTS is not immune against y-outliers and, thus, requires larger values of k to achieve a similar performance than LTA. For 68 instances, LLTA performs better than LTS and for k < 4 , LLTA is always better. Note that these results are heavily affected by dataset # 3 (Fig. 11a). LLTA performs comparable to LTA in most datasets. Exceptions are datasets # 2 (Fig. 10a), # 7 (Fig. 15a) and # 10 (Fig. 18a), where LTA is consistently better than LLTA and dataset # 9 (Fig. 17a), where LLTA outperforms LTA. It is remarkable that both LTA and LLTA have a quite similar performance, given that LTA has more degrees of freedom (because it can choose the leverage points among all data points compared to LLTA which is restricted to the index set of leverage points O ). Even more surprising is that LLTA shows a (strictly) better tMAE compared to LTA for 12 2 regression functions! Note that the computed regression functions are evaluated with respect to the tMAE (measuring the better half of the Fig. 9 Coleman data set residuals) and not with respect to the objective function (measuring all residuals except the k outliers). The difference in evaluation metrics and objective function also explains why the tMAE curves are not monotone decreasing in k (while we still observe a decreasing trend with an increase in k). This nonmonotone behavior can be seen, for example, in dataset # 1 (Fig. 9a). Figures 9b,10b,11b,12b,13b,14b,15b,16b,17b,18b and 19b show the number of visited branch-and-bound nodes, necessary to solve the corresponding instances. The trend here is clear: While LLTA shows a linear growth in k, both LTA and LTS follow an exponential curve. This behavior is by design, as the primary motivation to introduce LLTA is the significant reduction of the computational burden. In many instances, the LLTA can be solved in the root node.
The computational time is highly related to the number of visited branch-andbound nodes. Therefore,Figs. 9c,10c,11c,12c,13c,14c,15c,16c,17c,18c and 19c show a similar trend than the number of visited branch-and-bound nodes. In addition, we observe that the computational efforts to solve LTS are significantly greater than solving LTA. A similar percentage increase in runtime is observable for LTA compared to LLTA. Note that LLTA can solve any instance in at most 0.66 s, except the one instance for dataset #10 and k = 0 (Fig. 18c). In average, LLTA is 699.58 faster compared to LTA and 797.08 faster compared to LTS. However, the average speedup is mainly driven by dataset #3 where LLTA is 7543 faster than LTA Instance # 3 (Fig. 11a) is of particular interest due to its size. With 75 rows, this dataset contains about three times more rows than any other dataset (cf. Table 1). When inspecting Fig. 11b, c, we see that for k ≥ 11 , none of the instances for LTA and LTS can be solved to optimality within the time limit. This explains why we do not observe an exponential grows in the number of visited branch-and-bound nodes for large k for this dataset.

Summary and outlook
We introduce the Leveraged Least Trimmed Absolute Deviations (LLTA) which is based on the Least Trimmed Absolute Deviations (LTA). We make use of two observations. First, LTA is by design immune against y-outliers. Second, the leverage points can be computed beforehand in contrast to the y-outliers, because the y-outliers depend on the constructed regression function. As such, LLTA combines the advantages of LTA while considering only leverage points as potential x-outliers. This has the consequence that the proposed regression model LLTA is immune against both leverage points and y-outliers. At the same time, the computational burden is much lower compared to LTA.  . 19 Modified data on wood specific gravity Our computational results on known benchmark instances show that LLTA has a comparable performance of the computed regression models compared to LTA. At the same time, LLTA can be solved much faster compared to LTA. The computed regression models by LLTA tend to outperform the ones computed by the wellknown least trimmed squares (LTS). For small k, this effect is drastic. In addition, LLTA can be solved several orders of magnitude faster than LTS.