1.4 代数多重网格方法的理论分析

1.4.1 从一个差分方程到一个线性代数问题

一个差分方程，如果需要利用计算机进行求解，首先需要将其离散化。例如，求解一个简单的差分方程

满足边界条件为

其中， L 是一个差分算子； 是一个单一变量或者多个变量的函数（ （ x ， y ，…， t ））。使用一个简单的一维空间的例子， 。差分后可得

等效于求解

1.4.2 单一网格上求解线性方程的问题

定义线性系统 Ax = b ，残差 r _k = b-Ax _k ，预调整表达式 P ≈ A 。

预调整表达式 P 要尽量接近 A ，同时允许快速递归。Jacobi选择 A 的三角矩阵作为 P 是一个极端，速度快但是与 A 之间相差较大。另一个极端是 P = A 。将矩阵 A 分开得到 Ax = b 的新的表达式

相应地得到递归的另一种表达式

从任意给定的 x ₀ 开始，然后从 Px ₁ =（ P-A ） x ₀ + b 得到 x ₁ 。继续使用同样的矩阵得到 x ₂ 。所以通常通过 P = LU 分解来得到 P 的三角因子。有时 P 本身就是三角的，或者 L 、 U 有时近似于 A 的三角因子。 P 如果满足以下两个条件，则可以得到有效解。

（1）新的 x _{k +1} 必须能迅速计算出来；

（2） e _k = x-x _k 必须尽快地收敛到0。

用式（1.19）中减去式（1.18），可得

等效于

可以发现 e _k 的计算与 b 无关。每一步只需要将 e _k 乘以一个矩阵 M 即可。 x _k 收敛到 x 的速度完全依赖于矩阵 M 。矩阵 M 收敛必须满足

矩阵中绝对值最大的特征值也被称为谱半径 ρ （ M ）= 。矩阵 M 收敛要求 ρ （ M ）＜1。收敛率跟最大的特征值有关。当一个矩阵比较大时， ρ （ M ）=0.9甚至于 ρ （ M ）=0.99时是最好的。

当初始错误 e ₀ 刚好等于矩阵的特征向量时，下一个错误值 e ₁ = Me ₀ = λe ₀ 。随后的每一步的错误都要乘以一个 λ 。所以我们需要 ＜1。正常情况 e ₀ 是所有特征向量的线性组合结果。当乘以 M 时，每个特征向量都乘以自己的特征值。 k 步之后，每个向量都乘以 λ ^k ，其中最大的是（ ρ （ M ）） ^k 。

如果不使用预调整表达式，则 M = I-A 。为了收敛， A 的所有的特征值都必须在以1为中心的单位圆的内部。Jacobi矩阵 M _J = I-D ^-1 A = D ^-1 （ L + U ）的特征值为

，所以肯定收敛，但是收敛速度比较慢。

j =1时特征值最大，其频谱半径为

从式（1.24）中可以发现，低频时收敛较慢。而加权 Jacobi 方法相对较好， M _w = I-wD ^-1 A =[（1 -w ） I + wM _J ]，比较结果如图1.16所示。

从图1.16中的曲线可以看出，对于低频分量 衰减都不会很小，特别当 k =1时，误差分量相对最大；而对于高频分量 ，则可以通过 w 进行选择而控制衰减因子。当 w =2/3时，一次迭代就可以使得高频误差分量衰减到初始误差分量的1/3。这些结果表明，低频误差分量是迭代收敛慢的原因，而高频误差分量则在迭代中衰减得很快。换句话说，随着迭代次数的增加，误差的光滑性也进一步增强。对于这种现象，其他的迭代方法也是存在的。

1.4.3 多网格上求解线性方程的解决方案

多网格方法是在离散过程中利用网格属性影响收敛速度的特性的一种迭代求解方法。从傅里叶分析中可以看出，网格宽度（与空间频率成反比）是主导收敛的一个重要因素。

重新在粗网格上就选择初始的猜测值，然后在粗网格上选取其他任意一个点。粗网格上只有细网格上一半的网格点数。细网格上的傅里叶模式可以重新写成

这意味着在 Ω ^h 上的误差分量变成了 Ω ^{2 h} 上的误差分量，如图1.16所示。从图1.16上可以看出，细网格上的低频分量在粗网格上来看比在细网格上要振荡。从中可以看出，将细网格上的光滑误差分量限制到粗网格上，然后在粗网格上进行迭代处理就可以消除部分误差分量，基于这种思想我们可以一层一层地做下去直到消除所有的误差分量为止，最后还需要将求解问题一层一层地返回，直到最初的细网格上。

以一个二层网络为例来说明。对误差 e 使用式（1.22）进行平滑，如图1.17所示。根据 ，可以得到如图1.18所示的结果。

从图1.19中也可以发现当层次越高时，求解越迅速越精确。因为此时数据量小，进行松弛时的时空耗费小，而且收敛快。而在最细层中初始值和精确值之间存在较大的误差，而且毫无规律性。但是提取的粗网格层也是有限的，当提取的粗网格层过粗，或者原始数据波动特别迅速时，此时使用多重网格方法就存在较大的问题。当层次越高，其大致形状跟原始数据之间存在较大差异，插值回去时会出现较大误差，这时就需要进行重新考虑了。

1.4.4 粗网格在图分类方法中的应用

基于特征值方法的图分类方法的本质是在求解一个最优化问题。将图分成两类，表达形式如式（1.26）。

式（1.26）等效于

其中， a _ij = -w _ij ， a _ij =∑ _j w _ij ， b _ij = δ _ij 。得到的 A 就是图的拉普拉斯矩阵。

根据拉格朗日定理可知， x 的求解等效于求解式（1.28）的特征值。

当矩阵比较大的时候，这时需要考虑对该图在粗网格层上来进行处理。通过粗网格划分，组成一个新的由粗网格上的节点构成的矩阵 A _c 。

现在，优化式（1.29）等效于求解

1.4.5 利用拉普拉斯矩阵提取粗网格

拉普拉斯矩阵在目前的流形学习中也有很重要的地位，使用拉普拉斯作为流形学习的基本思想是在高维空间中距离很近的点在低维空间中投影的距离也应该离得很近，通过对图像进行拉普拉斯投影到高维空间，然后通过求解特征向量的方式将其投影回低维空间，可以得到图像中像素聚类的一些关系。因此，这里使用拉普拉斯矩阵来进行相应的粗网格的提取。

拉普拉斯矩阵的优点主要有：（1）拉普拉斯矩阵是个“0和”矩阵，最小的特征值为0，其他的特征值都大于0；（2）拉普拉斯矩阵的谱比邻接矩阵的谱更能反映图的图论性质，因此拉普拉斯矩阵的谱分析理论目前研究较多，较为成熟。

一个实对称矩阵可以表示成

其中 =1。

从式（1.31）中可以看出，一个矩阵可由它的特征向量完全表示，而每个向量所对应的特征值，就代表了矩阵在这个向量上的贡献率。

矩阵的能量可以表示为

根据式（1.32）的关系，对特征值进行排序，保留其中较大的特征值，也就可以保留对应的一些较大非对角线元素。

针对原始图像构建邻接矩阵，对邻接矩阵求解拉普拉斯矩阵。求解拉普拉斯矩阵的特征值和特征向量，根据一定的准则选择几个特征值，将这几个特征值对应的节点选取出来作为粗网格的节点。这种方法与主成分分析（PCA）方法的过程和结果是比较类似的。在PCA中，通过构建几个分量的协方差矩阵，求解特征方程得到特征值，然后按照寻找几个较大的特征值就可以得到分量在一个低维空间的近似的分布，可以达到近似和降维的目的。

在邻接矩阵中需要选择从大到小的特征值进行排序，提取粗网格节点。而在拉普拉斯矩阵中需要提取从小到大的特征值进行排序，这里给出证明。需要证明的问题表述为：度相等的拉普拉斯矩阵和邻接矩阵的特征值之间存在 λ ′= d-λ ，其中 d 为度矩阵的特征值， λ 为拉普拉斯矩阵的特征值， λ ′为邻接矩阵的特征值。

证明

L 矩阵为拉普拉斯矩阵，可以分离成度矩阵 D 和邻接矩阵 W 的差，其中 w _ii =0， w _ij 为邻接矩阵。

可以得到

进一步可得

其中 D ′中 = D _nn -λ 。 D _nn =…= D ₀₀ = d 时，可以将其转换成

这表示度相等的矩阵的拉普拉斯矩阵和邻接矩阵的特征值之间存在

对于其他度不相等的情况，通过理论很难去证实。因为目前大部分谱分析的研究都集中在研究每个特征值的上界和下界，而没有对两个具体的特征值之间的关系进行系统的研究。因为考虑到上面的研究结果，在此对两个特征值之间的关系进行理论和实践的研究。通过选择随机矩阵的方式来进行描述。通过随机产生邻接矩阵，构建度矩阵和拉普拉斯矩阵，分别计算特征值，根据大量随机矩阵的方式去进行分析。图1.20和图1.21分别为度矩阵相同（矩阵规模为2500）和度矩阵不同（矩阵规模为50）时各个特征值之间的关系图。图中， λ （ L ）为从大到小，而 λ （ W ）从小到大以及 λ （ D ）从小到大。图1.20中存在 λ （ L ）+ λ _n-k （ W ）= λ _k （ D ）的关系，只有几个节点除外，主要是因为边界效应（ d ₁ = d _n = d _i -1， i =2～ n ）的原因，这也验证了式（1.38）的推断。

通过图 1.20 与图 1.21 的对比分析，图 1.20 描述的是当度相同的时候，存在 λ _n -k （ L ）+ λ _k （ W ）= λ _k （ D ）的关系。也就是 λ _n-k （ L ）+ λ _k （ W ）和 λ _k （ D ）的差值能够从一定程度上反映图中不同度的分布情况，两者差值为0的临界值对度的分布进行进一步分析具有一定的指导意义。我们定义这个临界值为 k ₁ ，在 k ₁ 之下 λ _n-k （ L ）+ λ _k （ W ）≤ λ _k （ D ） （其中 λ （ L ）为从小到大排列，而 λ （ W ） 从大到小排列）。从大量的实验分析可看出， k ₁ 之上就是 λ _k （ L ）+ λ _n-k （ W ）＞ λ _k （ D ），但也有一些例外，也就是中间存在一些较小的波动，因此定义另一个临界值 k ₂ ≥ k ₁ +1，在 k ₂ 之上存在 λ _k （ L ）+ λ _n-k （ W ）＞ λ _k （ D ）。实际上这个波动的产生跟矩阵的内在特征也有很大的关系。

定理 1-3： 如果 A 和 A + E 是 n × n 的实对称矩阵，其中矩阵的特征值的排序如 λ _n ≤ λ _{n -1} …≤ ₁ λ ，有

证明：必定存在一个点 k ₁ ，使得 λ _{n-k +1} （ L ）+ λ _k （ W ）≤ λ _k （ D ）。

因为拉普拉斯矩阵和邻接矩阵的关系为 L + W = D ，根据定理1-3，有

因为 λ _n ≤ λ _{n -1} …≤ ₁ λ ，因此必定存在一个点 k ₁ ，使得

根据这个问题定义一个新的临界值 k ₁ 和 k ₂ 。使用 k ₁ 作为选择拉普拉斯特征值个数的依据。使用下面的例子来验证选择的正确性。如图1.22所示，计算得到 k ₁ =4，因此根据特征值排序，得到特征值较小的四个节点，如图1.22（b）所示。使用粗网格提取的邻接图的结果如图1.23所示。可见其中粗网格点都被提取出来了。尤其是图1.23（b）中B节点，从图1.22（b）中发现其度为3，跟E节点的度一样，但是因为它在原始图中的度比E节点的高，因此还是保留在粗网格之中了。