基因识别及隐马模型

显然 iδT1((α))就是[m解{a码x}问 t题(的)t解。] eoi1
学习算法
HMM的学习算法有许多种，包括 EM （期望最大化）算法，不同形式的梯度 下降法，模拟退火算法等。下面我们主 要介绍 EM 算法。
EM算法的用途
EM算法（The expectation maximization algorithm）是根 据不完整数据作最大似然估计 的一般方法。这种不完整可以 是缺失了某些数据，也可以是 存在某些无法观测的隐藏变量。
即每次叠代后非完整数据的似然性都会 增加，而似然性显然是有界的，因此叠 代必然收敛。
隐马模型的参数估计
利用EM算法进行隐马模型的参数估计， 又称为Baum-Welch算法。它要解决的主 要问题，是求
W arg max P(OW )
W
隐马模型复习
隐马模型：
– 存在隐序列H，由下列参数生成：
0 arg max L( X )

为简化计算，实际工作中常使用的是对 数似然估计。即：
0 arg max ln L( X )

具体计算公式为：
d ln L( ) 0 d
如果参数不止一个，相应令偏导数为0即 可。
例题1 正态分布的最大似然估计
设x1 ,x2 ,… xn是取自正态总体 N(μ ,σ 2) 的简单随机子样，μ 与σ 2是未知参数，求
P(O | w)
通过计算上述概率的极大值，可以得到 位置 i 处隐序列最可能的状态。 但在解码问题中，我们更关心的是最可 能的隐序列。它是一个整体，不能分解 开一位一位算。解决这个问题，需要使 用Viterbi算法。
Viterbi算法
定义 i ( ) max P(o1oi , hi | w) Hi 其中Hi 为一切长为i、以α结束的隐序列。 显然δi(α)对应的隐序列就是最可能生成 前i个字母明序列的隐序列。递推公式为：
对上式两边同时乘以P(y|x,θt)，并对所有
y求和得到
log P(x ) P( y x, t )log P(x, y ) P( y x, t )log P( y x, )
y
y
Q( t ) P( y x, t ) log P( y x, )
y
因此有：
log P(x ) log P(x t ) Q( t ) Q( t t )
当且仅当θt=θ，或P(y|x,θt)=P(y|x,θ)时等 式成立。
若取θt+1=argmaxQ(θ|θt)， 显然有 Q(θt+1|θt)>Q(θt|θt)，因此
log P(x t1) log P(x t )
解：每粒种子发芽与否可视为两点分布： 发芽， 则 X=1, 其概率为 p 不发芽，则 X=0, 其概率为1-p 由似然函数的构造，有： L(p)=P(X=x1|p)·P(X=x2|p)…P(X=xn|p)
由于共有m粒发芽，(n-m)粒不发芽， ∴ L(p)=pm(1-p)n-m
dL( p) mpm1 (1 p)nm (n m) (1) p m (1 p)nm1 dp p m1 (1 p)nm1[m(1 p) (n m) p]
n
(xi ) 0
i 1
ˆ

1 n
n i 1
xi

x
代入（2），得：
1
2
n
(xi x)2
i 1
n
ˆ 2

1 n
n i 1
(xi

x)2

S
2 n
即：μ 和σ 2的极大似然估计分别
为
x
和
S
2 n
。
例题2 二项分布的最大似然估计
取n粒种子作发芽试验，其中有m粒发芽， 求发芽率p的最大似然估计。
似然函数的概念
我们从某个未知分布得到了一组观察值 X={x1,x2,…xn}。此未知分布是由一组参数 Θ={θ1,θ2, … θm}决定。定义似然函数为：
n
由于对数L(函X数) 是 单i1 调P(递xi 增)的 P，( X因此) 不会改变
极值点的位置；而且它可以把连乘变成连加， 从而大大简化计算。因此实际工作中常使用 的是对数似然函数。即：
μ 和σ 的极大似然估计。
解：由于
f (x, , 2 )
1
2
exp

1
2
2
(x

)2

故有似然函数
L(, 2 )
(
1
2
)n
exp
1
2
2

n i1
(xi

)2

取对数，有：
ln L(, 2 ) n ln(2
2
2
)

1
2
2
n
(xi
i 1
)2
∴ 似然方程为：
ln L(, 2 )



ln
L(,
2
)
2
0 n
2
1
2
2
(2)
n i1
( xi

)

0

2 2
2

1
2 4

n i1
( xi
)2Hale Waihona Puke Baidu

0
(1) (2)
由（1）解得：
可能性问题：给定模型参数，当观察到 一个明序列时，这一明序列确实由给定 模型生成的概率有多大？
解码问题：给定模型参数，当观察到一 个明序列时，这一明序列所对应的最可 能的隐序列是什么？
学习问题：观察到足够多明序列时，如 何估计转移概率和生成概率？
基本算法
可能性问题：已知参数为w的隐马模型和 长为T的序列O＝{o1o2…ot…oT}；求在给 定模型下观测到序列O的概率。
{ ,T ,{ }}
– 观测序列O称为明序列，由下列参数生成：
{H , Ea ,{a}}
前向算法
定义在模型w下，隐序列第i位置状态为α， 明序列为 o1…oi时的概率为：
f i ( ) P(hi , o1 oi | w)
则有递推公式：
f i1( )
f i ( )teoi1
令上式等于0，由于 pm1(1 p)nm1 0，有：
m mpˆ npˆ mpˆ 0
pˆ m
n
即：发芽率p的极大似然估计为
pˆ m n
EM算法的一般概念
X：观测得到的数据，它是不完全的； Z=(X,Y)，是完全的；其中Y是缺失的或隐藏
的。它的联合分布为：
(i) arg max Q(, (i1) )

在EM算法中，上述E步骤和M步骤不断 重复，每一次重复似然值都会增加，直 到收敛到似然函数的局部极大值。
EM算法收敛性的证明
由于 P(x, y) P(y x,)P(x) 对数似然函数可以改写成： log P(x) log P(x, y) log P(y x,)
bi ( )
{ }
初值为：
f 1( ) eo1
则有：
P(O | w) f T () { }
后向算法
定义在模型w，隐序列第i位置状态为α条 件下，明序列为 oi＋1…oT时的概率为：
bi ( ) P(oi1 oT | hi , w) 则有递推公式：
例：基因组编码区的隐马模型
设基因组由两种功能区域组成，即编码 区和非编码区。分别由字母c，n代表。 转移矩阵为同种字母延伸或变为另一种 字母的概率。初始状态概率为第一个字 母出现c或n的概率。明序列由A，C，G， T四个字母组成，生成概率分别为编码区 和非编码区四个字母出现的概率。
隐马氏模型的三种典型问题
解：对于任意一条隐序列H，有：
T
P(O, H w) (t eoi ) i 1
其中α为隐序列H第i-1位置的状态， β为i 位置的状态。当i＝1时，tαβ应由πβ代替。
因此，所求概率为：
T
P(O w) P(O, H w)
(t eoi )
H
H i1
其中H为一切可能的隐序列。
初值为：
bT ( ) 1
则有：
P(O | w) b1( ) eo1 { }
使用前向或后向算法，计算P(O|w)的复 杂度约为O(N2)级。
联合应用前向和后向算法，可以容易地 计算隐序列位置i处于状态α的概率：
pi () P(hi | O, w) f i ()bi ()
Q(, (i1) ) E[ln P( X ,Y ) X , (i1) ]
ln P( X ,Y ) f ( y X , (i1) )dy yY
函数Q(θ,θ(i-1) )就是EM算法中，E步骤所 要计算的值。在上述函数中，θ(i-1)是上一 步计算的结果，是确定的值。而θ是接着 的M步骤中要调整的参数，调整它以便 使似然函数的值增加。即：
P(Z ) P(X ,Y ) P(Y X ,) P(X )
对于这个联合分布函数，我们可以定义它的似 然函数： L( Z) L( X ,Y) P(X ,Y )
称为完整数据的似然函数。而X的似然函数
L( X ) P(X )
称为非完整数据的似然函数。
注意完整数据似然函数的自变量中，X和 θ为常数，而Y是随机变量，且服从于某种 由X和θ所决定的分布。参数X是观测数据， 是完全确定的；参数θ则是我们需要估计 的，在计算过程中它会不断调整。设θ(i-1) 为当前我们使用的θ估计值，则我们可以 定义完整数据对数似然的期望：
由于H的数量随T的增加以指数增长， 这种算法实际并不可行。
前向算法
定义在模型w下，隐序列第i位置状态为α， 明序列为 o1…oi时的概率为：
f i ( ) P(hi , o1 oi | w)
则有递推公式：
f i1()
f i ( )teoi1
{ }
初值为：
n
ln L( X ) ln P(xi ) i 1
最大似然估计
目标：根据已知观测数据估计总体参数。 似然函数可视为以未知参数为自变量的函
数，它的统计意义又是在未知参数下观测 到该组数据的概率，因此我们很自然地取 使似然函数达到最大值的未知参数为估计 值。这就是最大似然估计。即：
f 1( ) eo1
则有：
P(O | w) f T () { }
后向算法
定义在模型w，隐序列第i位置状态为α条 件下，明序列为 oi＋1…oT时的概率为：
bi ( ) P(oi1 oT | hi , w)
则有递推公式：
bi ( ) bi1( )t eoi1 { }
生物信息学中的典型数学问题
基因识别 与隐马氏模型
基因结构及重要信号
基因识别模型
问题：从基因组序列出发，识别其中的 基因（基因组注释）。
输入：基因组序列 输出：各种基因元件（包括启动子，外
显子，内含子等）在基因组上的位置。 如果需要，也可输出翻译的蛋白质序列， 预测结果的可靠程度等。
基本概念
随机过程：一族无穷多个、相互有关联的
随机变量。记为：X (t), t T
由于参数 t 经常代表时间，故称为随机过 程。T常为自然数，整数或区间。当参数 取值为整数时，也称为随机序列。
马尔可夫过程：取值为整数的随机过程， 若 t = i 时刻的取值只与时刻 i-1 取值有关， 则称为马尔可夫过程。
隐马氏模型（HMM）
模型的数学描述：
– 存在一个隐序列H，它是不可观测的，且由以 下参数生成： { ,T ,{ }}
其中πα为初始状态出现概率；Tαβ为转移概率， 即tαβ = P(hi= β | hi-1= α)； α, β属于{σ}； {σ}为 字符集，即隐序列由哪些字符组成。
– 观测的结果称为明序列O，它由隐序列按照生 成概率e αa生成。其中eαa =P(a | α)；{a}为明序 列字符集。
log P(x ) log P(x t ) Q( t ) Q( t t ) P( y x, t ) log P( y x, t )
y
P( y x, )
上式中最后一项是P(y|x,θt)与P(y|x,θ)之间 的Kullback-Leibler距离，它总是大于0 的，因此