马尔可夫模型

1、从状态到马尔可夫链

利用计算机科学中有限状态机(FSA)模型，对于线性组合罚分系统给出一个解，具体说，将序列比对描述为三个不同状态之间不断转换的过程。并定义如下，作为动态规划求解的比对函数。

1. M(i,j) 表示在 Xi 比对到 Yj，也就是两个残基对在一起的时候，第一条序列X从第1位到第 i 位、第二条序列Y从第1位到第j位最好的比对分数。
2. X(i,j) 和 Y(i,j) 则分别表示在 Xi 或 Yj残基比对到空位时，序列X从第1位到第 i 位、序列Y从第1位到第 j 位最好的比对分数。

可以用动态规划矩阵来将三组转换关系表示位三个平面之间的单元格填充，只是对于每一个平面而言，回溯关系即可能来自于本平面，也可能来自于其他平面，或者说，另外一个状态。

马尔可夫链： 一个基于概率的随机过程模型，用来刻画一组之间存在关联的随机事件。具体说，用来描述一组离散状态之间在不同时刻的转移关系，值得注意的是，这里的状态转换关系不需要是唯一确定的，只需要可以由一个概率分布描述即可。唯一的要求是，t时刻状态的概率分布，由且只由之前有限的m个时刻状态的概率分布确定，称之为m阶马尔可夫链。事实上，我们通常可以考虑其最简化的情形：1阶马尔可夫链，也就是当前的状态与且只与其前一个状态有关。

事实上我们可以第任意比对计算出概率，也就说，我们通过引入马尔可夫链，给出序列比对的概率解释。

2、隐马尔可夫模型

很容易发现，上述的马尔可夫模型不足以真正完成序列比对。因为现有的状态模型只是区分了空位状态X，Y以及Match状态M。而没有考虑具体的残基。 因此我们需要进一步引入隐马尔可夫模型。所谓隐马尔可夫模型，是指在状态的基础上，增加了符号(Token)的概念。每个状态都可以以不同的概率产生一组可以观察到的符号。也就是说，除了状态转移概率之外，隐马尔可夫模型进一步引入“生成概率”的概念，每个状态都有自己的生成概率分布，可以按照不同的概率产生一组可以被观测到的符号。与马尔可夫模型不同的是，状态路径是看不到的。 这也是“隐”的含义。相反我们需要根据观测到的符号，来推测对应的状态。
例如，当我们观测到的字符串是aabc时，对应的状态路径可能是1,1,2,3，也可能是1,2,3,3、1,3,3,3甚至1,1,1,1、2,2,2,2、3,3,3,3等。其实每个状态都可能产生a，b，c，但问题在于，不同的状态产生的a,b,c的概率不同。
如图：

对于每个可能的状态路径，都可以按照这个方法来计算其产生aabc字符串的可能性，其中概率最大的那个状态路径，也是最可能产生这个字串的状态路径。

马尔可夫链引入空位状态X，Y以及Match状态M，可以很好的处理gap与gap extending 问题，但是没有考虑到具体残基，因此引入隐马尔可夫模型来补充这一点。我们用生成概率来处理残基，M状态生成的符号是所有可能的残基替代，其生成概率写作为$Pab$,而X和Y状态生成的符号则是所有可能插入的残基，其生成概率写成$Q_a$,这样我们就可以很方便的同时考虑状态和具体的残基了。

进而将序列比对问题重新描述为一个针对特定马尔可夫模型与符号串寻找最可能状态路径的问题。具体来说，我们可以应用动态规划的思想来分布求解。我们定义：

1. PM(i,j)表示在Xi比对到Yj，也就是两个残基对在了一起时候，第一条序列X从第1位到第i位、第二条序列Y从第1位到第j位最大的概率。
2. PX(i,j)和PY(i,j)分别表示在Xi或Yj残基比对到空位时，序列X从第1位到第i位、序列第1位到第j位最大概率。
   这样我们就可以根据状态转换图来定义每步的迭代函数。

并进一步应用动态规划的back tracking来得到最终比对。具体操作上，只需要之前每一步迭代时的最大值变为求和，再将最后的求最大值变换为求和就可以了。

隐马尔可夫模型通过符号观测序列来反推这个特点，它的应用领域并不局限于序列比对，在现代生物信息研究中，隐马尔可夫模型更多的是作为预测器来使用。

3、用隐马尔可夫模型建立预测模型

我们可以对每个可能的状态路径，计算产生观测符号串的可能性，而其中概率最大的那个状态路径也就是最可能这个字串的路径。
接下来使用一个例子来说明。做一个简单的基因预测，给定一段基因组DNA序列，来预测其中的编码区。按照上述的隐马尔可夫模型，先要区分不能直接观测的隐藏状态和可以直接观测的显符号，在这个例子中，很容易看出，给定的基因组DNA序列是可以观测到的符号串，而编码/非编码则是不能直接观测的隐状态。因此，可以画出编码和非编码区域的状态转换图。

假定我们经过对训练集的分析，分别填好了转移概率举证和生成概率矩阵。

我们现在需要根据这些数据，来对一个未知的给定基因组序列反推出最可能的状态路径，也就是概率最大的那个状态路径，因此我们还是和之前一样利用动态规划算法。由于我们在计算中需要做大量乘法，很容易出现下溢可能，因此对数据做了对数处理。

现在用一个测序的序列正式开始！序列为CGAAAAAATCG。

选出最终概率最大的值，进行回溯，拿到最终的回溯路径：

得到最终的结果：

注：详细内容可参考：北京大学生物信息学