Minimizing the maximal residual in mathematica

With this:

d1 = 10; d2 = 4;
mat = RandomReal[{-1, 1}, {d1, d2}];
vec = RandomReal[{-1, 1}, d1];

LeastSquares[mat, vec] returns the x, that minimizes Plus @@ ((mat.x-vec)^2)

What is the best way to make mathematica return the x that minimizes Max @@ ((mat.x-vec)^2)

I ended up with this, thank you Daniel 😉

With[{L = Length[First[mat]]},
LinearProgramming[Prepend[ConstantArray[0, L], 1],
Prepend[#, 1] & /@ Riffle[mat, -mat], Riffle[vec, -vec],
Prepend[ConstantArray[-\[Infinity], L], 0]]]




Nice. I might “borrow” (okay, steal) this next time I want to use LinearProgramming directly and avoid NMinimize overhead. I’d give it an upvote but I notice I already did that.
– Daniel Lichtblau
Feb 28 ’14 at 19:14


2 Answers


It is equivalent to minimize the absolute values. This can be set up as an explicit linear programming problem. The advantage over the approach of @bobthechemist (which is good, and I voted up) is that the problem can then be shipped to special case LP code.

vars = Array[x, d2];
linearexprs = mat.vars – vec;
constraints =
Join[Thread[max >= linearexprs], Thread[max >= -linearexprs]];
NMinimize[{max, constraints}, Join[{max}, vars]]

(* Out[790]= {0.634141, {max -> 0.634141, x[1] -> 0.679669,
x[2] -> 0.244346, x[3] -> -0.539059, x[4] -> -1.09431}} *)



Amazing, thank you, I tried to set this up with LinearProgramming[], but I think it’s not possible, do you think I’m wrong?
– Coolwater
Feb 26 ’14 at 19:35



It can be set up to use LinearProgramming directly, but that’s usually more work than I’m willing to do. Maybe someone out there has a program to convert an LP in NMinimize format to one directly usable by LinearProgramming? Of course then you’d be essentially replicating work that NMinimize does in coming up with that conversion itself.
– Daniel Lichtblau
Feb 26 ’14 at 20:15

Since LeastSquares can be written as

NMinimize[Plus @@ ((mat.{x1, x2, x3, x4} – vec)^2), {x1, x2, x3, x4}]

Then you can use a similar approach to minimizing your desired function:

NMinimize[Max @@ ((mat.{x1, x2, x3, x4} – vec)^2), {x1, x2, x3, x4}]

Although whether or not this is the “best” way is likely up for debate.