(Conjugate) Gradient Descent
Published:
共轭梯度方法，既是数值算法中的重要内容，也是优化的一类重要算法。本文介绍从梯度方法到共轭梯度方法的分析，着重介绍收敛性对比分析和证明过程。
Line Search Method
Wolfe Condition
线搜索算法在一个射线上寻找满足某些条件的点进行优化，其中$\alpha$为步长或者称为学习率，而$p_k$为一个下降方向，例如一个最简单的选取是最速下降方法，
\[\phi(\alpha) = f(x_k + \alpha p_k)\]步长的选取需要满足某些条件，
\[\begin{align} f(x_k+ \alpha_k p_k) &\le f(x_k) + c_1 \alpha_k \nabla f_k^T p_k,\text{Decrease Condition} \\ f(x_k+\alpha_k p_k)^T p_k &\ge c_2 \nabla f_k^T p_k , \text{Curvature Condition} \\ \text{With } c_1 &\in (0,1),c_2 \in (c_1,1), \nabla f_k^T p_k <0 \end{align}\]可以证明，一定存在满足上述条件的$\alpha$, 首先当$f(x)$有下界的时候，由于直线$f(x_k) + c_1 \alpha_k \nabla f_k^T p_k$无下界，且根据一阶条件，可以知道两者一定存在交点，再根据中值定理可以得到结果一定满足，
\[\begin{align} f(x_k + \alpha' p_k) &= f(x_k) + \alpha' c_1 \nabla p_kf_k^T,\exists \alpha' \\ f(x_k + \alpha' p_k) - f(x_k) &= \alpha' c_1 \nabla p_kf_k^T \\ \nabla f(x_k +\alpha '' p_k)^T \alpha' p_k &= c_1 \alpha' c_1 \nabla p_kf_k^T,\exists \alpha'' \end{align}\]只需要取中值定理构造出的位置$\alpha’’$, 其满足下面的条件，
\[\begin{align} f(x_k+ \alpha_k'' p_k) &\le f(x_k) + c_1 \alpha_k'' \nabla f_k^T p_k \\ f(x_k+\alpha_k'' p_k)^T p_k &\ge c_2 \nabla f_k^T p_k \\ \end{align}\]当被优化的函数$f(x)$满足L-Smooth条件的时候，我们尝试证明线搜索算法的收敛性，
\[\begin{align} (\nabla f_{k+1} - \nabla f_k)^T p_k &\ge (c_2-1) \nabla f_k^T p_k , \text{By Curvature Condition}\\ (\nabla f_{k+1}- \nabla f_k)^T p_k & \le \alpha_kL \Vert p_k \Vert_2^2 ,\text{By Lipschitz Assumption} \\ \alpha_k &\ge \frac{(c_2-1) \nabla f_k^T p_k}{L \Vert p_k \Vert_2^2} \\ f_{k+1} &\le f_k + c_1 \nabla f_k^T p_k \frac{(c_2-1) \nabla f_k^T p_k}{L \Vert p_k \Vert_2^2} \text{,By Decease Condition} \\ f_{k+1} &\le f_k + \frac{c_1(c_2-1) \Vert \nabla f_k \Vert _2^2 \cos \theta_k^2}{L} \text{,Define } \theta_k \\ f_{k+1} & \le f_0 - c \sum_{j=0}^k \cos \theta_j^2 \Vert \nabla f_k \Vert_2^2,\text{Let }c = \frac{c_1(1-c_2)}{L} \\ \sum_{j=0}^k \cos \theta_j^2 \Vert \nabla f_k \Vert_2^2 &\le \infty ,\text{By } f \text{ is Lower Bounded} \\ \lim_{n \rightarrow \infty} \cos \theta_j^2 \Vert \nabla f_k \Vert_2^2 &= 0 \end{align}\]对于最速下降方法，也即每次寻找梯度方向作为$p_k$，由下式可见该算法显然收敛到驻点的位置，而对于其他方法即便不是最速下降方向，利用上述条件也可以给出收敛性的证明，
\[\begin{align} \cos \theta_k =1, \lim_{n \rightarrow \infty} \Vert \nabla f_k \Vert_2^2 &= 0 \end{align}\]Gradient Descend
基于梯度的方法是最速下降方法的拓展，适用于任意的凸函数，当不可导的时候应该将梯度拓展为次梯度或者使用近端梯度下降算法，本节主要关注于最简单的梯度下降算法。
首先正式地给出凸函数的几个等价定义，
\[\begin{align} \alpha f(x) +(1-\alpha)f(y) &\ge f(\alpha x+(1-\alpha)y) \\ f(y) &\ge f(x) + \nabla f(x)^T (y-x) \\ (\nabla f(y) - \nabla f(x))^T (y -x) &\ge 0 \end{align}\]凸函数的全局极小值等价于局部极小值，且局部极小值等价于梯度为0的驻点。
另一个常见的定义是强凸函数，
\[\begin{align} g(x) &= f(x) - \frac{\mu}{2} \Vert x \Vert_2^2 \text{ is convex} \\ f(y) &\ge f(x) + \nabla f(x)^T (y-x) + \frac{\mu}{2} \Vert y-x \Vert_2^2 \\ (\nabla f(y) - \nabla f(x))^T(y-x) &\ge \mu \Vert y-x \Vert_2^2 \end{align}\]以及Lipschitz连续函数，
\[\begin{align} \Vert \nabla f(x) - \nabla f(y) \Vert_2 &\le L \Vert x - y\Vert_2 \\ g(x) &= \frac{L}{2} \Vert x \Vert_2^2 - f(x) \text{ is convex} \\ f(y) &\le f(x) + \nabla f(x)^T (y-x) + \frac{L}{2} \Vert y-x \Vert_2^2 \\ (\nabla f(y) - \nabla f(x))^T (y -x) &\le L \Vert y-x \Vert_2^2 \end{align}\]对于上述两个性质的函数，另外一个重要的性质是其梯度的范数可以被控制住，
\[\begin{align} f(x_{\star}) &= \min_y f(y) \\ & \ge \min_y f(x) + \nabla f(x)^T (y-x) + \frac{\mu}{2} \Vert y-x \Vert_2^2 \\ &= f(x) - \frac{1}{2\mu} \Vert \nabla f(x) \Vert_2^2,\text{Let } y = x - \frac{1}{\mu} \nabla f(x) \\ f(x_{\star}) &= \min_y f(y) \\ &\le \min_y f(x) + \nabla f(x)^T (y-x) + \frac{L}{2} \Vert y-x \Vert_2^2 \\ &= f(x) - \frac{1}{2L} \Vert \nabla f(x) \Vert_2^2,\text{Let } y = x - \frac{1}{L} \nabla f(x) \\ \end{align}\]合起来可以得到，
\[\begin{align} \frac{1}{2L} \Vert \nabla f(x) \Vert_2^2 \le f(x) - f(x_{\star}) \le \frac{1}{2\mu} \Vert \nabla f(x) \Vert_2^2 \end{align}\]另外一个重要的性质称为矫顽力(coercivity），此处略过推导过程，
\[\begin{align} (\nabla f(x) - \nabla f(y))^T(x -y) &\ge \frac{1}{L} \Vert \nabla f(x) - \nabla f(y) \Vert_2^2 \\ (\nabla f(x) - \nabla f(y))^T(x -y) &\le \frac{1}{\mu} \Vert \nabla f(x) - \nabla f(y) \Vert_2^2 \\ \end{align}\]对于同时具有上述两个性质的函数，还可以得到，
\[\begin{align} h(x) &=f(x) - \frac{\mu}{2} \Vert x \Vert_2^2 \text{ is } L- \mu \text{ smooth and convex} \\ 0 &\le (\nabla h(x) - \nabla h(y))^T (x - y ) \\ &= (\nabla f(x) - \nabla f(y))^T(x-y) - \mu \Vert x - y\Vert_2^2 \\ &\le (L- \mu) \Vert x- y\Vert_2^2 \end{align}\]更进一步可以验证，
\[\begin{align} (\nabla h(x) - \nabla h(y))^T(x-y) &\ge \frac{1}{L - \mu} \Vert \nabla h(x) - \nabla h(y) \Vert_2^2 \\ ((\nabla f(x) - \nabla f(y)) - (\mu x- \mu y))^T (x- y ) &\ge \Vert(\nabla f(x) - \nabla f(y)) - (\mu x- \mu y) \Vert_2^2 \\ (L+ \mu) (\nabla f(x) - \nabla f(y))^T (x-y) & \ge L \mu \Vert x- y\Vert_2^2 + \Vert \nabla f(x) - \nabla f(y) \Vert_2^2 \end{align}\]对于梯度下降方法，在满足优化函数Lipschitz连续的情况下，可以证明其收敛性，但却是一个次线性收敛的结果，在本节中我们简单地采取学习率为$\frac{1}{L}$,
\[\begin{align} x_{k+1} &= x_k - \frac{1}{L} \nabla f(x_k) \\ f(x_{k+1}) - f(x_\star) &= f(x_k - \frac{1}{L} \nabla f(x_k)) - f(x_{\star} ) \\ &\le f(x_k ) - \frac{1}{2L} \nabla f(x_k)^T \nabla f(x_k) -f(x_{\star} ), \text{By Lipschitz } \\ &\le - \nabla f(x_k)^T(x_{\star} - x_k) - \frac{1}{2L} \nabla f(x_k)^T \nabla f(x_k),\text{By Convex } \\ &= \frac{L}{2} ( \Vert x_{\star} - x_k \Vert_2^2 - \Vert x_{\star} - x_k - \frac{1}{L} \nabla f(x_k) \Vert_2^2) \\ &= \frac{L}{2} ( \Vert x_{\star} - x_k \Vert_2^2 - \Vert x_{\star} - x_{k+1} \Vert_2^2) \\ f(x_k) - f(x_{\star}) &\le \frac{L}{2k}\sum_{i=1}^k ( \Vert x_{\star} - x_k \Vert_2^2 - \Vert x_{\star} - x_{k+1} \Vert_2^2) \\ &= \frac{L}{2k}( \Vert x_{\star} - x_0 \Vert_2^2 - \Vert x_{\star} - x_{k} \Vert_2^2) \\ &\le \frac{L}{2k} \Vert x_{\star} - x_0 \Vert_2^2 \end{align}\]在此基础上可以发现其达到$\epsilon$-最优需要的迭代次数为，
\[N \ge O(\frac{1}{\epsilon})\]而如果附加上优化函数为强凸函数的假设，梯度算法可以做到线性收敛，
\[\begin{align} x_{k+1} &= x_k - \frac{1}{L} \nabla f(x_k) \\ f(x_{k+1}) - f(x_\star) &= f(x_k - \frac{1}{L} \nabla f(x_k)) - f(x_{\star} ) \\ &\le f(x_k ) -f(x_{\star} ) - \frac{1}{2L} \Vert \nabla f(x_k) \Vert_2^2\\ &\le f(x_k ) -f(x_{\star} ) - \frac{\mu}{L} (f(x_k ) -f(x_{\star} )) \\ &= (1- \frac{\mu}{L})(f(x_k ) -f(x_{\star} )) \end{align}\]在此基础上可以发现其达到$\epsilon$-最优需要的迭代次数为，
\[N \ge O(\log \frac{1}{\epsilon})\]而对于距离最优点的距离，需要用到矫顽力的性质，并且需要选取更小的步长，$\alpha_k = \frac{2}{\mu+L}$
\[\begin{align} \Vert x_{k+1} - x_{\star} \Vert_2^2 &= \Vert x_{k} - \frac{2}{L+ \mu} \nabla f(x_k) - x_{\star} \Vert_2^2 \\ &= \Vert x_k - x_{\star} \Vert_2^2 -\frac{4}{L+\mu}(x_{k}-x_{\star})^T \nabla f(x_k) + \frac{4}{(L+\mu)^2}\Vert \nabla f(x_k) \Vert_2^2 \\ & \le \Vert x_k - x_{\star} \Vert_2^2 + \frac{4}{(L+\mu)^2}\Vert \nabla f(x_k) \Vert_2^2 - \frac{4}{(L+ \mu)^2}( L \mu \Vert x_k - x_{\star}\Vert _2^2 + \Vert \nabla f(x_k) \Vert_2^2) ,\text{By Coercivity } \\ &= (1- \frac{4 L \mu}{(L+\mu)^2}) \Vert x_k - x_{\star} \Vert_2^2 \end{align}\]可以发现，如果选取了合适的步长，梯度下降算法在连续强凸函数上面可以做到线性收敛。
Cojugate Gradient Method
Steeping Gradient Descend
共轭梯度方法最开始源于二次正定优化问题，或者解正定线性系统，
\[\begin{align} &Ax = b, A >0 \\ \text{Or } &\min_x f(x) = \min_x\frac{1}{2} x^T Ax - b^T x \end{align}\]最简单的想法是使用最速下降方法，每次寻找的方向为梯度的方向，之后在该方向上寻找使得函数值最小的步长，
\(\begin{align} p_k &= -\nabla f(x_k)= r_k = b - Ax_k \\ \alpha _k &= \min_{\alpha} f(x_k+ \alpha r_k) \\ \end{align}\) 利用梯度求解最优的步长，
\[\begin{align} \nabla f(x_k+ \alpha_k r_k) &= 0 \\ r_k^T (A(x_k + \alpha_k r_k) - b) &= 0 \\ \end{align}\]可以解得，
\[\begin{align} \alpha_k &= \frac{r_k^T(b-Ax_k)}{r_k^TA r_k} = \frac{r_k^T r_k}{r_k^T A r_k} \end{align}\]有意思的是可以发现最速下降方法每次的搜索方向实际上为正交方向，
\[\begin{align} r_{k}^T r_{k+1} &= r_k^T (b - A x_{k+1} )\\ &= r_k^T (b - A(x_k +\alpha_k r_k)) \\ &= r_k^T(r_k -A \alpha_k r_k) \\ &= r_k^T r_k - r_k^TA r_k \alpha_k \\ &= 0, \text{With } \alpha_k = \frac{r_k^T r_k}{r_k^T A r_k} \end{align}\]Cojugate Gradient Descent
我们已经看到，在最速下降方法中，每一步选取的下降方向为正交方向，而共轭梯度下降方法中下降方向的选取中用到了矩阵$A$的信息，共轭方向被定义为，
\[p_i^T A p_j = 0 ,\forall i \ne j\]共轭方向虽然不是正交的，但在矩阵$A$为正定矩阵的前提下其为线性无关的方向，可以通过反证法，
\[\begin{align} \text{Assume that } &\exists \alpha_i \ne 0, \sum_{j=1}^N \alpha_ip_j = 0 \\ \Rightarrow & \sum_{j=1}^N \alpha_i Ap_j = 0 \\ \Rightarrow &\alpha_i p_i^T A p_i = 0 \\ \Rightarrow &\alpha_i = 0 \\ \end{align}\]下降的步长同样可以通过选取在该方向上的最优值得到，
\[\begin{align} \nabla f(x_k +\alpha_k p_k) &= 0 \\ p_k^T (A(x_k - \alpha_kp_k)-b) &= 0\\ \alpha_k &= \frac{p_k^T r_k}{p_k^T A p_k} \end{align}\]对于所有的下降方向为基进行展开，可以神奇地发现共轭梯度下降方法可以保证在有限步内收敛，通过对比在基展开下的系数和步长，可以发现共轭方法中寻找到的步长正好是基的展开系数，
\[\begin{align} x_{\star}- x_0 &= \sum_{i=1}^N \beta_i p_i \\ \beta_k & = \frac{p_k^TA(x_{\star}-x_0)}{p_k^T A p_k} \\ &= \frac{p_k^TA(x_{\star}-x_k +x_k-x_0)}{p_k^T A p_k} \\ &= \frac{p_k^TA(x_{\star}-x_k)}{p_k^T A p_k}, \text{With } p_k^T A(X_k -x_0) = 0 \\ &= \frac{p_k^T r_k}{p_k^T A p_k} = \alpha_k \end{align}\]类似地，利用归纳法也可以证明共轭梯度方法实际上每步在由$[p_0,…,p_k]$张成的子空间上最小化目标函数，本质上共轭梯度下降方法是一种Krylov子空间迭代法，利用归纳法还可以证明$p_k,r_k$分别是该子空间的一组共轭基和正交基，证明暂略。从这个角度来说，共轭梯度下降方法显然是优于最速下降方法的，因为最速下降方法仅仅保证了相邻两次下降的残差是正交的，但共轭梯度下降方法通过下降方向的合理选取，保证了任意两次下降的方向都是相互正交的，并且每一次下降都做到了该方向上的极致。
\[\begin{align} \text{span}(r_0,...,r_k) = \text{span}(p_0,...,p_k) = \text{span}(r_0,Ar_0,...,A^kr_0) \end{align}\]Convergence
本节关注于算法的收敛性，从收敛性的角度证明共轭梯度方法的优势所在。
下面来看最速下降方法在该问题上的收敛性，其证明首先需要用到向量的$A$范数，容易验证当$A$为正定矩阵的时候，其满足范数的定义，
\[\Vert x \Vert_A = \sqrt{x^T Ax}\]对于$A$范数与矩阵函数的关系，有一个常用的引理，
\[\Vert p(A) x\Vert_A \le \max_i \vert p(\lambda_i) \vert \Vert x \Vert_A\]对向量$x$在$A$的特征向量上做展开之后放缩即可，
\[\begin{align} \Vert p(A) x\Vert_A^2 &= \Vert \sum_{i=1}^n \beta_i \phi_i p(A) \Vert_A^2 \\ &= \Vert\sum_{i=1}^n \beta_i \phi_i p(\lambda_i)\Vert_A^2 \\ &= \sum_{i=1}^n \lambda_i \beta_i^2 p^2(\lambda_i) \\ &\le \max_{i} p^2(\lambda_i) \sum_{i=1}^n \lambda_i \beta_i^2 \\ &= \max_i p^2(\lambda_i) x^T Ax \end{align}\]和上述引理配合使用的是Chebshev多项式，在所有满足$p_k(0)=1$和$k$阶多项式中，
\[\begin{align} \min_{p_k} \max_{a \le \lambda \le b} \vert p_k(\lambda ) \vert = \frac{T_k(\frac{b+a-2x}{b-a})}{T_k(\frac{b+a}{b-a})} \end{align}\]证明用到了ChebShev多项式的交错性质，首先上述多项式显然是满足条件的一个解，如果存在一个更优的多项式，可以发现该多项式与原多项式之差在$\lambda=0$ 上有一个零点，且在区间$[a,b]$上根据ChebShev多项式的性质也有$n$个零点，则该差值只能为零次多项式，那么通过反证法就说明了上述问题的解。关于多项式的逼近，可以参见 函数逼近与拟合
下面我们再证明一个小引理，说明了$x_k,r_k$的收敛性之间的关系，因此我们不仅仅关心点是否收敛，还关心最终的残差或者梯度是否在很小的误差范围内，
\[\begin{align} f(x) - f(x_{\star}) &= (\frac{1}{2}x^T Ax - b^T x) -(\frac{1}{2}x_{\star}^T A x_{\star} - b^T x_{\star})\\ &= (\frac{1}{2}x^T Ax - x_{\star}^T A x) -(\frac{1}{2}x_{\star}^T A x_{\star} - x_{\star}^T A x_{\star}) \\ &=\frac{1}{2} (x - x_{\star})^T A(x- x_{\star}) \\ &= \frac{1}{2} \Vert x - x_{\star} \Vert_A^2 \end{align}\]我们知道在线搜索方法中每次都寻找到了该方向上的最优值，因此我们可以知道，
\[\begin{align} \frac{1}{2} \Vert x_{k+1} - x_{\star} \Vert_A^2 &= f(x_{k+1}) - f(x_{\star}) \\ &=f(x_{k} + \alpha_k r_k) - f(x_{\star}) \\ &= \min_{\alpha} f(x_{k} + \alpha r_k) - f(x_{\star}) \\ &= \min_{\alpha} \frac{1}{2} \Vert x_k+ \alpha r_k - x_{\star} \Vert_A^2 \\ &= \min_{\alpha} \frac{1}{2} \Vert (I-\alpha A)(x_k - x_{\star})\Vert_A^2 \end{align}\]利用$A$范数与矩阵多项式函数的引理，并且利用Chebshev多项式的最佳一致逼近性质可以得到，
\[\begin{align} \Vert x_{k+1} - x_{\star} \Vert_A &\le \min_{\alpha}\max_{\lambda_i}\vert (1-\alpha \lambda_i) \vert \Vert (x_k - x_{\star})\Vert_A \\ &\le \min_{\alpha}\max_{\lambda_1 \le \lambda \le \lambda_n} \vert 1-\alpha \lambda \vert \Vert (x_k - x_{\star})\Vert_A \\ &= \min_{p_1(\lambda)} \max_{\lambda_1 \le \lambda \le \lambda_n} \vert p_1(\lambda ) \vert \Vert (x_k - x_{\star})\Vert_A \\ &=\frac{1}{T_1(\frac{\lambda_n + \lambda_1}{\lambda_n- \lambda_1})} \Vert (x_k - x_{\star})\Vert_A \\ &= \frac{\lambda_n - \lambda_1}{\lambda_n + \lambda_1} \Vert (x_k - x_{\star})\Vert_A \\ &= \frac{\kappa-1}{\kappa+1} \Vert (x_k - x_{\star})\Vert_A ,\text{With }\kappa = \frac{\lambda_n}{\lambda_1} \end{align}\]据此可以得到其收敛性，为了达到$\epsilon$-最优的解，
\[\begin{align} f(x_{N}) - f(x_{\star}) &\le f(x_0) - f(x_{\star}) (\frac{\kappa-1}{\kappa+1})^{2N} \le \epsilon (f(x_0) - f(x_{\star})) \\ (\frac{\kappa-1}{\kappa+1})^{2N} & \le \epsilon \\ 2N \log (1- \frac{2}{\kappa+1}) & \le \log \epsilon \\ N &\le -(\kappa+1) \log \epsilon, \kappa \rightarrow \infty \\ N &\ge O[(\kappa+1) \log \frac{1}{\epsilon}] \end{align}\]经过一番功夫我们得到了最速下降方法的收敛率，下面我们对比共轭梯度方法的收敛率，
利用Krylov子空间的结论，以及Chebshev多项式的性质，其中$p_k(0)=1$，
\[\begin{align} \Vert x_k - x_{\star} \Vert _A &= \min_{p_k} \Vert p_k(A) A^{-1} r_0 \Vert_A \\ &= \min_{p_k} \Vert p_k(A)(x_0 - x_{\star}) \Vert_A \\ &= \min_{p_k} \max_{\lambda_i} \vert p_k(\lambda_i) \vert \Vert (x_0 - x_{\star}) \Vert_A \\ & \le \min_{p_k} \max_{\lambda_1 \le \lambda \le \lambda_n} \vert p_k(\lambda) \vert \Vert (x_0 - x_{\star}) \Vert_A \\ &= \frac{1}{T_k(\frac{\lambda_n+\lambda_1}{\lambda_n- \lambda_1})} \Vert (x_0 - x_{\star}) \Vert_A \\ & \le \frac{1}{\frac{1}{2}[(x+\sqrt{x^2-1})^k+(x-\sqrt{x^2-1})^k]} \Vert (x_0 - x_{\star}) \Vert_A ,\text{With }x = \frac{\lambda_n+\lambda_1}{\lambda_n- \lambda_1} \\ & \le (\frac{\sqrt{\kappa}-1}{\sqrt{\kappa}+1})^k \Vert (x_0 - x_{\star}) \Vert_A \end{align}\]完全类似地我们可以证明对于共轭梯度方法，为了达到$\epsilon$-最优解，
\[N \ge O[(\sqrt{\kappa}+1) \log \frac{1}{\epsilon}]\]因此从收敛的角度考虑，共轭梯度下降方法也是一个更优的方法，实际上面的推导过程仅仅得到了一个较为宽松的界，采用其他技术可以得到更紧的界，但上述的证明已经足以对比两种方法。