Find best fitted approximated parameters

Question

I want to find the best fitted approximated parameters for the following case:

$$\begin{align*} 2a + 3b &= 5\\ a + b &= 4\\ a + b &= 3.9\\ 2a + 3b &= 5.1 \end{align*}$$

When I use

Solve[{2 a + 3 b == 5, a + b == 4}, {a, b}]`

I can get exact values of a and b.

If I add two more equations to obtain the best approximate values of a and b for the 4 equations, Solve doesn't work.

Could you please give me your comments?

Answer 1

Clear["Global`*"]

Your equations are of the form f = a*x + b*y

f[x_, y_] := a*x + b*y

Your equations represent the {x, y, f} data set

data = {{2, 3, 5}, {1, 1, 4}, {1, 1, 3.9}, {2, 3, 5.1}};

Use FindFit to fit your data to the function f

param = FindFit[data, f[x, y], {a, b}, {x, y}]

(* {a -> 6.8, b -> -2.85} *)

The fitted function is

ff[x_, y_] = f[x, y] /. param

(* 6.8 x - 2.85 y *)

Plotting the result

Show[
 Plot3D[
  ff[x, y], {x, 0.95, 2.05}, {y, 0.95, 3.05},
  PlotStyle -> Opacity[0.5]],
 Graphics3D[{Red, AbsolutePointSize[5],
   Point[data]}],
 AxesLabel -> (Style[#, 14, Bold] & /@ {x, y, f})]

Answer 2

For this case, you can reformulate the system in matrix-vector format and then use LeastSquares[]:

mat = {{2, 3}, {1, 1}, {1, 1}, {2, 3}}; vec = {5, 4, 39/10, 51/10};
LeastSquares[mat, vec]
   {34/5, -57/20}

N[%]
   {6.8, -2.85}

which is of interest if one wants a solution that is optimal in the $2$-norm sense.

---

Sometimes, however, one might want to find solutions that are optimal in other norms, e.g. the $1$-norm or $\infty$-norm. One can of course use NMinimize[] for that:

NMinimize[Norm[mat.{x, y} - vec, 1], {x, y}]
   {0.19999999999999973, {x -> 6.690832044871543, y -> -2.77644284441384}}

NMinimize[Norm[mat.{x, y} - vec, ∞], {x, y}]
   {0.05000000029087026, {x -> 6.800000000759247, y -> -2.8500000004683765}}

but it is better to use specialized methods in these cases, which rely on their equivalence to a linear programming problem (see e.g. this paper):

L1Solve[A_?MatrixQ, b_?VectorQ, opts___] := Module[{m, n, id},
  {m, n} = Dimensions[A]; id = IdentityMatrix[m, Head[A]]; 
  Take[LinearProgramming[PadLeft[ConstantArray[1, 2 m], 2 m + n], 
                         Join[A, -id, id, 2], PadRight[List /@ b, {Automatic, 2}], 
                         PadRight[ConstantArray[-∞, n], 2 m + n], opts], n]]

LInfinitySolve[A_?MatrixQ, b_?VectorQ, opts___] := Module[{n = Last[Dimensions[A]]},
  Rest[LinearProgramming[UnitVector[n + 1, 1], 
                         PadLeft[Join[-A, A], {Automatic, n + 1}, 1], Join[-b, b], 
                         ConstantArray[-∞, n + 1], opts]]]

L1Solve[mat, vec]
   {33/5, -27/10}

N[%]
   {6.6, -2.7}

LInfinitySolve[mat, vec]
   {34/5, -57/20}

N[%]
   {6.8, -2.85}

Somewhat coincidentally, the $\infty$-norm solution coincides with the $2$-norm solution in this case.

(On another note, the documentation used to have notes on finding $1$-norm and $\infty$-norm solutions.)