掰开揉碎推导Normal Equation
Normal Equation是一种基础的最小二乘方法，本文将从线性代数的角度来分析Normal Equation（而不是从矩阵求导 matrix derivative 的角度）。
很多作者（特别是智商比较高的）在推导公式的时候有意无意的忽略了思考过程，只留下漂亮的步骤。这让很多读者（比如说我）跟不上节奏，最后一头雾水。本文将从求解“貌似无解”的方程组入手，再讲讲投影（Projection）的使用，最后进入到Normal Equation的应用。我的目的是让和我一样蠢的孩子对\overrightarrow\theta = (A^TA)^{-1}A^T\overrightarrow{y} 这个重要公式有一个Big Picture——即使忘记了也可以重头推出。
更新记录：
更新1 增加了对A^TA的使用解释（偏导数证明）
一、求解不可解的方程组
先看一个最最简单的例子——
例1.0 如图，在R^2空间中有两个向量，求一个常数\theta 使两个向量满足\theta \overrightarrow{a} =\overrightarrow{b} 。
这个方程明显不可解，因为\overrightarrow{b}与\overrightarrow{a} 不共线，无法通过对\overrightarrow{a}数乘得到\overrightarrow{b}。
再看下一个比较简单的例子——
例2.0 在R^3空间中的平面P有一组基\overrightarrow{a_1}和\overrightarrow{a_2}，如图所示，求出常数\theta_1与\theta_2使向量\overrightarrow{b}满足条件\theta_1\overrightarrow{a_1}+\theta_2\overrightarrow{a_2} =\overrightarrow{b} 。
这个方程也明显不可解，因为\overrightarrow{b}不在平面P上，而\overrightarrow{a_1}与\overrightarrow{a_2}的线性组合只能得到平面上的向量。
以上两个问题非常的典型，因为在解决实际问题的时候，我们很难得到Perfect Solution，我们只能尽力而为的争取Best Solution。以上两个例子明显没有做到perfect(连基本的方向都错了)，那么如何找到best solution呢？
二、投影的应用
思路很简单：我们只要找到一个\theta^*使\overrightarrow{a}方向上的向量\overrightarrow{p} =\theta^* \overrightarrow{a} 距离\overrightarrow{b}最近。
回到最简单的例子
如图，在R^2空间中有两个向量，求一个常数\theta 使两个向量满足\theta \overrightarrow{a} =\overrightarrow{b} 。
现在应该如何寻找\theta \overrightarrow{a} =\overrightarrow{b} 的解呢？
最好的方法就是抛弃\overrightarrow{b} 向量中垂直\overrightarrow{a}的分量，只要计算\theta使\theta \overrightarrow{a}等于向量\overrightarrow{b} 在\overrightarrow{a}
方向的分量（即\overrightarrow{b} 在\overrightarrow{a} 上的投影(Proj)\overrightarrow{p} ），同时我们把向量\overrightarrow{b} 垂直\overrightarrow{a} 方向的分量称为\overrightarrow{e} (error)。
原来的问题\theta \overrightarrow{a} =\overrightarrow{b} 变成了求解\theta^* \overrightarrow{a} =\overrightarrow{p} （\theta^*是\theta的估计量）
因为\overrightarrow{p} 与\overrightarrow{e} 合成了\overrightarrow{b} 向量（\overrightarrow{e} +\overrightarrow{p}= \overrightarrow{b} ），而且\overrightarrow{e} 垂直于\overrightarrow{a} (\overrightarrow{e}\bot \overrightarrow{a})，所以我们得出了一个非常重要的结论（敲黑板）！！！核心啊！！！
\overrightarrow{a}^T(\overrightarrow{b}-\theta^*\overrightarrow{a} ) = 0
这个方程的核心就是写成向量内积形式的\overrightarrow{e} 与\overrightarrow{a} 的垂直关系，只不过\overrightarrow{e} 被拆开书写。其实这个方程也可以写作\overrightarrow{a}\times (\overrightarrow{b}-\theta^*\overrightarrow{a} ) = 0 ，但是写作转置向量\overrightarrow{a}^T的形式可以让这个方程更自然的拓展到高维。好了，我们继续改写方程……
\theta^*\overrightarrow{a}^T \overrightarrow{a} = \overrightarrow{a}^T\overrightarrow{b}
\theta^* = \frac{\overrightarrow{a}^T\overrightarrow{b}}{\overrightarrow{a}^T\overrightarrow{a}}
在这一步我们就得到了best的\theta，但考虑到这并不perfect，所以我们称之为\theta^*。
P.S.如果想用投影矩阵P来简化从\overrightarrow{a}转换到\overrightarrow{p} 的过程，可以把\theta^*的结果带入到\overrightarrow{p} =\theta^* \overrightarrow{a}中。我们发现投影矩阵P在形式上就等于乘数\theta^* = \frac{\overrightarrow{a}^T\overrightarrow{b}}{\overrightarrow{a}^T\overrightarrow{a}}，即\overrightarrow{p} 满足\overrightarrow{p} =P\overrightarrow{a}。
现在我们再看看怎么在R^3中解决不可解方程。
例2.0 在R^3空间中的平面P有一组基\overrightarrow{a_1}和\overrightarrow{a_2}，如图所示，求出常数\theta_1与\theta_2使向量\overrightarrow{b}满足条件\theta_1\overrightarrow{a_1}+\theta_2\overrightarrow{a_2} =\overrightarrow{b} 。
平面P有基向量\overrightarrow{a_1}和\overrightarrow{a_2}，故P可以表示成基的线性组合\theta_1\overrightarrow{a_1}+\theta_2\overrightarrow{a_2} ，即 \begin{bmatrix} a_1 a_2 \end{bmatrix} \begin{pmatrix} \theta_1 \\ \theta_2 \end{pmatrix}
令基向量组成的矩阵A=\begin{bmatrix} | \ \ |\\ a_1 a_2\\ | \ \ |\\ \end{bmatrix}，参数组成的向量\overrightarrow{\theta} = \theta_{1 \dots n} \ (n = 2)，与平面垂直的误差向量\overrightarrow{e} = \overrightarrow{b}-A \overrightarrow{\theta^*} 。（这里插一句话，最小二乘法的核心就是找出一个\theta^*就是让\left|| A\theta-b \right||^2最小化）
我们发现在R^2中的问题\theta \overrightarrow{a} =\overrightarrow{b} 在这里拓展成为了A\overrightarrow{\theta }=\overrightarrow{b} 。
相应的，\theta^* \overrightarrow{a} =\overrightarrow{p} 问题在这里拓展成了A\overrightarrow{\theta^*}=\overrightarrow{p}，其中\overrightarrow{p}=\overrightarrow{\theta_1}\overrightarrow{a_1}+\overrightarrow{\theta_2}\overrightarrow{a_1}。
还是一样的套路，我们还是从垂直关系入手——因为P \bot \overrightarrow{e} ，而且P \in \overrightarrow{\theta_1} \overrightarrow{a_1} +\overrightarrow{\theta_1} \overrightarrow{a_2} ，所以有以下方程组——
\begin{cases} \overrightarrow{a_1}^T(\overrightarrow{b}-A\overrightarrow{\theta^*})=0\\ \overrightarrow{a_2}^T(\overrightarrow{b}-A\overrightarrow{\theta^*})=0 \end{cases}整理成矩阵的形式——
\begin{bmatrix} -\overrightarrow{a_1}^T-\\ -\overrightarrow{a_2}^T-\\ \end{bmatrix} (\overrightarrow{b}-A\overrightarrow{\theta^*})=0（敲黑板！！！敲黑板！！！）
A^T(\overrightarrow{b}-A\overrightarrow{\theta^*})=0
写到这里回头看看R^2情景下的核心公式\overrightarrow{a}^T(\overrightarrow{b}-\theta^*\overrightarrow{a} ) = 0 ，可以这家伙换一套马甲又出现了！！！看来方程A^T(\overrightarrow{b}-A\overrightarrow{\theta^*})=0是一种高维的拓展。我们可以把R^2中的\overrightarrow{a}看成一个只有一列的矩阵。
我们继续整理这个公式——
A^T\overrightarrow{b} = A^TA\overrightarrow{\theta}\overrightarrow{\theta} = (A^TA)^{-1}A^T\overrightarrow{b}
写到这里我们就没什么可以干的了。
有人可能想说——明明还可以继续化简啊！！！
\begin{equation} \begin{split} \overrightarrow{\theta} &= (A^TA)^{-1}A^T\overrightarrow{b}\\ & = A^{-1}(A^T)^{-1}A^T\overrightarrow{b}\\ & = A^{-1}\overrightarrow{b} \end{split} \end{equation}但实际的情况中，我们不能保证矩阵A总是方阵（square），但是A^TA总是可以保证是方阵。因为只有方阵才有逆矩阵，所以我们只能保证有(A^TA)^{-1}，而不能保证有A^{-1}。
所以我们只能回到\overrightarrow{\theta} = (A^TA)^{-1}A^T\overrightarrow{b}这里。如果你有读过Andrew Ng著名的公开课CS229的Lecture Notes，你一定记得他用矩阵求导得出的Normal Equation——
你会发现除了\overrightarrow{y} 和\overrightarrow{b} 不一样以外，我们已经把Normal Equation（\overrightarrow{\theta} = (A^TA)^{-1}A^T\overrightarrow{b}）推出来了……我居然在下一部分还没有开始讲就把内容说完了，场面一度非常尴尬啊。可见从投影推出Normal Equation是一件多么自然的事情啊~~~我都不知道哪里切开。
说到这里先总结一下投影的几个意义（敲黑板）！！！
A\overrightarrow{\theta}的所有可能结果都在一个固定的区域中，在线性代数中我们称这个区域为列空间(column space),列空间顾名思义就是矩阵各列的所有线性组合\overrightarrow{\theta_1}\overrightarrow{a_1}+\overrightarrow{\theta_2}\overrightarrow{a_2}+\dots+\overrightarrow{\theta_n}\overrightarrow{a_n}。在1-D的情况下列空间就是一条线，在2-D的情况下列空间就是一个平面。但是我们的数据哪里会这么恰好的落在矩阵的列空间里呢？天底下哪有这样的好事啊！！！
特别是在数据量特别大的情况下，矩阵
A会成为一个
n\gg m的超级高大的
n\times m矩阵（如下图）。在这种等式数量远大于未知数数量的情况中，我们很难满足每一个等式的约束。
A = \begin{bmatrix} 1 & 1\\ 1 & 2\\ . & .\\ . & .\\ . & .\\ 1 & 3\\ \end{bmatrix}但是目标不再在空间里并不代表不能求出解，只能说没有perfect solution（语出Gilbert Strang），但是我们努力一下还是可以做到最好的(best solution)。我们用投影向量\overrightarrow{p}来寻找最合适的\overrightarrow{\theta^*}。\overrightarrow{\theta^*}就是并不存在的完美解\overrightarrow{\theta}的估计值。
三、Normal Equation应用
既然Normal Equation在上文都推导完了，这里我们就随便带几个数据来玩玩咯。
练手案例 找一条直线来拟合点 (1,1)、(2,2)、(3,2)
我们如果用一条直线来拟合的话，设h(\theta) = \theta_0 + \theta_1x_1 ，我们先得到以下值——
\overrightarrow{\theta} = \theta_{0 \dots n} \ (n = 1)
A = \begin{bmatrix} 1 & 1\\ 1 & 2\\ 1 & 3\\ \end{bmatrix}A^T = \begin{bmatrix} 1 & 1 & 1 \\ 1 & 2 & 3 \\ \end{bmatrix}A^TA = \begin{bmatrix} 3 & 6 \\ 6 & 14\\ \end{bmatrix}
\overrightarrow{h(\theta)}= (1,2,2)^T我们发现A\overrightarrow{\theta} =\overrightarrow{h(\theta)} 很遗憾的没有解，于是我们左右各乘上A^T，祭出了投影大招——A^TA\overrightarrow{\theta^*}=A^T\overrightarrow{h(\theta^*)} 。
再把这个方程变换成Normal Equation:\overrightarrow{\theta^*} = (A^TA)^{-1}A^T\overrightarrow{h(\theta^*)}
带入数值在Matlab中小跑一下就得到了结果\overrightarrow{\theta^*} = \begin{pmatrix} \frac 23\\ \\ \frac 12 \\ \end{pmatrix}
即直线h(x) = \frac 23 + \frac 12 x 是上述三个点的拟合结果。
四、其他想说的话
1.关于A^TA\overrightarrow{\theta^*}=A^T\overrightarrow{h(\theta^*)} 的暴力使用
在前一步可以不用判断是否可解，可以直接使用A^TA\overrightarrow{\theta^*}=A^T\overrightarrow{h(\theta^*)} 。事实上，在最小二乘时遇到长方形矩阵A^T，我们就可以用上A^TA替代A^T计算。这是是一种路子很野的但是很简单实用的经验规则，可以简单实验如下——
用直线h(\theta) = \theta_0 + \theta_1x_1 拟合三个点 (1,1)、(2,2)、(3,2)时，自然希望真实值和估计值的误差\left \| h(\theta^*)-h(x)\right \|^2越小越好。
\begin{equation} \begin{split} error=\left \| h(\theta^*)-h(x)\right \|^2 &= (\theta_0 +\theta_1 -1 )^2+(\theta_0 +2\theta_1-2 )^2+(\theta_0 +3\theta_1 -2 )^2\\ \end{split} \end{equation}分别对\theta_0和\theta_1求偏导数等于的零的值——
\begin{equation} \begin{split} \partial error/\partial\theta_0&= 2(\theta_0 +\theta_1 -1 )+2(\theta_0 +2\theta_1-2 )+2(\theta_0 +3\theta_1 -2 )\\ & = 6\theta_0+12\theta_1-10\\ & = 0\\ \end{split} \end{equation}
\begin{equation} \begin{split} \partial error/\partial\theta_0&= 2(\theta_0 +\theta_1 -1 )+2(\theta_0 +2\theta_1-2 )\times 2+2(\theta_0 +3\theta_1 -2 )\times 3\\ & = 12\theta_0+28\theta_1-22\\ & = 0\\ \end{split} \end{equation}整理以上公式我们得到了方程组——
\begin{cases} 3\theta_0+6\theta_1 = 5\\ 6\theta_0+14\theta_1 = -11\\ \end{cases}再整理一下，把这个方程写成矩阵乘法的形式——
\begin{bmatrix} 3&6\\ 6&14 \end{bmatrix} \begin{pmatrix} \theta_0\\ \theta_1 \end{pmatrix}= \begin{pmatrix} 5\\ -11\end{pmatrix}在最后一步整理以后我们发现刚才千辛万苦算出来的\begin{bmatrix} 3 & 6 \\ 6 & 14\\ \end{bmatrix}就是上文的A^TA啊！！！
说明这个经验方法是可以信得过的！！！
2.关于化简的问题
因为投影的性质非常美妙，如果矩阵A是各行线性无关的方阵(square)，说明存在A^{-1}，则Normal Equation会变成如下形式——
\begin{equation} \begin{split} \overrightarrow{\theta} &= (A^TA)^{-1}A^T\overrightarrow{b}\\ &=A^{-1}{A^T}^{-1}A^T\overrightarrow{b}\\ &=A^{-1} \overrightarrow{b}\\ \end{split} \end{equation}说明如果存在一个perfect solution，该解不会受到影响。
3.多次投影有影响吗？
已经在空间中的向量乘上投影矩阵P仍然等于本身，二次投影不会有任何副作用！也就是说P^2 = P。证明如下——
\begin{equation} \begin{split} P^2 &= (A(A^TA)^{-1}A^T)(A(A^TA)^{-1}A^T)\\ &=A(A^TA)^{-1}(A^TA)(A^TA)^{-1}A^T\\ &=A(A^TA)^{-1}A^T\\ &=P \end{split} \end{equation}五、参考资料
1.Gilbert Strang Introduction to Linear Algebra 4.2 Projection 4.3 Least Squares Approximations
2.Andrew Ng CS229 Lecture Note 1 Supervised learning/The normal equations
六、最后的话
列空间没展开讲不知道有没有必要。
笔力不够好，想象中应该写的更简单易懂的，但是没有达到效果，会再更新。
欢迎拍砖！！！