龚光鲁,钱敏平著应用随机过程教程及其在算法与智能计算中的应用 清华大学出版社,2003 第9章以图象信息为背景的随机场,迭代 Markov系统 随机场的一个直观背景,是以图像为基本事件的随机元素.另一方面,迭代 Markov系 统是取值为压缩映射随机映射的迭代,其目的是得到其不变元素,以用来描述一个图像 1有限格点上的 Markov随机场与图象 1.1有限格点上的 Markov随机场 定义9.1(相邻系统与子团) 设G为平面有限格点集合:G={():1≤i≤M1l≤j≤N2}.对于G中的任意两 点x≠y,可以定义它们相邻,或者不相邻.这样定义了S中的一种相邻概念以后,就得到 了G上的一个相邻系统.注意,两点是否相邻完全可以是抽象的,并不需要在直观上比邻 (有时为了看起来有对称性,也可以认为N1与1恒同,N2也与1恒同,即认为G为环面 格点.) 在本书中,我们把x,y相邻记为x~y,记与x的相邻的点的全体为a(x).G的子集 C称为G的∂-子团,如果C中的任意两个元素都相邻.含x的子团记为 一个相邻系统由G中点对的相邻概念所决定,因此也由全体子团所决定,所以人们常 用∂来代表一个相邻系统,或用全体子团{C}来代表一个相邻系统 例9.2(1)若∂(x)=φ,则只有单点集是G的∂-子团 (2)若O(x)=G\{x,则G的所有子集都是∂一子团 (3)一阶相邻(最紧相邻,简称紧邻):*的相邻点。可以如下图 此时G的全部-子团形式为 (4)二阶相邻 的相邻点o可以如下图
229 龚光鲁，钱敏平著 应用随机过程教程及其在算法与智能计算中的应用 清华大学出版社, 2003 第 9 章 以图象信息为背景的随机场, 迭代 Markov 系统 随机场的一个直观背景, 是以图像为基本事件的随机元素．另一方面, 迭代 Markov 系 统是取值为压缩映射随机映射的迭代, 其目的是得到其不变元素, 以用来描述一个图像. 1 有限格点上的 Markov 随机场与图象 1. 1 有限格点上的 Markov 随机场 定义９.１（相邻系统与子团） 设G 为平面有限格点集合: {( , ):1 ,1 } 1 N2 G = i j £ i £ N £ j £ . 对于 G 中的任意两 点 x ¹ y , 可以定义它们相邻, 或者不相邻. 这样定义了S 中的一种相邻概念以后, 就得到 了G 上的一个相邻系统. 注意, 两点是否相邻完全可以是抽象的, 并不需要在直观上比邻 (有时为了看起来有对称性, 也可以认为 N1与 1 恒同, N2也与 1 恒同, 即认为 G 为环面 格点.) 在本书中, 我们把 x, y 相邻记为 x y ¶ ~ , 记与 x 的相邻的点的全体为¶(x) .G 的子集 C 称 为 G 的 ¶ - 子 团 , 如 果 C 中的任意两个元素都相邻 . 含 x 的子团 记 为 C(x) = U¶(x) . 一个相邻系统由 G 中点对的相邻概念所决定, 因此也由全体子团所决定, 所以人们常 用¶ 来代表一个相邻系统, 或用全体子团{C} 来代表一个相邻系统. 例９．２ (1) 若 ¶(x) = f , 则只有单点集是G 的¶ -子团. (2) 若 ¶(x) = G \ {x},则G 的所有子集都是¶ -子团.. (3) 一阶相邻(最紧相邻, 简称紧邻): * 的相邻点o 可以如下图 o o o o * 此时G 的全部¶ -子团形式为: f , * , * - * , * * | (4)二阶相邻: * 的相邻点o 可以如下图
此时G的全部∂-子团形式为 及其转动, 及其转动,|× 高阶相邻等可以如下图所示 54345 42124 631*136 其中数字代表分阶为相邻的阶.S的一个相邻系统也就是以G的点为顶点的一个双向图 定义9.3(随机场) 以平面矩形格点集G中的格点为参数,取值于{0…,L-1}(也可以是实直线R,或 个离散的集合)的随机变量族,称为取值于{0,…,L-1}的G一随机场,简称为随机场 按照统计物理中的术语,我们称{0,…,L-1}为自旋空间 于是,随机场的一个样本,就是在格点集G的每个格点上”坐着”{0,…,L-1}中某 个数.随机场的样本的一个直观背景是图像.黑白图象是连续变化的一个灰度图.而在实际 图象处理中的黑白灰度都是经过离散化采样的.设黑白灰度可分成L个(一般L=256)不同 的等级.{0,…,L-1}中的一个数值可以解释为黑白的灰度.最粗放的是黑白轮廓图,在每 个采样点上只有”黑”,“白”两个象素之一(即L=1的情形).一般G为256×256个采样 点(即M1=N2=256)那么,一幅黑白灰度图就可以看成在G的每一个格点上的一个灰 度(象素, Pixel),即给每一个格点赋予0,…,L-1中的一个值.所以,G一随机场的样本, 即做一次试验可能得到的结果,可以解释为G上的一个黑白灰度图 在这个意义上,随机场就也可以解释为随机图,或者说,它是在图像空间取值的一个随 机元素(简称为随机元).这是随机场的第2种解释.从相对灰度(每个格点上的灰度除以总 灰度(灰度的和)的角度也可以说,随机场的分布是格点集G上的一个随机的概率分布 定义9.4(组态)从数学的角度抽象,随机场的一个样本就是S={0,…,L-1}°的 个元素,记为a,B…,称为一个组态( onfiguration),即a=(a2:x∈G)
230 o o o o o o o o * 此时G 的全部¶ -子团形式为: f , * , * - * , * * | , * * \ 及其转动, * - * * | \ 及其转动, * - * ´ * - * | | . 高阶相邻等可以如下图所示 6 5 4 3 4 5 4 2 1 2 4 6 3 1 1 3 6 4 2 1 2 4 5 4 3 4 5 6 * 其中数字代表分阶为相邻的阶. S 的一个相邻系统也就是以G 的点为顶点的一个双向图. 定义９.３（随机场） 以平面矩形格点集G 中的格点为参数, 取值于{0,L,L -1} (也可以是实直线 R , 或 一个离散的集合)的随机变量族, 称为取值于{0,L,L -1} 的 G －随机场, 简称为随机场. 按照统计物理中的术语, 我们称{0,L,L -1} 为自旋空间. 于是, 随机场的一个样本, 就是在格点集 G 的每个格点上 ”坐着” {0,L,L -1} 中某一 个数. 随机场的样本的一个直观背景是图像．黑白图象是连续变化的一个灰度图. 而在实际 图象处理中的黑白灰度都是经过离散化采样的. 设黑白灰度可分成 L 个 (一般 L=256)不同 的等级. {0,L,L -1} 中的一个数值可以解释为黑白的灰度. 最粗放的是黑白轮廓图，在每 个采样点上只有 ”黑”, “白”两个象素之一 (即L = 1的情形)． 一般G 为 256´ 256 个采样 点 (即N1 = N2 = 256 ). 那么, 一幅黑白灰度图就可以看成在G 的每一个格点上的一个灰 度(象素,Pixel), 即给每一个格点赋予0,L,L - 1中的一个值. 所以, G －随机场的样本, 即做一次试验可能得到的结果, 可以解释为G 上的一个黑白灰度图. 在这个意义上, 随机场就也可以解释为随机图, 或者说, 它是在图像空间取值的一个随 机元素(简称为随机元). 这是随机场的第 2 种解释. 从相对灰度(每个格点上的灰度除以总 灰度(灰度的和))的角度也可以说, 随机场的分布是格点集G 上的一个随机的概率分布. 定义９.４（组态） 从数学的角度抽象, 随机场的一个样本就是 G S ={0, ,L -1} D L 的一 个元素 , 记 为 a, bL , 称为一个 组 态 (configuration), 即 ( : x G) a = ax Î
ax∈{0.…,L-1}称为组态a的x坐标,再记a=(2:x∈J(cG),称为组态a的J 坐标. 所以,随机场也可以看成是一个随机组态,或者看成为在组态空间S={0…,L-1}° 上取值的随机元.这是随机场的第3种解释.从分布的角度也可以说,随机场的分布是组态 空间上的一个概率分布 由此可见,一个灰度图就可以近似地视为一个组态 定义9.5(- Markov随机场) 随机场{x:x∈G}称为∂- Markov场,如果对于任意格点x,y∈G及任意组态 a∈S={…,L},恒有 P(x=a1|5,=ay,一切y≠x)=P(x=a215,=a,一切y~x).(9.1) a- Markov场的含义是:在任意格点x的余格点位置上取已知值的条件下,随机场在格 点x处取值的概率规律只与格点x的∂相邻点有关 由定义可以直接证明下述定理 定理9.6随机场{x:x∈G}是∂- Markov场的等价条件为:对任意格点集 AcG及任意组态α∈S={0,…,L-1}°,恒有 P(5x=ax,x∈A|5=a,yA)=P(5x=x,x∈A|5y=a,y∈(A) 其中a(4)={4中点的邻点全体} 1.2相邻系统的Gbs分布与Gbs随机场(∂-邻位势 Gibbs场) 定义9.7(格点集G的子集与组态的相互作用) 对于给定的一个子格点集A(cG),组态上的一个与之相关的函数U4(x) (a∈{0,…,L-1})称为A对于组态a的相互作用,如果它满足:对于任意两个组态 a,B∈{0,…L-1}°,只要5,(a)=5(B)(Wx∈A),就有U4a)=U(B)(U4(a) 代表组态a在A处的”坐标”对a的作用,所以要求对在A处”坐标”相同的任意组态的相 互作用一样) 定义9·8(位势函数)在补充定义U。=0后,相互作用 qUA: ACG) 的全体称为位势.对于给定的相邻系统∂,如果对于任意A,只要它不是∂一子团,就有 U4=0,那么{c:C为G的∂-子团}就称为∂-邻位势.如果∂-相邻正好就是直观上 的紧邻,则这样的的∂一邻位势简称紧邻位势. 31
231 Î{0, ,L -1} a x L 称为组态a 的 x 坐标. 再记 ( : x J ( G)) aJ = ax Î Ì , 称为组态a 的 J 坐标. 所以, 随机场也可以看成是一个随机组态, 或者看成为在组态空间 G S = {0,L, L -1} 上取值的随机元. 这是随机场的第 3 种解释. 从分布的角度也可以说, 随机场的分布是组态 空间上的一个概率分布. 由此可见，一个灰度图就可以近似地视为一个组态． 定义９.５ （¶ -Markov 随机场） 随机场{ : x G} x x Î 称为 ¶ - Markov 场, 如果对于任意格点 x, y ÎG 及任意组态 G a Î S = {1,L, L} , 恒有 P( | , y x) x x = ax x y = ay 一切 ¹ P( | , y ~ x) x x y y ¶ = x = a x = a 一切 . (9. 1) ¶ -Markov 场的含义是：在任意格点x 的余格点位置上取已知值的条件下，随机场在格 点 x 处取值的概率规律只与格点 x 的¶ 相邻点有关． 由定义可以直接证明下述定理 定理９．６ 随机场{ : x G} x x Î 是 ¶ - Markov 场的等价条件为: 对任意格点集 A Ì G 及任意组态 G a Î S = {0,L,L - 1} , 恒有 P( , x A | , y A) x x = ax Î x y = ay Ï P( , x A | , y (A)) x x y y = x = a Î x = a Î ¶ , 其中 (A) {A中点的邻点全体} D ¶ = . 1．２ 相邻系统的 Gibbs 分布与 Gibbs 随机场 (¶ -邻位势 Gibbs 场) 定义９．７（格点集G 的子集与组态的相互作用） 对 于给定 的一个子 格点集 A(Ì G) , 组态上的一个 与之相关的 函数 (a) UA ( G a Î{0,L,L - 1} ) 称为 A 对于组态a 的相互作用, 如果它满足: 对于任意两个组态 G a,b Î{0,L, L -1} , 只要 x (a) x (b) x = x ("x Î A ), 就有 (a) (b) UA =UA ( (a) UA 代表组态a 在 A 处的 ”坐标”对a 的作用, 所以要求对在 A 处”坐标”相同的任意组态的相 互作用一样). 定义９．８（位势函数） 在补充定义 = 0 Uf 后, 相互作用 {U : A G} A Ì 的全体称为位势. 对于给定的相邻系统¶ , 如果对于任意 A ，只要它不是¶ -子团，就有 UA = 0 ，那么 {U :C为G的¶ -子团} C 就称为¶ -邻位势. 如果¶ -相邻正好就是直观上 的紧邻，则这样的的¶ -邻位势简称紧邻位势
定义9.9(能量函数)对于∂一邻位势,组态的函数 称为∂一邻位势的能量函数 定义9.10(∂一邻位势的 Gibbs分布与Gibs随机场) 在组态空间S=(0,…,L-1}°上定义的概率分布兀=丌(a)(a∈S),称为一邻位 势(在能量函数下)的Gbbs分布(这里是离散分布),其中 Hu为a-邻位势的能量函数,而Z为一个归一化常数,称为 Gibbs分布的配分常数 Z 以Gibs分布为分布的随机场称为Gbs场,即 P(5 就是说,Gibs分布π是随机场在所有的格点处取值的联合分布 在第6,7章中的 Ising模型中的 Gibbs分布,就是这里的 Gibbs场的概率分布的特殊情 形,即L=2,自旋空间为{-1,}的情形,其对应的位势是紧邻的(而把外场看成格点对自 己的作用的和). Ising模型中 Gibbs分布中有一个倒温度参数,这是为了调节用的.这里 无妨省去 再则,由于在不同的邻域中的科学家有他们习惯用的记号,所以,在本节中的记号 6,7章中用的有些不同,请读者见谅由于记号改动所带来的不便.我们列出其对比 组态配称分布倒温度随机场 第6,7章5,na, 本节a,Bm 在第6,7章中,由于 Ising模型较为简单,我们只提出了 Gibbs分布的概念,而并未提出随 机场的概念.事实上, Gibbs分布是 Gibbs随机场的不变分布(有的文献中 Gibbs分布与 Gibs随机场的术语是混用的).因此在第6,7章中,我们如很多书上那样用ξ,n表示组态 而在本节中,我们强调的是随机场,它用ξ,n表示更为合适,所以我们就与很多书上不一样, 改用a,B表示组态 定理9.11( Hammersley- Clifford基本定理 对于给定的相邻系统∂,∂- Markov随机场与∂一邻位势 Gibbs场是一样的.】 这个定理的意义在于,用∂- Markov随机场来描述图象是一个直观认识,而用∂一邻位 势 Gibbs场则可归结为其能量函数,后者有时可以用简单的参数来表达(例如相互作用具有 平移不变性的时候).那时就可用几个参数来描述图象.从而可以达到信息压缩的目的.于 32
232 定义９.９（能量函数） 对于¶ -邻位势, 组态的函数 å ¶ D = C HU UC 子团 (9. 2) 称为¶ -邻位势的能量函数. 定义９.１０（¶ -邻位势的 Gibbs 分布与 Gibbs 随机场） 在组态空间 G S = {0,L , L - 1} 上定义的概率分布p = p (a),(a Î S ), 称为 ¶ - 邻位 势(在能量函数下)的 Gibbs 分布(这里是离散分布), 其中 1 ( ) ( ) a p a HU e Z - = , HU 为¶ -邻位势的能量函数, 而 Z 为一个归一化常数, 称为 Gibbs 分布的配分常数: åÎ - = S HU Z e a (a ) ． (9. 3) 以 Gibbs 分布为分布的随机场称为 Gibbs 场，即 P(x =a ,"x ÎG) = p (a) x x ． 就是说，Gibbs 分布p是随机场在所有的格点处取值的联合分布． 在第 6,7 章中的 Ising 模型中的 Gibbs 分布, 就是这里的 Gibbs 场的概率分布的特殊情 形, 即 L = 2 , 自旋空间为{-1,1} 的情形，其对应的位势是紧邻的(而把外场看成格点对自 己的作用的和). Ising 模型中 Gibbs 分布中有一个倒温度参数, 这是为了调节用的. 这里 无妨省去. 再则, 由于在不同的邻域中的科学家有他们习惯用的记号, 所以, 在本节中的记号与第 6,7 章中用的有些不同, 请读者见谅由于记号改动所带来的不便. 我们列出其对比: a b p x h a p b , 6,7 , , 本节 第 章 组态 配称分布 倒温度 t Xt x 随机场 . 在第 6,7 章中, 由于 Ising 模型较为简单, 我们只提出了 Gibbs 分布的概念, 而并未提出随 机场的概念. 事实上, Gibbs 分布是 Gibbs 随机场的不变分布（ 有的文献中 Gibbs 分布与 Gibbs 随机场的术语是混用的）．因此在第 6,7 章中, 我们如很多书上那样用x ,h 表示组态. 而在本节中, 我们强调的是随机场, 它用x ,h 表示更为合适, 所以我们就与很多书上不一样, 改用a, b 表示组态. 定理 9.1１ (Hammersley-Clifford 基本定理) 对于给定的相邻系统¶ , ¶ -Markov 随机场与¶ -邻位势 Gibbs 场是一样的. 】 这个定理的意义在于, 用¶ -Markov 随机场来描述图象是一个直观认识, 而用¶ -邻位 势 Gibbs 场则可归结为其能量函数, 后者有时可以用简单的参数来表达(例如相互作用具有 平移不变性的时候). 那时就可用几个参数来描述图象. 从而可以达到信息压缩的目的. 于
是,此定理将直观认识转化为参数表示,以便于数学处理 1.3图象处理的随机过程方法的思路原则概述 正如前面所述,一幅黑白图象是某个组态空间中的一个组态,在图像传输过程中,需要对于灰度图的 信息进行处理与压缩.对于接收到的信息还需要恢复( Retrieving)为原图像.常见的方法有 第1种方法(用能量函数的 Gibbs分布对图像建模) 按黑白图像的相对灰度(即每个格点的灰度除以图像的总灰度和)大小,把图像看成格点集G上的组态 空间S={01…,L}°中的一个组态,或者看成G上的一个概率分布.由于自然景物等的图象往往具有 一些局部相联系(空间不变性)的结构,这些局部结构可以抽象地归结为组态间按某种相邻系统的能量函数 从而可以把一个图像用一个能量函数取到最小值的组态,也就是用 Gibbs分布达到最大值处的组态,从而 可以通过能量函数或 Gibbs分布来对图像建模.这个方法的要害处是:如何由图像给出其能量函数.能量 函数表示图像能起一定的信息压缩作用.而反过来由能量函数恢复图像的问题,则可以用模拟退火或其它 Markov Monte Carlo方法运作 第2种方法(用低阶 Markov场的样本对图像建模) 将黑白图像看成随机场的一个样本.由定理9.11知道任意随机场都可以看成 Markov场.而图像的 局部结构又说明,用较低阶相邻的 Markov场的样本来描述图象更为合理. 第3种方法(用迭代随机映射的不动点对图像建模) 对于一个黑白轮廓图像,构造与之相系的有限个压缩映射,使得由此压缩映射组给出的一个迭代映射 的不动点近似地是此轮廓图.而对于一个黑白灰度图,构造与之相系的一个随机地取有限个压缩映射的随 机映射,使此图像近似地是此随机映射的不变分布.(这在第3节中将较多地介绍).这种方法实际上是处 理图像的分形方法的另一角度的发展 一般地,彩色图可以看成红,黄,兰三张灰度图的叠 收到了不清晰的图像,也就是受到了干扰的图像,就需要对图像进行滤波,称为图像的清洗.其中常 见的一种清洗方法为 Baves方法,即把未清洗的图像看成先验分布,而寻找转移机制,使清洗后的恢复图象 是此先验分布的 Bayes分布(后验分布)图像受到污染的典型情形,例如,被测量的图象象素ax受到检测 与录制系统的干扰,使实际测到的象素为 B2=g(∑ha,)+" 其中g为一个已知函数,b,为已知的通道系数,而噪声x的统计规律为已知这时要恢复图象象素a 就需要进行滤波 此外,图象边缘信息也可能受到干扰而产生变形,也需要用边缘检测来恢复图象的信息.因此有时应 在图象的先验信息中加进“边-随机场”的先验分布 有时在处理图像前可能还需要对图像进行有部分相重叠的分割( segmentation) Gibbs分布的样本的 Gibbs采样法( Gibbs Sampler) 要对 Gibbs场作一些统计估计(例如参数估计),就先要得到 Gibbs分布的配分常数.但 是由于组态空间S={1…,L}°的元素数通常非常庞大,例如可以有100,这就使直接计算
233 是, 此定理将直观认识转化为参数表示, 以便于数学处理. １．３ 图象处理的随机过程方法的思路原则概述 正如前面所述，一幅黑白图象是某个组态空间中的一个组态．在图像传输过程中, 需要对于灰度图的 信息进行处理与压缩. 对于接收到的信息还需要恢复（Retrieving）为原图像．常见的方法有: 第 1 种方法 (用能量函数的 Gibbs 分布对图像建模) 按黑白图像的相对灰度( 即每个格点的灰度除以图像的总灰度和)大小, 把图像看成格点集 G 上的组态 空间 G S = {0,1,L, L} 中的一个组态，或者看成G 上的一个概率分布. 由于自然景物等的图象往往具有 一些局部相联系(空间不变性)的结构, 这些局部结构可以抽象地归结为组态间按某种相邻系统的能量函数. 从而可以把一个图像用一个能量函数取到最小值的组态，也就是用 Gibbs 分布达到最大值处的组态，从而 可以通过能量函数或 Gibbs 分布来对图像建模. 这个方法的要害处是：如何由图像给出其能量函数．能量 函数表示图像能起一定的信息压缩作用．而反过来由能量函数恢复图像的问题，则可以用模拟退火或其它 Markov Monte Carlo 方法运作． 第 2 种方法 (用低阶 Markov 场的样本对图像建模) 将黑白图像看成随机场的一个样本. 由定理９.１１知道任意随机场都可以看成 Markov 场. 而图像的 局部结构又说明,用较低阶相邻的 Markov 场的样本来描述图象更为合理. 第 3 种方法 (用迭代随机映射的不动点对图像建模) 对于一个黑白轮廓图像，构造与之相系的有限个压缩映射，使得由此压缩映射组给出的一个迭代映射 的不动点近似地是此轮廓图．而对于一个黑白灰度图，构造与之相系的一个随机地取有限个压缩映射的随 机映射，使此图像近似地是此随机映射的不变分布．（这在第３节中将较多地介绍）．这种方法实际上是处 理图像的分形方法的另一角度的发展． 一般地，彩色图可以看成红，黄，兰三张灰度图的叠置． 收到了不清晰的图像, 也就是受到了干扰的图像, 就需要对图像进行滤波, 称为图像的清洗. 其中常 见的一种清洗方法为 Bayes 方法, 即把未清洗的图像看成先验分布, 而寻找转移机制, 使清洗后的恢复图象 是此先验分布的 Bayes 分布(后验分布). 图像受到污染的典型情形, 例如, 被测量的图象象素ax 受到检测 与录制系统的干扰, 使实际测到的象素为 åζ = + ( ) ( ) y x x y x wx b g h a , 其中 g 为一个已知函数, y h 为已知的通道系数, 而噪声 wx 的统计规律为已知. 这时要恢复图象象素ax 就需要进行滤波. 此外, 图象边缘信息也可能受到干扰而产生变形, 也需要用边缘检测来恢复图象的信息. 因此有时应 在图象的先验信息中加进 “边-随机场”的先验分布. 有时在处理图像前可能还需要对图像进行有部分相重叠的分割 (segmentation). １. ４ Gibbs 分布的样本的 Gibbs 采样法(Gibbs Sampler) 要对 Gibbs 场作一些统计估计（例如参数估计）, 就先要得到 Gibbs 分布的配分常数. 但 是由于组态空间 G S = {1,L, L} 的元素数通常非常庞大, 例如可以有 6 10 10 , 这就使直接计算
配分常数Z时的求和项太多,从而常会出现大量被求和的项的值小于计算机误差的精度, 以致使计算在实际上失去可能.用 Markov链 Monte carlo方法(MCMC)可以解决这个困难 即先近似地得到 Gibbs分布,再进行参数估计.所以MMC是一个重要的计算手段 MCMC方法就是构造以 Gibbs分布为可逆分布的 Markov链的转移矩阵.因为这里的随机 场的 Gibbs分布远比 Ising模型的 Gibbs分布复杂,我们只构造离散时间的 Markov转移矩 阵,而且还希望如 Glauber动力学那样,每次只在一个格点上改变(称为按 Gibbs方式转移) 最后还要保证由这样得到的转移矩阵的各行都收敛于 Gibbs分布 为了使构造更为透明,我们先对任给的格点子集J定义如下的“部分转移” B)={Z en(v)(a,B在/外坐标相同) (其它) 其中β,as表示:把组态a中在J处的坐标换成组态B中的相应坐标,而Z(a)是一个归 化常数(“局部的配分常数) ∑e (9.5) 这种“部分转移”表示组态间只在位置J外的坐标相同时才可能转移 在此基础上我们定义转移矩阵如下:记G中的全部格点的一个给定的排序为 {x1…,xa,(G|=N1N2).任给组态a,B∈S,定义 P,=(P,(a, B))a,Bes,P= 则P是一个转移矩阵.注意由于P中的转移只在一个格点x,上进行,由P构造的 Markov 链的一次转移,可以视为G|次连续转移的结果,其中每次转移只在一个格点上进行,而且 要做遍所有的格点 下面的遍历收敛定理是 Markov链的 Monte carlo方法的理论基础.而 Markov链 Monte Carlo方法主要用于获得 Gibbs分布的样本,以便进一步得到Gibs场的各种统计平均量. 对于任意格点x∈G,任意组态a∈S,记 Hu(r,aGuri)=m,=min n, Hu(nraGyr) 那么 (Hu(,:)-m,I Pir(a, B -Hu(B! -m, e max( Hu(B)Hu(y)BGx=yGx)A G 其中δx为能量函数HU在格点x上的振幅,Δ为Hu的振幅,其确切定义为
234 配分常数 Z 时的求和项太多，从而常会出现大量被求和的项的值小于计算机误差的精度， 以致使计算在实际上失去可能．用 Markov 链 Monte Carlo 方法 (MCMC)可以解决这个困难， 即先近似地得到 Gibbs 分布, 再进行参数估计. 所以 MCMC 是一个重要的计算手段. MCMC 方法就是构造以 Gibbs 分布为可逆分布的 Markov 链的转移矩阵. 因为这里的随机 场的 Gibbs 分布远比 Ising 模型的 Gibbs 分布复杂, 我们只构造离散时间的 Markov 转移矩 阵, 而且还希望如Glauber动力学那样, 每次只在一个格点上改变（称为按 Gibbs方式转移）, 最后还要保证由这样得到的转移矩阵的各行都收敛于 Gibbs 分布. 为了使构造更为透明, 我们先对任给的格点子集 J 定义 如下的“部分转移” ï î ï í ì = - 0 ( ) ( , ) ( ) 1 ( , ) ( ) \ 其它 e 在J外坐标相同 p Z HU J S J J J a b a b a b a , (9. 4) 其中 b JaS \ J 表示: 把组态a 中在 J 处的坐标换成组态 b 中的相应坐标, 而 (a) ZJ 是一个归 一化常数(“局部的”配分常数): å - = J HU J S J J Z e b b a a ( ) \ ( ) . (9. 5) 这种 “部分转移”表示组态间只在位置 J 外的坐标相同时才可能转移. 在此基础上我们定义转移矩阵如下：记 G 中的全部格点的一个给定的排序为 { , , ),(| | ) 1 | | G N1N2 x x L G = . 任给组态a,b Î S , 定义 P J J S p = a b a ,bÎ ( ( , )) , P = Õ= | | 1 G k P ( xk } , (9. 6) 则 P 是一个转移矩阵. 注意由于 P ( xk } 中的转移只在一个格点 i x 上进行, 由 P 构造的Markov 链的一次转移,可以视为| G |次连续转移的结果, 其中每次转移只在一个格点上进行, 而且 要做遍所有的格点. 下面的遍历收敛定理是 Markov 链的 Monte Carlo 方法的理论基础. 而 Markov 链 Monte Carlo 方法主要用于获得 Gibbs 分布的样本, 以便进一步得到 Gibbs 场的各种统计平均量. 对于任意格点 x ÎG ,任意组态a Î S , 记 ( ) min ( } HU xa G\{ x} mx x HU hxaG\{ x} g h D D = = . (9. 7) 那么 ( , ) p{x} a b åÎ - - - - = y G H m H m U y G y x U x G x x e e [ ( ) ] [ ( ) ] \{ } \{ } b a b a - -D D - - = ³ = ³ e G e G G e x HU HU G x G x | | 1 | | 1 | | max( ( ) ( ): ) \{ } \{ } d b g b g , 其中 x d 为能量函数 HU 在格点 x 上的振幅, D 为 HU 的振幅, 其确切定义为
δ,= max a(x2=a I Hu(a)-Hu(B) 由此推出 C(Px)≤l-G| min aB Pirl(a,B)≤1 从而得到 (P)=C(Px1)…C(Px)≤(1 (9.10) 定理9.12(遍历收敛定理)以上定义的转移矩阵P的 Dobrushin收缩系数满 C(P)≤1-e(△为H的振幅) 记以P为转移矩阵的 Markov链为n,则它以Gibs分布兀为可逆分布,而且有 证明显见r(a)p(a,B)=r(B)P(B,a)·由P=(p(,B)a,Bes的定义立得 丌(a)p(a,B)=(B)p(B,a).因此兀为可逆分布.又由C(P)0∑g(x)=1例如均匀分布)的随机数作随机选取来代替.这时转移矩阵P应 作相应的改动,代替P1x(a,B)的是,以g(x)的概率对只在x坐标不同的组态间进行转移,即令 =(P),P=(p(ax,B) 其中 ,B)=8(x)(x,B)(若o=Bm:) (其它情形) 这时我们仍有 1 1.5 Gibbs分布的模拟退火 本段讨论求得达到组态空间上能量函数最小值的组态的概率方法,就是用能量函数的 Gibs分布作模拟退火的方法.考虑依赖于温度Tn的Gibs分布,以它们作为不变分布构造 各步为非时齐Gibs转移矩阵,由此确定的一个非时齐 Markov链{n}.我们可以用它作模 拟退火,使得当n充分大以后,此非时齐的 Markov链n以大概率落在能量函数最小的组态 5
235 max | ( ) ( )| \{ } \{ } d x a b HU a HU b G x G x = = - , x G x d Î D D =max . (9. 8) 由此推出 -D C P £ - G p £ - e ( {x} ) 1 | | min a ,b {x} (a,b) 1 . (9. 9) 从而得到 | | { } { } ( ) ( ) ( ) (1 ) 1 G x x C P C P C P e M -D = L £ - | | 1 G e -D £ - . (9. 10) 定理 9. １２ （遍历收敛定理） 以上定义的转移矩阵P 的 Dobrushin 收缩系数满 足 | | ( ) 1 G C P e -D £ - ( D 为 HU 的振幅). 记以 P 为转移矩阵的 Markov 链为 n x , 则它以 Gibbs 分布p 为可逆分布, 而且有 P p n n T ¾ ¾®1 ®¥ . 证明 显见 p (a) (a,b ) p (b ) (b,a) pJ = pJ . 由 P S p = a b a,bÎ ( ( , )) 的定义立得 p (a) p(a, b ) = p (b ) p(b ,a) . 因此p 为可逆分布. 又由 C(P) 0,å ( ) = 1 xÎG g x g x 例如均匀分布) 的随机数作随机选取来代替. 这时转移矩阵 P 应 作相应的改动，代替 ( , ) { } a b xk p 的是, 以 g (x) 的概率对只在 x 坐标不同的组态间进行转移, 即令 P ( D = ~ P | | ) G , ~ P = p a b a,bÎS ~ ( ( , )) , 其中 î í ì = = 0 ( ) ( ) ( , ) ( ) ( , ) { } \{ } \{ } ~ 其它情形 p x 若 G x G x g x p a b a b a b . 这时我们仍有 P p n n T ¾ ¾®1 ®¥ . １. ５ Gibbs 分布的模拟退火 本段讨论求得达到组态空间上能量函数最小值的组态的概率方法, 就是用能量函数的 Gibbs 分布作模拟退火的方法. 考虑依赖于温度Tn的 Gibbs 分布, 以它们作为不变分布构造 各步为非时齐 Gibbs 转移矩阵，由此确定的一个非时齐 Markov 链{ }n x . 我们可以用它作模 拟退火, 使得当n 充分大以后，此非时齐的 Markov 链 n x 以大概率落在能量函数最小的组态
所组成的集合中 更具体地说,对于如下的依赖于温度Tn的O一邻位势Gibs分布 -Hula) -Hu(a) 相应地定义 pj(a, B)=z.(a) (a,B在/外坐标相同) (9.11) (其它) -Hu(Busy) (a)= ∑ 再令 1(=(P1(a,B)s,P=∏P,(9.12) TEG 则P(是一个转移矩阵.我们将以它为非时齐转移矩阵列所确定的非时齐 Markov链记为 {n} 定理9.13( Gibbs分布的模拟退火的收敛性定理) 如果温度T满足如下的 Geman- Geman条件:对0<T√0,存在n0使n≥n0后有 Tn≥ (9.13) In 其中Δ是能量函数Hu(a)的振幅 A=maxes Hu(a)-m H() 那么在任意初始分布下,以P为第n步转移矩阵列的非时齐的 Markov链{n},在时 刻n的绝对分布n有极限 n→ 其中丌)为集合{a∈S:H(x)=min}上的离散均匀分布.从而还有 lm n P(Hy(sn)=minges hy(a))=1 证明我们验证 Dobrushin- Isaacson- Madsen定理的条件成立.条件(A.1)的验证与第 8章3.2段的模拟退火定理中的证明完全一样.下面验证 Dobrushin条件(A.2)’.这里 定理9.12中的不等式变为 36
236 所组成的集合中. 更具体地说, 对于如下的依赖于温度Tn的¶ -邻位势 Gibbs 分布 ( ) 1 ( ) 1 ( ) a p a U n H T n n e Z - = , åÎ - = S H T n U n Z e a (a ) 1 , 相应地定义 ï î ï í ì = - 0 ( ) ( , ) ( ) 1 ( , ) ( ) 1 , ( ) \ 其它 e 在J外坐标相同 p Z U J S J n H T n J n J a b a b a b a , (9. 11) å - = J U J S J n H T n J Z e b b a a ( ) 1 , \ ( ) . 再令 P ( ( ) { } D = n x S n = p x a b a,bÎ ( ) { } ( ( , )) , P = (n) ÕxÎG P ( ) { } n x , (9. 12) 则 P (n) 是一个转移矩阵. 我们将以它为非时齐转移矩阵列所确定的非时齐 Markov 链记为 { }n x . 定理９.１３（Gibbs 分布的模拟退火的收敛性定理） 如果温度Tn满足如下的 Geman-Geman 条件: 对0 < ¯ 0 Tn , 存在n0 使n ³ n0后有 n G Tn ln | | D ³ , (9. 13) 其中D 是能量函数 (a) HU 的振幅: D max (a) min (a) a a D = ÎS HU - ÎS H , 那么在任意初始分布 m0下, 以P (n) 为第n 步转移矩阵列的非时齐的 Markov 链{ }n x , 在时 刻n 的绝对分布 mn 有极限 ® (¥) mn p , 其中 (¥ ) p 为集合{a Î : (a) = min} S HU 上的离散均匀分布. 从而还有 lim n®¥ P(HU (xn ) = min aÎS HU (a)) = 1. 证明 我们验证 Dobrushin-Isaacson-Madsen 定理的条件成立. 条件(A.1)的验证与第 8 章 3. 2 段的模拟退火定理中的证明完全一样. 下面验证 Dobrushin 条件(A. 2)’. 这里 定理９.１２中的不等式变为
其中△为能量函数HU的振幅.于是对于PP+1)…P,我们有 C(PP)…P)≤C(PC(P+)…C(P)≤∏(1-e) 由无穷乘积的理论,当n→∞时,上面右方趋于0的充要条件为级数∑e是发散的 但是由 Geam- Geam条件,我们有 所以 Dobrushin条件(A.2)成立 [注1] Gibbs场的模拟退火在图象的清污(图象的滤波)中的应用 设图象为组态a,而在传输过程这受到了污染,使得接受方观测到的图象为组态B.需要用数学方法 进行滤波,以便恢复图象的清晰度,我们把在观测为组态B的条件下,状态a的分布,记为P(a|B) (就是后验分布 Posteriori distribution,即 Baves分布).一般它也可用写成Gbbs分布的形式 P(aB) Z(B) 其中H(a|B)称为条件能量函数.典型的例子如 H(aB)=0∑(a(x)-a(y)2+∑a(x)-B(x)2 其中第一项强调了图形的光滑性,第二项为保真度.滤波的目的就是要找组态α,使条件能量函数达到最 H(a|β)= min h(a|B 这个最佳组态O可以通过模拟退火方法求得 又例如,在边缘检测中,β就是离散化点上的象素,而α就是“边随机链”的一个样值,与滤波问 的不同的是,此时的后验图形只是先验图形的一部分,即边的位置.在复杂的图形处理技术中,例如,X断 层摄影术的处理,情况就远为复杂.适当的利用分割技术可以得到较好的恢复效果 [注2]图象的处理有许多数学方法, Gibbs分布方法只是其中的一种.另一种更为常用的方法是,把 图象的一些特征参数作为观察量(这同时也达到了图象信息压缩的目的.用它们来估计 Markov链的状态 就是重建图象.这就是隐 Markov方法(M),一般地说,如果特征参数取得合适,HM方法是十分有效的. 例如,纹理图象分析方法( Texture Analysis)提取特征参数.其它的方法还有:下节将介绍的随机函数迭
237 - D £ - Tn G n C P e | | ( ) ( ) 1 , 其中 D 为能量函数 HU 的振幅. 于是对于 (j) (j 1) (n) P P LP + , 我们有 £ + ( ) ( j) ( j 1) (n ) C P P LP £ + ( ) ( ) ( ) ( j) ( j 1) (n) C P C P LC P (1 ) | | | | n T G G k j e D - = Õ - . 由无穷乘积的理论, 当n ® ¥时, 上面右方趋于 0 的充要条件为级数å D - n T G n e | | 是发散的. 但是由 Gemam-Gemam 条件, 我们有 å D - n T G n e | | å - ³ n n e ln = å = ¥ n n 1 . 所以 Dobrushin 条件(A. 2)’成立. ［注 1］ Gibbs 场的模拟退火在图象的清污(图象的滤波)中的应用 设图象为组态a , 而在传输过程这受到了污染, 使得接受方观测到的图象为组态 b . 需要用数学方法 进行滤波, 以便恢复图象的清晰度. 我们把在观测为组态 b 的条件下, 状态a 的分布, 记为P(a | b ) (就是后验分布 Posteriori distribution, 即 Bayes 分布). 一般它也可用写成 Gibbs 分布的形式: ( | ) ( ) 1 ( | ) a b b a b H e Z P - = , (9. 14) 其中 H (a | b )称为条件能量函数. 典型的例子如 H (a | b ) å å ζ = - + - ( ) 2 2 2 ( ( ) ( )) 2 1 ( ( ) ( )) x y x x y a x b x s q a a , 其中第一项强调了图形的光滑性, 第二项为保真度. 滤波的目的就是要找组态 ^ a , 使条件能量函数达到最 小: ( | ) min ( | ) ^ H a b = a H a b . 这个最佳组态 ^ a 可以通过模拟退火方法求得. 又例如, 在边缘检测中, b 就是离散化点上的象素, 而a 就是 “边随机链”的一个样值, 与滤波问题 的不同的是, 此时的后验图形只是先验图形的一部分, 即边的位置. 在复杂的图形处理技术中, 例如, X 断 层摄影术的处理, 情况就远为复杂. 适当的利用分割技术可以得到较好的恢复效果． [注 2] 图象的处理有许多数学方法, Gibbs 分布方法只是其中的一种. 另一种更为常用的方法是, 把 图象的一些特征参数作为观察量(这同时也达到了图象信息压缩的目的. 用它们来估计 Markov 链的状态 , 就是重建图象. 这就是隐 Markov 方法(HMM), 一般地说, 如果特征参数取得合适, HMM 方法是十分有效的. 例如, 纹理图象分析方法(Texture Analysis)提取特征参数. 其它的方法还有: 下节将介绍的随机函数迭
代系统方法.除此以外,还有小波分析等非概率方法, [注3]当格点集G为无限集时,组态的数目就不再是可数的了.此时就应归入连续状态的 Markov 场的范畴.在统计物理中应用的模型,就是这种模型.有时自旋空间S也可以是连续的,例如,圆周或其 它集合.这时就可能存在多个不变分布,称为有相变,因为每个分布代表了统计物理中质点的一个宏观 相”.统计物理中也正是对这种相变的存在性更感兴趣. [注4]若每次发展的修改的格点集合J不是单点集,且各次的修改集可以不同,那么,此时的 Markov 链未必有遍历定理,这时并不能保证此 Markov链5n的时间平均收敛到不变分布,因此,即使n很大,也 不能认为5n的样本近似为不变分布的样本 2时间离散状态连续的 Markov链 概率空间再访 设g2={0:O-基本事件},7是样本空间的某些子集组成的集合类(事件体),满足 对于任意AAn∈?(n21),恒有UAn∩A1,A=92-A∈只 (这样的事件体?代表了可以通过理论推算,或实际测量能够得到概率的全体随机事件).在 概率理论中?称为Ω上的一个σ-代数.(Ω,7)称为概率空间.于是概率P就是定义在 σ-代数?上的一个非负函数 定义9.14( Borel集与 Borel函数) 设是R的某些子集组成的一个-代数,且满足以下条件 (1)R4中的任意开集G∈ (2)任意一个有R的某些子集组成的σ-代数另,只要它含有一切开集,则就有 那么因称为R“上 Borel集合类,其中的集合称为Bore|集.即 Borel集合类是以开集为元 素的最小σ-代数.(以前在本书中的所谓“常见的”集合,均指 Borel集).Rl上的 Borel 集称为一维 Borel集 再设w为R上的一个满足下述条件的函数类 1)对极限封闭:若∫n∈w且fn(x)→f(x),则有∫∈ (2)任意连续函数∫∈w, (3)若R上一个函数类∠只要对极限封闭,而且包含全体连续函数,则有wC∠
238 代系统方法. 除此以外, 还有小波分析等非概率方法. [注 3] 当格点集 G 为无限集时, 组态的数目就不再是可数的了. 此时就应归入连续状态的 Markov 场的范畴. 在统计物理中应用的模型, 就是这种模型. 有时自旋空间 S 也可以是连续的, 例如, 圆周或其 它集合. 这时就可能存在多个不变分布, 称为有相变, 因为每个分布代表了统计物理中质点的一个宏观 “相”. 统计物理中也正是对这种相变的存在性更感兴趣. [注 4] 若每次发展的修改的格点集合 J 不是单点集, 且各次的修改集可以不同, 那么, 此时的 Markov 链未必有遍历定理, 这时并不能保证此 Markov 链 n x 的时间平均收敛到不变分布, 因此, 即使n 很大, 也 不能认为 n x 的样本近似为不变分布的样本. 2 时间离散状态连续的 Markov 链 2. 1 概率空间再访 设W ={w :w - 基本事件}，F 是样本空间W 的某些子集组成的集合类(事件体), 满足 对于任意 A, An ÎF (n ³ 1)，恒有 U I ¥ = ¥ = = - Î 1 1 , , n n n An A A W A D F. (这样的事件体F 代表了可以通过理论推算, 或实际测量能够得到概率的全体随机事件). 在 概率理论中 F 称为 W 上的一个s - 代数. (W , F ) 称为概率空间. 于是概率P 就是定义在 s - 代数 F 上的一个非负函数. 定义９.１４（Borel 集与 Borel 函数） 设 B 是 d R 的某些子集组成的一个s - 代数，且满足以下条件: (1) d R 中的任意开集GÎB , (2) 任意一个有 d R 的某些子集组成的s - 代数 G, 只要它含有一切开集, 则就有 B Ì G． 那么 B 称为 d R 上 Borel 集合类, 其中的集合称为 Borel 集. 即 Borel 集合类是以开集为元 素的最小s - 代数. (以前在本书中的所谓 “常见的”集合, 均指 Borel 集). 1 R 上的 Borel 集称为一维 Borel 集． 再设 M 为 d R 上的一个满足下述条件的函数类: (1) 对极限封闭: 若 f n ÎM 且 f (x) f (x) n ® , 则有 f ÎM, (2) 任意连续函数 f ÎM, (3) 若 d R 上一个函数类 L, 只要对极限封闭, 而且包含全体连续函数, 则有 MÌ L．
