TWI420007B - 基因測序序列的組合系統及方法 - Google Patents

Description

基因測序序列的組合系統及方法

本發明係關於一種基因測序序列的組合系統及方法，尤指一種基因測序(DNA sequencing)的資料分析方法。從頭開始組合(De novo assembly)基因測序(DNA sequencing)產生的核苷酸鹼基序列的短字串。將可能有錯誤的短核苷酸鹼基序列字串拼接成正確的長基因序列。

按基因測序序列的組合技術，是一種應用在基因序列的拼圖技術。基因測序所產生的序列，如以下的例子(26421212個序列)：

序列1:TCCTGTATATTCTAAACTTAGAGATTGTTCAT；

序列2:CATAAACATCTTTATAAAATACTAATAGAAAG；

序列3:AAAGGAGAGAACGTCGTCGTTTTCGTCGAAGT；

序列4:ACAACCCTAACTCTTTTTTTTTTGGCTATTGT；

…

序列26421209:

TCTTCCGCCGTCGCAACTTTACCCAACGCCGC；

序列26421210:

ACCGCAAAAGCAAGATGATTCATTGTGTATCC；

序列26421211:

CTGGATCACAGCATCCACACGCACAAATATC；

序列26421212:

CCAATGGATTCTTTCTTTACTAACAATATCGA。

上述的基因測序序列的組合問題和普通拼圖(Jigsaw Puzzle)的拼接不同之處主要有：

(1)基因測序序列的組合是一維的字串拼接。基因測序資料的碎片是只包含四種(A,G,C,T)核苷酸鹼基的字串，而拼接的時候是把不同的序列碎片依照其一致的部份重疊起來。

(2)基因測序序列的組合是巨大數量的碎片拼接。常見圖像拼圖的碎片數可能是300片或500片。1000片的拼圖就很難拼接。而基因測序的組合所要拼接的序列數量是巨大的，往往有1,000,000到1,000,000,000，甚至隨著技術進步可以產生更多的序列碎片。

要拼接數目如此龐大的測序序列，往往需要非常大量的記憶體來記住過程中所產生的重疊群產物。而且，因為基因測序的過程可能會產生序列的錯誤。因此，如何判斷序列中的錯誤也是一項重要的待解決課題。

傳統的基因測序從頭開始組合方法包含以下三種：

(1)重疊-排列-一致法(Overlap-Layout-Consensus)；

(2)De Bruijn圖(De Bruijn graph)；及

(3)貪婪延伸演算法(Greedy extension algorithm)。

以下以5個基因序列的拼接例子，分別對以上三種習知方法做進一步說明：

r1 CCCTTCCAAC；

r2 ATTTAATCCC；

r3 TTAATCCCTT；

r4 TTCCAACAGC；及

r5 AACAGCCGCC

(1)重疊-排列-一致法是最傳統的方式。包括三個階段：

階段一：將兩兩的基因測序序列重疊看看，找出其距離。如r2和r1可以對齊中間3個鹼基，記成d(r2,r1)=-3，並如下表。

如r1和r2可以對齊中間8個鹼基，記成d(r2,r3)=-8，並如下表。

階段二：將兩兩的基因測序序列重疊看看，找出其距離。以基因序列間的距離來建立所有基因序列的有向圖，如圖1所示。

階段三：在有向圖中找出有一致關係的排列順序(如以最小擴充樹法)，如圖2所示。圖2中，如以最小擴充樹法可以得到序列拼接的順序是由左至右重疊r1,r2,r3,r4,r5。這五個基因序列的重疊稱為一個重疊群(contig)，由這個重疊群中每個一致的鹼基可以得到最後組合序列：ATTTAATCCCTTCCAACAGCCGCC，如圖3所示。

(2) De Bruijn圖(De Bruijn graph)：係把序列以每k個組成一個節點，如圖4所示。將不同的基因序列而具相同的節點予以合併，可以得到如圖5所示之結果。De Bruijn圖是把圖5中相鄰的節點合併成一個更大的節點。因為圖中示例只形成一個序列，所以把相鄰的節點合併最後會形成單一的節點，如圖6所示。如果，合併的結果形成複數的節點，則最後找尋一筆畫的路徑(Eulerian path)來做為最後可以合併的序列。應用De Bruijn圖之專利技術有美國第7,071,324號、第7,034,143號、第6,865,491號、第6,689,563號及第5,683,881號專利案。

(3)貪婪延伸演算法(Greedy extension algorithm)：係選取一個基因序列如r1: CCCTTCCAAC，看看其字尾(postfix)是不是別的字的字首(prefix)。如是，便將其重疊上去。如TTCCAAC是r1的字尾，是r4: TTCCAACAGC的字首。所以合併r1及r2變成重疊群(1,4): CCCTTCCAACAGC。重疊群(1,4)的字尾AACAGCC是r5: AACAGCCGCC的字首，所以把重疊群(1,4)及r5合併成重疊群(1,4,5): CCCTTCCAACAGCGCC。重疊群(1,4,5)沒有任何字尾是別的序列的字首，所以就停止。換從序列如r2: ATTTAATCCC開始接。r2的字尾TTAATCCC是r3: TTAATCCCTT的字首，所以合併r2及r3變成重疊群(2,3): ATTTAATCCCTT。最後再把重疊群(2,3)及重疊群(1,4,5)合併，可以得到ATTTAATCCCTTCCAACAGCCGCC。

以上三種傳統的方式都需要不斷進行合併(merge)的動作。把比較短的基因序列之重疊群合併成比較長的基因序列之重疊群。合併的過程需要大量的記憶體來存放拼接過程的暫時結果。然而，當資料量很大時，往往需要很大的記憶體來存放拼接過程的結果。甚至要多達數百Giga的記憶體才能進行。因此，當基因序列資料大時，往往受限於記憶體的限制無法完成拼接的動作。而且，當基因序列中有鹼基是錯誤時，往往就無法被組合。

再者，關於基因序列之重組或分析的技術有很多，例如中華民國第I326431號專利案，美國第7,809,509號專案，以及如附件一之參考文獻[1]至[10]所發表之技術內容。然而，目前所見的在先技藝，尚未發現有如本發明之領先技術者。

本發明之目的，在於提供一種組合基因測序序列的系統及方法，用以解決傳統方法所產生之二個問題：(1)合併的動作需要大量的記憶體；及(2)容許序列中有錯誤的鹼基也可以進行組合。

為解決上述問題1，本發明之技術手段係提供一種雙向延伸組合的方法來拼接各別的基因序列以形成一目標基因序列。這個技術是發展來將各別基因序列同時向一個待接基因序列的二側延伸接續下去。因為是向待接基因序列二側延伸，因此，我們可以任意選取一個待接之基因序列(為一基因序列或由數個基因序列接續組合而成的一重疊群)開始進行其二側的雙向延伸接續其他候選基因序列的動作。最後可以找出位在同一個重疊群上的其他基因序列，並將它們組合成一目標基因序列(即由更多基因序列接續組合而成的一更長的重疊群)。

為解決上述問題2，本發明之技術手段係提供一個容錯序列的篩選機制來找出正確的序列。由於基因測序產生的序列可能會有以下二種錯誤：

(1)鹼基配對失誤(mismatch)之錯誤：原基因序列中某個鹼基被錯誤配對成其他鹼基。例如：ACATTAAGCCTT是原本的基因序列，經過基因測序處理所產生的序列為AGATTAAGCCTT。也就是第二個鹼基C被錯誤配對成G。

(2)鹼基插入或刪除(insert/deletion)之錯誤：經過基因測序處理所產生的基因序列比原基因序列中多出額外的鹼基，或減少某個鹼基。例如ACATTAAGCCTT是原本的基因序列，第4~5鹼基是連續的鹼基T。如經過基因測序處理所產生的序列為ACATAAGCCTT，比原本的基因序列在相同位置減少一個T，此情形稱之為鹼基刪除的錯誤。反之，如經過基因測序處理所產生的基因序列為ACATTTAAGCCTT比原本的基因序列在相同位置多一個T，此情形稱之為鹼基插入的錯誤。

本發明之容錯序列篩選器，是在所有的候選基因序列中找出可以接續基因序列的正確基因序列，而容許基因序列中有錯誤的鹼基。

為讓本發明之目的及其他特徵能更加清楚，以下茲舉出一些較佳實施例，並配合所附圖式圖7到圖16，作詳細說明。在這些實施例的說明中，為了簡明解釋原理，所以在不同實例使用不同的序列長度。以及不同的索引鍵長度。

I.本發明之系統

如圖7所示，本發明之基因序列組合系統10係用以將一個基因序列的集合拼接成目標的基因序列，其一種具體實施例係包括有一輸入介面11、一索引器12、一雙向延伸組合器13、一容錯序列篩選器14、一重疊群建構器15及一輸出介面16。茲將前述各元件詳述如下。

輸入介面11係從儲存在資料庫或記憶體中之檔案110讀入複數個基因序列(其可以是由基因測序系統所產生的複數個基因序列)，用以給予輸入的基因序列編號，並建立基因序列的左索引結構及右索引結構。本發明一種具體實施例中，其輸入介面11會讀入基因序列二次。第一次輸入基因序列，取其基因序列前後各N個鹼基做為索引鍵值資料111，並將索引鍵值資料111置入索引器12中存放。索引鍵值可以字串或轉換成數值表示。輸入介面11有一個序列使用記錄陣列，來記錄序列是否已被排入重疊群中。

索引器12，其用以儲存該複數個基因序列之索引值資料111，該索引值供找出可能可以接續在一個待接基因序列之二側的候選基因序列。其可以是置於記憶體中的一個索引陣列，或者是置於硬碟的索引檔，也可以是一個置放在遠端的資料庫，作用是輸入經輸入介面11並索引的基因序列(即短鹼基序列130)，及輸出與索引對應之多重候選基因序列122。

雙向延伸組合器13，用以將經由容錯序列篩選器所決定之該選定基因序列接續在待接基因序列至少一側而延伸成一個更長的基因序列，直至決定該目標基因的鹼基序列為止。本發明實施例中，係以雙向延伸組合器13取出待接基因序列(或目前已組合重疊群序列)二側各M個鹼基長度分別做為一延伸測試視窗21/22，本發明係以一個基因序列的長度做為延伸測試視窗之長度，該二延伸測試視窗21/22分別供自該索引器中找尋出可以附加在該延伸測試視窗的該基因序列以做為該候選基因序列。本發明具體實施例中，雙向延伸組合器13會從長度為1開始位移延伸測試視窗，將位移後的延伸測試視窗中的基因序列分成新的索引鍵131及容錯比對區域132，如圖11所示。其中，新的索引鍵131用以向索引器12查詢可能的候選基因序列。而容錯比對區域132提供給容錯序列篩選器14，用以比對出正確而可供延伸接續的選定基因序列。

容錯序列篩選器14，用以決定候選基因序列為可接續在待接基因序列之二側的選定基因序列。其係根據由雙向延伸組合器13輸入的容錯比對區域132及多重候選基因序列122，請配合參看圖13、14所示的篩選過程，留下帶有位置且正確可供延伸接續之選定基因序列141給重疊群建構器15。

重疊群建構器15將選定基因序列141依其位置重疊排列，建構出重疊群(contig)151。透過輸出介面16將此重疊群151輸出到檔案161中。

Ⅱ.本發明之方法

請配合參看圖7至16所示，本發明之基因序列組合方法的一種具體實施例，係包括有以下所述之步驟。

步驟S201：由輸入介面11輸入複數個基因序列，給予輸入的所有基因序列一個編號，並且建立此基因序列的左索引結構及右索引結構，並儲存在索引器12。

步驟S211：輸入介面11從序列使用記錄陣列中找出一個未使用的序列，先和其鄰近的序列比對，確定每個鹼基的正確性後，此序列做為雙方延伸組合器13進行雙向延伸的啟始待接基因序列112。因為單一的基因序列可能會有錯誤，因此可以使用數個連續相鄰的複數個基因序列來做為啟始待接基因序列片段。尋找連續相鄰的複數個基因序列，是用位移的索引鍵找尋彼此鹼基都一致的基因序列來先重疊成啟始待接基因序列片段。其中，如果相鄰的基因序列彼此鹼基不一致，就不能做為雙向延伸的啟始待接基因序列片段的二側。

步驟S221及S222是左右對稱的運算程序，在此以向右的實施例做說明。雙向延伸組合器13取出待接基因序列20二側各M個鹼基長度分別做為左延伸測試視窗21及右延伸測試視窗22，並從長度1開始位移測試視窗21/22，將位移後的測試視窗21/22中的基因序列分成新的索引鍵131及容錯比對區域132，以新的索引鍵131向索引器12查詢可能的候選基因序列122，並將容錯比對區域132提供給容錯序列篩選器14，用以比對出正確而可供延伸接續的選定基因序列。

如圖9所示，本發明以待接基因序列20(或為目前完成群組之基因序列片段)的二側，分別做為右延伸測試視窗及左延伸測試視窗。如圖10所示，本發明滑動左延伸測試視窗21及右延伸測試視窗22，用以找出可以接續在目前已知待接基因序列20左右二側的候選基因序列122。

步驟S231及S232是左右對稱的運算程序，圖12係以產生向右延伸之候選基因序列為例，目前比對視窗23中的比對參考序列型式為CACAGCAGTAAGTTTCCAATATATGGT。此序列中，CACAGCA做為索引鍵，而GTAAGTTTCCAATATATGGT是用以進行容錯比對的區域。從索引器12中找出所有左側索引鍵為CACAGCA的基因序列。這些基因序列也分成索引鍵及容錯比對的區域。比較延伸視窗23及候選基因序列122的比對參考序列型式，計算出其不同鹼基的數目，如果不同鹼基的數目小於一個閥值T，則該基因序列被選為可能延伸之候選基因序列。

步驟S241及S242是左右對稱的運算程序，由前一個步驟產生的可能延伸之候選基因序列，必須進一步測試是否有測序錯誤。本發明之方法是把所有被找出的候選基因序列依其可能延伸的位置重疊排列，計算每個位置其ACGT鹼基所佔的比率，即統計不同候選基因序列排列後相同位置的鹼基，以判斷是否是測序產生的序列錯誤或者該候選基因序列並不是接在此位置的序列。對單一的序列而言，如果其某個位置的鹼基和其他基因序列的相同位置之鹼基不同，會有二種情形：第一種情形是此基因序列是正確的候選基因序列，但是發生鹼基配對失誤的測序錯誤；第二種情形是此一基因序列並不是可以接在此位置的候選基因序列。圖13及14分別圖示說明此二種情形。在圖13中，序列r1,r2,r3,r4,r6各有1~2個鹼基和其他候選基因序列不一致。然而其重疊時，各別位置的錯誤鹼基沒有超過一定百分比，如1/5。此時，錯誤鹼基被視為鹼基配對失誤的測序錯誤。此外，在圖14中，序列r1,r2,r3,r4,r5,r6,各有1~2個鹼基和其他候選序列不一致。當其重疊時，r1,r2,r3有一個相同位置的錯誤鹼基超過一定的比率，如1/5。該位置的鹼基都是A，相較於其他基因序列在此位置的鹼基都是T，因此r1,r2,r3等基因序列被判定為不是接在此位置的候選基因序列。步驟S241及S242也偵測是否發生鹼基插入或刪除的測序錯誤。在圖15及16顯示鹼基插入或刪除的測序錯誤偵測。

步驟S251及S252是左右對稱的判斷程序步驟，如果前一個步驟產生一些的候選基因序列，可以附加到已知待接基因序列之右側，則重新進行步驟S221。如果前一個步驟產生一些的候選基因序列，可以附加到已知待接基因序列之左側，則重新進行步驟S222。

步驟S261，當待接基因序列二側都無法繼續附加新的基因序列，則把所有找到的可延伸之選定基因序列依其位置重疊成重疊群(contig)。並輸出重疊群每個位置最判定的鹼基以成為組合的目標序列。

如圖9及10所示，係為本發明之雙向延伸組合器進行序列組合的實施例圖。此實施例說明本發明找尋可以拼接在一起的基因序列群的主要方法。由一個小的啟始待接基因序列向二端延伸，找出可以接在適當位置的基因序列。

如圖11所示，係為本發明之容錯序列篩選器的簡化實施例圖。此一實施例說明容錯序列篩選器和延伸測試視窗的關係。延伸測試視窗是啟始待接基因序列二側的比對序列。雙向延伸組合器會位移此延伸測試視窗，並將延伸測試視窗內的基因序列分成索引鍵及容錯比對區域。

如圖12所示，係本發明找出可以用以延伸基因序列的候選基因序列之容錯比對方法。

如圖13及14所示，係顯示對候選基因序列進行篩選，偵測是否發生鹼基配對失誤的測序錯誤。圖13顯示發生鹼基配對失誤的測序錯誤情形，圖14顯示非鹼基配對失誤的測序錯誤情形。

圖15及16顯示對候選基因序列進行篩選，偵測是否發生鹼基插入或刪除的測序錯誤。鹼基插入或刪除之錯誤的偵測，係將原來比對的序列型式轉換成差別序列型式進行比對。在延伸測試視窗的比對參考型式ref會被轉換成差別序列型式dref，方法是掃描基因序列。連續相同的鹼基被視為單一鹼基。例如比對參考型式ref=GTAAGTTTCCAATATATGGT，其差別序列型式dref=GTAGTCATATATGT，也就是在ref中的連續二個AA鹼基在dref中只表示成單一個A鹼基。同理，在ref中的連續二個TTT鹼基在dref中只表示成單一個T鹼基。在進行候選基因序列篩選時，候選基因序列r1的比對參考型式GTAAAGTTTCCAATATATGGT，其差別序列型式dr1=GTAGTCATATATGT。比對二個差別序列型式(dref,dr1)是一致的，因此r1的比對參考型式會被取代成er1=GTAAGTTTCCAATATATGGT。如此，r1被視為可以接在此位置的可延伸之選定基因序列。

雖然本發明已以較佳實施例揭露如上，然其並非用以限定本發明，任何熟悉此項技藝者，在不脫本發明之精神和範圍內，當可做些許更動與潤飾，因此本發明之保護範圍當視後附之申請專利範圍所界定為準。

10．．．基因序列的組合系統

11．．．輸入介面

110,161．．．檔案

111．．．索引鍵值資料

112．．．啟始待接基因序列

12．．．索引器

122．．．候選基因序列

13．．．雙向延伸組合器

130．．．短鹼基序列

131．．．索引鍵

132．．．容錯比對區域

14．．．容錯序列篩選器

141．．．選定基因序列

15．．．重疊群建構器

151．．．重疊群

16．．．輸出介面

21．．．左延伸測試視窗

22．．．右延伸測試視窗

圖1為習知有向圖；

圖2為習知有向圖中找出有一致關係的排列順序之示意圖；

圖3為習知以重疊-排列-一致的方式組合序列的示意圖；

圖4為習知De Bruijn示意圖；

圖5為習知De Bruijn圖中相鄰節點合併一大節點示意圖；

圖6為習知De Bruijn圖合併而成的序列示意圖；

圖7為本發明之基因序列組合系統的一種實施例示意圖；

圖8為本發明之基因序列組合方法的一種實施例流程圖；

圖9為本發明雙向延伸組合器具有左右延伸測試視窗之示意圖；

圖10為本發明雙向延伸組合器進行序列組合的簡化實施例圖；

圖11為本發明容錯序列篩選器的簡化實施例圖；

圖12為本發明找出可以用以延伸序列的候選序列方法示意圖；

圖13為本發明對候選序列進行篩選及偵測是否發生鹼基配對失誤的測序錯誤之一種示意圖；

圖14為本發明對候選序列進行篩選及偵測是否發生鹼基配對失誤的測序錯誤之另一種示意圖；

圖15為本發明對候選序列進行篩選及偵測是否發生鹼基插入或刪除的測序錯誤之一種示意圖；及

圖16為本發明對候選序列進行篩選及偵測是否發生鹼基插入或刪除的測序錯誤之另一種示意圖。

附件一：參考文獻。

10‧‧‧基因序列的組合系統

11‧‧‧輸入介面

110,161‧‧‧檔案

111‧‧‧索引鍵值資料

112‧‧‧啟始待接基因序列

12‧‧‧索引器

122‧‧‧候選基因序列

13‧‧‧雙向延伸組合器

131‧‧‧索引鍵

132‧‧‧容錯比對區域

14‧‧‧容錯序列篩選器

141‧‧‧選定基因序列

15‧‧‧重疊群建構器

151‧‧‧重疊群

16‧‧‧輸出介面

Claims (12)

1. 一種基因序列組合系統，其用以拼接經一基因測序系統所產生之複數個基因序列，以決定一目標基因的鹼基序列，該系統包括：一索引器，其用以儲存該複數個基因序列之索引值資料，該索引值供找出可能可以接續在一個待接基因序列之二側的候選基因序列；一容錯序列篩選器，其用以決定該候選基因序列為可接續在該待接基因序列之二側的選定基因序列；及一雙向延伸組合器，其用以將經由該容錯序列篩選器所決定之該選定基因序列接續在該待接基因序列至少一側而延伸成一個更長的基因序列，直至決定該目標基因的鹼基序列為止。
2. 如請求項1所述之基因序列組合系統，其中該雙向延伸組合器包括有分別對應於該待接基因序列二側的二延伸測試視窗，該二延伸測試視窗分別供自該索引器中找尋出可以附加在該延伸測試視窗的該基因序列以做為該候選基因序列。
3. 如請求項1所述之基因序列組合系統，其中該容錯序列篩選器從該候選基因序列中偵測是否有基因測序的序列錯誤。
4. 如請求項3所述之基因序列組合系統，其中該容錯序列篩選器所偵測的該序列錯誤，包括鹼基配對失誤之錯誤及鹼基插入或刪除之錯誤。
5. 如請求項4所述之基因序列組合系統，其中該鹼基配對失誤之錯誤的偵測，係將複數個該候選基因序列依其可能延伸的位置排列，統計不同序列排列後相同位置的鹼基，以判斷是否是測序產生的序列錯誤或者該候選基因序列並不是接在此位置的序列。
6. 如請求項4所述之基因序列組合系統，其中該鹼基插入或刪除之錯誤的偵測，係將原來比對的序列型式轉換成差別序列型式進行比對。
7. 如請求項6所述之基因序列組合系統，其中該差別序列型式是將基因序列中相同的連續鹼基表示成單一的鹼基。
8. 如請求項4所述之基因序列組合系統，其中該容錯序列篩選器進行該鹼基插入或刪除之錯誤的偵測時，當成功比對該選定基因序列及延伸測試視窗的差別序列型式，則該選定基因序列所插入的多餘鹼基或刪除的缺少鹼基會被取代成正確的鹼基數目。
9. 如請求項1所述之基因序列組合系統，其中該索引器是以該基因序列的部份片段做為索引鍵值，用以將所有的該基因序列分類存放，於組合時，再依其索引鍵值供取回該基因序列做使用。
10. 一種基因序列組合方法，其用以拼接經一基因測序系統所產生之複數個基因序列，以決定一目標基因的鹼基序列，其包括：步驟(A)提供如請求項1所述之系統；步驟(B)輸入該複數個基因序列並建立索引，將其索引值資料儲存於該索引器；步驟(C)從未使用的該序列中產生一待接基因序列；步驟(D)以該雙向延伸組合器取出該待接基因序列二側各一預定個鹼基長度分別做為一左延伸測試視窗及一右延伸測試視窗，並以一預定長度分別左移該左延伸測試視窗及右移該右延伸測試視窗，將每次位移後的該左延伸測試視窗及該右延伸測試視窗中的基因序列分成新的一索引鍵及一容錯比對區域，以該新的索引鍵向該索引器查詢可能的候選基因序列；步驟(E)以該容錯序列篩選器決定該候選基因序列為可接續在該待接基因序列之二側的選定基因序列；步驟(F)以該雙向延伸組合器將經由該容錯序列篩選器所決定之該選定基因序列接續在該待接基因序列至少一側而延伸成一個更長的基因序列，當延伸成功時，再重複步驟(D)至(F)，直至該待接基因序列二側都無法繼續接續為止；及將所有找到的可延伸之該選定基因序列依其位置重疊成重疊群，並輸出該重疊群每個位置最有可能的鹼基而成該目標基因的鹼基序列。
11. 如請求項10所述之基因序列組合方法，其中步驟(B)係經由一輸入介面輸入該複數個基因序列，給予每一該基因序列一編號，並且建立該基因序列的左索引結構及右索引結構。
12. 如請求項10所述之基因序列組合方法，其中步驟(C)係使用連續相鄰的複數個基因序列來做為啟始待接基因序列片段，而該連續相鄰的複數個基因序列是用位移的索引鍵找尋彼此鹼基都一致的基因序列並重疊成該啟始待接基因序列片段。