哈夫曼树|Huffman Tree

在数据通信中,如果允许对不同的字符采用不等长的二进制码进行表示,那么这种编码方式称为 可变长度编码(variable-length code)。该编码方式对使用频率较高的字符赋以短编码,从而使得字符的平均编码长度减小,实现了压缩数据的效果。

哈夫曼编码 就是一种被广泛运用的数据压缩编码,通常可以节省 20%~90% 的空间。

如果没有一个编码是另一个编码的前缀,则称这样的编码是 前缀码(prefix code)
例如 二进制串 001011101,如果根据前缀码的唯一性(没有二义性)的思路进行解码,那么可以得到0,0,101,1101,这个解码过程是没有歧义的。

前缀码被证明是可以保证达到最优数据压缩率的一种编码,而哈夫曼编码就是基于前缀码实现的。

字符解码

在上面的例子中,我们并没有给出解码的过程。这里再展开阐述。
前缀码的解码依赖于二叉树,即我们可以通过构造二叉树的方式轻松对前缀码进行解码。

构造的二叉树应具有如下特点:

  1. 叶子结点用于表示所给的的字符;
  2. 字符的二进制码通过根结点到对应的叶子结点的路径表示,“0”表示转向左子结点,“1”表示转向右子结点。

注:该二叉树并不是BST,其自己子结点并未有序排列,内部结点也不包含key.

如下图所示,按照这个解码逻辑以及给出的树,我们可以轻松解码刚刚的001011101aabe
解码示例

文件的最优编码方案总是对应一棵 满二叉树。即任何树内的非叶子结点都有左右两个孩子。

编码与形式化

我们为二叉树定义了一个衡量其是否最优的指标:带权路径长度

WPL=i=1nwiliWPL=\sum_{i=1}^n w_il_i

其中, 表示的是一棵树中结点被赋予某种意义的数值, 路径长度 是从树根到任意结点经过的边数,而 WPL 就是我们定义的带权路径长度,由叶子结点与其上的权的乘积组成。wiw_i就是第ii 个叶子结点所带权值weightweightlil_i 就是第ii 个叶子结点到根结点的路径长度lengthlength.

如果一棵二叉树使得 WPL 最小,那么这棵树我们就称为 最优二叉树

根据前缀码的特点,以及二叉树的特性,我们引入 WPL 就可以得到类似的指标:平均传输位数
二叉树的权用字符的出现频率表示,路径长度即是二进制码位数。

如果一棵二叉树的平均传输位数最小,那么其对应的编码方式我们就称为 最优前缀码


哈夫曼为了构造最优二叉树,利用了贪心算法的思想。

假定CC 是一个含有nn 个字符的字符集合。cC\forall c\in C,用c.freqc.freq 表示CC 中的字符cc 出现的频率
于是,目标函数便是最优化 平均传输位数(WPL):

WPL=cCc.freq×dT(c)WPL=\sum_{c\in C}c.freq\times d_T(c)

其中,dT(c)d_T(c) 表示字符cc 在二叉树TT 中的深度。

算法伪代码如下:

  • 输入:C={c1,c2,...,cn}C=\{c_1,c_2,...,c_n\}
  • 输出:QQ 构造完成的二叉树

Algorithm:   HUFFMAN(C)1.  nC2.  QC  //Q为按频率 C.freq 递增的优先队列3.  for  i  =1  to  n1  do4.  zAllocate-Node()  //生成一个新结点z5.  z.leftxExtract-min(Q)6.  z.rightyExtract-min(Q)7.  z.freqx.freq+y.freq7.  Insert(Q,z)8.  return  Q\begin{aligned} &\text{Algorithm: }\;\text{HUFFMAN}(C)\\\\ 1.&\;n\leftarrow|C|\\ 2.&\;Q\leftarrow C\;//Q\text{为按频率 }C.freq\text{ 递增的优先队列}\\ 3.&\;\mathbf{for}\;i\;=1\;\mathbf{to}\;n-1\;\mathbf{do}\\ 4.&\;\qquad z\leftarrow \text{Allocate-Node}()\;//\text{生成一个新结点}z\\ 5.&\;\qquad z.left\leftarrow x\leftarrow\text{Extract-min}(Q)\\ 6.&\;\qquad z.right\leftarrow y\leftarrow\text{Extract-min}(Q)\\ 7.&\;\qquad z.freq\leftarrow x.freq+y.freq\\ 7.&\;\qquad \text{Insert}(Q,z)\\ 8.&\;\mathbf{return}\;Q \end{aligned}

哈夫曼树的构造算法就是:自底向上地,将C|C| 个叶子结点通过关于c.freqc.freq最小优先队列QQ 进行两两合并。合并得到的新对象,记其频率为构成该对象的原来两个对象的频率之和,放回队列中。
一共进行了C1|C|-1 次合并后,得到结果。

算法流程示例

如果优先队列QQ 采用 最小二叉堆 实现,则每一次循环取出当前最小元的过程需要花费O(logn)O(\log n),循环次数O(n1)O(n-1),事先装载队列需要O(n)O(n).
综合来看,该算法的时间复杂度为O(nlogn)O(n\log n),空间复杂度O(n)O(n).

正确性证明

哈夫曼树的构造方法被证明完全能得到最优解,即其构造的二叉树的WPL一定最小。因此,我们也把最优二叉树叫做哈夫曼树

下面给出其正确性的证明。

因为哈夫曼算法采用的是贪心策略,因此我们要证明其正确性需要从 贪心选择性质最优子结构 两方面进行论证。

贪心选择性质

贪心选择性质 引理1

如果x,yC\exists x,y \in C,且x.freqx.freqy.freqy.freq 最小,则存在一个CC 的最优前缀码使得x,yx,y 的二进制码位数相等,且只有最后一位的数值不等。

引理1的证明

证明TT 是关于CC 的任意一个最优前缀码对应的最优二叉树。令a,ba,bTT 中深度最大的两个叶结点(aba\neq b
不妨设a.freqb.freq  &&  x.freqy.freqa.freq \leq b.freq\;\&\&\; x.freq\leq y.freq.
由于a,ba,b 虽然深度最深,但其频率不一定最小,且一定不会比x,yx,y 还小。所以有:

x.freqa.freq      y.freqb.freqx.freq \leq a.freq\;\;\;\\y.freq \leq b.freq

  • x.freq=b.freqx.freq=b.freq 时,由于x.freqy.freqx.freq\leq y.freq ,所以x,bx,b 都是频率最低的,而已知yy 也是频率最低的两个字符中的一个,因此有y=by=by=xy=x ,即y.freq=x.freqy.freq=x.freq,最后有:

x.freqa.freqb.freq=x.freq=y.freqa.freq=b.freq=x.freq=y.freq\begin{aligned}x.freq \leq a.freq \leq b.freq =x.freq=y.freq\\\\\Rightarrow a.freq=b.freq=x.freq=y.freq\end{aligned}

这表明,如果将a,b,x,ya,b,x,y 的位置交换,不改变 WPL 的值,问题得证。

  • x.freqb.freqx.freq\neq b.freq 时,即xbx\neq b 时,我们将x,ax,a 进行交换,得到区别于TT 的新树TT',将y,by,b 进行交换,得到区别于TT' 的新树TT''.
    (注:如果x=b,yax=b,y\neq a 那么上面的这个交换方法得到的TT'' 中,x,yx,y 不是兄弟结点,即不会出现只有最后一位不等的情况)

于是:

WPL(T)WPL(T)=cCc.freq×dT(c)cCc.freq×dT(c)=(a.freqx.freq)[dT(a)dT(x)]0\begin{aligned}WPL(T)-WPL(T')&=\sum_{c\in C}c.freq\times d_T(c)-\sum_{c\in C}c.freq\times d_{T^*}(c)\\&=(a.freq-x.freq)[d_T(a)-d_T(x)]\\&\geq 0\end{aligned}

类似地,有WPL(T)WPL(T)0WPL(T')-WPL(T'')\geq 0,从而得到WPL(T)WPL(T)WPL(T)\geq WPL(T'').
说明,交换后得到的新树TT'' 的WPL 小于或等于 最优二叉树TT。而最优二叉树的WPL是理论最小的,所以只有二者相等。即交换之后的新树仍为最优二叉树,问题得证。

最优子结构

最优子结构性质 引理2

C=(C{x,y})z,z.freq=x.freq+y.freqC'=(C-\{x,y\})\cup {z},\quad z.freq=x.freq+y.freq,且TT'CC' 的最优二叉树,则将TT' 中的结点zz 替换为以x,yx,y 为左右子结点的内部节点从而生成的树TTCC 的最优二叉树。

引理2的证明

证明 根据题意,有cC{x,y},dT(c)=dT(c)\forall c\in C-\{x,y\}, d_T(c)=d_{T'}(c) ,并且还有dT(x)=dT(y)=dT(x)+1d_T(x)=d_T(y)=d_{T'}(x)+1.
于是就有:

WPL(T)=WPL(T)+x.freq+y.freqWPL(T)=WPL(T)x.freqy.freq\begin{aligned}WPL(T)&=WPL(T')+x.freq+y.freq\\WPL(T')&=WPL(T)-x.freq-y.freq\end{aligned}

利用反证法
如果TT 不是CC 的最优二叉树,则存在TT^* 使得WPL(T)<WPL(T)WPL(T^*)\lt WPL(T).
定理1 可知,TT^* 中存在频率最小的兄弟结点x,yx,y 。现将它们及其父结点替换为新的结点zz,并且令z.freq=x.freq+y.freqz.freq=x.freq+y.freq,得到新的树TT^{*'}
从而,有:

WPL(T)=WPL(T)x.freqy.freq<WPL(T)x.freqy.freq=WPL(T)\begin{aligned}WPL(T^{*'})&=WPL(T^*)-x.freq-y.freq\\&\lt WPL(T)-x.freq-y.freq=WPL(T')\end{aligned}

这说明,存在TT^{*'} 使得CC' 的 WPL 最小,与已知条件矛盾。


哈夫曼算法的正确性引理1引理2 可知,哈夫曼算法能够得到最优解.

编程实现

采用 C++ 语言进行编程,此处利用了 STL 标准库中的优先队列(头文件<queue>

因为哈夫曼树是二叉树的一种,因此考虑采用结构体的方式进行存储:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
typedef struct HT{
    int freq; //权重
    char key; //字符内容
    HT *left; //左子结点 
    HT *right; //右子结点
}HuffmanTree;

struct cmp{
    bool operator()(HuffmanTree *a, HuffmanTree *b){
        return a->freq > b->freq;//小的优先级高 
    }
}; //设置比较函数用于优先队列根据freq排序

typedef struct {
    char ch;
    int freq;
}CharSet; //需要编码的字符结构体

// 优先队列
typedef priority_queue<HuffmanTree*,vector<HuffmanTree*>,cmp> HTqueue;

然后是根据伪代码实现具体细则:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
void CreateHuffmanTree(HTqueue &Q, CharSet C[], int n){
    // 将字符集 C 导入 Q 并生成哈夫曼树
    // n 为字符个数
    for(int i = 0; i < n; i++){
        HuffmanTree *tmp = new HuffmanTree();
        tmp->key = C[i].ch;
        tmp->freq = C[i].freq;
        tmp->left = NULL;
        tmp->right = NULL;
        Q.push(tmp);
    }
    
    for(int i = 0; i < n-1; i++){
        HuffmanTree *tmp = new HuffmanTree();
        HuffmanTree *left = Q.top();  Q.pop();
        HuffmanTree *right = Q.top();  Q.pop();
        tmp->left = left;
        tmp->right = right;
        tmp->freq = left->freq + right->freq;
        tmp->key = '\0';
        Q.push(tmp);
    }
}

编写函数实现哈夫曼编码:

1
2
3
4
5
6
7
8
9
10
11
void OutputHuffmanCode(HuffmanTree *CodeTree, string code=""){
    //先序遍历,并且利用code存储编码
    if(CodeTree->key != '\0'){
            cout << CodeTree->key << ": (" << CodeTree->freq << ")";
            cout << code << endl;
    }
    if(CodeTree->left != NULL && CodeTree->right != NULL){
        OutputHuffmanCode(CodeTree->left,code+"0");
        OutputHuffmanCode(CodeTree->right,code+"1");
    }    
}

实例演示
考虑字符集<a:45,b:13,c:12,d:16,e:9,f:5><a:45,b:13,c:12,d:16,e:9,f:5>
手工绘制得到如下哈夫曼结果:
100a: 4555253014c: 12b: 13f: 5e: 9d: 16a:0b:101c:100d:111e:1101f:1100

下面通过调用算法查看.

主函数:

1
2
3
4
5
6
7
8
int main(void){
    HTqueue Q;
    CharSet C[6]={{'a',45},{'b',13},{'c',12},{'d',16},{'e',9},{'f',5}};
    CreateHuffmanTree(Q, C, 6);
    OutputHuffmanCode(Q.top());
    cout << "END" << endl;
    while(1);
}

输出结果:

1
2
3
4
5
6
7
a: (45)0
c: (12)100
b: (13)101
f: (5)1100
e: (9)1101
d: (16)111
END

Prefect!

最优二叉检索树|Optimal Binary Search Tree

哈夫曼树,也被称作最优二叉树,是通过建立一棵二叉树,对指定频率的字符进行编码,从而使得字符的平均编码长度达到最低,实现了数据的压缩。

而最优二叉检索树则是为了实现对于给定频率的字符,建立BST使得搜索效率最高。

例如想要实现对某英文文本中一个单词的快速查找并翻译,我们可以建立一棵 二叉搜索树(BST) 存储英文单词,而对于英文中出现频率很高的 the 这一单词,我们希望它在 BST 中更靠近根结点,这样可以很快检索到它。这种使得所有搜索过程中,访问结点个数最少(也就是搜索效率最高)的 BST,就被称为 最优二叉检索树(OBST)。

形式化表达

给定nn 个不同关键字(key)的已排序序列K=<k1,k2,...,kn>K=<k_1,k_2,...,k_n>,有k1<k2<<knk_1\lt k_2\lt\cdots\lt k_n.
每个关键词附有一个属性ki.pk_i.p 表示该关键字的搜索频率,即一个待搜索的关键字xxkik_i 的概率.
此外,还存在有n+1n+1 个伪关键字d0,d1,...,dnd_0,d_1,...,d_n 坐落在KK空隙中,i0,in\forall i\neq0,i\neq ndid_i 表示在(ki,ki+1)(k_i,k_{i+1}) 内的所有关键字,而d0:(,k1),  dn:(kn,+d_0:(-\infty,k_1),\;d_n:(k_n,+\infty).
伪关键字也有属性di.pd_i.p 表示它所表示的那个范围的搜索频率,即xxdid_i 表示范围内的概率.

自然有:

i=1nki.p+i=0ndi.p=1\sum_{i=1}^nk_i.p+\sum_{i=0}^nd_i.p=1

对上述关键字序列建立任意一棵 BST,显然代表空隙的did_i 始终作为叶子结点出现。一个示例如下图所示:

k2k1k4k3k5d0d1d2d3d4d5

通过 BST 搜索某个结点xx 时,需要访问的结点数就是dT(x)+1d_T(x)+1,我们就将访问结点数作为一个关键字的搜索代价。(其中,dT(x)d_T(x) 表示字符xx 在二叉搜索树TT 中的深度)

那么与WPLWPL 类似,一个最优二叉检索树的期望代价

E[Search cost in T]=i=1n[dT(ki)+1]×ki.p+i=0n[dT(di)+1]×di.p=1+i=1ndT(ki)×ki.p+i=0ndT(di)×di.p\begin{aligned} E[\text{Search cost in }T]&=\sum_{i=1}^n[d_T(k_i)+1]\times k_i.p+\sum_{i=0}^n[d_T(d_i)+1]\times d_i.p\\ &=1+\sum_{i=1}^nd_T(k_i)\times k_i.p+\sum_{i=0}^nd_T(d_i)\times d_i.p \end{aligned}

求解最优二叉检索树就是找到这样一棵二叉搜索树TT^* ,使得E(T)E(T^*) 最小.

T=argmaxTE(T)T^*=\arg\max_{T}E(T)

问题分析

构造一棵任意的 BST 是简单的。我们只需对KK 进行全排列,以排列后的序列第一个元素作为根结点,依照 BST 的构建方法逐步添加结点,直到所有的kik_i 加入树中,然后再填补空隙did_i 作为叶子结点即可。这样能得到所有可能的 BST(含重复),即便如此,搜索所有可能的解空间需要消耗的时间也十分巨大。

《算法导论》一书中指出,nn 个结点的二叉树(虽然并不是BST)数量在Ω(4n/n3/2)\Omega(4^n/n^{3/2}) ,所以暴力穷举求解是十分困难的,存在指数爆炸问题。因此下面将介绍一种使用动态规划求解 OBST 的算法。

状态转移方程

现在我们从后往前思考,对于一棵以krk_r 为根结点的 BST,其搜索代价必然包含krk_r 本身的代价,以及由krk_r 所划分的左右子树T1,  T2T_1,\;T_2 的搜索代价.

其中,T1T_1 是由全体小于krk_r 的关键字子序列K1:=K<k1..kr1>K_1:=K<k_1..k_{r-1}> 构成,当然还包括伪关键字(空隙)的子序列D1=D<d0..dr1>D_1=D<d_0..d_{r-1}>
T2T_2 是由全体大于krk_r 的关键字子序列K2:=K<k1..kr1>K_2:=K<k_1..k_{r-1}> 及其对应的空隙子序列D2D_2 构成。

krK1:<k1,k2,...,kr-1>K2:<kr+1,kr+2,...,kn>TT1T2

并且,以krk_r 作为根结点时的二叉树TT 比之其左右子树T2,  T2T_2,\;T_2 深度多了正好一个单位
从而,有:

E(T)=[E(T1)+w(1,k1)]+[E(T2)+w(k+1,n)]+kr.p=E(T1)+E(T2)+[w(1,k1)+w(k+1,n)+kr.p]=E(T1)+E(T2)+w(1,n)\begin{aligned} E(T)&=[E(T_1)+w(1,k-1)]+[E(T_2)+w(k+1,n)]+k_r.p\\ &=E(T_1)+E(T_2)+[w(1,k-1)+w(k+1,n)+k_r.p]\\ &=E(T_1)+E(T_2)+w(1,n) \end{aligned}

其中,w(1,k1)w(1,k-1) 是左子树深度增加一位之后,比起E(T1)E(T_1)多出来的搜索代价,右子树同理。不难给出对于任意新增代价w(i,j)w(i,j) 的表达式:

w(i,j)=l=ijkl.p+l=i1jdl.pw(i,j)=\sum_{l=i}^jk_l.p+\sum_{l=i-1}^jd_l.p

事实上,T1,  T2T_1,\;T_2 正是由我们需要解决的问题中划分出来的更小规模的子问题。我们只需要搜索KK 中,以哪一个结点作为根结点可以使得E(T)E(T) 最大,那么问题就可以解决。

更一般化地,如果需要从子序列<ki,ki+1,...,kj><k_i,k_{i+1},...,k_j> 中求最小期望代价,那么可以将序列的起始点和终止点作为 状态,在该状态下的最优值记作dp(i,j)dp(i,j) 则有:

dp(i,j)={minikj{dp(i,k1)+dp(k+1,j)+w(i,j)},jidi1.p,j=i1dp(i,j)=\begin{cases} \min\limits_{i\leq k\leq j}\{dp(i,k-1)+dp(k+1,j)+w(i,j)\},&j\geq i\\ d_{i-1}.p,&j=i-1 \end{cases}

该状态方程指出,当划分出空子树时,即j=i1j=i-1 时,其空隙的频率作为叶子结点出现,该频率就是唯一的搜索代价。
jij\geq i 时,就从iji\to j 中搜索出以哪个结点kk 为根结点时的代价最小。

值得注意的是,我们的转移方程虽然建立起来了,但是我们还需关注自底向上的数据是如何建立的。

例如,想要计算dp[1,5]dp[1,5],根据转移方程:
左边需要遍历dp[1,0],dp[1,1],d[1,2],dp[1,3],dp[1,4]dp[1,0],dp[1,1],d[1,2],dp[1,3],dp[1,4]
右边需要遍历dp[2,5],dp[3,5],dp[4,5],dp[5,5],dp[6,5]dp[2,5],dp[3,5],dp[4,5],dp[5,5],dp[6,5].

而这些都需要我们一开始计算并保存过。下面给出将动态规划表旋转后的结果,该旋转使得数组的对角线元素在水平线上。

654321543214321321211ij dpTable0516

图中的数字并不是实际结果,而是一个标记。可见只有在先算完图中标记为11 的数据之后,才能去计算标记为22 的数据,以此类推。
如果记ll 为上表中的 标记值-1,则:

  1. l=0l=0 时,i=1,2,3,4,5,6.  j=0,1,2,3,4,5i=1,2,3,4,5,6.\;j=0,1,2,3,4,5 (此时j=i1j=i-1 ,所以应提前处理)
  2. l=1l=1 时,i=1,2,3,4,5.  j=1,2,3,4,5i=1,2,3,4,5.\;j=1,2,3,4,5
  3. l=2l=2 时,i=1,2,3,4.  j=2,3,4,5i=1,2,3,4.\;j=2,3,4,5
  4. l=3l=3 时,i=1,2,3.  j=3,4,5i=1,2,3.\;j=3,4,5
  5. l=4l=4 时,i=1,2.  j=4,5i=1,2.\;j=4,5
  6. l=5l=5 时,i=1.  j=5i=1.\;j=5dp[1,5]dp[1,5] 即为所求)

根据以上自底向上的计算顺序,可总结出算法的循环体中应该有:
i:1nl+1,  j:il+1i:1\to n-l+1,\;j:i-l+1

伪代码与复杂度

在算法中,为了简单起见,避免ki.p,  di.pk_i.p,\;d_i.p 这种多余的符号,我们规定输入关键字序列KK 的个数nn,其频率p[1..n]p[1..n],空隙出现的频率q[1..n+1]q[1..n+1]

利用二维数组dp[1..n+1,0..n]dp[1..n+1,0..n] 作为动态规划表(即备忘录)。根据转移方程可知,只在ji1j\geq i-1 时,该数组的内容才会被用到。

dp[..]dp[..] 的行列数均是n+1n+1 是因为存在伪关键字单独作为结点的情况,如dp[n+1,n]dp[n+1,n] 表示以knk_n 为根结点时,其右子树为空树,但补上空隙dnd_n 作为子结点。同理,dp[1,0]dp[1,0] 表示k1k_1 为根结点时,左子树为空树时,补上d0d_0

dp[..]dp[..] 类似,我们利用二维数组w[1..n+1,0..n]w[1..n+1,0..n] 存储新增代价,避免后续的重复计算。新增代价的迭代更新也很简单:

w[i,j]={w[i,j1]+p[j]+q[j],jidi1.p,j=i1w[i,j]=\begin{cases} w[i,j-1]+p[j]+q[j],&j\geq i\\ d_{i-1}.p,&j=i-1 \end{cases}

为记录结果(追溯解),我们利用二维数组root[1..n,1..n]root[1..n,1..n] 存储在下标子域root[i,j]root[i,j] 时使得代价最小时选取的根结点rr. 显然,只有在1ijn1\leq i\leq j\leq n 时,该数组的内容才有意义。

  • 输入:关键字序列KK 的个数nn,其频率p[1..n]p[1..n],空隙出现的频率q[1..n+1]q[1..n+1]
  • 输出:OBST的最小代价,根结点的选取表root[1..n,1..n]root[1..n,1..n]

Algorithm:   Optimal-BST(n,p,q)1.  初始化二维数组 dp,  w,  root2.  for  i  =1  to  n+1  do3.  dp[i,i1]q[i1]4.  w[i,i1]q[i1]5.  for  l  =1  to  n  do6.  for  i  =1  to  nl+1  do7.  ji+l18.  dp[i,j]+9.  w[i,j]w[i,j1]+p[j]+q[j]10.  for  r  =i  to  j  do11.  tmpdp[i,r1]+dp[r+1,j]+w[i,j]12.  if  dp[i,j]>tmp  then13.  dp[i,j]tmp14.  root[i,j]r15.  return  dp[1,n]  and  root\begin{aligned} &\text{Algorithm: }\;\text{Optimal-BST}(n,p,q)\\\\ 1.&\;\text{初始化二维数组 }dp,\;w,\;root\\ 2.&\;\mathbf{for}\;i\;=1\;\mathbf{to}\;n+1\;\mathbf{do}\\ 3.&\;\qquad dp[i,i-1]\leftarrow q[i-1]\\ 4.&\;\qquad w[i,i-1]\leftarrow q[i-1]\\ 5.&\;\mathbf{for}\;l\;=1\;\mathbf{to}\;n\;\mathbf{do}\\ 6.&\;\qquad\mathbf{for}\;i\;=1\;\mathbf{to}\;n-l+1\;\mathbf{do}\\ 7.&\;\qquad\qquad j\leftarrow i+l-1\\ 8.&\;\qquad\qquad dp[i,j]\leftarrow +\infty\\ 9.&\;\qquad\qquad w[i,j]\leftarrow w[i,j-1]+p[j]+q[j]\\ 10.&\;\qquad\qquad\mathbf{for}\;r\;=i\;\mathbf{to}\;j\;\mathbf{do}\\ 11.&\;\qquad\qquad\qquad tmp\leftarrow dp[i,r-1]+dp[r+1,j]+w[i,j]\\ 12.&\;\qquad\qquad\qquad\mathbf{if}\;dp[i,j]\gt tmp\;\mathbf{then}\\ 13.&\;\qquad\qquad\qquad\qquad dp[i,j]\leftarrow tmp\\ 14.&\;\qquad\qquad\qquad\qquad root[i,j]\leftarrow r\\ 15.&\;\mathbf{return}\;dp[1,n]\;\mathbf{and}\;root \end{aligned}

显然,算法嵌套了三层循环,每层执行的操作在O(1)O(1),因此该算法的时间复杂度为Θ(n3)\Theta(n^3),而利用了3个近似为n×nn\times n 的的二维数组作为备忘和结果,所以空间复杂度为O(n2)O(n^2).

实例与编程

假设给定概率分布如下:

ii012345
p[i]p[i]0.150.100.050.100.20
q[i]q[i]0.050.100.050.050.050.10

下给出 c 系列的编程实现:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
#include <iostream>
#include <algorithm>
#include <string.h>
#include <vector>
#include <queue>
#include <cstdlib>
using namespace std;
#define M 100

double dp[M][M];
double w[M][M];
int root[M][M];

double Optimal_BST(double p[], double q[], int n){
// 初始化
memset(dp,0,sizeof(dp));
memset(root,0,sizeof(root));
memset(w,0,sizeof(w));

// 单独处理最底层
for(int i = 1; i <= n+1; i++){
dp[i][i-1] = q[i-1];
w[i][i-1] = q[i-1];
}
// 从 l=1 开始逐层填表
for(int l = 1; l <= n; l++){
for(int i = 1; i <= n-l+1; i++){
int j = i+l-1;
dp[i][j] = INT32_MAX;
w[i][j] = w[i][j-1]+p[j]+q[j];
for(int r = i; r <= j; r++){
double tmp = dp[i][r-1]+dp[r+1][j]+w[i][j];
if(dp[i][j] > tmp){
dp[i][j] = tmp;
root[i][j] = r;
}
}
}
}
return dp[1][n];
}
int main(void){
double p[6] = {0, 0.15, 0.10, 0.05, 0.10, 0.20};
double q[6] = {0.05, 0.10, 0.05, 0.05, 0.05, 0.10};

cout << "min E = " << Optimal_BST(p, q, 5) << endl << endl;
for(int l = 0; l <= 5; l++){
for(int i = 1; i <= 5-l+1; i++){
int j = i+l-1;
cout << "(" << i << "," << j << "): dp = "; printf("%.2f\t",dp[i][j]);
cout << "w = "; printf("%.2f\t",w[i][j]);
cout << "root = " << root[i][j] << endl;
}
}
return 0;
//while(1);
}

得到如下旋转过后的各数组结果:

实例的求解结果

从而得到对于该问题的最小代价为2.752.75.

最小生成树|Minimum Spanning Tree,MST

对于一个 无向带权连通图G=(V,E,W)G=(V,E,W) ,我们需要寻找GG 的一个无环子集TT 使得TGT\subseteq GTT 包含GG 的所有顶点的同时又具有最小权重,即w(T)=(u,v)Tw(u,v)w(T)=\sum\limits_{(u,v)\in T}w(u,v) 的值最小。

TT 无环且连通所有顶点可知,TT 必然是一棵树,我们就称TT 是 图GG最小生成树Minimum Spanning Tree,MST。

生成树的性质

生成树具有如下性质
G(V,E,W)G(V,E,W)nn 阶连通图,即V=n|V|=n,那么:

  1. TTGG 的生成树 当且仅当TT 无环且有n1n-1 条边;
  2. 若某条边eE(G)e\in E(G),且GG 的一个生成树TT 中不含有ee,即eE(T)e\notin E(T),那么E(T){e}E(T)\cup \{e\} 必含有一个回路,记为CC;(注:用E(X)E(X) 表示某个图XX边集
  3. 紧接第二点,如果去掉回路CC 中任意一条边e,eee',\quad e'\neq e,最终得到的TT' 仍是生成树。

如下图所示:

生成树的性质

根据该性质,我们可以初步得出一种找到 MST 的策略:
首先得到一棵生成树TT ,然后对其加一条非树边ee 以形成回路CC.
如果能在CC 中去掉一条边ee',其中有w(e)<w(e)w(e')\lt w(e),那么得到的新树TT' 就有:w(T)w(T)=w(e)w(e)<0w(T')-w(T)=w(e')-w(e)\lt0,即是w(T)<w(T)w(T')\lt w(T).

从而可以根据该策略不断改进,最终得到 MST.

于是,求出一棵生成树的问题也就可以按照如下的通用算法进行描述:

Algorithm:   Generic-MST(G,E,w)1.  A2.  while  A does not from a spanning tree  do3.  find an edge (u,v) that is safe for A4.  AA{(u,v)}5.  return  A\begin{aligned} &\text{Algorithm: }\;\text{Generic-MST}(G,E,w)\\\\ 1.&\;A\leftarrow\varnothing\\ 2.&\;\mathbf{while}\;A\text{ does not from a spanning tree}\;\mathbf{do}\\ 3.&\;\qquad \text{find an edge }(u,v)\text{ that is safe for }A\\ 4.&\;\qquad A\leftarrow A\cup\{(u,v)\}\\ 5.&\;\mathbf{return}\;A \end{aligned}

此处的AA 即是 MST 的边集.
此处的safesafe 表示的是这样一条边(u,v)(u,v) 使得A{(u,v)}A\cup\{(u,v)\} 仍然是 MST 边的子集,把这样的边称为 安全边.

下面将具体介绍两种算法用于求解最小生成树,并且二者均是利用了贪心策略的思想,但具体如何找到安全边的方法则不一样。

Prim算法

Prim 算法的核心思想是将图的顶点VV 划分为两个子集S,VSS,V-S.
如果边e=(u,v)e=(u,v) 中,uS,vVSu\in S,v\in V-S 则称ee 横跨VV 的一个切割(S,VS)(S,V-S).

而 Prim算法 要求不断找到横跨某个切割的最小权重边,这样的边称为 轻量级边。可以证明这样的轻量级边就是一条安全边

Prim算法 将当前的轻量级边划入子集AA ,对应的端点加入SS 中,从而构造出下一个切割,再寻找下一条轻量级边,直到所有结点都划入SS ,最终的结果即为 MST.

说人话就是,初始化最终 MST 的点集和边集为空。
先选定一个点,找与该点距离最近的边,然后将此边的连接的另一个端点加入点集,此边加入边集。
接下来,又继续找距离当前这个点集最近的边,依次加入,直到所有顶点都加入点集。

算法伪码

  • 输入:连通图GG,边权ww,最小生成树的根结点rr
  • 输出:最小生成树的边集AA
  • 其他:
    • v.keyv.key 表示树中与vv 相连的边权重,若不存在则赋为\infty
    • v.πv.\pi 表示vv 在树中的父结点,有:e=(v,v.π)e=(v,v.\pi) 的权值是v.keyv.key
    • G.VG.V 表示图GG 的点集(数学中表示为V(G)V(G)
    • G.Adj[v]G.Adj[v] 表示图GG 中与结点vv 相邻的结点的集合(数学中表示为NG(v)N_G(v)

Algorithm:   MST-Prim(G,w,r)//r是起始点,最终作为树的根结点1.  for  each uG.V  do2.    u.key+3.    u.πNIL4.  r.key05.  QG.V  //Q为按权重 V.key 递增的优先队列6.  A7.  while  Q  do8.  uExtract-min(Q)9.  AA{(u,u.π)}10.  for  each vG.Adj[u]  do11.  if  vQ and w(u,v)<v.key  then12.  v.πu13.  v.keyw(u,v)14.  return  A\begin{aligned} &\text{Algorithm: }\;\text{MST-Prim}(G,w,r)\\ &//r\text{是起始点,最终作为树的根结点}\\\\ 1.&\;\mathbf{for}\;\text{each }u\in G.V\;\mathbf{do}\\ 2.&\;\qquad\;u.key\leftarrow+\infty\\ 3.&\;\qquad\;u.\pi\leftarrow NIL\\ 4.&\;r.key\leftarrow0\\ 5.&\;Q\leftarrow G.V\;//Q\text{为按权重 }V.key\text{ 递增的优先队列}\\ 6.&\;A\leftarrow\varnothing\\ 7.&\;\mathbf{while}\;Q\neq\varnothing\;\mathbf{do}\\ 8.&\;\qquad u\leftarrow\text{Extract-min}(Q)\\ 9.&\;\qquad A\leftarrow A\cup\{(u,u.\pi)\}\\ 10.&\;\qquad\mathbf{for}\;\text{each }v\in G.Adj[u]\;\mathbf{do}\\ 11.&\;\qquad\qquad\mathbf{if}\;v\in Q\text{ and }w(u,v)\lt v.key\;\mathbf{then}\\ 12.&\;\qquad\qquad\qquad v.\pi\leftarrow u\\ 13.&\;\qquad\qquad\qquad v.key\leftarrow w(u,v)\\ 14.&\;\mathbf{return}\;A \end{aligned}

其中,第10到13行是对队列中的结点的属性进行更新,使得每次循环都能保证队列中的结点vv 到MST子集AA 的权重v.keyv.key 最小,此时保存好在AA 中与vv 相连的结点,即是其父结点v.πv.\pi

如果优先队列QQ 采用 最小二叉堆 实现,则每一次循环取出当前最小元的过程(利用最小二叉堆实现时)需要花费O(logV)O(\log |V|)while 循环次数O(V)O(|V|)for 循环总次数O(E)O(|E|),事先装载队列需要O(V)O(|V|).

综合来看,该算法的时间复杂度为O(VlogV+ElogV)=O(ElogV)O(|V|\log |V|+|E|\log |V|)=O(|E|\log |V|),空间复杂度O(V)O(|V|).

如果改用斐波那契堆,则时间复杂度还可以进一步优化为O(E+VlogV)O(|E|+|V|\log|V|).

正确性证明

提出命题k<n=V\forall k\lt n=|V|,存在一棵 MST 包含有算法前kk 步选择边而纳入的子集AA

  1. k=1k=1 时,假设最小生成树TT^* 不包含由选择的根结点rr 组成的某条轻量级边(r,i)(r,i),即(r,i)T(r,i)\notin T^*.
    T=(T{(r,j)}){(r,i)}T=(T^*-\{(r,j)\})\cup\{(r,i)\} .
    因为(r,i)(r,i) 是一条轻量级边,即w(r,i)=min(r,t)Ew(r,t)w(r,j)w(r,i)=\min\limits_{(r,t)\in E}w(r,t)\leq w(r,j),所以替换后得到的生成树仍然是最优解

  2. 假设命题成立,证明命题对k+1k+1 的情况也成立。

设算法当前得到的子集AA 中的顶点集为SS,那么第k+1k+1 步将会从集合VSV-S 中选取下一个点ik+1i_{k+1} 纳入SS,且保证ik+1i_{k+1}SS 有最小权边。如图所示:

12S45iV-Sk+1ilee'

假设该最小权边记为e=(ik+1,il)e=(i_{k+1},i_l),则有:eTe\in T,若不然,则将ee 加入到TT^* 中形成回路,这条回路中一定包含有横跨(S,VS)(S,V-S) 的另一条边ee',将其从回路中删除得到新的生成树:

T=(T{e}){e}T=(T^*-\{e'\})\cup\{e\}

因为ee 是最小权边,所以w(e)w(e),w(T)w(T)w(e)\leq w(e'),w(T)\leq w(T^*),所以TT 仍然是最小生成树,即算法的第k+1k+1 步仍然能够得到MST的子集。

当算法执行完毕后,结果A=TA=T 是图GG 的一个最小生成树。

编程实现

采用 C++ 进行编程实现.

首先创建给出结果的边集结构体(Edge),以及由 NodeEdge 数组组合而成的图(Graph).

此外还有为Prim算法设计的优先队列,因为C++中没有特定的Decrease-key() 方法,因此直接采用 vector +遍历的方法模拟优先队列。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
typedef struct NODE{
    int tag;
    int key; //最小权重
    int pi; //父结点
}Node;

typedef struct EDGE{
    int u,v;
}Edge;

typedef struct GRAPH{
    int num; //图结点个数
    vector<int> nodes;
    vector<Edge> edges;
}Graph;

然后是具体的Prim算法实现:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
Node PopMinKey(vector<Node> &Q){
//找到key最小的结点从队列中删除,并返回该结点
    int min = INT32_MAX;
    Node u;
    int index;
    for(int i = 0; i < Q.size(); i++){
        if(min > Q[i].key){
            min = Q[i].key;
            u = Q[i];
            index = i;
        }
    }
    Q.erase(Q.begin()+index);
    return u;
}


void UpdateQueue(vector<Node> &Q, int target, int sourse, int weight){
//对队列 Q 中的结点 key值 进行更新
    for(int i = 0; i < Q.size(); i++){
        if(Q[i].tag == target){
            Q[i].key = weight;
            Q[i].pi = sourse;
            break;
        }
    }
}


vector<Edge> MST_Prim(Graph G, int w[], int r = 0){

    vector<Node> Q;
    for(int i = 0; i < G.num; i++){
        Q.push_back({G.nodes[i], INT32_MAX, -1});
    }
    Q[0].key = 0; Q[0].pi = 1; //初始化队列

    vector<Edge> A;

    while(!Q.empty()){
        Node u = PopMinKey(Q);
        A.push_back({u.tag, u.pi});
        for(int i = 0; i < G.edges.size(); i++){
            if(G.edges[i].u == u.tag){
                UpdateQueue(Q, G.edges[i].v, u.tag, w[i]);
            }
            if(G.edges[i].v == u.tag){
                UpdateQueue(Q, G.edges[i].u, u.tag, w[i]);
            }
        }
    }
    return A;
}

实例演示
考虑如下图所示的实例,其正确答案已给出。

31245663656421553124566365642155312456636564215531245663656421312456365642131245635421

编写主函数:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
int main(void){
    Graph G;
    int u,v;
    G.num = 6;
    for(int i = 0; i < 6; i++){
        G.nodes.push_back(i+1);
    }
    for(int i = 0; i < 10; i++){
        cin >> u >> v;
        G.edges.push_back({u,v});
    }
    int w[10]= {6,1,5,5,3,5,6,4,2,6};
    vector<Edge> MST = MST_Prim(G,w);
    for(int i = 1; i < MST.size(); i++){
        cout << "(" << MST[i].u << ", " << MST[i].v << ")" << endl;
    }
    while(1);
}

结果如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
输入(边集):
1 2
1 3
1 4
2 3
2 5
3 4
3 5
3 6
4 6
5 6

输出(边集):
(3, 1)
(6, 3)
(4, 6)
(2, 3)
(5, 2)

Prefect!

Kruskal算法

Kruskal算法的核心是与 Prim算法 不同。
Kruskal算法并不要求每一步所选最小权边都是相互连通的。换言之,Kruskal算法中第kk 步产生的子集AA 不一定是一棵子树,而有可能是森林。当然,即使是森林也不可能出现环。

不难发现,Prim算法因为每一步仅选取一个点及其轻量级边纳入AA,所以可以保证AA 中的边一定不构成回路,且是一棵树,而Kruskal算法则需要有一个判断是否构成回路的策略。这个策略可以利用并查集这种数据结构解决。

综上,Kruskal算法的思想就是,先对边集进行升序排列,然后考察当前最小边纳入AA 之后是否构成回路,如果构成回路则跳过并寻找下一条边。

算法伪码

  • 输入:连通图GG,边权ww
  • 输出:最小生成树的边集AA

Algorithm:   MST-Kruskal(G,w)1.  for  each uG.V  do2.    Make-Set(u)3.  QG.E  //Q为按权重 w(e) 递增的优先队列4.  A5.  repeat6.  (u,v)Extract-min(Q)7.  if  Find-Set(u)Find-Set(u)  then8.  AA{(u,u.π)}9.  Union(u,v)10.  until  A=V111.  return  A\begin{aligned} &\text{Algorithm: }\;\text{MST-Kruskal}(G,w)\\\\ 1.&\;\mathbf{for}\;\text{each }u\in G.V\;\mathbf{do}\\ 2.&\;\qquad\;\text{Make-Set}(u)\\ 3.&\;Q\leftarrow G.E\;//Q\text{为按权重 }w(e)\text{ 递增的优先队列}\\ 4.&\;A\leftarrow\varnothing\\ 5.&\;\mathbf{repeat}\\ 6.&\;\qquad (u,v)\leftarrow\text{Extract-min}(Q)\\ 7.&\;\qquad\mathbf{if}\;\text{Find-Set}(u)\neq\text{Find-Set}(u)\;\mathbf{then}\\ 8.&\;\qquad\qquad A\leftarrow A\cup\{(u,u.\pi)\}\\ 9.&\;\qquad\qquad\text{Union}(u,v)\\ 10.&\;\mathbf{until}\;|A|=|V|-1\\ 11.&\;\mathbf{return}\;A \end{aligned}

其中:

  • Make-Set(u)\text{Make-Set}(u) 用于创建集合,第 1到2 行用于为每一个结点单独创建一个集合(树)
  • Find-Set(u)\text{Find-Set}(u) 返回用于表示元素uu 所在的集合的元素,第 7 行 用于判断加入(u,v)(u,v) 后是否构成环,即合并两颗树后是否仍是树
  • Union(u,v)\text{Union}(u,v) 用于将u,vu,v 所在集合进行合并

分析得知,一共循环O(E)O(|E|) 次,每个循环下,生成集合、取出最小元、合并集合都在O(logE)O(\log |E|) 下完成。因此该算法的时间复杂度为O(ElogE)O(|E|\log |E|),如果注意到EV2|E|\leq|V|^2,则还可以表示为O(ElogV)O(|E|\log|V|)
空间复杂度为O(V)O(|V|)

正确性证明

利用数学归纳法.
证明 对于任意结点数为V=n|V|=n 的无向图GG,贪心算法都能得到一棵最小生成树.

  1. n=2n=2 时,算法能够直接从连接图GG 的两个端点中找到一条最小权边,显然能够得到最小生成树;
  2. 假设 对于任意nn 阶图,算法都能得到MST,下面说明对于n+1n+1 阶图也能得出.

对于任意n+1n+1 阶图GG ,将其最小权边e=(i,j)e=(i,j) 短接,如下图所示,将图化为nn 阶图GG'。根据归纳假设,算法能够给出GG' 的一棵 MST 记为TT'

132456短接13254-6GG'13254-6T'132456恢复TMST

T=T{e}T=T'\cup\{e\} 就是GG 的一棵 MST,若不然,Te\exists T^*\ni eGG 的 MST.

注:就算eTe\notin T^* 也可以通过加边再减边的方式得到包含ee 的 MST.
(这部分更详细的解释已在Prim算法的正确性证明时给出)

此时,有w(T)<w(T)w(T{e})<w(T{e})=w(T)w(T*)\lt w(T)\Rightarrow w(T^*-\{e\})\lt w(T-\{e\})=w(T')
这与TT'GG' 的最小生成树矛盾.

综上所述,Kruskal 算法能够给出任意nn 阶图的MST.

编程实现

采用 C++ 进行编程实现.

首先声明必要的结构体以及以边权大小排序的优先队列:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
typedef struct EDGE{
    int u,v;
    int w;
}Edge;

typedef struct GRAPH{
    int num; //图结点个数
    vector<int> nodes
    vector<Edge> edges;
}Graph;

struct cmp{
    bool operator()(Edge a, Edge b){
        return a.w > b.w; //小的优先级高 
    }
};

typedef priority_queue<Edge,vector<Edge>,cmp> EdgeQue;

利用 vector 自建并查集及其相关的合并查找函数:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
typedef struct SET{
    int root;
    int value;
}Set;

vector<Set> MakeSet(vector<int> a){
    vector<Set> set(a.size());
    for(int i = 0; i < a.size(); i++){
        set[i].root = a[i];
        set[i].value = a[i];
    }
    return set;
}

void SetUnion(vector<Set> &set, int a, int b){
    int i = max(set[a].root, set[b].root);
    int j = min(set[a].root, set[b].root);
    for(int k = 0; k < set.size(); k++){
        if(set[k].root == i){
            set[k].root = j;
        }
    }
}

最后是 Kruskal算法的实现:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
vector<Edge> MST_Kruskal(Graph G){
    vector<Edge> A;
    vector<Set> set = MakeSet(G.nodes);
    EdgeQue Q;
    for(int i = 0; i < G.edges.size(); i++){
        Q.push({G.edges[i].u, G.edges[i].v, G.edges[i].w});
    }
    do{
        Edge e = Q.top(); Q.pop();
        if(set[e.u].root != set[e.v].root){
            A.push_back({e.u, e.v});
            SetUnion(set, e.u, e.v);
        }
    } while (A.size() < G.num-1);

    return A;
}

利用和 Prim进行演示时相同的实例,编写主函数:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
int main(void){
    Graph G;
    int u,v;
    G.num = 6;
    int w[10]= {6,1,5,5,3,5,6,4,2,6};
    for(int i = 0; i < 6; i++){
        G.nodes.push_back(i+1);
    }
    for(int i = 0; i < 10; i++){
        cin >> u >> v;
        G.edges.push_back({u,v,w[i]});
    }
    vector<Edge> MST = MST_Kruskal(G);
    for(int i = 0; i < MST.size(); i++){
        cout << "(" << MST[i].u << ", " << MST[i].v << ")" << endl;
    }
    while(1);
}

结果如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
输入(边集):
1 2
1 3
1 4
2 3
2 5
3 4
3 5
3 6
4 6
5 6

输出(边集):
(1, 3)
(4, 6)
(2, 5)
(3, 6)
(2, 3)

Prefect!

参考

  1. 算法导论(原书第三版)》
  2. 算法设计与分析|中国大学MOOC
  3. 哈夫曼树(赫夫曼树、最优树)详解
  4. C++优先队列自定义排序总结|CSDN