Cataloged from PDF version of article.The robust linear regression problem using Huber's piecewise-quadratic M-estimator
function is considered. Without exception, computational algorithms for this problem have been
primal in nature. In this note, a dual formulation of this problem is derived using Lagrangean
duality. It is shown that the dual problem is a strictly convex separable quadratic minimization
problem with linear equality and box constraints. Furthermore, the primal solution (Huber's
M-estimate) is obtained as the optimal values of the Lagrange multipliers associated with the dual
problem. As a result, Huber's M-estimate can be computed using off-the-shelf optimization software

Pinar, M. C.

Applied Mathematics Letters

AbstractThe robust linear regression problem using Huber's piecewise-quadratic M-estimator function is considered. Without exception, computational algorithms for this problem have been primal in nature. In this note, a dual formulation of this problem is derived using Lagrangean duality. It is shown that the dual problem is a strictly convex separable quadratic minimization problem with linear equality and box constraints. Furthermore, the primal solution (Huber's M-estimate) is obtained as the optimal values of the Lagrange multipliers associated with the dual problem. As a result, Huber's M-estimate can be computed using off-the-shelf optimization software

Pinar, M.Ç.

Elsevier - Publisher Connector 

Duality in robust linear regression using Huber's M-estimator 

Bilkent University Institutional Repository

~ )  Pergamon Appl. Math. Lett. Vol. 10, No. 4, pp. 65-70, 1997 Copyright©1997 Elsevier Science Ltd Printed in Great Britain. All rights reserved 0893-9659/97 $17.00 + 0.00 P Ih  S0893-9659(97)00061-X Dual i ty  in Robust  L inear Regress ion  Us ing  Huber ' s  M-Est imator  M. (~. P INAR Industrial Engineering Department, Bilkent University 06533 Bilkent, Ankara, Turkey must af ap©bilkent, edu. tr (Received October 1996; accepted December 1996) Communicated by R. H. Byrd Abst ract - -The  robust linear regression problem using Huber's piecewise-quadratic M-estimator function is considered. Without exception, computational gorithms for this problem have been primal in nature. In this note, a dual formulation of this problem is derived using Lagrangean duality. It is shown that the dual problem is a strictly convex separable quadratic minimiza- tion problem with linear equality and box constraints. Furthermore, the primal solution (Huber's M-estimate) is obtained as the optimal values of the Lagrange multipliers associated with the dual problem. As a result, Huber's M-estimate can be computed using off-the-shelf optimization software. Keywords- -Lagrangean duality, Huber's M-estimator, Robust regression, Quadratic program- ming. 1. INTRODUCTION There has been considerable interest in the theory and algor i thms for robust  est imat ion in the past  two decades. In part icular ,  Huber 's  M-est imator  [1] has received a great  deal of at tent ion from both theoret ical  and computat ional  points of view. Robust  est imat ion is concerned with identifying "outliers" among data  points and giving them less weight. Huber 's  M-est imator  is essential ly the least squares est imator,  which uses the g l -norm for points that  are considered outl iers with respect to a certain threshold. Hence, the Huber  criterion is less sensit ive to the presence of outl iers. More precisely, Huber 's  M-est imate is a minimizer x* E ~n of the function where f(x) = i----1 (1) 12 ~-=.t , if Itl < % p(t) = "Y 1 (2) I t l -  2-r, if Itl _~ % with a tuning constant ~ > O, and a scaling factor a that  depends on the data  to be est imated.  The residual r~(x) is defined as r~(x) -- b~ - a~x, (3) Special thanks are due to R. H. Byrd for his careful handling of the paper. One referee pointed out the QP form of the primal added to the revised version. Another eferee clarified several issues concerning robust estimators and data analysis, and suggested directions for future research. I am grateful to both. 65 66 M. ~. PINAR for all i = 1, . . .  ,m with r -- b - ATx .  To view this minimization problem in a more familiar format, define a "sign vector" ST(X) = [S~I(X),...,S~m(X)] (4) with and where -1,  if ri(x) _< -7 ,  sv~(x) = 0, if Iri(x)[ < 7, 1, if ri(x) >_ 7, (5) Ws = diag (w l , . . . ,  win), (6) wi = 1 - s~. (7) Now, assuming a unit a for the moment, Huber's M-estimation problem can be expressed as the following minimization problem. PROBLEM [P]. 1 T T[ 1 ] minimize F(x )  - ~Tr W, r  + s~ r - ~Ts~ , (8) where the argument x of r is dropped for notational convenience. Clearly, F measures the "small" residuals ([ri(x)[ < 7) by their squares, while the "large" residuals are measured by the gl func- tion. Thus, F is a piecewise quadratic function, and it is once continuously differentiable in Nn. The contribution of the present note is to introduce a dual approach to this computationally intensive problem. A simple derivation is used to give a dual problem which turns out to be a problem familiar to the numerical optimization community: linearly constrained separable convex quadratic programming. To the best of the author's knowledge, this duality relationship has not been noticed before. This simple result reduces the M-estimation problem to one that is easily solved using off-the-shelf optimization software. The interest in Huber's M-estimator derives from its perceived utility as a robust estimation procedure. In this context, 7 is an important quantity as it controls the spread of the residuals. This suggests that it should be related to the scaling factor a. Some algorithms [2,3] estimate aonly once at the beginning using some rules of thumb, and keep it fixed throughout the compu- tation. Another approach is to estimate a iteratively by treating it as an independent variable. This is used by Huber in [4] who suggests a to be computed from an auxiliary equation involving a, 7, and x. Shanno and Rocke [5] and Ekblom [6] also use this approach in their respective papers. An alternative way to connect hese two parameters and to iteratively estimate a is to set a = 1 for any values of 7 and a, and to replace 7 with 7a. Hence, there is no loss of infor- mation in using 7 only, while 7 and a are both allowed to vary during the computation. Clark and Osborne [7] use this observation in their algorithm that computes both a and V. This is a continuation method, so that at each stage they have the M-estimate for the current 7. Clark and Osborne also describe a partitioning algorithm to resolve degenerate situations in the continua- tion algorithm. The continuation method is essentially based on tracing a curve of M-estimates as a function of 7. Their approach consists of two stages. 1. Construction of a solution for a particular 7. The easiest cases are 7 = 0 (~1 estimate) and 7 = oo (least squares). These are both used in [7]. 2. Continuation with respect o 7 to solve the data analysis problem. A potential contribution of the present paper in the above context is that it suggests a new method for obtaining Huber's M-estimate for a particular value of 7. Then, continuation with respect o 7 of a solution obtained by the dual method can be pursued using methods of para- metric quadratic programming. In this connection, the rule of thumb choices that have been Huber's M-Estimator 67 suggested for ~ in the literature [2-4] can be used to see if they provide better starting points for the continuation method than least squares and ei estimates. This is a topic for future research. When a and 9, are fixed, most algorithms that have proved successful for the computation ofHuber's M-estimate have been iterative in nature. Since F is only once continuously differentiable, research concentrated ondeveloping successful applications of Newton's method to this problem, and on studying ideas such as the iteratively reweighted least squares (IRLS). Among these, Huber and Dutter's method [8] and Huber's [4] apply Newton's method to the nonlinear equation system: which represent the optimality conditions. The IRLS algorithm is attributed to Beaton and Tukey [9]. This algorithm has been discussed in several other papers as well [2,3,10-12]. A brief review of these algorithms i given in [7]. Madsen and Nielsen gave finite modified Newton algorithms for the minimization of the Huber function in [13]. These algorithms capitalize on the piecewise quadratic nature of the function. The idea is that if the sign vector associated with a minimizer x*, s* say, were known, then the minimizer could be computed using one step of Newton's method. Since sign vectors correspond to a subdivision of ~n into a finite number of subregions, the methods of Madsen and Nielsen reduce to a search for the correct subregion. 2. A DUAL  PROBLEM In this section, a dual problem to [P] is derived using Lagrangean duality. The interested reader is directed to the book by Rockafellar [14] for a detailed exposition of Lagrangean duality. We use a single parameter 7 to mean 3'a as discussed in the previous ection. Consider the problem [P] in a slightly different form: minF(x) =_ @ (b - AT  x)  . (10) Let u = b - ATx ,  and rewrite the problem as: min @(u) s.t. b - ATx  = u. Associating the multipliers y E ~m with the equality constraints b-ATx  = u, we get the following Lagrangean problem: max min {@(u) + yT (b - AT  x -- u) } .  (11) X,~ This is equivalent to: max [yTb + min ((I)(u)- yTu} ÷ m~n (_yT A Tx}]. y u Observing that x is a free variable, the term min {_yT AT x} yields the constraint It remains to simplify the term Ay=O.  (12) (13) (i4) min - ysu}. (15) 68 M. ~. PINAR Simple calculus hows that min{O(u)--yTu)} = --~TY Y, U --OO, Hence, the dual problem is the following. PROBLEM [D]. 1 T max bmy-  ~TY Y i f -1  < y _< 1, otherwise. (16) s.t. Ay = 0 - l<y<l .  An alternative way to arrive at the above dual is to pose the primal problem [P] as a quadratic programming problem: 1 m ~(  ) max + q ' -  i=.l i=l  s.t. -p  - q < b - ATx  < p -{- q p<_Te q>_O. where e denotes a vector with all components unity. This is not surprising since the dual of a quadratic program is again a quadratic program. However, this alternative derivation is sub- stantially longer since it requires ome transformations on the resulting dual to be cast in the form [D]. So, it is not included into the present paper. The remarkable fact about the duality result is that the dual we have derived is a quadratic programming problem with a strictly concave separable objective function, linear equality, and box constraints. This is important in two respects. First, this dual problem is an intensely researched, numerically well-solved problem. An excellent software system is available from Stanford University for the solution of quadratic programming problems [15,16]. Furthermore, numerical procedures for solving this problem are part of almost any subroutine library. Such procedures are also available through matrix manipulation packages such as Matlab TM and Oc- tave [17]. Second, it has been shown to be polynomially solvable [18], and efforts to turn related algorithms into reliable software on this front are also under way. In particular, quadratic pro- gramming versions of the software systems CPLEX TM and LoQo are available [19]. Hence, any advances made in the fast and accurate solution of the quadratic programming problem will benefit Huber's M-estimation problem. Note that any numerical procedure which yields a primal-dual solution to [D] gives a mini- mizer x* of the primal problem, which is really of interest here. To see this, it is instructive to derive a dual problem to [D] using Lagrangean duality. Associating multipliers v E !l~" with the constraints Ay = O, we get the following Lagrangean problem: min max ~ b T y - 1 T } v -*<_~<* [. ~TY Y +vT( - -Ay+O)  . (17) Now, rearranging terms and using r =- b - ATv,  simple calculus hows that: max y~r~ - = -1<~,<1 ~7Y~ 27r"  (18) 57, > 7. But, this is precisely Huber's M-estimator function. Therefore, the optimal values of the dual multipliers associated with the equality constraints in [D] give precisely Huber's M-estimate. Huber's M-Estimator 69 Further insight into the relationship between [P] and [D] is gained through Theorem 1 below which links the opt imal  solutions of [P] and [D]. This theorem shows that  the opt imal  solution to [D] is obtained as the first derivative of the function F with respect o r at any pr imal opt imal  point. Before stat ing this result, we quote the following property that  was proved in [20]. LEMMA 2.1. Let  x be a minimizer  of  F for some value of  ~ > O, and let s = sx(x) with Ws defined accordingly. Also, let r = b - AVx.  Then, ri is constant for all i such that si = O. Furthermore,  sw is constant for any x that  minimizes F.  By strict concavity, it is clear that  the optimal solution to [D] is unique. THEOREM 2.1. Let  x be a minimizer  o fF  for some value of  "l > O, and let s = s~(x) with Ws defined accordingly. Then, the unique optimal solution of [D] is given as: w.r (x )  y - - -  + s. (19) PROOF. Associate multipliers x with the equality constraints and nonnegative multipliers a,/3 E N'~ with the box constraints. From [6, Theorem 28.3], it is well known that  for y to be an opt imal  solution for [D], and for (x, a,/3) to be a Lagrange multiplier vector, it is necessary and sufficient that  (y, x, a,/3) be a saddlepoint of the Lagrangean of [D] as defined in [6, p. 280]. This condition holds if and only if the components of (y, x, a,/3) satisfy the following conditions: ~y-  b+ ATx  + a- /3  = 0, (20) Ay = 0, (21) ai(yi  - 1) = 0, for all i = 1 . . . . .  m, and (22) /3i(-y~ - 1) -- 0, for all i = 1 . . . . .  m. (23) From the first condition (20), we obtain b - A-rx o~ - /3  u - (24) Let r - b - A-r x and consider three cases. CASE 1. If  --1 < Yi < 1, clearly, ai =/3i = 0. This implies that  [ri(x)[ < ~ with s~i(x) = 0, i.e., = CASE 2. If  Yi = 1, this implies that /3 i  = 0. Hence, r i (x)  = ~ + ai. Therefore, r i (x)  >_ "y with s~i(x) = 1, i.e., yi = s~i(x) = 1. CASE 3. If  Yi = --1, one has ai  = 0. Hence, similarly to Case 2 above, r i (x)  <_ -~ with s~i(x) = -1 ,  i.e., Yi = s~(x)  = -1 .  Therefore, an alternative xpression for y is given as: + (25) The result now follows from the previous lemma. | 70 M. ~. PINAR REFERENCES I. P. Huber, Robust Statistics, John Wiley, New York, (1981). 2. J.B. Birch, Some convergence properties of iterated reweighted least squares in the location model, Comm. Statist. B9, 359-369 (1980). 3. P.W. Holland and R.E. Welsch, Robust regression using iteratively reweighted least squares, Comm. Statist. A6, 813-827 (1977). 4. P. Huber, Robust regression: Asymptotics, conjectures and Monte Carlo, Annals of Statistics 1, 799-821 (1973). 5. D.F. Shanno and D.M. Rocke, Numerical methods for robust regression: Linear models, SIAM J. on Scientific and Statistical Computing 7, 86-97 (1986). 6. H. Ekblom, A new algorithm for the Huber estimator in linear models, BIT 28, 123-132 (1988). 7. D.I. Clark and M.R. Osborne, Finite algorithms for Huber's M-Estimator, SIAM J. on Scientific and Sta- tistical Computing 7, 72-85 (1986). 8. P. Huber and R. Dutter, Numerical solution of robust regression problems, In COMPSTAT I974 Proc. Symposium on Computational Statistics, (Edited by G. Brushmann), pp. 165-172, Physike Verlag, Berlin, (1974). 9. A.E. Beaton and J.W. Tukey, The fitting of power series, meaning polynomials, illustrated on band- spectroscopic data, Technometrics 16, 147-185 (1974). 10. J.B. Birch, Effects of the starting value and stopping rule on robust estimates obtained by iterated weighted least squares, Comm. Statist. Bg, 133-140 (1980). 11. R.H. Byrd and D.A. Pyne, Some results on the convergence of the iteratively reweighted least squares algorithm for robust regression, In Proc. Statistical Computing Section, pp. 87-90, American Statistical Association, (1979). 12. S.A. Ruziusky and E.T. Olsen, L1 and Leo minimization via a variant of Karmarkar's algorithm, IEEE Transactions on Acoustics, Speech and Signal Processing 37, 245-253 (1989). 13. K. Madsen and H.B. Nielsen, Finite algorithms for robust linear regression, BIT 30, 682-699 (1990). 14. tLT. Rockafellar, Convex Analysis, Princeton University Press, Princeton, N J, (1970). 15. P.E. Gill, W. Murray, M.A. Saunders and M.H. Wright, User's guide for SOL/QSPOL (Version 3.2), Tech- nical Report SOL 84-6, Systems Optimization Laboratory, Department of Operations Research, Stanford University, Stanford, CA, (1984). 16. P.E. Gill, S.J. Harnmarling, W. Murray, M.A. Saunders and M.H. Wright, User's guide for LSSOL (Version 1.0): A Fortran package for constrained linear least squares and convex quadratic programming, Technical Report SOL 86-1, Systems Optimization Laboratory, Department of Operations Research, Stanford Univer- sity, Stanford, CA, (1986). 17. J.W. Eaton, OCTAVE: A High-Level Interactive Language for Numerical Computation, (1994). 18. Y. Ye and E. Tse, An extension of Karmarkar's projective algorithm for convex quadratic programming, Mathematical Programming 44, 157-179 (1989). 19. R.J. Vanderbei and T.J. Carpenter, Symmetric indefinite systems for interior point methods, Mathematical Programming 58, 1-32 (1993). 20. K. Madsen, H.B. Nielsen and M.(~. Pmar, New characterizations of £1 solutions to overdetermined systems of linear equations, Operations Research Letters 16, 159-166 (1994). 21. D.I. Clark, The mathematical structure of Huber's M-estimator, SIAM J. on Scientific and Statistical Computing 6, 209-219 (1985). 

Duality in robust linear regression using Huber's M-estimator

The robust linear regression problem using Huber's piecewise-quadratic M-estimator function is considered. Without exception, computational algorithms for this problem have been primal in nature. In this note, a dual formulation of this problem is derived using Lagrangean duality. It is shown that the dual problem is a strictly convex separable quadratic minimization problem with linear equality and box constraints. Furthermore, the primal solution (Huber's M-estimate) is obtained as the optimal values of the Lagrange multipliers associated with the dual problem. As a result, Huber's M-estimate can be computed using off-the-shelf optimization software

~ )  Pergamon 
Appl. Math. Lett. Vol. 10, No. 4, pp. 65-70, 1997 
Copyright©1997 Elsevier Science Ltd 
Printed in Great Britain. All rights reserved 
0893-9659/97 $17.00 + 0.00 
P Ih  S0893-9659(97)00061-X 
Dual i ty  in Robust  L inear Regress ion  
Us ing  Huber ' s  M-Est imator  
M. (~. P INAR 
Industrial Engineering Department, Bilkent University 
06533 Bilkent, Ankara, Turkey 
must af ap©bilkent, edu. tr 
(Received October 1996; accepted December 1996) 
Communicated by R. H. Byrd 
Abst ract - -The  robust linear regression problem using Huber's piecewise-quadratic M-estimator 
function is considered. Without exception, computational gorithms for this problem have been 
primal in nature. In this note, a dual formulation of this problem is derived using Lagrangean 
duality. It is shown that the dual problem is a strictly convex separable quadratic minimiza- 
tion problem with linear equality and box constraints. Furthermore, the primal solution (Huber's 
M-estimate) is obtained as the optimal values of the Lagrange multipliers associated with the dual 
problem. As a result, Huber's M-estimate can be computed using off-the-shelf optimization software. 
Keywords- -Lagrangean duality, Huber's M-estimator, Robust regression, Quadratic program- 
ming. 
1. INTRODUCTION 
There has been considerable interest in the theory and algor i thms for robust  est imat ion in the 
past  two decades. In part icular ,  Huber 's  M-est imator  [1] has received a great  deal of at tent ion 
from both theoret ical  and computat ional  points of view. Robust  est imat ion is concerned with 
identifying "outliers" among data  points and giving them less weight. Huber 's  M-est imator  is 
essential ly the least squares est imator,  which uses the g l -norm for points that  are considered 
outl iers with respect to a certain threshold. Hence, the Huber  criterion is less sensit ive to the 
presence of outl iers. 
More precisely, Huber 's  M-est imate is a minimizer x* E ~n of the function 
where 
f(x) = 
i----1 
(1) 
12 
~-=.t , if Itl < % 
p(t) = "Y 1 (2) 
I t l -  2-r, if Itl _~ % 
with a tuning constant ~ > O, and a scaling factor a that  depends on the data  to be est imated.  
The residual r~(x) is defined as 
r~(x) -- b~ - a~x, (3) 
Special thanks are due to R. H. Byrd for his careful handling of the paper. One referee pointed out the QP form 
of the primal added to the revised version. Another eferee clarified several issues concerning robust estimators 
and data analysis, and suggested directions for future research. I am grateful to both. 
65 
66 M. ~. PINAR 
for all i = 1, . . .  ,m with r -- b - ATx .  To view this minimization problem in a more familiar 
format, define a "sign vector" 
ST(X) = [S~I(X),...,S~m(X)] (4) 
with 
and 
where 
-1,  if ri(x) _< -7 ,  
sv~(x) = 0, if Iri(x)[ < 7, 
1, if ri(x) >_ 7, 
(5) 
Ws = diag (w l , . . . ,  win), (6) 
wi = 1 - s~. (7) 
Now, assuming a unit a for the moment, Huber's M-estimation problem can be expressed as 
the following minimization problem. 
PROBLEM [P]. 
1 T T[ 1 ] 
minimize F(x )  - ~Tr W, r  + s~ r - ~Ts~ , (8) 
where the argument x of r is dropped for notational convenience. Clearly, F measures the "small" 
residuals ([ri(x)[ < 7) by their squares, while the "large" residuals are measured by the gl func- 
tion. Thus, F is a piecewise quadratic function, and it is once continuously differentiable in Nn. 
The contribution of the present note is to introduce a dual approach to this computationally 
intensive problem. A simple derivation is used to give a dual problem which turns out to be a 
problem familiar to the numerical optimization community: linearly constrained separable convex 
quadratic programming. To the best of the author's knowledge, this duality relationship has not 
been noticed before. This simple result reduces the M-estimation problem to one that is easily 
solved using off-the-shelf optimization software. 
The interest in Huber's M-estimator derives from its perceived utility as a robust estimation 
procedure. In this context, 7 is an important quantity as it controls the spread of the residuals. 
This suggests that it should be related to the scaling factor a. Some algorithms [2,3] estimate a
only once at the beginning using some rules of thumb, and keep it fixed throughout the compu- 
tation. Another approach is to estimate a iteratively by treating it as an independent variable. 
This is used by Huber in [4] who suggests a to be computed from an auxiliary equation involving 
a, 7, and x. Shanno and Rocke [5] and Ekblom [6] also use this approach in their respective 
papers. An alternative way to connect hese two parameters and to iteratively estimate a is to 
set a = 1 for any values of 7 and a, and to replace 7 with 7a. Hence, there is no loss of infor- 
mation in using 7 only, while 7 and a are both allowed to vary during the computation. Clark 
and Osborne [7] use this observation in their algorithm that computes both a and V. This is a 
continuation method, so that at each stage they have the M-estimate for the current 7. Clark and 
Osborne also describe a partitioning algorithm to resolve degenerate situations in the continua- 
tion algorithm. The continuation method is essentially based on tracing a curve of M-estimates 
as a function of 7. Their approach consists of two stages. 
1. Construction of a solution for a particular 7. The easiest cases are 7 = 0 (~1 estimate) and 
7 = oo (least squares). These are both used in [7]. 
2. Continuation with respect o 7 to solve the data analysis problem. 
A potential contribution of the present paper in the above context is that it suggests a new 
method for obtaining Huber's M-estimate for a particular value of 7. Then, continuation with 
respect o 7 of a solution obtained by the dual method can be pursued using methods of para- 
metric quadratic programming. In this connection, the rule of thumb choices that have been 
Huber's M-Estimator 67 
suggested for ~ in the literature [2-4] can be used to see if they provide better starting points for 
the continuation method than least squares and ei estimates. This is a topic for future research. 
When a and 9, are fixed, most algorithms that have proved successful for the computation of
Huber's M-estimate have been iterative in nature. Since F is only once continuously differentiable, 
research concentrated ondeveloping successful applications of Newton's method to this problem, 
and on studying ideas such as the iteratively reweighted least squares (IRLS). Among these, 
Huber and Dutter's method [8] and Huber's [4] apply Newton's method to the nonlinear equation 
system: 
which represent the optimality conditions. The IRLS algorithm is attributed to Beaton and 
Tukey [9]. This algorithm has been discussed in several other papers as well [2,3,10-12]. A brief 
review of these algorithms i given in [7]. 
Madsen and Nielsen gave finite modified Newton algorithms for the minimization of the Huber 
function in [13]. These algorithms capitalize on the piecewise quadratic nature of the function. 
The idea is that if the sign vector associated with a minimizer x*, s* say, were known, then the 
minimizer could be computed using one step of Newton's method. Since sign vectors correspond 
to a subdivision of ~n into a finite number of subregions, the methods of Madsen and Nielsen 
reduce to a search for the correct subregion. 
2. A DUAL  PROBLEM 
In this section, a dual problem to [P] is derived using Lagrangean duality. The interested 
reader is directed to the book by Rockafellar [14] for a detailed exposition of Lagrangean duality. 
We use a single parameter 7 to mean 3'a as discussed in the previous ection. 
Consider the problem [P] in a slightly different form: 
minF(x) =_ @ (b - AT  x)  . (10) 
Let u = b - ATx ,  and rewrite the problem as: 
min @(u) 
s.t. b - ATx  = u. 
Associating the multipliers y E ~m with the equality constraints b-ATx  = u, we get the following 
Lagrangean problem: 
max min {@(u) + yT (b - AT  x -- u) } .  (11) 
X,~ 
This is equivalent to: 
max [yTb + min ((I)(u)- yTu} ÷ m~n (_yT A Tx}]. 
y u 
Observing that x is a free variable, the term 
min {_yT AT x} 
yields the constraint 
It remains to simplify the term 
Ay=O.  
(12) 
(13) 
(i4) 
min - ysu}. (15) 
68 M. ~. PINAR 
Simple calculus hows that 
min{O(u)--yTu)} = --~TY Y, 
U --OO, 
Hence, the dual problem is the following. 
PROBLEM [D]. 
1 T 
max bmy-  ~TY Y 
i f -1  < y _< 1, 
otherwise. 
(16) 
s.t. Ay = 0 
- l<y<l .  
An alternative way to arrive at the above dual is to pose the primal problem [P] as a quadratic 
programming problem: 
1 m ~(  ) 
max + q ' -  
i=.l i=l  
s.t. -p  - q < b - ATx  < p -{- q 
p<_Te 
q>_O. 
where e denotes a vector with all components unity. This is not surprising since the dual of a 
quadratic program is again a quadratic program. However, this alternative derivation is sub- 
stantially longer since it requires ome transformations on the resulting dual to be cast in the 
form [D]. So, it is not included into the present paper. 
The remarkable fact about the duality result is that the dual we have derived is a quadratic 
programming problem with a strictly concave separable objective function, linear equality, and 
box constraints. This is important in two respects. First, this dual problem is an intensely 
researched, numerically well-solved problem. An excellent software system is available from 
Stanford University for the solution of quadratic programming problems [15,16]. Furthermore, 
numerical procedures for solving this problem are part of almost any subroutine library. Such 
procedures are also available through matrix manipulation packages such as Matlab TM and Oc- 
tave [17]. Second, it has been shown to be polynomially solvable [18], and efforts to turn related 
algorithms into reliable software on this front are also under way. In particular, quadratic pro- 
gramming versions of the software systems CPLEX TM and LoQo are available [19]. Hence, any 
advances made in the fast and accurate solution of the quadratic programming problem will 
benefit Huber's M-estimation problem. 
Note that any numerical procedure which yields a primal-dual solution to [D] gives a mini- 
mizer x* of the primal problem, which is really of interest here. To see this, it is instructive to 
derive a dual problem to [D] using Lagrangean duality. 
Associating multipliers v E !l~" with the constraints Ay = O, we get the following Lagrangean 
problem: 
min max ~ b T y - 1 T } v -*<_~<* [. ~TY Y +vT( - -Ay+O)  . (17) 
Now, rearranging terms and using r =- b - ATv,  simple calculus hows that: 
max y~r~ - = -1<~,<1 ~7Y~ 27r"  (18) 
57, > 7. 
But, this is precisely Huber's M-estimator function. Therefore, the optimal values of the dual 
multipliers associated with the equality constraints in [D] give precisely Huber's M-estimate. 
Huber's M-Estimator 69 
Further insight into the relationship between [P] and [D] is gained through Theorem 1 below 
which links the opt imal  solutions of [P] and [D]. This theorem shows that  the opt imal  solution 
to [D] is obtained as the first derivative of the function F with respect o r at any pr imal opt imal  
point. Before stat ing this result, we quote the following property that  was proved in [20]. 
LEMMA 2.1. Let  x be a minimizer  of  F for some value of  ~ > O, and let s = sx(x) with Ws 
defined accordingly. Also, let r = b - AVx.  Then, ri is constant for all i such that si = O. 
Furthermore,  sw is constant for any x that  minimizes F.  
By strict concavity, it is clear that  the optimal solution to [D] is unique. 
THEOREM 2.1. Let  x be a minimizer  o fF  for some value of  "l > O, and let s = s~(x) with Ws 
defined accordingly. Then, the unique optimal solution of [D] is given as: 
w.r (x )  
y - - -  + s. (19) 
PROOF. Associate multipliers x with the equality constraints and nonnegative multipliers 
a,/3 E N'~ with the box constraints. From [6, Theorem 28.3], it is well known that  for y to 
be an opt imal  solution for [D], and for (x, a,/3) to be a Lagrange multiplier vector, it is necessary 
and sufficient that  (y, x, a,/3) be a saddlepoint of the Lagrangean of [D] as defined in [6, p. 280]. 
This condition holds if and only if the components of (y, x, a,/3) satisfy the following conditions: 
~y-  b+ ATx  + a- /3  = 0, (20) 
Ay = 0, (21) 
ai(yi  - 1) = 0, for all i = 1 . . . . .  m, and (22) 
/3i(-y~ - 1) -- 0, for all i = 1 . . . . .  m. (23) 
From the first condition (20), we obtain 
b - A-rx o~ - /3  
u - (24) 
Let r - b - A-r x and consider three cases. 
CASE 1. If  --1 < Yi < 1, clearly, ai =/3i = 0. This implies that  [ri(x)[ < ~ with s~i(x) = 0, i.e., 
= 
CASE 2. If  Yi = 1, this implies that /3 i  = 0. Hence, r i (x)  = ~ + ai. Therefore, r i (x)  >_ "y with 
s~i(x) = 1, i.e., yi = s~i(x) = 1. 
CASE 3. If  Yi = --1, one has ai  = 0. Hence, similarly to Case 2 above, r i (x)  <_ -~ with 
s~i(x) = -1 ,  i.e., Yi = s~(x)  = -1 .  
Therefore, an alternative xpression for y is given as: 
+ (25) 
The result now follows from the previous lemma. | 
70 M. ~. PINAR 
REFERENCES 
I. P. Huber, Robust Statistics, John Wiley, New York, (1981). 
2. J.B. Birch, Some convergence properties of iterated reweighted least squares in the location model, Comm. 
Statist. B9, 359-369 (1980). 
3. P.W. Holland and R.E. Welsch, Robust regression using iteratively reweighted least squares, Comm. Statist. 
A6, 813-827 (1977). 
4. P. Huber, Robust regression: Asymptotics, conjectures and Monte Carlo, Annals of Statistics 1, 799-821 
(1973). 
5. D.F. Shanno and D.M. Rocke, Numerical methods for robust regression: Linear models, SIAM J. on Scientific 
and Statistical Computing 7, 86-97 (1986). 
6. H. Ekblom, A new algorithm for the Huber estimator in linear models, BIT 28, 123-132 (1988). 
7. D.I. Clark and M.R. Osborne, Finite algorithms for Huber's M-Estimator, SIAM J. on Scientific and Sta- 
tistical Computing 7, 72-85 (1986). 
8. P. Huber and R. Dutter, Numerical solution of robust regression problems, In COMPSTAT I974 Proc. 
Symposium on Computational Statistics, (Edited by G. Brushmann), pp. 165-172, Physike Verlag, Berlin, 
(1974). 
9. A.E. Beaton and J.W. Tukey, The fitting of power series, meaning polynomials, illustrated on band- 
spectroscopic data, Technometrics 16, 147-185 (1974). 
10. J.B. Birch, Effects of the starting value and stopping rule on robust estimates obtained by iterated weighted 
least squares, Comm. Statist. Bg, 133-140 (1980). 
11. R.H. Byrd and D.A. Pyne, Some results on the convergence of the iteratively reweighted least squares 
algorithm for robust regression, In Proc. Statistical Computing Section, pp. 87-90, American Statistical 
Association, (1979). 
12. S.A. Ruziusky and E.T. Olsen, L1 and Leo minimization via a variant of Karmarkar's algorithm, IEEE 
Transactions on Acoustics, Speech and Signal Processing 37, 245-253 (1989). 
13. K. Madsen and H.B. Nielsen, Finite algorithms for robust linear regression, BIT 30, 682-699 (1990). 
14. tLT. Rockafellar, Convex Analysis, Princeton University Press, Princeton, N J, (1970). 
15. P.E. Gill, W. Murray, M.A. Saunders and M.H. Wright, User's guide for SOL/QSPOL (Version 3.2), Tech- 
nical Report SOL 84-6, Systems Optimization Laboratory, Department of Operations Research, Stanford 
University, Stanford, CA, (1984). 
16. P.E. Gill, S.J. Harnmarling, W. Murray, M.A. Saunders and M.H. Wright, User's guide for LSSOL (Version 
1.0): A Fortran package for constrained linear least squares and convex quadratic programming, Technical 
Report SOL 86-1, Systems Optimization Laboratory, Department of Operations Research, Stanford Univer- 
sity, Stanford, CA, (1986). 
17. J.W. Eaton, OCTAVE: A High-Level Interactive Language for Numerical Computation, (1994). 
18. Y. Ye and E. Tse, An extension of Karmarkar's projective algorithm for convex quadratic programming, 
Mathematical Programming 44, 157-179 (1989). 
19. R.J. Vanderbei and T.J. Carpenter, Symmetric indefinite systems for interior point methods, Mathematical 
Programming 58, 1-32 (1993). 
20. K. Madsen, H.B. Nielsen and M.(~. Pmar, New characterizations of £1 solutions to overdetermined systems 
of linear equations, Operations Research Letters 16, 159-166 (1994). 
21. D.I. Clark, The mathematical structure of Huber's M-estimator, SIAM J. on Scientific and Statistical 
Computing 6, 209-219 (1985). 


file:///data/core-remote/dit/data/elsevier/pdf/159/aHR0cDovL2FwaS5lbHNldmllci5jb20vY29udGVudC9hcnRpY2xlL3BpaS9zMDg5Mzk2NTk5NzAwMDYxeA==.pdf

Duality in robust linear regression using Huber's M-estimator

Abstract

Similar works

Full text

Available Versions

Elsevier - Publisher Connector

Bilkent University Institutional Repository

Bilkent University Institutional Repository

Elsevier - Publisher Connector