如何理解正规方程

在《如何理解线性回归》（建议先看）一文中，我们提到了使用梯度下降来求解最小二乘问题，但在线性回归中，还有一种不需要迭代的方法来求解最小二乘问题，这就是正规方程。这是基于矩阵求导来及计算的（不会LaTeX看起来有点糟心）
还是上次的数据

x = [[1, 1, 1, 1, 1], [1, 2, 3, 4, 5]]
y = [4.8, 7.2, 8.8, 11.2, 12.8]
plt.scatter(x[1], y)
plt.xlim(0,6)
plt.ylim(0,14)
plt.show()

如何理解正规方程

和线性回归前半部分的理解一样，我们的目标是最小化J，只是这里使用的J是J=1/2*∑(y-y_h)^2（我们用y表示真实值，y_h表示预测值），标准的写法是这样的，之前我们在用梯度下降时J没有1/2也没有关系。那么现在我们就要用矩阵的形式来表示J，并取它的导数，当导数等于零时，也就可以得到最终的正规方程，下面是J的矩阵表示

J = 1/2 * ∑(y - y_h)^2 
  = 1/2 * tr[(Xθ - Y)T (Xθ - Y)]

tr是什么，tr表示矩阵的迹，也就是矩阵主对角线元素的和，表示为trA=∑aii，其中有些重要的定理

tr(AB) = tr(BA)
tr(A+B) = tr(A) + tr(B)
d(tr(AB))/d(A) = BT
d(tr(ATB))/d(A) = B
tr(A) = tr(AT)
tr(a) = a（a ∈ R）
d(tr(ABATC))/d(A) = CAB+CTABT

就这个问题，我们可以有如下的表示

$xi= \left( \begin{matrix} 1\\ x \end{matrix} \right)\tag{1}这里的x就是上面的[1, 2, 3, 4, 5]$ $θ= \left( \begin{matrix} b\\ a \end{matrix} \right)\tag{2}这是y=a*x+b的参数$ $X= \left( \begin{matrix} x0T\\ x1T\\ x2T\\ x3T\\ x4T \end{matrix} \right)\tag{3}xiT就是xi的转置$ $Y= \left( \begin{matrix} y0\\ y1\\ y2\\ y3\\ y4\\ \end{matrix} \right)\tag{4}这里的yi就是上面的[4.8, 7.2, 8.8, 11.2, 12.8]$

上面的x是一个2x1的向量，这样xTθ就是一个1x2的向量乘一个2x1，就得到一个答案，就是a*x+b，这个很眼熟吧。所以XTθ-Y就包含了所有的y-y_h=a*x+b-y_h，转置相乘再去迹就可以得到所有式子的平方和。这个式子比较好理解
麻烦的就是下面d(J)/d(θ)。首先暴力拆开(Xθ-Y)T(Xθ-Y)

  (Xθ - Y)T (Xθ - Y) 
= (θTXT - YT)(Xθ - Y)
= θTXTXθ - θTXTY - YTXθ + YTY

所以这时要解d(J)/d(θ)，我们算一下维度看，θTXTXθ是(1x2 x 2x5 x 5x2 x 2x1)=(1x1)，θTXTY是(1x2 x 2x5 x 5x1)=(1x1)，还有YTXθ是(1x5 x 5x2 x 2x1)=(1x1)，YTY是(1x5 x 5x1)=(1x1)，都是些1x1的阵。苦逼没学过矩阵的求导，只有马上去学了，至少可以看出最后的YTY与θ半毛钱关系都没有，所以就可以忽略了，得到下面的式子

d(J)/d(θ) = 1/2 * d(tr(θTXTXθ - θTXTY - YTXθ + YTY)) / d(θ)
          = 1/2 * d(tr(θTXTXθ - θTXTY - YTXθ)) / d(θ)

然后我们发现θTXTY是YTXθ的转置，也就是(YTXθ)T=θTXTY。由迹的一些定理（上面的2和5），我们可以这样处理一下

d(J)/d(θ) = 1/2 * d(tr(θTXTXθ - θTXTY - YTXθ + YTY)) / d(θ)
          = 1/2 * d(tr(θTXTXθ - θTXTY - YTXθ)) / d(θ)
          = 1/2 * (d(tr(θTXTXθ) / d(θ) - d(tr(θTXTY) / d(θ) - d(tr(YTXθ) / d(θ))
          = 1/2 * (d(tr(θTXTXθ) / d(θ) - 2*d(tr(θTXTY) / d(θ))

然后现在需要求导了，这里有一些有用的求导的定理

d(a)/d(x) = 0
d(x)/d(x) = I （I是单位矩阵）
d(Ax)/d(x) = AT
d(xTA)/d(x) = A
d(ATBA)/d(A) = (B+BT)A

现在就可以继续向下做了，由求导的5和迹的4可以得到如下的结果

d(J)/d(θ) = 1/2 * d(tr(θTXTXθ - θTXTY - YTXθ + YTY)) / d(θ)
          = 1/2 * d(tr(θTXTXθ - θTXTY - YTXθ)) / d(θ)
          = 1/2 * (d(tr(θTXTXθ) / d(θ) - d(tr(θTXTY) / d(θ) - d(tr(YTXθ) / d(θ))
          = 1/2 * (d(tr(θTXTXθ) / d(θ) - 2*d(tr(θTXTY) / d(θ))
          = 1/2 * (2*XTXθ - 2*XTY)
          = XTXθ - XTY

这是让d(J)/d(θ)=0，就是导数等于零，我们知道在J最小时拟合最佳，而J的导数为零时对应J的一个极小值，这里就是它的最小值，得到

 d(J)/d(θ) = 0
XTXθ - XTY = 0
      XTXθ = XTY
         θ = ((XTX)-1)XTY

这样就得到了正规方程θ=((XTX)-1)XTY，这里我们算算θ的维度，((XTX)-1)XTY得到(2x5 x 5x2 x 2x5 x 5x1)=(2x1)，可以看到这里的θ就是一个2x1的向量，里面包含了所有参数的值，这里只有a和b，所以这里的θ[0]=b，而θ[1]=a
我们回到开头的那个问题，看看数据和数据的维度

th = np.zeros((2, 1))  # 这是θ
X = np.array(x).T
Y = np.array(y).reshape(5, 1)
print(th.shape)
print(X.shape)
print(Y.shape)

这里的输出是

(2, 1)
(5, 2)
(5, 1)

现在我们构造这个正规方程并且计算，在python-numpy中，np.dot()是矩阵相乘，np.linalg.inv()是求逆矩阵，现在输出看看

th = np.dot(np.dot(np.linalg.inv(np.dot(X.T, X)), X.T), Y)
print('a='+str(th[1])+' b='+str(th[0]))

结果为，也就是y=2*x+2.96

a=[2.] b=[2.96]

我们画个图看看

plt.scatter(x[1], y)
x2 = np.linspace(0, 6, 100)
fy = th[1]*x2+th[0]
plt.plot(x2, fy)
plt.xlim(0,6)
plt.ylim(0,14)
plt.show()

如何理解正规方程

这个数据是按照y=2*x+3构造的，这个结果和梯度下降的结果是差不多的，这就是矩阵求导最小二乘的解法
梯度下降和正规方程都适用于求最小二乘的一种方法，相同点在于它们都求了导数，目标都是最小化J；不同在于梯度下降需要设定学习率和需要多次迭代，而正规方程不需要。一般使用时，n<10000时且为线性模型我们可以用正规方程，当n>10000时由于逆矩阵的运算复杂度高，效率很低；而且正规方程只适用于线性模型，梯度下降法则可以适用各种模型，是更常用的做法

如何理解正规方程

相关推荐