Least Squares Parabola

Solving an over-determined linear system A x = b is generally impossible, but finding an x0 so that A x0 is as close as possible to b leads to the concept of a least squares solution. One can see that such an x0 must satisfy the normal equations
A_tr A x0 = A_tr b
where A_tr is the transpose of the matrix A. Generically this system will have a unique solution. We use this worksheet to study this in the case of fitting a set of points to a parabola.

Problem : Find the best fitting parabola y = ax^2 +bx+c through the following points: (1,2),(1,3),(2,4),(3,5),(4,6)

> pointlist := [[1,2],[1,3],[2,4],[3,5],[4,6]];

Solution: The over-determined system here is
1a +1b+c=2
1a +1b+c=3
4a +2b+c=4
9a +3b+c=5
16a+4b +c=6
In matrix form Ax= B:

> with (linalg):

Warning, the protected names norm and trace have been redefined and unprotected

> A:= array([[1,1,1],[1,1,1],[4,2,1],[9,3,1],[16,4,1]]);

> x:= vector ([a, b, c]);

> B:= vector([2,3,4,5,6]);

> evalm( A &* x = B);

The normal equations Nx=m.

> N:= evalm ( transpose(A) &* A );

> m:= evalm ( transpose(A) &* B);

> x0:= evalm (inverse(N) &* m);

> parabola := x -> x0[1] * x^2 + x0[2] * x + x0[3];

> plot ({pointlist, parabola});

>

>