yunlongwang的个人博客分享 http://blog.sciencenet.cn/u/yunlongwang

博文

数学模型也可以像玩具一样好玩(记一个漂亮的马尔科夫链)

已有 8716 次阅读 2014-3-27 13:13 |系统分类:科普集锦|关键词:学者| 马尔科夫链, 数学模型。

上周,楼主开始尝试着写博文。承蒙科学网上各位老师的关照和指导,楼主介绍贝叶斯定理的科普文得到了好多的点击量。这对我真是一个莫大的鼓励,非常感谢。我今天就冒昧的再写一篇关于马尔科夫链的小文章。


楼主本学期平均每周都要带一次概率论习题课,一次模电实验。呵呵,两门课,系里穷呀,拿一个研究生当两个用,还美其名曰“给你的是两个‘半助教’”。 言归正传,每次上课呢,我似乎都要强调一些老生常谈,比如说拜托來以前把电路搭好了,上课请不要迟到,拜托仔细看电路图,别忘了接电源线哦,亲,巴拉巴拉巴拉。。。

 

我就想啊,每次我强调一次,我就影响的几个人,这几个人下次可能就不会犯低级错误;我再强调一次,又多几个人。那么我需要强调多少次,才能让每一个同学都知道我想表达的意思呢?

 

类似的问题在我们生活中还是能见到的,比如啊,我想写博文,做科普,我每写一篇博文,每个用户都有一定的几率接触到“贝叶斯”这个概念, 那我需要平均需要写多少篇博文才能让大家都知道“贝叶斯”?生物我不太懂,但是我揣测着,有没有可能每次喂点药,培养皿里每个细胞就以一定概率变性,那平均喂多少次药这个培养皿就只有变性的细胞了呢。昨天和隔壁实验室的沈喆兄讨论射频ID的问题:  每次我把带着扫描装置小车推过货架,每个物品都有一定几率被扫描到,那平均需要多少次才能得到每个货物的信息?

 

总之,就是有限数量个节点,每执行一次操作,每个节点都有一定几率发生不可逆的变化,那么平均需要多少次操作才能完成整个系统的变化?有一个类似但有本质区别的问题是:在执行M次操作之后,发生变化的节点数目是多少?由于第二个问题仅仅是几何分布加二项式分布,我就不写了。

 

好啦,背景介绍结束,楼主开始建模了。


第一步,我们做一点假设,假设有N=3各节点,每个节点都是“不健康的”。每次操作之后,每个节点从“非健康节点”变成“健康节点”或“非健康节点”的概率分别是p1-p,而“健康节点”不可能变成“非健康节点”。假设各个节点的变化都是独立的。

 

第二步,定义随机变量X为系统的状态,表示“健康节点”在系统中的数目, $X\in\{0,1,2,3\}$ 。那么我们可以把该系统用下图中的马尔科夫链【1】【2】来表示:


 

Fig.1 描述本模型的马尔科夫链

 

其中 $p_{ij}$ 表示从任意t时刻到t+1时刻的过程中,系统从i状态转换到j状态的转移概率。比如 $p_{01}$ 表示再执行此操作之前,系统有0个“健康节点”,而操作之后,系统变成有1个“健康节点”的概率,注意哦,这是仅仅一次操作。此处 $p_{01}=C_3^1p^1(1-p)^2$ ,不难推出,当 $i\leq j$ ,所有 $p_{ij}" style="font-family:宋体;font-size:14.44444465637207px;$ 都服从二项式分布;当 $i>j$ ,所有 $p_{ij}=0$ 。为了避免数学公式给阅读带来的不快,我就不写通式了。您就知道转移概率都很容易求就行了。最后,定义矩阵 $\mathbf{P}$ 为该系统的转移概率矩阵。

 

第三步,我们想要知道的,是从0状态进化到3状态的平均时间。因为3状态是一个吸收状态(就是一旦到达3状态,就不可能再变到别的状态了)。上谷歌搜absorbing Markov chain。可以直接得到这类问题的解法,可惜是英文的【4】。为了方便参加数学建模比赛的大学生朋友, 楼主这就不是很严谨的给出平均时间的解法。定义 $t_{i3}" style="font-family:宋体;font-size:14.44444465637207px;$ 为从i状态到3状态所花费的时间,定义 $\mathbb{E}t_{i3}$ $t_{i3}" style="font-family:宋体;font-size:14.44444465637207px;$ 的数学期望,则我们有以下方程组(感谢S.C.同学):

$\begin{align} \mathbb{E}t_{23}=&1+p_{22}\mathbb{E}t_{23}\nonumber\\ \mathbb{E}t_{13}=&1+p_{12}\mathbb{E}t_{23}+p_{11}\mathbb{E}t_{13}\nonumber\\ \mathbb{E}t_{03}=&1+p_{02}\mathbb{E}t_{23}+p_{01}\mathbb{E}t_{13}+p_{00}\mathbb{E}t_{03}\nonumber \end{align}$


让我解释一下这几个式子,比如第一个,就表示因为从2状态只能到达3状态或者2状态,到达3状态就不需要再花费时间,到达2状态就还需要 $\mathbb{E}t_{23}$ 这么多的平均时间。第二个式子类似,1状态可以到达1状态,2状态,或者3状态。好啦好啦,有结果啦,三个方程三个未知数,可以解出来我们需要的 $\mathbb{E}t_{03}$ 。 用矩阵来表示:

$\begin{align} \mathbb{E}\mathbf{t}=&\mathbf{1}+\mathbf{P}^{(sub)}\mathbb{E}\mathbf{t}\nonumber\\ \mathbb{E}\mathbf{t}=&(\mathbf{I}-\mathbf{P}^{(sub)})^{-1}\mathbf{1}\nonumber \end{align}$


其中 $\mathbf{1}$ 是一个全部元素为1的列向量, $\mathbf{I}$ 3乘以3的单位矩阵, $\mathbf{P}^{(sub)}$ 是刨掉 $\mathbf{P}$ 矩阵最后一行和最后一列后的分块矩阵, $\mathbf{t}=[t_{03},t_{13},t_{23}]^{\top}$ 是由 $t_{i3}" style="font-family:宋体;font-size:14.44444465637207px;$ 构成的列向量。我为什么要把方程总结成矩阵形式呢,因为本例中N=3,如果有更多的节点,那么用上面的那个矩阵,就可以很方便的拓展啦。

 

辛苦辛苦,您读到这里想必已经有点累了吧,我这就进行总结。

 

总结:1. 很遗憾,我想不出有什么简单的办法来求 $t_{i3}" style="font-family:宋体;font-size:14.44444465637207px;$ 的分布。用递归的方法,我想是可以得到解析形式的。但考虑到解法的复杂性,其已经超出了科普的范畴。有兴趣的老师请帮我想一想,要是有结果,请您给我留言告诉我。

 

2. 现实问题远远比数学模型复杂,楼主在这里加入了很多隐含的假设,以至于数学模型如此简单。但是我理解的数学建模是一个抓主要矛盾的过程,有些细枝末节不妨用假设将其忽略,否则模型有时候过于复杂反而得不到有意义的结果。当然,具体问题具体分析。

 

3. 不一定对啊,我感觉建模有时候还是要建成自己能理解能解决的问题。本例中,楼主用马尔科夫链,还算是比较熟悉吧。

 

4. 我觉得这个马尔科夫链比较有意思的是最后解平均时间的几个方程,虽然在熟悉随机过程的朋友看来,根本不值一提,我还是斗胆把它展示出来。因为我想有些大学生朋友可能会用得到,比如参加数学建模比赛的电子系学生。

 

5. 我没有完全分析这个问题,至少还可以发散出类似“平均经过多少时间整个系统可以有一半变成健康的?”这种问题。


6. 有时候觉得贝叶斯定理,马尔科夫链这些小模型就跟玩具一样,似乎生活中随处可见呀。

 

后记: 我本科的时候参加过数学建模比赛,当时深感中文的科普资料不多,百度百科的专业词条有时候并不准确。比如均方误差现在还是错的3】,我很认真的给它改过,编辑没让我通过,说我夹杂广告。而我实际上只是有一条参考文献是一本书。我后来写电邮和他们争辩,他们不理我,于是我不开心,于是我再也不写百度百科了。

 

后来我上研究生,我看wiki上有些词条日文的解释都比中文的多,而我们有十倍于他们的人口。想象我是一名大一的化学专业学生,我想了解“马尔科夫链”,似乎除了看书,我很难在网络上学到这一概念,因为写科普的老师不多。但是,我看书看不懂。当然,我英语要是够好,我可以阅读老美写的很多很多的类似博文的“brief introduction on 什么什么”。

 

我要非常非常感谢S.C.同学,我从上初一就开始问你数学题,真是太谢谢你了。

 

1】马小洪 老师的博客:http://blog.sciencenet.cn/blog-560320-714691.html

2】百度百科词条:http://baike.baidu.com/view/340221.htm

3wiki上关于均方误差的词条:

http://zh.wikipedia.org/wiki/%E5%9D%87%E6%96%B9%E5%B7%AE 

【4】 http://en.wikipedia.org/wiki/Absorbing_Markov_chain




https://m.sciencenet.cn/blog-624263-779648.html

上一篇:在国外大学监考有感
下一篇:大学老师组织高中课外活动

10 杨正瓴 徐传胜 陆泽橼 蒋迅 王春艳 李天成 强涛 zjz631026 ccgoodluck rosejump

该博文允许注册用户评论 请点击登录 评论 (15 个评论)

数据加载中...
扫一扫,分享此博文

Archiver|手机版|科学网 ( 京ICP备07017567号-12 )

GMT+8, 2024-5-19 23:48

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部