total least squares derivation with matrices Taken from a computer vision book: "to minimize the sum of the perpendicular distances between points and lines, we need to minimize $$ \sum_i (ax_i + by_i +c)^2$$ subject to $a^2 +b^2 =1$. Now using a Lagrangian multiplier $\lambda$, we have a solution if $$ \left( \begin{array}{ccc}
\overline{x^2} & \overline{xy} & \overline{x} \\
\overline{xy} & \overline{y} & \overline{y} \\
\overline{x} & \overline{y} & 1 \end{array} \right)\left[ \begin{array}{c}
   a \\
   b \\
   c
  \end{array} \right]  = \lambda \left[ \begin{array}{cc}
   2a\\
   2b \\
   0
  \end{array} \right]$$
How is the book getting these matrices? 
Also, the notion is that $\overline{u} = \frac{\sum u_i}{k}$. (Yeah, I don't know what $k$ stands for. I can only assume this is an average.)
It goes onto say that $c = -a\overline{x} - b\overline{y}$, and that we can substitute this back to get the eigenvalue problem $$\left[ \begin{array}{cc}
   \overline{x^2} -\overline{x}~\overline{x} &  \overline{xy} -  \overline{x}\overline{y}\\
    \overline{xy} -  \overline{x}\overline{y} &  \overline{y^2} -  \overline{y} ~\overline{y} \\
  \end{array} \right] \left[\begin{array}{cc}
   a\\
   b
  \end{array} \right] = \mu \left[ \begin{array}{cc}
   a\\
   b
  \end{array} \right].$$
I don't see what they substituted into, and how the answer is derived.
 A: We want to minimize
$$
S=\sum_i(ax_i+by_i+c)^2\tag{1}
$$
Varying $(1)$ yields
$$
\frac12\delta S=
\begin{bmatrix}\delta a&\delta b&\delta c\end{bmatrix}
\left(\sum_i\begin{bmatrix}x_i\\y_i\\1\end{bmatrix}\begin{bmatrix}x_i&y_i&1\end{bmatrix}\right)
\begin{bmatrix}a\\b\\c\end{bmatrix}\tag{2}
$$
Since, $a^2+b^2=1$, we restrict the variations to those so that
$$
0=\begin{bmatrix}\delta a&\delta b&\delta c\end{bmatrix}\begin{bmatrix}a\\b\\0\end{bmatrix}\tag{3}
$$
We want to find $a,b,c$ that cancel $(2)$ under the constraints in $(3)$. Orthogonality implies that we need
$$
\left(\sum_i\begin{bmatrix}x_i\\y_i\\1\end{bmatrix}\begin{bmatrix}x_i&y_i&1\end{bmatrix}\right)
\begin{bmatrix}a\\b\\c\end{bmatrix}
=\lambda\begin{bmatrix}a\\b\\0\end{bmatrix}\tag{4}
$$
Equation $(4)$ is the first equation in the question.
The bottom row of $(4)$ says
$$
\left(\sum_i\begin{bmatrix}x_i&y_i&1\end{bmatrix}\right)
\begin{bmatrix}a\\b\\c\end{bmatrix}=0\tag{5}
$$
The top two rows of $(4)$ say
$$
\left(\sum_i\begin{bmatrix}x_i\\y_i\end{bmatrix}\begin{bmatrix}x_i&y_i&1\end{bmatrix}\right)
\begin{bmatrix}a\\b\\c\end{bmatrix}
=\lambda\begin{bmatrix}a\\b\end{bmatrix}\tag{6}
$$
we can rewrite $(6)$ as
$$
\left(\sum_i\begin{bmatrix}x_i\\y_i\end{bmatrix}\begin{bmatrix}x_i&y_i\end{bmatrix}\right)
\begin{bmatrix}a\\b\end{bmatrix}
+\left(\sum_i\begin{bmatrix}x_i\\y_i\end{bmatrix}\right)c
=\lambda\begin{bmatrix}a\\b\end{bmatrix}\tag{7}
$$
and $(5)$ as
$$
\left(\sum_i\begin{bmatrix}x_i&y_i\end{bmatrix}\right)
\begin{bmatrix}a\\b\end{bmatrix}+\left(\sum_i1\right)c=0\tag{8}
$$
If we set $n=\sum\limits_i1$ and substitute $c$ computed from $(8)$ into $(7)$, we get
$$
\left(\sum_i\begin{bmatrix}x_i\\y_i\end{bmatrix}\begin{bmatrix}x_i&y_i\end{bmatrix}-\frac1n\sum_i\begin{bmatrix}x_i\\y_i\end{bmatrix}\sum_i\begin{bmatrix}x_i&y_i\end{bmatrix}\right)
\begin{bmatrix}a\\b\end{bmatrix}
=\lambda\begin{bmatrix}a\\b\end{bmatrix}\tag{9}
$$
Equation $(9)$ is the last equation in the question.
A: Take the function
$$
F(a,b,c,\lambda)=\sum_i (ax_i+by_i+c)^2-\lambda(a^2+b^2-1)
$$
so that
\begin{align}
\frac{\partial F}{\partial a}&=2\sum_i (ax_i+by_i+c)x_i-2\lambda a=0\\
\frac{\partial F}{\partial b}&=2\sum_i (ax_i+by_i+c)y_i-2\lambda b=0\\
\frac{\partial F}{\partial c}&=2\sum_i (ax_i+by_i+c)=0\\
\frac{\partial F}{\partial \lambda}&=-(a^2+b^2-1)=0
\end{align}
or better
\begin{align}
&a\sum_i x_i^2 +b\sum_i x_iy_i+c\sum_i x_i=\lambda a\\
&a\sum_i x_iy_i+b\sum_i y_i^2 +c\sum_i y_i=\lambda b\\
&a\sum_i x_i   +b\sum_i y_i   +nc=0\\
&a^2+b^2=1
\end{align}
