2. 贵州大学 计算机科学系, 贵州 贵阳 550025
2. Department of Computer Science, Guizhou University, Guiyang 550025, China
约束可满足问题(constraint satisfaction problem,简称CSP)是人工智能中的一个重要研究领域[1-3].一个CSP问题是由n个变量和m个约束组成.每个变量从其定义域D中取值.从变量集合中随机选取k个不同的变量构成一个约束,并且有一个相应的不协调集合Q来限制这k个变量的取值为不可满足的赋值[4].从可能的$\left( \begin{matrix} n \\ k \\ \end{matrix} \right)$个约束中随机选取m个构成CSP的一个随机实例.求解这个随机实例就是至少找到一个解,也就是这n个变量的一组赋值,使得这m个约束同时被满足;或者证明无解.可满足性判定问题(satisfiability,简称SAT)是典型的CSP问题,即给定一个合取范式公式(conjunctive normal form,简称CNF)F,SAT判定问题指是否存在一组指派使得F为真.一般情况下,随机CSP是NP完全的.在随机CSP问题中,约束个数m与变量个数n之比值a=m/n是一个重要参数,它不仅对实例的可满足性产生影响,还对求解实例的难易程度产生影响.
· 一方面,大量的理论分析和实验结果表明,NP完全问题存在可满足性相变现象,即存在一个阈值ad,几乎所有的随机实例在a<ad时是高概率可满足的(satisfiable,简称SAT);相反地,几乎所有的随机实例在a>ad时是高概率不可满足的(unsatisfiable,简称UNSAT)[5-8] .称阈值ad为相变点.
· 另一方面,大多数优秀的算法能够求解的实例与a密切相关.在相变点ad附近,实例的求解难度最大.如比较典型的随机3SAT问题,基于纯文字规则的算法只能有效求解a<1.6区域的随机3SAT实例(如DPLL算法),基于正文字或负文字出现次数最多而设计的算法可有效求解a<3.52区域的随机3SAT实例(如Zchaff算法),WALKSAT算法可有效求解a稍高于3.95区域的随机3SAT实例[9, 10].
最近几年,统计物理中的方法得到了广泛应用,基于该方法设计的求解可满足性问题的信息传播算法,在接近相变点时仍然有非常出色的表现.如,置信传播(belief propagation,简称BP)算法可有效求解a<3.95区域的随机3SAT实例;调查传播(survey propagation,简称SP)算法能够有效求解a<4.26区域的随机3SAT实例 [11-17],几乎能够有效求解接近相变点的实例.这弥补了以往其他随机算法在接近相变点时就会失效的缺憾.
信息传播算法不仅对固定定义域的随机CSP实例非常有效(如随机kSAT),还可用于求解具有增长定义域的随机CSP实例.许可等人在NP完全问题相变现象的研究中,提出了具有精确相变和难解实例的RB(k,n,a,rc,p)模型,该模型是一个典型的具有增长定义域的随机CSP实例产生模型,被广泛应用于构造难解性问题[18].随后, 文献[19-21]分别从理论和实验上说明在相变点的阈值附近,求解RB(k,n,a,rc,p)模型实例的难度随着问题规模的增加以指数形式增长.文献[22]表明,基于变量熵的BP算法在求解RB(k,n,a,rc,p)模型实例时仍然非常有效,几乎能够求解接近相变点的实例.
尽管信息传播算法求解CSP实例时非常有效,但仍有两个问题有待理论上深入研究:
1) 信息传播算法是否收敛.
2) 如果信息传播算法收敛,得到的结果是否是变量边缘概率分布的有效近似;如果信息传播算法不收敛,这种情况下算法失效.
BP算法之所以能够有效求解CSP实例,缘于该算法能够在CSP实例上正确收敛,进而得到变量的有效边缘概率分布,用以固定部分变量的赋值.然而在某些实例上,信息传播算法不是总有效,常表现为不收敛.我们已经知道:对于因子图为树结构和带有一个环的图结构的CSP实例,BP算法能够有效收敛[23, 24].但对于因子图带有多个环的CSP实例,BP算法不是总收敛.最近几年,一些学者初步研究了BP算法的收敛性,取得了一些结果.例如,文献[25]表明,对于任意结构的高斯图模型实例,BP算法能够正确收敛.Tatikonda等人分析了基于计算树的Gibbs测度序列的唯一性与BP算法的收敛性之间的关系,给出了BP算法收敛的有效条件[26].Heskes利用Bethe自由能量最小值的唯一性,给出了BP算法有唯一固定点的有效条件[27].Ihler直接从BP算法的信息更新函数出发,获得了传递错误信息的界,并利用该界给出了BP算法收敛的有效条件[28].Mooij等人利用商空间中信息更新函数的压缩映射原理导出了BP算法收敛的有效条件,指出,如果更新函数的系数矩阵的谱半径严格小于1,则BP算法收敛[29].Shi等人得到了传递错误信息的一个更紧致的界,利用这个界研究了BP算法的动态行为,并给出了算法收敛的有效条件[30].文献[31]中分析了求解最小代价网络流(minimum cost network flow,简称MCF)的BP算法的收敛性,研究结果表明,在亚多项式时间内,BP算法收敛到MCF问题的最优解,并且这个最优解具有唯一性.与此同时,文献中还提供了一个简化了的BP算法,能 够对MCF问题给出一个多项式时间随机近似方案.文献[32]中给出了BP算法收敛的一个概率条件,研究结果表明,BP算法能够给出最大权重匹配问题和最小代价网络流问题的最优解,且算法迭代次数高概率地有多项式界.文献[33]中分析了BP算法的复杂性,提出一个SBP(stochastic BP)算法,该算法的消息更新复杂性较BP算法的复杂性低了一个指数级别.
上述有关置信传播算法的收敛性证明方法大多都忽略了产生问题的模型,而是把注意力集中在具体问题本身,增加了这一问题的研究难度.事实上,基于实例的产生模型来研究置信传播算法的收敛性是非常有意义的.研究这一问题的优点至少包含以下3点:
(1) 为了研究一般性复杂问题,往往需要构造该问题的实例产生模型,而该模型一般为一个随机模型(如,CSP问题的实例产生模型B,D,RB,随机3-SAT实例产生模型等),利用随机模型的特征来分析问题所对应的因子图的某些结构性质的概率界.我们知道,置信传播算法是基于实例因子图的信息迭代算法,因子图的结构影响算法的收敛性.因此,利用因子图中某些结构性质的概率界来估算置信传播算法迭代步数的上界,以提高算法运行的效率,这是非常有意义的.
(2) 实例产生模型一般具有规则结构,规则图中的许多良好的性质对于研究实例因子图的结构性质具有潜在的作用.
(3) 基于这种规则结构的实例,置信传播算法的收敛性研究则更为具体和直观,有助于分析置信传播算法中的消息传递方式和消息迭代的稳定性.
RB模型是一个具有规则因子图结构的随机实例产生模型,目前得到了国内外学者的重视,是国际上应用最为广泛的难解实例产生模型,该模型的难解实例被应用于各种国际算法竞赛(如,ACM竞赛等)[34].基于此,我们研究基于RB模型的置信传播算法的收敛性.文献[22]中提出了一种基于变量熵的置信传播算法来求解RB模型实例,实验结果表明,该算法能够有效求解接近相变点的RB可满足实例.然而,未能从理论上解释这种现象.
RB模型是一个具有增长定义域的实例产生模型,已有的关于算法的收敛性分析都是基于具有固定定义域的情形.文献[29]中利用向量空间中函数的压缩映射原理给出了置信传播算法收敛的一般条件,但该判定条件的计算过程比较复杂,且适合于因子图结构较为简单的CSP实例.本文借鉴文献[29]中的证明方法,充分利用RB实例产生模型中不相容赋值的概率分布,分析基于变量熵的置信传播算法的收敛性,给出了该算法收敛的一个概率条件.具体来讲,在RB(k,n,a,rc,p)模型中,取k=2,$\alpha >\frac{1}{k}$,rc>0均为常数,且满足$k{{\text{e}}^{-\tfrac{\alpha }{{{r}_{c}}}}}\ge 1$,概率p控制每个约束的不协调赋值的个数.我们证明,如果p∈(0,n-2a),则基于变量熵的BP算法在RB(k,n,a,rc,p)模型产生的随机二元可满足实例集上高概率收敛.最后,在RB(k,n,a,rc,p)模型上选取了几组不同的数据进行数值模拟,实验表明该结论有效.更进一步地,当问题规模n增大时,在RB(k,n,a,rc,p)模型的可满足区域,实验收敛区间趋于一个固定范围,而理论收敛区间逐渐变窄.原因在于,RB(k,n,a,rc,p)模型是一个具有增长定义域的随机CSP实例产生模型,不协调赋值的数目与参数 p及问题规模n有关.
1 问题描述设有一个变量集合X={x1,…,xn},每个变量从其定义域D={1,…,dn}中取值.从变量集合中随机选取k(k≥2)个不同的变量$\{x_{i}^{1},...,x_{i}^{k}\}$构成约束Ci.对于每个约束Ci,都有一个相应的不协调赋值集合${{Q}_{{{C}_{i}}}}\subset {{D}^{k}}$来限制这k个变量的取值.若这k个变量的取值属于集合${{Q}_{{{C}_{i}}}},$则称此约束是不可满足的;反之,则称此约束是可满足的.由以下两步生成一个RB(k,n,a,rc,p)模型实例:
1) 从可能的$\left( \begin{matrix} n \\ k \\ \end{matrix} \right)$个约束中随机选取m=rcnlnn个约束构成集合$\{{{C}_{i}}\}_{i=1}^{m}$,其中,rc>0是常数;
2) 对于每个约束Ci,从可能的$d_{n}^{k}$个赋值中随机不重复地选取$|{{Q}_{{{C}_{i}}}}|=pd_{n}^{k}$个赋值作为不协调赋值集合,其中,0<p<1称为约束紧度,dn=na,a>0为常数.
这m个约束就组成RB(k,n,a,rc,p)模型的一个实例.如果存在这n个变量的一组赋值,使得所有的约束同时被满足,则这组赋值就是该实例的一个解.
在文献[18]中,用一阶矩和二阶矩方法证明了RB(k,n,a,rc,p)模型存在精确的可满足性相变现象.具体来讲,用Pr{SAT}表示RB(k,n,a,rc,p)模型实例是可满足的概率,则有如下结论成立.
定理1 [18]. 设${{p}_{c}}=1-{{\text{e}}^{-\tfrac{\alpha }{{{r}_{c}}}}}$.如果k,a,rc满足$k{{\text{e}}^{-\tfrac{\alpha }{{{r}_{c}}}}}\ge 1$,且$\alpha >\frac{1}{k}$,rc>0均为常数,则
$\underset{n\to \infty }{\mathop{\lim }}\,\Pr \{SAT\}=\left\{ \begin{array}{*{35}{l}} 0, & \text{if }p>{{p}_{c}} \\ 1, & \text{if }p<{{p}_{c}} \\ \end{array} \right.$ | (1) |
由此可见,当变量规模n充分大以后,RB(k,n,a,rc,p)模型的随机实例在pc处发生了可满足性的突变.也就是说:当p<pc时,RB(k,n,a,rc,p)模型产生的随机实例以接近概率1可满足;当p>pc时,RB(k,n,a,rc,p)模型产生的随机实例以接近概率1不可满足.当k≥2时,RB(k,n,a,rc,p)模型实例的判定问题是NP完全的.因此,为了理解方便,本文只考虑k=2的RB(k,n,a,rc,p)模型实例F.该方法容易推广到k>2的多元实例的情形.
二元(k=2) RB(k,n,a,rc,p)模型实例的因子图表示形式如图 1所示.因子图是一个二部图,图中的结点有两类:一类是变量结点(图中用圆环表示,记作i,j,k,…);另一类是约束结点(图中用方框表示,记作a,b,c,d,…).
在因子图的每条边(i,a)上,定义变量发送给约束的消息ui→a(si)以及约束发送给变量的消息ha→i(si).其中,ui→a(si)表示在没有约束a的情况下,变量i取值为si∈D的概率,ha→i(si)表示约束a要求变量i取值为si∈D的概率.V(i)表示包含变量i的约束结点集合,V(a)表示约束a中的变量结点集合.
求解二元RB(k,n,a,rc,p)模型实例的BP算法的消息迭代方程为
$u_{j\to a}^{(t)}({{s}_{j}})=\frac{1}{{{Z}^{j\to a}}}\prod\limits_{b\in V(j)\backslash a}{\eta _{b\to j}^{(t)}({{s}_{j}})}$ | (2) |
$\eta _{a\to i}^{(t+1)}({{s}_{i}})=\frac{1}{{{Z}^{a\to i}}}\sum\limits_{{{s}_{j}}\in D,j\in V(a)\backslash i}{{{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}u_{j\to a}^{(t)}({{s}_{j}})}$ | (3) |
其中,Zj→a和Za→i为归一化因子,${{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}=\left\{ \begin{array}{*{35}{l}} 0,\text{ if }({{s}_{i}},{{s}_{j}})\in {{Q}_{a}} \\ 1,\text{ o}\text{.w}\text{.}\ \\ \end{array} \right.$为标识函数.因此,概率为
$\Pr \{{{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}=0\}=p,\text{ }\Pr \{{{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}=1\}=1-p.$ |
将公式(2) 和公式(3) 合并,得到如下迭代方程:
$\eta _{a\to i}^{(t+1)}({{s}_{i}})=\frac{1}{{{Z}^{a\to i}}}\sum\limits_{{{s}_{j}}\in D,j\in V(a)\backslash i}{{{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}\prod\limits_{b\in V(j)\backslash a}{\eta _{b\to j}^{(t)}({{s}_{j}})}}$ | (4) |
当公式(4) 收敛时,可以根据固定点$\eta _{a\to i}^{*}({{s}_{i}})$得到变量i取值si的边缘概率分布p(si)的近似值bi(si).具体如下:
$p({{s}_{i}})\approx {{b}_{i}}({{s}_{i}})=\frac{1}{{{Z}^{i}}}\prod\limits_{a\in V(i)}{\eta _{a\to i}^{*}({{s}_{i}})}$ | (5) |
其中,${{Z}^{i}}=\sum\nolimits_{{{s}_{i}}\in D}{\prod\nolimits_{a\in V(i)}{\eta _{a\to i}^{*}({{s}_{i}})}}$为归一化因子.对于每个变量i,计算熵:
${{S}_{i}}=-\sum\nolimits_{{{s}_{i}}\in D}{{{b}_{i}}({{s}_{i}})\log ({{b}_{i}}({{s}_{i}}))}$ | (6) |
取出熵值Si最小的变量i,将变量i的取值固定到分量$s_{i}^{*}={{\max }_{{{s}_{i}}}}{{b}_{{{i}^{*}}}}({{s}_{i}})$上.
文献[22]中给出了求解RB(k,n,a,rc,p)模型实例的基于变量熵的BP算法,在该算法中,将变量分为3种类型.
· A类型:已经被固定赋值的变量.
· B类型:与A类型的变量出现在同一个约束中的变量.
· C类型:其他变量.
具体算法如下:
输入:二元RB模型的一个实例,最大迭代次数tmax和精度e.
输出:实例的解或者算法不收敛.
Step 1. t←0,对于因子图中的每条边(a,i),随机初始化信息$\eta _{a\to i}^{(0)}({{s}_{i}})\in [0,1].$
Step 2. t←t+1,根据公式(4) 计算$\eta _{a\to i}^{(t)}({{s}_{i}}).$
Step 3. 如果$|\eta _{a\to i}^{(t)}({{s}_{i}})-\eta _{a\to i}^{(t-1)}({{s}_{i}})|<\varepsilon $,则令$\eta _{a\to i}^{*}({{s}_{i}})=\eta _{a\to i}^{(t)}({{s}_{i}})$,并转到step 5;如果t>tmax,则输出不收敛.
Step 4. 返回到Step 2.
Step 5. 利用公式(5) 计算每个变量的边缘概率分布bi(si).
Step 6. 对于每个变量i,计算熵${{S}_{i}}=-\sum\limits_{{{s}_{i}}\in D}{{{b}_{i}}({{s}_{i}})\log ({{b}_{i}}({{s}_{i}}))},$挑选熵最小的变量${{i}^{*}}=\underset{i\in V}{\mathop{\arg }}\,\min {{S}_{i}},$将i的取值固定到分量$s_{i}^{*}=\underset{{{s}_{i}}\in D}{\mathop{\arg }}\,\max {{b}_{{{i}^{*}}}}({{s}_{i}}).$
Step 7. t←0,对于因子图中的每条边(a,i),按如下方式初始化信息$\eta _{a\to i}^{(0)}({{s}_{i}}):$
(1) 对于A类型的变量,不做处理.
(2) 对于B类型的变量,如果i的取值si与固定变量j的取值$s_{j}^{*}$有$({{s}_{i}},s_{j}^{*})\in {{Q}_{a}}$,则令$\eta _{a\to i}^{(0)}({{s}_{i}})=0$;否则,$\eta _{a\to i}^{(0)}({{s}_{i}})=1.$
(3) 对于C类型的变量,随机初始化信息$\eta _{a\to i}^{(0)}({{s}_{i}})\in [0,1].$
Step 8. t←t+1:
(1) 对于A类型的变量,不做处理.
(2) 对于B类型的变量,$\eta _{a\to i}^{(t)}({{s}_{i}})=\eta _{a\to i}^{(t-1)}({{s}_{i}}).$
(3) 对于C类型的变量,对于因子图中的每条边(a,i),根据公式(4) 计算$\eta _{a\to i}^{(t)}({{s}_{i}}).$
Step 9.
(1) 对于A类型的变量,不做处理.
(2) 对于B类型和C类型的变量,如果$|\eta _{a\to i}^{(t)}({{s}_{i}})-\eta _{a\to i}^{(t-1)}({{s}_{i}})|<\varepsilon ,$则令$\eta _{a\to i}^{*}({{s}_{i}})=\eta _{a\to i}^{(t)}({{s}_{i}}),$并转到Step 11;
如果t>tmax,则输出不收敛.
Step 10. 返回到Step 8.
Step 11. 对于自由变量(没有被固定的变量),利用公式(5) 计算每个自由变量的边缘概率分布bi(si).
Step 12. 对于每个自由变量i,计算熵${{S}_{i}}=-\sum\limits_{{{s}_{i}}\in D}{{{b}_{i}}({{s}_{i}})\log ({{b}_{i}}({{s}_{i}}))},$挑选熵最小的变量${{i}^{*}}=\underset{i\in V}{\mathop{\arg }}\,\min {{S}_{i}},$将i*的取值固定到分量$s_{i}^{*}=\underset{{{s}_{i}}\in D}{\mathop{\arg }}\,\max {{b}_{{{i}^{*}}}}({{s}_{i}}).$
Step 13. 返回到Step 7,直到所有变量都被赋值.
Step 14. 对于每个约束a,如果变量的取值都有$(s_{i}^{*},s_{j}^{*})\notin {{Q}_{a}},$则输出这组赋值作为解;否则,输出算法不收敛.
注意到:在算法的迭代过程中,消息更新对3种类型的变量分别做不同的处理.对于A类型变量的消息不做处理,也即,关于A类型的变量在算法迭代中消息保持不变;对于B类型变量的消息,在相继两次迭代中也保持不变;对于C类型变量的消息,是利用公式(4) 进行更新.从某种意义上说,在基于变量熵的BP算法的消息迭代过程中,A类型变量和B类型变量可能加速或延迟算法收敛的时间,但并不对算法是否收敛产生影响.因此,算法是否收敛完全取决于消息更新方程公式(4) .实验结果表明,该算法求解RB(k,n,a,rc,p)模型实例时非常有效,几乎能够求解接近相变点的实例[22].然而对于算法的性能,至今缺少系统的理论解释.BP算法是一类消息迭代算法,算法的收敛性直接影响算法的性能.对算法收敛性分析,有助于算法性能的进一步优化,对实际应用中信息传播算法的设计具有指导意义.本文分析了基于变量熵的BP算法求解RB(k,n,a,rc,p)模型实例的收敛性,给出了算法收敛的一个概率条件.
2 主要结论及证明基于变量熵的置信传播算法求解RB(k,n,a,rc,p)模型实例时不是总有效,常表现为不收敛.事实上,对于二元RB(k,n,a,rc,p)模型实例,有如下结论:
定理2. 在二元RB(k,n,a,rc,p)模型中,如果p∈(0,n-2a),则基于变量熵的BP算法在该模型产生的随机实例集上高概率收敛.
在二元RB(k,n,a,rc,p)模型中,每个变量的取值为${{s}_{i}}\in D=\{{{s}_{1}},{{s}_{2}},...,{{s}_{{{d}_{n}}}}\}.$相应地,BP算法中的每个消息可表示为向量形式$\eta _{a\to i}^{(t)}=(\eta _{a\to i}^{(t)}({{s}_{1}}),\eta _{a\to i}^{(t)}({{s}_{2}}),...,\eta _{a\to i}^{(t)}({{s}_{{{d}_{n}}}})).$
为了方便处理,我们对消息更新方程公式(4) 取对数,即$\lambda _{a\to i}^{(t+1)}=\log \eta _{a\to i}^{(t+1)},$则公式(4) 表示为
$\lambda _{a\to i}^{(t+1)}({{s}_{i}})=\log \left( \sum\limits_{{{s}_{j}}\in D,j\in V(a)\backslash i}{{{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}{{h}^{a\backslash i}}({{s}_{j}})} \right)$ | (7) |
${{h}^{a\backslash i}}({{s}_{j}})=\exp \left( \sum\limits_{b\in V(j)\backslash a}{\lambda _{b\to j}^{(t)}({{s}_{j}})} \right)$ | (8) |
每个消息$\lambda _{a\to i}^{(t+1)}$是Va→i=$\mathbb{R}$|D|中的向量.对于s∈D,用$\lambda _{a\to i}^{(t+1)}(s)$表示$\lambda _{a\to i}^{(t+1)}$的某个分量的值.定义全局向量空间$V=\underset{i\in a\in F}{\mathop{\oplus }}\,\,\,{{V}^{a\to i}},$其中,⊕表示向量的方向和.定义函数f:V→V,假设算法中的消息更新并行执行.自然,可以用函数f的迭代过程来刻画算法中的消息更新过程.其中,消息更新方程公式(7) 可表示为函数f的分量.
设(V,||×||)是一个有限尺度的实向量赋范空间,由范数||×||诱导出V上的一个距离:d(v,w)=||v-w||.若存在0≤k<1,对于任意x,y∈V,有d(f(x),f(y))≤kd(x,y)成立,则称映射f:V→V是一个关于d的压缩映射.由代数知识,如果f:V→V是一个关于完备度量空间(V,d)的压缩映射,那么f有唯一固定点x∞∈V,且与初始信息无关.也即,对于任意x∈V,有序列x,f(x),f2(x),…,x∞,称f收敛到x∞.
引理 1. 设函数f:V→V,d是V上的度量.对于某个N∈∞,假设fN是关于d的压缩函数,那么f收敛到唯一固定点x∞,且与初始信息无关.
证明:对于任意的x∈V,考虑如下N-1个序列:
$\begin{align} & x,{{f}^{N}}(x),{{f}^{2N}}(x),... \\ & f(x),{{f}^{N+1}}(x),{{f}^{2N+1}}(x),... \\ & \ \ \vdots \\ & {{f}^{N-1}}(x),{{f}^{2N-1}}(x),{{f}^{3N-1}}(x),... \\ \end{align}$ |
因为fN是关于d的压缩函数,所以上述每一个序列都收敛到x∞.因此,序列x,f(x),f2(x),…一定收敛到x∞,所以函数f收敛到唯一固定点x∞,且与初始信息无关.
根据中值定理,如果函数f:V→V可微,那么对于任意x,y∈V,有公式(9) 成立:
$||f(y)-f(x)||\ \ \le \ \ ||y-x||\cdot \underset{z\in [x,y]}{\mathop{\sup }}\,||{f}'(z)||$ | (9) |
其中,[x,y]={cx+(1-c)y:c∈[0,1]}.
因此,如果$\underset{v\in V}{\mathop{\sup }}\,||{f}'(v)||<1,$那么f是压缩函数,收敛于固定点x∞∈V,且与初始信息无关.
更进一步地,如果$\underset{v\in V}{\mathop{\sup }}\,||{f}'(v)||<1,$则BP算法收敛,且与初始消息无关.因此,关键问题是构造范数||f'(v)||.
赋范空间(V,||×||)可诱导一个关于映射A:V→V的矩阵范数:$||A||=\underset{\begin{smallmatrix} \,v\in V \\ ||v||\le 1 \end{smallmatrix}}{\mathop{\sup }}\,||Av||$.
l1-范数定义为$||A|{{|}_{1}}=\underset{j\in \{1,2,...,N\}}{\mathop{\max }}\,\sum\limits_{i=1}^{N}{|{{A}_{ij}}|}$,l∞-范数定义为$||A|{{|}_{\infty }}=\underset{i\in \{1,2,...,N\}}{\mathop{\max }}\,\sum\limits_{j=1}^{N}{|{{A}_{ij}}|}.$
对于子向量空间Va→i,设其上的范数表示为||×||a→i.向量空间V上的范数定义为
$||v||=\sum\limits_{a\to i}{||{{v}^{a\to i}}|{{|}_{a\to i}}}$ | (10) |
设l∈V,f'(l):V→V为f:V→V在l处的导函数.
在l1-范数下,则$||{f}'(\lambda )|{{|}_{1}}\ =\underset{b\to j}{\mathop{\max }}\,\sum\limits_{a\to i}{||{f}'{{(\lambda )}_{a\to i,b\to j}}||_{a\to i}^{b\to j}}$,而$||{f}'{{(\lambda )}_{a\to i,b\to j}}||_{a\to i}^{b\to j}=\underset{\begin{smallmatrix} v\in {{V}^{b\to j}} \\ ||v|{{|}_{b\to j}}\le 1 \end{smallmatrix}}{\mathop{\sup }}\,||{f}'{{(\lambda )}_{a\to i,b\to j}}v|{{|}_{a\to i}}.$
对于二元RB(k,n,a,rc,p)模型实例的每一个约束a∈F,假设$\forall i\in a,\forall {{s}_{i}}\in D,\exists {{s}_{j}}\in D:\ {{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}>0,$那么方程公式(7) 的导数为
$\frac{\partial \lambda _{a\to i}^{(t+1)}({{s}_{i}})}{\partial \lambda _{b\to i}^{(t)}({{s}_{j}})}={{1}_{V(j)\backslash a}}(b){{1}_{V(a)\backslash i}}(j)\times \frac{{{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}{{h}^{a\backslash i}}({{s}_{j}})}{\sum\nolimits_{{{{{s}'}}_{j}}\in D}{{{\delta }_{({{s}_{i}},{{{{s}'}}_{j}})\in {{Q}_{a}}}}{{h}^{a\backslash i}}({{{{s}'}}_{j}})}}$ | (11) |
其中,标识函数${{1}_{X}}(x)=\left\{ \begin{array}{*{35}{l}} 1, & \text{if }x\in X \\ 0, & \text{o}\text{.w}\text{.} \\ \end{array} \right..$因此,
${{\left\| {f}'(\lambda ) \right\|}_{1}}=\underset{b\to j}{\mathop{\max }}\,\sum\limits_{a\to i}{{{1}_{V(j)\backslash a}}(b){{1}_{V(a)\backslash i}}(j){{B}_{a\to i,b\to j}}({{h}^{a\backslash i}}(\lambda ))}$ | (12) |
其中,${{B}_{a\to i,b\to j}}({{h}^{a\backslash i}}(\lambda ))=\left\| \frac{{{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}{{h}^{a\backslash i}}(\lambda )}{\sum\nolimits_{{{{{s}'}}_{j}}\in D}{{{\delta }_{({{s}_{i}},{{{{s}'}}_{j}})\in {{Q}_{a}}}}{{h}^{a\backslash i}}(\lambda )}} \right\|_{a\to i}^{b\to j}.$因此,只要满足$\underset{\lambda \in V}{\mathop{\sup }}\,||{f}'(\lambda )|{{|}_{1}}<1,$BP算法就收敛.
定义矩阵A中的元素为
${{A}_{a\to i,b\to j}}=\underset{{{h}^{a\backslash i}}>0}{\mathop{\sup }}\,{{B}_{a\to i,b\to j}}({{h}^{a\backslash i}}(\lambda ))$ | (13) |
假设||×||a→i是Va→i=$\mathbb{R}$|D|上的l∞范数,由代数知识,Aa→i,b→j可表示为如下形式[29]:
${{A}_{a\to i,b\to j}}=N(\delta ,i,j)=\underset{\begin{smallmatrix} \ \ {{s}_{i}}\ne {{{{s}'}}_{i}} \\ \ \ {{s}_{i}},{{{{s}'}}_{i}}\in D \end{smallmatrix}}{\mathop{\sup }}\,\underset{\begin{smallmatrix} \ {{s}_{j}}\ne {{{{s}'}}_{j}} \\ \ {{s}_{j}},{{{{s}'}}_{j}}\in D \end{smallmatrix}}{\mathop{\sup }}\,\frac{\sqrt{{{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}{{\delta }_{({{{{s}'}}_{i}},{{{{s}'}}_{j}})\in {{Q}_{a}}}}}-\sqrt{{{\delta }_{({{{{s}'}}_{i}},{{s}_{j}})\in {{Q}_{a}}}}{{\delta }_{({{s}_{i}},{{{{s}'}}_{j}})\in {{Q}_{a}}}}}}{\sqrt{{{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}{{\delta }_{({{{{s}'}}_{i}},{{{{s}'}}_{j}})\in {{Q}_{a}}}}}+\sqrt{{{\delta }_{({{{{s}'}}_{i}},{{s}_{j}})\in {{Q}_{a}}}}{{\delta }_{({{s}_{i}},{{{{s}'}}_{j}})\in {{Q}_{a}}}}}}$ | (14) |
根据函数的压缩映射原理及公式(9) 、公式(12) ~公式(14) ,有如下推论 :
推论 1. 设$\forall i\in a,\forall {{s}_{i}}\in D,\exists {{s}_{j}}\in D:\ {{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}>0.$若$\underset{b\to j}{\mathop{\max }}\,\sum\limits_{a\to i}{{{1}_{V(j)\backslash a}}(b){{1}_{V(a)\backslash i}}(j)N(\delta ,i,j)}<1,$则BP算法收敛到唯一固定点,且与初始信息无关.
定义矩阵M的元素为Ma→i,b→j=1V(j)\a(b)1V(a)\i(j)N(d,i,j),l1,…,ln是矩阵M的特征值,l(M)={l1,…,ln}为矩阵M的谱,r(M)=sup{|li|:li∈l(M)}为矩阵M的谱半径.
定理3. 设函数f:V→V可微,进一步假设|f'(x)|≤M,M是一个不依赖于x的非负矩阵.如果r(M)<1,那么对于任意x∈V',f收敛到唯一固定点x∞.也即,x,f(x),f2(x),…,x∞(注意,|f'(x)|≤M表示由f'(x)构成的雅可比矩阵中的每个元素的绝对值不大于M中对应元素的值).
证明:对于任意的x∈V',假设n=1,2,…,我们有:
$||({{f}^{n}}{)}'(x)||=\left\| \prod\limits_{i=1}^{n}{{f}'({{f}^{i-1}}(x))} \right\|={{\left\| \left| \prod\limits_{i=1}^{n}{{f}'({{f}^{i-1}}(x))} \right| \right\|}_{1}}\le {{\left\| \prod\limits_{i=1}^{n}{|{f}'({{f}^{i-1}}(x))}| \right\|}_{1}}\le {{\left\| \prod\limits_{i=1}^{n}{M} \right\|}_{1}}=||{{M}^{n}}|{{|}_{1}}$ | (15) |
由Gelfand引理 [35],$\underset{n\to \infty }{\mathop{\lim }}\,||{{M}^{n}}|{{|}_{1}}^{1/n}=\rho (M)$.我们选取e>0,使得r(M)+e<1,对于某个N∈∞,||MN||1=(r(M)+e)N<1.因此,||(fN)'(x)||<1.所以,||(fN)(x)||是压缩函数.
由引理 1可知,函数f(x)收敛到唯一固定点x∞,且与初始信息无关.
因此,根据定理3及公式(11) ,有如下结论:
推论 2. 设$\forall i\in a,\forall {{s}_{i}}\in D,\exists {{s}_{j}}\in D:\ {{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{a}}}}>0.$如果矩阵M的谱半径r(M)严格小于1,则BP算法收敛到唯一固定点,且与初始信息无关.
对于一个具体的二元RB(k,n,a,rc,p)模型的随机实例F,依据RB(k,n,a,rc,p)模型实例的生成方法,每个约束I∈F有n2ap个不协调赋值和n2a(1-p)个协调赋值.也即,有n2ap个${{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{I}}}}=0$的赋值和 n2a(1-p)个${{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{I}}}}=1$的赋值.下面给出$\forall i\in I,\forall {{s}_{i}}\in D,\exists {{s}_{j}}\in D:\ {{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{I}}}}>0$高概率成立的p的取值界,有如下引理 :
引理 2. 设$p=1-\frac{\text{e}}{{{n}^{c}}}$,c为常数.则c=a是$\forall i\in I,\forall {{s}_{i}}\in D,\exists {{s}_{j}}\in D:\ {{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{I}}}}>0$高概率成立的临界点.
证明:设约束I的一个赋值为(si,sj)∈DxD,所有的赋值构成赋值矩阵S,其元素为sij=(si,sj).位于矩阵S的不同行、不同列元素构成一个排列${{s}_{1{{i}_{1}}}},{{s}_{2{{i}_{2}}}},...,{{s}_{{{d}_{n}}{{i}_{{{d}_{n}}}}}},$这样的排列共有dn!个.如果在n2a(1-p)个协调赋值中至少含有一个排列${{s}_{1{{i}_{1}}}},{{s}_{2{{i}_{2}}}},...,{{s}_{{{d}_{n}}{{i}_{{{d}_{n}}}}}}$中的元素,则都有$\forall i\in I,\forall {{s}_{i}}\in D,\exists {{s}_{j}}\in D:\ {{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{I}}}}>0$成立.设随机变量x表示在n2a(1-p)个协调赋值中含有元素${{s}_{1{{i}_{1}}}},{{s}_{2{{i}_{2}}}},...,{{s}_{{{d}_{n}}{{i}_{{{d}_{n}}}}}}$的事件发生的数目,则
$\begin{align} & \underset{n\to \infty }{\mathop{\lim }}\,E(x)=\underset{n\to \infty }{\mathop{\lim }}\,\frac{{{n}^{\alpha }}!\cdot \left( \begin{matrix} {{n}^{2\alpha }}-{{n}^{\alpha }} \\ (1-p){{n}^{2\alpha }}-{{n}^{\alpha }} \\ \end{matrix} \right)}{\left( \begin{matrix} {{n}^{2\alpha }} \\ (1-p){{n}^{2\alpha }} \\ \end{matrix} \right)} \\ & \quad \quad \,\quad \ \,=\underset{n\to \infty }{\mathop{\lim }}\,\frac{{{n}^{\alpha }}!\cdot ((1-p){{n}^{2\alpha }})((1-p){{n}^{2\alpha }}-1)...((1-p){{n}^{2\alpha }}-{{n}^{\alpha }}+1)}{{{n}^{2\alpha }}({{n}^{2\alpha }}-1)...({{n}^{2\alpha }}-{{n}^{\alpha }}+1)} \\ & \quad \quad \,\quad \ \,=\underset{n\to \infty }{\mathop{\lim }}\,\frac{{{n}^{\alpha }}!\cdot ((1-p))\left( (1-p)-\frac{1}{{{n}^{2\alpha }}} \right)...\left( (1-p)-\frac{{{n}^{\alpha }}-1}{{{n}^{2\alpha }}} \right)}{\left( 1-\frac{1}{{{n}^{2\alpha }}} \right)\left( 1-\frac{2}{{{n}^{2\alpha }}} \right)...\left( 1-\frac{{{n}^{\alpha }}-1}{{{n}^{2\alpha }}} \right)} \\ & \quad \quad \,\quad \,\,=\underset{n\to \infty }{\mathop{\lim }}\,{{n}^{\alpha }}!\cdot {{(1-p)}^{{{n}^{\alpha }}}}\quad \\ & \quad \quad \,\quad \,\,\approx \underset{n\to \infty }{\mathop{\lim }}\,\frac{{{n}^{\alpha \left( {{n}^{\alpha }}+\frac{1}{2} \right)}}}{{{\text{e}}^{{{n}^{\alpha }}}}}\cdot {{(1-p)}^{{{n}^{\alpha }}+1}} \\ & \quad \quad \,\quad \,\,\le \underset{n\to \infty }{\mathop{\lim }}\,\frac{{{n}^{\alpha ({{n}^{\alpha }}+1)}}}{{{\text{e}}^{{{n}^{\alpha }}}}}\cdot {{(1-p)}^{{{n}^{\alpha }}+1}}\text{ (}\alpha >\frac{1}{2}\text{)} \\ \end{align}$ | (16) |
因此,当c>a时,$\underset{n\to \infty }{\mathop{\lim }}\,E(x)=0$.由一阶矩方法可知,Pr(x≥1) =0.
不妨设v1=(si,sj),v2=(sk,sl)是约束I的两个赋值,Pr(<v1,v2>)表示赋值v1和v2同时满足约束I的概率.
· 如果赋值v1和v2相同,则$\Pr (\langle {{v}_{1}},{{v}_{2}}\rangle )={\left( \begin{matrix} d_{n}^{k}-1 \\ |{{Q}_{I}}| \\ \end{matrix} \right)}/{\left( \begin{matrix} d_{n}^{k} \\ |{{Q}_{I}}| \\ \end{matrix} \right)}\;=1-p;$
· 如果赋值v1和v2不同,则$\Pr (\langle {{v}_{1}},{{v}_{2}}\rangle )={\left( \begin{matrix} d_{n}^{k}-2 \\ |{{Q}_{I}}| \\ \end{matrix} \right)}/{\left( \begin{matrix} d_{n}^{k} \\ |{{Q}_{I}}| \\ \end{matrix} \right)}\;={{(1-p)}^{2}}+O\left( \frac{1}{{{n}^{2\alpha }}} \right).$
设xi,xj是n2a(1-p)个协调赋值中分别含有元素${{s}_{1{{i}_{1}}}},{{s}_{2{{i}_{2}}}},...,{{s}_{{{d}_{n}}{{i}_{{{d}_{n}}}}}}$和${{s}_{1{{j}_{1}}}},{{s}_{2{{j}_{2}}}},...,{{s}_{{{d}_{n}}{{j}_{{{d}_{n}}}}}}$的两组赋值,所以,
$\begin{align} & E({{x}^{2}})=\sum\limits_{i=j}{E({{x}_{i}},{{x}_{j}})}+\sum\limits_{i\ne j}{E({{x}_{i}},{{x}_{j}})} \\ & \text{ }=E(x)+({{({{n}^{\alpha }}!)}^{2}}-{{n}^{\alpha }}!){{\left[ {{(1-p)}^{2}}+O\left( \frac{1}{{{n}^{2\alpha }}} \right) \right]}^{{{n}^{\alpha }}}} \\ & \text{ }=E(x)+({{({{n}^{\alpha }}!)}^{2}}-{{n}^{\alpha }}!){{(1-p)}^{2{{n}^{\alpha }}}}+O\left( \frac{1}{{{n}^{2\alpha }}} \right) \\ \end{align}$ | (17) |
进而有:
$\underset{n\to \infty }{\mathop{\lim }}\,\frac{E({{x}^{2}})}{{{E}^{2}}(x)}=\underset{n\to \infty }{\mathop{\lim }}\,\frac{E(x)+({{({{n}^{\alpha }}!)}^{2}}-{{n}^{\alpha }}!){{(1-p)}^{2{{n}^{\alpha }}}}}{{{E}^{2}}(x)}=\underset{n\to \infty }{\mathop{\lim }}\,\left( \frac{1}{E(x)}+1-\frac{1}{{{n}^{\alpha }}!} \right)$ | (18) |
当c<a时,$\underset{n\to \infty }{\mathop{\lim }}\,E\left( x \right)=\infty $.因此,
$\underset{n\to \infty }{\mathop{\lim }}\,\frac{E({{x}^{2}})}{{{E}^{2}}(x)}\le \underset{n\to \infty }{\mathop{\lim }}\,\left[ \frac{1}{E(x)}+1-\frac{1}{{{n}^{\alpha }}!} \right]=1+o(1)$ | (19) |
由二阶矩方法,Pr(x=0) =0.
所以,c=a是$\forall i\in I,\forall {{s}_{i}}\in D,\exists {{s}_{j}}\in D:\ {{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{I}}}}>0$高概率成立的临界点.
引理 2表明,当p<1-en-a时,则$\forall i\in I,\forall {{s}_{i}}\in D,\exists {{s}_{j}}\in D:\ {{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{I}}}}>0$高概率成立.注意到,在RB(k,n,a,rc,p)模型中的可满足相变点为${{p}_{c}}=1-{{\text{e}}^{-\tfrac{\alpha }{{{r}_{c}}}}},$而$1-{{\text{e}}^{-\tfrac{\alpha }{{{r}_{c}}}}}<1-e{{n}^{-\alpha }}.$根据定理1和引理 2,在RB(k,n,a,rc,p)模型的可满足区域,总有$\forall i\in I,\forall {{s}_{i}}\in D,\exists {{s}_{j}}\in D:\ {{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{I}}}}>0$高概率成立.因此,推论 2可用作为判定BP算法在RB(k,n,a,rc,p)模型的可满足实例上是否高概率收敛的充分条件.我们有如下推论 :
推论 3. 设${{p}_{c}}=1-{{\text{e}}^{-\tfrac{\alpha }{{{r}_{c}}}}},$其中,$\alpha >\frac{1}{k}$,rc>0均为常数,且满足$k{{\text{e}}^{-\tfrac{\alpha }{{{r}_{c}}}}}\ge 1$.当p<pc时,对于二元RB(k,n,a,rc,p)模型
实例F,如果矩阵M的谱半径r(M)严格小于1,则BP算法高概率收敛,且与初始信息无关(注意,这里的高概率指引理 2中的高概率).
注意到,矩阵M中的元素由N(d,i,j)决定,为此给出N(d,i,j)=0高概率成立的条件.
引理 3. 设p=n-c,c为常数.当c>2a时,N(d,i,j)=0高概率成立.
证明:设x表示约束I不协调赋值的数目的随机变量.由马尔可夫不等式得:Pr{x≥1}≤E[x].而
$E\left[ x \right]={{n}^{2a}}\times p={{n}^{2a}}\times {{n}^{-c}}={{n}^{2a-c}}$ | (20) |
当c>2a时,$\underset{n\to \infty }{\mathop{\lim }}\,E(x)=0$.所以,Pr(x≥1) =0,进一步有Pr(x<1) =1.因此,根据公式(14) 有Pr{N(d,i,j)=0}=1.
在矩阵分析中,我们知道:关于非负方阵M的谱半径r(M),如果存在一个矩阵范数||×||使得||M||<1,那么$\underset{n\to \infty }{\mathop{\lim }}\,{{M}^{n}}=0;$如果$\underset{n\to \infty }{\mathop{\lim }}\,{{M}^{n}}=0$当且仅当r(M)<1.下面给出BP算法在二元RB(k,n,a,rc,p)模型上高概率收敛的条件.
定理 4. 设p=n-c,c为常数.当c>2a时,BP算法在二元RB(k,n,a,rc,p)模型实例集上高概率收敛.
证明:由引理 2可知:当p1<1-en-a时,高概率地有$\forall i\in I,\forall {{s}_{i}}\in D,\exists {{s}_{j}}\in D:\ {{\delta }_{({{s}_{i}},{{s}_{j}})\in {{Q}_{I}}}}>0$成立.而引理 3中,当p2< n-2a时,高概率地有$||M|{{|}_{\infty }}=\underset{i\in \{1,2,...,N\}}{\mathop{\max }}\,\sum\limits_{j=1}^{N}{|{{M}_{ij}}|}<1.$所以,高概率地有r(M)<1成立.容易看出,p2≤p1.因此,当p<n-2a时,BP算法在二元RB(k,n,a,rc,p)模型实例集上高概率收敛.至此,定理2得证.
定理2表明,如果p∈(0,n-2a),则基于变量熵的BP算法在二元RB(k,n,a,rc,p)模型产生的随机实例集上高概率收敛.
上述有关证明过程及结论容易扩展到多元RB模型实例上.
3 数值实验及分析为了对本文定理2中的结论作出实验验证,数值实验中选用RB模型实例.在RB模型中,我们选择了两组不同的参数,分别为a=0.8,rc=3和a=0.6,rc=1,用来生成随机CSP实例.对于每组参数取5个不同规模n的CSP随机实例,见表 1和表 2(表 1中,取a=0.8,rc=3,对于不同的问题规模n,相应的dn和m的值;表 2中,取a=0.6,rc=1,对于不同的问题规模n,相应的dn和m的值).在算法中,取tmax=103,e=10-4.实验运行的软、硬件环境为:Intel Core i7-2600 CPU,主频3.40GHZ,内存4.00GB,操作系统WIN7,采用VC++6.0集成开发工具,使用C和MATLAB语言编程.
在RB(k,n,a,rc,p)模型中,可满足相变点${{p}_{c}}=1-{{\text{e}}^{-\tfrac{\alpha }{{{r}_{c}}}}}.$对于第1组数据,pc=0.234.对于第2组数据,pc=0.451.数值实验中,每个数据点由100个RB(k,n,a,rc,p)模型的随机实例构成.第1组数据和第2组数据的运行结果分别如图 2和图 3所示.
图 2(a)是基于变量熵的BP算法在第1组数据上成功收敛的概率随着参数p变化的情形.当p<0.145时,算法高概率收敛.当p>0.145时,算法从收敛突变到几乎不收敛.并且在p>0.195时,算法成功收敛的概率逐渐增加.我 们发现:在0.160<p<0.195区间内,算法几乎高概率不收敛,即算法失效.事实上,此区间恰好为RB(k,n,a,rc,p)模型的难解区间,实例的求解难度随着变量规模的增加成指数形式增长.
图 2(b)是基于变量熵的BP算法在第1组数据上能够有效求解实例的概率随着参数p变化的情形.当p<0.145时,算法几乎能够高概率地给出实例的解.当p>0.145时,算法不能给出实例的解:一方面,算法本身不收敛;另一方面,算法收敛但不能正确给出实例的解.
对比图 2(a)与图 2(b),我们发现,当p>0.195时,对于部分实例该算法收敛,但未能正确给出实例的解.
对于第2组数据,基于变量熵的BP算法的运行结果如图 3所示.当p<0.245时,算法高概率地收敛.当p>0.245时,算法从收敛突变到几乎不收敛.并且在p>0.420时,算法成功收敛的概率逐渐增加.在0.325<p<0.420 区间内,算法几乎高概率不收敛,即算法失效,如图 3(a)所示.
图 3(b)是算法在第2组数据上能够解的运行结果.当p<0.245时,算法几乎能够高概率地给出实例的解.当p>0.245时,算法不能给出实例的解.我们发现,当p>0.420时,算法在部分实例收敛,但未能给出实例的解.
为了进一步分析和验证本文结论,并与其他结论进行比较,我们选取RB模型中的参数见表 3.对于表 3中的第1组数据(n,a,dn,rc,m)=(20,0.4,3,4,240) ,根据定理2可计算出算法收敛的条件p的取值范围为(0,0.091) .同理,对于第2组数据(n,a,dn,rc,m)=(40,0.4,4,7,1033) 和第3组数据(n,a,dn,rc,m)=(60,0.4,5,10,2457) ,p的取值范围分别为(0,0.052) 和(0,0.038) ,这3组数据对应实例的相变点pc分别为0.095,0.056,0.039,具体参数值见表 4.不难看出,在定理2的收敛性判定条件中,p的取值上界几乎接近相变点pc.进一步说明,对于接近相变点的复杂因子图结构的二元RB模型实例,定理2的判定条件几乎有效.针对表 3中的第1组数据,选取p为0.078,由RB(k,n,a,rc,p)模型产生1 000个二元CSP实例进行数值实验,本文能够有效判定的概率为0.87,文献[29]中能够有效判定的概率为 0.91;同样,对于表 3中的第2组数据和第3组数据,本文能够有效判定的概率为1.00和0.93,而文献[29]中分别为0.84和0.79.表 4中给出了本文结论与文献[29]中结论的对比结果(每组参数生成1 000个RB模型实例).
当问题规模n增大时,理论收敛区间逐渐减小.原因在于,RB(k,n,a,rc,p)模型是一个具有增长定义域的随机CSP实例产生模型,不协调赋值的数目与参数p及问题规模n有关.
用难解实例测试算法性能及相关结论是非常必要的,我们使用国际上常用的由Lecoutre提供的benchmark案例进行测试[34].数值实验中选取了两组实例集frb30-15和frb40-19进行测试,结果表明,基于变量熵的置信传播算法在这些实例上基本不收敛.事实上,在这些实例集上,定理2中的条件已经不满足.也即,定理2中判定算法收敛的条件失效.当然,不能仅凭某些局部的实验数据就断定某种方法的优劣,每种方法都有其优越性和局 限性.
4 结束语本文分析了基于变量熵的置信传播算法在二元RB(k,n,a,rc,p)模型随机实例集上的收敛性,利用置信传播算法的信息更新函数在实向量空间中的压缩性原理,及更新函数的系数矩阵的谱半径等理论,证明了当p∈(0,n-2a)时,该算法在二元RB(k,n,a,rc,p)模型产生的随机CSP实例集上高概率收敛,给出了置信传播算法能够有效求解二元RB(k,n,a,rc,p)模型实例的理论解释.我们希望这些方法的引入可以提供一种新的思路,用来分析其他信息传播算法的收敛性及性能.进一步工作是:(1) 扩展基于变量熵的BP算法在RB(k,n,a,rc,p)模型实例集上的收敛区间;(2) 分析基于变量熵的BP算法在RB(k,n,a,rc,p)模型实例集上的有效性,也即,算法收敛后能否给出问题的有效解.
[1] | Martin D, Alan F, Michael M. A probabilistic analysis of randomly generated binary constraint satisfaction. Theory Computer Science, 2003, 290 (3) :1815–1828. [doi:10.1016/S0304-3975(02)00317-1] |
[2] | Creignou N, Eaude H. The SAT-UNSAT transition for random constraint satisfaction problems. Discrete Mathematics, 2009, 309 (8) :2085–2099. [doi:10.1016/j.disc.2008.04.025] |
[3] | Martin OC, Monasson R, Zecchina R. Statistical mechanics methods and phase transitions in optimization. Theory Computer Science, 2001, 265 (1-2) :3–67. [doi:10.1016/S0304-3975(01)00149-9] |
[4] | Tsang E. A glimpse of constraint satisfaction. Artificial Intelligence Review, 1999, 13 (3) :215–277. [doi:10.1023/A:1006558104682] |
[5] | Xu K, Li W. The SAT phase transition. Science in China, Series E, 1999, 42 (5) :494–501. [doi:10.1007/BF02917402] |
[6] | Fan Y, Shen J. On the phase transitions of random k-constraint satisfaction problems. Artificial Intelligence, 2011, 175 (3-4) :914–927. [doi:10.1016/j.artint.2010.11.004] |
[7] | Achlioptas D, Naro A, Peres Y. Rigorous location of phase transitions in hard optimization problems. Nature, 2005, 435 (7043) :759–764. [doi:10.1038/nature03602] |
[8] | Mertens S, Mèzard M, Zecchina R. Threshold values of random k-SAT from the cavity method. Random Structure Algorithms, 2006, 28 (3) :340–373. [doi:10.1002/rsa.20090] |
[9] | Moskewicz MW, Madigan CF, Zhao Y. Chaff:Engineering an efficient SAT solver. In:Proc. of the 38th Annual Design Automation Conf. New York:ACM Press, 2001. 530-535.[doi:10.1145/378239.379017] |
[10] | Aurell E, Gordon U, Kirkkpatrick S. Comparing beliefs, surveys, and random walks. Advances in Neural Information Processing Systems, 2004, 17 (1) :1–8. http://cn.bing.com/academic/profile?id=2139955802&encoded=0&v=paper_preview&mkt=zh-cn |
[11] | Mezard M, Parisi G. The cavity method at zero temperature. Journal of Statistical Physics, 2003, 111 (1-2) :1–34. [doi:10.1023/A:1022221005097] |
[12] | Mezard M, Zecchina R. Random k-satisifability problem:From an analytic solution to an efficient algorithm. Physical Review E, 2002, 66 (5) :No.056126. [doi:10.1103/PhysRevE.66.056126] |
[13] | Maneva E, Mossel E, Wainwright M. A new look at survey propagation and its generalizations. Journal of the ACM, 2007, 54 (4) :1089–1098. [doi:10.1145/1255443.1255445] |
[14] | Braunstein A, Mezard M, Zecchina R. Survey propagation:An algorithm for satisfiability. Random Structures and Algorithms, 2005, 27 (2) :201–226. [doi:10.1002/rsa.20057] |
[15] | Braunstein A, Zecchina R. Survey and belief propagation on random k-SAT. Lecture Notes in Computer Science, 2004, 2919 (1) :519–528. [doi:10.1007/978-3-540-24605-3_38] |
[16] | Yedidia JS, Freeman WT, Weiss Y. Understanding belief propagation and its generalizations. Artificial Intelligence, 2003, 8 (1) :239–269. http://cn.bing.com/academic/profile?id=2019599312&encoded=0&v=paper_preview&mkt=zh-cn |
[17] | Mezard M, Parisi G, Zecchina R. Analytic and algorithmic solution of random satisfiability problems. Science, 2002, 297 (5582) :812–815. [doi:10.1126/science.1073287] |
[18] | Xu K, Li W. Exact phase transitions in random constraint satisfaction problems. Journal of Artificial Intelligence Research, 2000, 12 (1) :93–103. http://cn.bing.com/academic/profile?id=2139300273&encoded=0&v=paper_preview&mkt=zh-cn |
[19] | Liu T, Lin XX, Wang CY, Su KL, Xu K. Large hinge width on sparse random hyper graphs. In:Proc. of the 22nd Int’l Joint Conf. on Artificial Intelligence. Barcelona:AAAI, 2011. 611-616.[doi:10.5591/978-1-57735-516-8/IJCAI11-109] |
[20] | Xu K, Li W. Many hard examples in exact phase transitions. Theory Computer Science, 2006, 355 (1) :291–302. [doi:10.1016/j.tcs.2006.01.001] |
[21] | Xu K, Boussemart F, Hemery F, Lecoutre C. Random constraint satisfaction:Easy generation of hard instances. Artificial Intelligence, 2007, 171 (1) :514–534. [doi:10.1016/j.artint.2007.04.001] |
[22] | Zhao CY, Zheng ZM. A belief propagation algorithm based on variable entropy for constraint satisfaction problems. Science China Information Sciences, 2012, 42 (9) :1170–1180(in Chinese with English abstract). http://www.cnki.com.cn/Article/CJFDTOTAL-PZKX201209010.htm |
[23] | Kschischang F, Frey B, Loeliger H. Factor graph and the sum-product algorithm. Information Theory, 2001, 47 (1) :498–519. [doi:10.1109/18.910572] |
[24] | Weiss Y. Correctness of local probability propagation in graphical models with loops. Neural Computation, 2000, 12 (1) :1–41. [doi:10.1162/089976600300015880] |
[25] | Weiss Y, Freeman WT. Correctness of belief propagation in Gaussian graphical models of arbitrary topology. Neural Computation, 2001, 13 (10) :2173–2200. [doi:10.1162/089976601750541769] |
[26] | Tatikonda S, Jordan MI. Loopy belief propagation and Gibbs measure. In:Proc. of the 18th Annual Conf. on Uncertainty in Artificial Intelligence. 2002. 493-500. |
[27] | Heskes T. On the uniqueness of loopy belief propagation fixed points. Neural Computation, 2004, 16 (11) :2379–2413. [doi:10.1162/0899766041941943] |
[28] | Ihler AT, Fisher JW, Willsky AS. Loopy belief propagation:Convergence and effects of message errors. Machine Learning Research, 2005, 6 (1) :905–936. http://cn.bing.com/academic/profile?id=2163364417&encoded=0&v=paper_preview&mkt=zh-cn |
[29] | Mooij JM, Kappen HJ. Sufficient conditions for convergence of the sum-product algorithm. IEEE Trans. on Information Theory, 2007,53(12):4422-4437.[doi:10.1109/TIT.2007.909166] |
[30] | Shi XQ, Schonfeld D, Tuninetti D. Message error analysis of loopy belief propagation for the sum-product algorithm. Computer and Information Science, 2010, 1009 (2) :1–30. [doi:10.1109/ICASSP.2010.5495075] |
[31] | Gamarnik D, Shah D, Wei Y. Belief propagation for min-cost network flow:Convergence and correctness. Operations Research, 2012, 60 (2) :410–428. [doi:10.1287/opre.1110.1025] |
[32] | Brunsch T, Cornelissen K, Manthey B, Roglin H. Smoothed analysis of BP for minimum cost flow and matching. Journal of Graph Algorithms and Applications, 2013, 17 (6) :647–670. [doi:10.7155/jgaa.00310] |
[33] | Norshams N, Wainwright MJ. Stochastic belief propagation:A low complexity alternative to the sum product algorithm. IEEE Trans. on Information Theory, 2013,59(4):1981-2000.[doi:10.1109/TIT.2012.2231464] |
[34] | Lecoute. http://www.cril.univ-artois.fr/~lecoutre/benchmarks.html |
[35] | Dieudonné J. Foundations of modern analysis. Vol.10-I, New York:Academic Press, 1969. |
[22] | 赵春艳, 郑志明. 一种基于变量熵求解约束满足问题的置信传播算法. 中国科学:信息科学, 2012,42 (9) :1170–1180. http://www.cnki.com.cn/Article/CJFDTOTAL-PZKX201209010.htm |