模體發(fā)現(xiàn)問題中OOPS模型的EM算法

黃影

摘 要 生物序列模體發(fā)現(xiàn)是從生物序列集合找出具有特定功能子序列的問題，它是生物信息學中最重要問題之一。OOPS模型是一種簡單常用的模體發(fā)現(xiàn)的序列模型，它表示每條序列中有一個模體出現(xiàn)。本文提出了一種適用于OOPS序列模型的EM算法，用于查找每條序列中所出現(xiàn)的模體。實驗結(jié)果表示，所提出的算法的時間性能和識別準確率均好于經(jīng)典的Project算法。

關(guān)鍵詞 模體發(fā)現(xiàn) 生物序列 EM算法

中圖分類號：TP181 文獻標識碼：A DOI：10.16400/j.cnki.kjdkx.2015.08.010

EM Algorithm of OOPS Model in Motif Discovery Problem

HUANG Ying

（School of Information Engineering， Xi'an University， Xi'an， Shaanxi 710065）

Abstract Biological sequence motif discovery from biological sequence set to find out having a problem specific functional sub-sequence， it is one of the most important issues in bioinformatics. OOPS model series model is a simple common motif finding， which represents each sequence has a motif appears. This paper presents a suitable OOPS series model EM algorithm， for finding that appears in each sequence motif. The results indicated that the time performance of the algorithm and the recognition accuracy rate of the proposed algorithm is better than the classical Project.

Key words motif discovery; biological sequence; EM algorithm

1背景

EM（expectation-maximization）算法是Dempster， Laird和Rubin（DLR）三個人在1977年正式提出的。主要是用于在不完全數(shù)據(jù)的情況下計算最大似然估計。在EM算法正式提出以來，人們對EM算法的性質(zhì)有更加深入的研究;并且在此基礎(chǔ)上，提出了很多改進的算法，在數(shù)理統(tǒng)計，數(shù)據(jù)挖掘，機器學習以及模式識別等領(lǐng)域有廣泛的應(yīng)用。

生物序列模體（motif）發(fā)現(xiàn)問題是從一個從相關(guān)生物序列集合找出具有某種特定功能的子序列的問題，它是生物信息學中最重要問題之一，也是最基本問題之一。針對motif在每條序列中出現(xiàn)或不出現(xiàn)或者出現(xiàn)多次，有三種對應(yīng)的模型，即：若每條序列中有且僅有一個motif出現(xiàn)，則被稱為OOPS模型，這是最簡單的模型;若每條序列中至多有一個motif出現(xiàn)，則被稱為ZOOPS模型;若每條序列允許多個motif出現(xiàn)或者沒有motif出現(xiàn)，則被稱為TCM模型，這是最復(fù)雜的模型。

本文實現(xiàn)基于最簡單的OOPS模型實現(xiàn)EM算法，來查找每條序列中所出現(xiàn)的motif。

2模體發(fā)現(xiàn)的EM算法

在最大化似然的問題中，給定被稱為數(shù)據(jù)參數(shù)似然值的函數(shù)（ ; ），或者就稱為似然函數(shù)：

（∣ ） = （∣ ） = （ ; ）

我們的目標就是找到使最大的 值，即

* = ? （ ; ）

EM算法的核心思想就是根據(jù)已有的數(shù)據(jù)來遞歸估計似然函數(shù)。

所有的EM算法都由兩個主要的步驟組成：

E-step

M-step

在E-step中，對于未知的潛在變量，使用當前參數(shù)的估計值和當前的觀察。在M-step中，重新計算出來一個新的參數(shù)的估計值。在每次迭代中，似然值都會增加，因此該過程總會達到一個漸近最大值;但EM可能陷入一個局部最優(yōu)。通過交替使用這兩個步驟，EM算法逐步改進模型的參數(shù)，使參數(shù)和訓(xùn)練樣本的似然概率逐漸增大，最后終止于一個極大點。

直觀地理解EM算法，它也可被看作為一個逐次逼近算法：事先并不知道模型的參數(shù)，可以隨機的選擇一套參數(shù)或者事先粗略地給定某個初始參數(shù)，確定出對應(yīng)于這組參數(shù)的最可能的狀態(tài)，計算每個訓(xùn)練樣本的可能結(jié)果的概率，在當前的狀態(tài)下再由樣本對參數(shù)修正，重新估計參數(shù)，并在新的參數(shù)下重新確定模型的狀態(tài)，這樣，通過多次的迭代，循環(huán)直至某個收斂條件滿足為止，就可以使得模型的參數(shù)逐漸逼近真實參數(shù)。

EM算法的主要目的是提供一個簡單的迭代算法計算后驗密度函數(shù)，它的最大優(yōu)點是簡單和穩(wěn)定，但容易陷入局部最優(yōu)。

如果給定一組已比對好的序列集，可以直接構(gòu)造一個表征motif的profile矩陣，即PWM;圖1是一個例子。

但對于未比對的序列，我們怎么構(gòu)造profile矩陣？因為在通常情況下，我們并不知道m(xù)otif是什么，這時我們就可以使用EM算法。

使用概率矩陣 = 表示motif概率，這里表示列的字符所出現(xiàn)的概率。用概率 = 來表示背景（即，序列中motif之外的部分）概率，用表示背景中的字符的出現(xiàn)概率。并定義矩陣，其元素表示motif出現(xiàn)在序列的位置的概率。例如，在圖2中，給定長度為6的3條DNA序列，這里motif的長度W=3，其對應(yīng)的、和Z：

圖1 給定一組已比對好的序列集，直接構(gòu)造表征motif的profile矩陣

圖2 給定長度為6的3條DNA序列，motif的長度W=3，對應(yīng)的、和

整個EM算法的框架如下：

EM algorithm

Given： length parameter W， training set of sequences

set initial values for P

do ? ?re-estimate Z from P ? ? ? ? ? （E-step）

re-estimate P from Z ? ? ? ? ? （M-step）

until change in P <

Return： P， Z

在E-step中， 由計算的估計即motif出現(xiàn)的最可能位置;在M-step中，由重新計算的估計即一個新的motif的PWM。

假定給定條長度為的序列： = （， ，…，）， = 1， …，，我們找出長度為的motif。則該序列在位置j出現(xiàn)motif的概率為：

在E-step：由計算

在E-step的第次迭代中，由下式計算獲得：

= （ = 1∣ ， ， ）

=

假定motif等可能的出現(xiàn)在序列的任意位置，則上式可化簡為：

=

在M-step：由計算

對于motif，對任意的，有

=

這里 = ?。

對于背景，上式的可用 = 計算得到。這里表示序列中堿基的總個數(shù)，是為避免計算中出現(xiàn)零概率而引入的擾動。

3實驗結(jié)果

我們在隨機生成的模擬數(shù)據(jù)上對我們的程序進行了測試。隨機生成t條滿足均勻分布的DNA序列，在每條序列中植入（l， ）motif，即motif長度為l，發(fā)生個變異。對所實現(xiàn)算法的正確性，我們使用性能系數(shù)進行評估：

performace coefficient =

這里表示所植入的motif在序列中的位點集， 表示由算法所找到的motif在序列中的位點集。

我們將我們的實現(xiàn)同project方法進行了比較，結(jié)果如表1所示：

表1 EM和project的pc值

表1中每項表示在十個樣本上的平均的性能系數(shù)，其中每組包含20條數(shù)據(jù)，每條序列長為600。每輪測試進行10次迭代。對于（12， 3）、（15， 4）、 （18， 5）、（21， 6）、（24， 7），我們的實現(xiàn)的性能系數(shù)分別為0.6474、0.8740、0.9810、0.9723、1，同相應(yīng)的project的方法相差不大，而它的性能系數(shù)分別為0.6576、0.8740、0.9905、0.9723、1。

同時，我們改變motif的長度，對兩個算法的執(zhí)行時間進行比較（單位：s），如表2所示：

表2 EM和project的運行時間

從表2中我們可以看出，我們的方法隨著所找motif長度的增加，時間增長比較慢;而相應(yīng)的project方法所用的時間增長很快。

參考文獻

[1] Arthur Dempster， Nan Laird， and Donald Rubin. Maximum likelihood from incomplete data via the EM algorithm.Journal of the Royal Statistical Society，Series B，1977.39（1）：1-38.