拉普拉斯矩阵和基尔霍夫定律

随着科技时代的发展，我们在日常生活和工作实践中使用了很多‘图’。‘图’形象直观，可以包含很多信息，节省了语言的阐述。图和数的结合成为解决很多实际问题的利器，而矩阵是联系图和计算机的有力工具，本文阐述在这方面至关重要的拉普拉斯矩阵，它的构造，性质和简单应用。

（一）图(G)，邻接矩阵（W），度矩阵(D)，拉普拉斯矩阵(L) 我们先看一张图：

n个节点A,B,C......由边连接起来，边表示节点与节点之间的联系，每一条边对应一个常数，如 W_{ab} 表示点A与点B联系的紧密程度，称为权重。权重用W表示，权重越大表示两个节点越亲密，没有联系的节点如A和D之间的权重为0。 这张n个节点的图，对应于一个n*n的矩阵，表示出图中所有的权重，这个矩阵就是邻接矩阵（W)。邻接矩阵的构造法则是：n个节点从上至下顺序排列，称为矩阵的A行，B行..... n个节点从左至右顺序排列，称为矩阵的A列，B列......， 行和列的每一个交点处填写相关的权重，如 W_{ab} , W_{ae} ......

与邻接矩阵并列的是度矩阵(D)。’度‘是与节点相关的一个常数，它等于与一个节点相连的所有边的权重之和。度用d表示，在图1 中 d_{a} = W_{ab} + W_{ae} ， d_{b} = W_{ab} + W_{bc} + W_{bd} ......。 度矩阵是一个n*n 的对角阵，各对角元素表示对应的节点的度。 拉普拉斯矩阵(L) 定义为度矩阵与邻接矩阵之差，即 L=D-W。它的每一行（列）对应于一个节点，表示出这个节点与其余节点之间的权重，并写出这个节点的度。 拉普拉斯矩阵（L)是一个非常好的矩阵，它具有一些优良的性质。 （1）L矩阵的每一行（列）的所含元素之和为0，由此推断它的行列式为0，L矩阵不可逆。 （2）L矩阵是实对称矩阵，所以它有n个实数特征值，其中有一个特征值为0，其对应的特征向量为 [1,1,1......1 ]^{T} 。这一点是容易验证的。 关于L矩阵，还有一些性质我们将在下一节中继续讨论。这一节最后还有一个问题要说明，在研究‘图’时，往往定义一个节点函数 f_{(i)} ，如图1中的 f_{(a)} , f_{(b)} ...... 这样可以给各个节点赋值。节点函数视情况可以有不同的含义，例如在电路图中把它定义为节点电势，可以用一个向量 [ f_{(a)} , f_{(b)} ...... ]^{T} 来表示节点函数。

（二）电路图和基尔霍夫定律 上节提到的‘图’在很多领域都有应用，图内的权重，节点函数等概念在不同的场合会有不同的含义。例如在电路图中，我们可以把‘边’看作一条条‘支路’，权重则是每一条支路的‘电导’，即这条支路电阻的倒数，单位是（欧姆 )^{-1} , \Omega^{-1} 或者S。 在图1中选择两个节点，例如A和B，接入一个电源，这时各个节点上都会产生电势，这些电势的值就是节点函数 f_{(i)} 的函数值。电势的单位是伏特，V。

如果两个节点的电势不相等，存在电势差，它们之间的支路里会产生电流，例如与 W_{ab} 对应的电流 I_{ab} = W_{ab} ( f_{(a)} - f_{(b)} )。电流的单位是安培，A。 我们考察一下拉普拉斯矩阵（L)和节点函数的乘积，发现其结果表示每一个节点处进出电流的总和，其中外出电流为正，进入电流为负。可以用一个向量 [ I_{a} , I_{b} ...... ]^{T} 表示各个节点的进出电流。

下面再考察 L矩阵与节点函数的另一个乘积，这里按习惯把字母符号的下标写为1，2, .... n。

这个式子是线性代数里的二次型，由于 L矩阵的特性（各行（列）的所有元素之和均为0），所以我们得到了右面的结果，推导的过程就不在这里展示了。对电路图而言，这个结果表示各条支路里电流所作电功率的总和，它是大于或等于0的。 由此得出 L矩阵的第三条性质：（3）由于二次型 f^{T} L f = \frac{1}{2} \sum_{i,j=1}^{n}{} W_{ij} ( f_{i} - f_{j} )^{2} \geq 0，所以 L矩阵是半正定矩阵，半正定矩阵的 n个特征值均大于或等于0。 （4）拉普拉斯矩阵（L)的瑞利熵 R(L, f) = ( f^{T} L f) /( f^{T} f) 的最大值就是 L矩阵的最大特征值；其最小值就是 L矩阵的最小特征值。这里不作证明。 （三）应用拉普拉斯矩阵（L)求解简单电路问题

这里的简单电路指的是一个电源和多个电阻组成的电路。这种电路计算原理并不复杂，但当电阻数量很多时，计算很繁琐，容易出错。若应用 L矩阵，则计算过程简捷了许多。以图2a为例，它包括12个节点，14条支路，各支路的电导（电阻的倒数）已经确定。现在任意选择两个节点如 H和 G接入一个电源，见图2b，用 L矩阵可以一步计算出各个节点的电势，进而计算出各条支路上的电流。

解题的过程是这样的： 首先写出该电路图的 L矩阵。然后考虑各个节点的进出电流，电源正极接 H极，先设定 I_{h} 等于10A（这个数字在后面是可以调整的），则电源负极G点的 I_{g} 等于 -10A，其余节点的进出电流和均为0，这样得到了电流向量 [0,0,0,0,0,0,-10,10,0,0,0,0 ]^{T} 。需要计算的是各节点的电势向量 [ f_{a} , f_{b} , f_{c} , f_{d} , f_{e} , f_{f} ,0, f_{h} , f_{i} , f_{j} , f_{k} , f_{l} ]^{T} ，其中 f_{g} 应该为 0。列出的线性方程组如下。

求解这样的线性方程组需要使用计算机。我使用的是 python numpy 软件，在输入程序之前对矩阵作了两个简化。 （1）已知 f_{g} =0，需要把这个 0和与之相关的 L矩阵G列元素划除，这样12元线性方程组降为11元。 （2）观察这个11元线性方程组的增广矩阵，它的每一列元素之和均为 0，我们把上面11行都加到第12行上，第12行变为全0，随即消去。这样方程的个数也降为11。 输入计算机的程序如下。注意矩阵括号内的数字不要输入。

计算的结果为 [4.75, 2.27, 1.86, 0.62, 6.07, 5.8, 6.54, 5.1, 3.6, 2.29, 1.31 ]^{T} ，这11个数和 f_{g} =0 填入图2b对应的节点即为各点的电势值。各支路的电流随即可以算出。这些数据是相对的，全部的电势和电流可以乘以同一个比例进行伸缩。 （四）拉普拉斯矩阵（L)应用于谱聚类 上面讨论的‘图’也可以看作 n个节点的一个集合 V，有时候需要把 V分切成 k个子集，也就是在图中选择若干条边将之切断，这样的操作称为聚类。以三分聚类为例，即把 V分切成X,Y,Z三个子集。 X \bigcup_{}^{} Y \bigcup_{}^{} Z =V X \bigcap_{}^{} Y=X \bigcap_{}^{} Z=Y \bigcap_{}^{} Z=O 聚类的要求是不同子集间相连的边（即被切断边）之‘权重和’尽可能低，而子集内‘边权重和’尽可能高。 聚类的方法可以是先确定一个分切方案，通过计算，不断优化，最后得到满意的结果。以图2的三分聚类为例，先找出几个节点组成 X子集并把它分切出去，这时被切断边之‘权重和’用符号 W_{(X,cX)} 表示，其中 cX表示 X的补集，cX=Y \bigcup_{}^{} Z。我们的要求是 W_{(X,cX)} / \left| X \right| 尽量小，式中 \left| X \right| 表示 X 子集中的节点的个数。 同理对子集 Y有 W_{(Y,cY)} / \left| Y\right| ; 对子集 Z有 W_{(Z,cZ)} / \left| Z \right| 。综合起来用一个指标 RatioCut(X,Y,Z)= \frac{1}{2} ( W_{(X,cX)} / \left| X \right| + W_{(Y,cY)} / \left| Y \right| + W_{(Z,cZ)} / \left| Z\right| ) 表示聚类的成效，这个和数越小越好。 欲求得 W_{(X,cX)} ，对 电路图可采取施加电势的方法，即接入一个1V的电源，X子集中的节点都接到电源的正极，电势为1；X补集中的节点都接到电源的负极，电势为0.。这时只有连接X与cX的边中流有电流，有电流的边要被切断，它们的‘权重和’就是 W_{(XcX)} 。

这里我们用到 L矩阵的性质 3， f^{T} L f = \frac{1}{2} \sum_{i,j=1}^{}{} W_{i,j} ( f_{i} - f_{j} )^{2} 。观察图2c中的电势差只有1或者0，所以 f^{T} L f 等于被切断边之权重和。由此得到等式 W_{(X,cX)} = h_{ix}^{T} L h_{ix} 这里的 h_{ix} 就是 f_{i} ，不过换了个名称叫指示向量。 h_{ix} 是一个n维向量，它的每一个元素对应一个节点，当节点属于 X子集时， h_{ix} =1;当节点属于X补集时， h_{ix} =0。同理对 Y有 h_{iy} 对Z有 h_{iz} 。现在可以用 L矩阵与指示向量计算 W_{(X,cX)} , W_{(Y,cY)} , W_{(Z,cZ)} 。

看图2d，实际计算时我们一次确定一个聚类的初步方案 ，得到 X,Y,Z的指示向量 h_{ix} , h_{iy} , h_{iz} ，把它们排成一个3*n 的矩阵 H^{T} ，称之为指示矩阵。用指示矩阵作乘法 H^{T} L H更方便，输入计算机的结果如下：

上面得到的3*3 矩阵，其对角元素分别为 W_{(X,cX)} , W_{(Y,cY)} , W_{(Z,cZ)} ，非对角元素可以忽视。这样聚类指标 RatioCut 就可以计算出来了。实际问题中往往节点很多，不会象图2的例子这么简单，计算机要对大量的聚类方案进行运算，根据聚类指标最小化的原则优化结果。文献里有不少这方面的运算程序 ，这里就不作深入探讨了。 最后考虑一下拉普拉斯矩阵（L)的性质（4），发现一个巧合，即其瑞利熵就等于我们上面所要求的聚类指标，（ h_{ix}^{T} L h_{ix} )/( h_{ix}^{T} h_{ix} )= W_{(X,cX)} / \left| x \right| 。而 L的瑞利熵又与 L的特征值有关。这样我们得到了一个标准，即聚类结果的指标应该接近 L矩阵的较小的特征值。 在图2中，L矩阵的特征值由 numpy 算出： \lambda =26.5, 0, 0.82, 1.46, 3.4, 4.47, 14.24, 13.54, 12.84, 7.97, 8.39, 10.35.。图2d的 W_{(X,cX)} / \left| X \right| =2.75, W_{(Y,cY)} / \left| Y \right| =3 W_{(Z,cZ)} / \left| Z\right| =2.25。这个结果是令人满意的。