鹼基序列比對系統及方法
2023-05-04 05:33:46
鹼基序列比對系統及方法
【專利摘要】本發明公開一種鹼基序列比對系統及方法。根據本發明一個實施例的鹼基序列比對系統，用於將包含第一序列及第二序列的一對鹼基序列比對到參考序列，包括：種子序列生成單元，從所述第一序列及所述第二序列中分別生成一個以上的片段，並由此構成第一種子序列集合及第二種子序列集合；映射值計算單元，將所述參考序列劃分為多個區間，並按所述多個區間分別計算包含於所述第一種子序列集合中的種子序列在對應區間內的第一映射值以及包含於所述第二種子序列集合中的種子序列在對應區間內的第二映射值；比對單元，選擇計算出的所述第一映射值及所述第二映射值均為基準值以上的第一區間，並在所述第一區間內搜尋所述第一序列及所述第二序列的映射位置。
【專利說明】鹼基序列比對系統及方法
【技術領域】
[0001]本發明的實施例涉及一種用於分析基因組的鹼基序列的技術。
【背景技術】
[0002]測序儀用於從原始鹼基序列中生成作為長度較短的鹼基序列的短片段(read)，此時將有一對短片段配對(pair)生成。這樣形成配對的短片段是在原始DNA中的預定距離內生成，且根據測序儀的種類的不同，形成為在參考序列中具有反向互補(reversecomplement)方向或相同方向。此時生成的兩個短片段之間的距離(insert size)及各短片段的長度為事先根據鹼基序列分析目的進行設定，且在相同實驗中生成的短片段都具有類似的值。在這些成對的短片段中把先生成的稱為5』短片段，後生成的稱為3』短片段，而且將5』短片段與3』短片段的方向為反向互補關係的稱為雙末端短片段(paired-end read),反之5』短片段與3』短片段具有相同方向的稱為配對短片段(mate-pair read)。
[0003]在對這種雙末端短片段或配對短片段進行比對(alignment)時要同時考慮以下三個條件。
[0004](I)各短片段與參考序列之間的鹼基序列同源性(homology)
[0005](2)兩個短片段的比對方向
[0006](3)兩個短片段的比對位置之間的距離
[0007]現有技術中的比對算法構成為根據條件(I)將兩個短片段分別比對到參考序列上之後，在兩個短片段的比對位置中選擇滿足上述條件(2)、(3)的位置。然而如果這樣進行雙末端短片段或配對短片段的比對，則為了首先獲得滿足上述條件(I)的各短片段的比對位置，對參考序列中不滿足上述條件(2)、(3)的位置也都要進行搜尋，因此存在不必要的計算量過多的問題。
【發明內容】
[0008]本發明實施例的目的在於提供一種能夠確保映射(mapping)準確度的同時通過改善映射過程的複雜度而提高處理速度的針對一對短片段的比對方案。
[0009]根據本發明一個實施例的鹼基序列比對系統，用於將包含第一序列及第二序列的一對鹼基序列比對到參考序列，包括:種子序列生成單元，從所述第一序列及所述第二序列中分別生成一個以上的片段，並由此構成第一種子序列集合及第二種子序列集合；映射值計算單元，將所述參考序列劃分為多個區間，並按所述多個區間分別計算包含於所述第一種子序列集合中的種子序列在對應區間內的映射值即第一映射值、以及包含於所述第二種子序列集合中的種子序列在對應區間內的映射值即第二映射值；比對單元，選擇計算出的所述第一映射值及所述第二映射值均為基準值以上的第一區間，並在所述第一區間內搜尋所述第一序列及所述第二序列的映射位置。
[0010]根據本發明另一實施例的鹼基序列比對系統，用於將包含第一序列及第二序列的一對鹼基序列比對到參考序列，包括:誤差估計單元，分別計算所述第一序列及所述第二序列的最小誤差估計值；比對單元，從所述第一序列與所述第二序列中選擇計算出的所述最小誤差估計值較小的序列，並計算該序列的相對所述參考序列的比對位置，且在以計算出的所述比對位置為基準而設定的可映射範圍內執行針對剩餘序列的全局比對。
[0011]根據本發明一個實施例的鹼基序列比對方法，用於在鹼基序列比對系統中將包含第一序列及第二序列的一對鹼基序列比對到參考序列，包括如下步驟:在種子序列生成單元中，從所述第一序列及所述第二序列中分別生成一個以上的片段，並由此構成第一種子序列集合及第二種子序列集合；在映射值計算單元中，將所述參考序列劃分為多個區間，並按所述多個區間分別計算包含於所述第一種子序列集合中的種子序列在對應區間內的映射值即第一映射值、以及包含於所述第二種子序列集合中的種子序列在對應區間內的映射值即第二映射值；在比對單元中，選擇計算出的所述第一映射值及所述第二映射值均為基準值以上的第一區間，並在所述第一區間內搜尋所述第一序列及所述第二序列的映射位置。
[0012]根據本發明另一實施例的鹼基序列比對方法，用於在鹼基序列比對系統中將包含第一序列及第二序列的一對鹼基序列比對到參考序列，包括如下步驟:在誤差估計單元中，分別計算所述第一序列及所述第二序列的最小誤差估計值；在比對單元中，從所述第一序列與所述第二序列中選擇計算出的所述最小誤差估計值較小的序列，並計算該序列的相對所述參考序列的比對位置；以及在所述比對單元中，在以計算出的所述比對位置為基準而設定的可映射範圍內執行針對剩餘序列的全局比對。
[0013]根據本發明的各實施例，在將雙末端短片段或配對短片段比對到參考序列時，預先選擇具有形成配對的可能性的區間，並在對應區間內執行針對所述雙末端短片段或配對短片段的比對，從而與現有方法相比可以顯著減少計算量。並且，還可以提供一種在對雙末端短片段或配對短片段進行比對時，不僅在特定鹼基被置換的情況下可以進行比對，而且在特定鹼基被插入或刪除而存在缺口( gap)狀不一致的情況下也可以進行比對的比對算法。
【專利附圖】
【附圖說明】
[0014]圖1為用於說明根據本發明一個實施例的鹼基序列比對方法100的圖。
[0015]圖2為用於舉例說明根據本發明一個實施例的鹼基序列比對方法100的步驟104中的最小誤差估計值(MEB)計算過程的圖。
[0016]圖3為用於詳細說明根據本發明一個實施例的鹼基序列比對方法100中的比對步驟114的順序圖。
[0017]圖4為用於詳細說明根據本發明一個實施例的鹼基序列比對方法100中的有效配對搜尋過程的順序圖。
[0018]圖5為示出根據本發明一個實施例的鹼基序列比對系統500的模塊圖。
[0019]圖6為示出根據本發明另一實施例的鹼基序列比對系統600的模塊圖。
[0020]符號說明:
[0021]500、600:喊基序列比對系統502:種子序列生成單兀
[0022]504:映射值計算單元506:比對單元
[0023]602:誤差估計單元604:比對單元【具體實施方式】
[0024]以下參照【專利附圖】
【附圖說明】本發明的【具體實施方式】。然而這僅僅是示例，本發明並不局限於此。
[0025]在對本發明進行說明時，如果遇到對有關本發明的公知技術的具體說明有可能不必要地幹擾本發明的主旨的情況，則省略其詳細說明。並且，後述的術語均為考慮本發明中的功能而進行定義的，其可能因使用者、運用人員的意圖或習慣等而有所不同。因此，要以整個說明書的內容為基礎對其進行定義。
[0026]本發明的技術思想由權利要求書確定，以下的實施例只是為了將本發明的技術思想有效地傳遞給本發明所屬【技術領域】中具有普通知識的人員而採用的一種手段。
[0027]在對本發明的實施例進行具體說明之前，首先對本發明中使用的術語進行如下說明。
[0028]首先，「短片段序列(read sequence)"(或者簡稱為「短片段(read)」)是指基因組測序儀(genome sequencer)中輸出的短鹼基序列數據。短片段的長度因基因組測序儀的種類而不同，通常構成為35?500bp (base pair)範圍的多種長度,在DNA鹼基的情況下，通常用A、C、G、T等四個字母表示。
[0029]在本發明的實施例中，基因組測序儀輸出一對配對(pair)的短片段。此時，將所述一對短片段中的第一個短片段稱為5』短片段而將第二個短片段稱為3』短片段，所述5』短片段與3』短片段的方向可形成為反向互補(reverse complement)關係(雙末端短片段)，或者形成為相同的方向(配對短片段)。例如，對於雙末端短片段而言，如果5』短片段為正向(forward)短片段,則3』短片段將是反向互補(reverse complement)短片段,與之相反，如果5』短片段為反向互補短片段，則3』短片段將是正向短片段。並且，對於配對短片段而言，如果5』短片段為正向短片段，則3』短片段也將是正向短片段，與之相反，如果5』短片段為反向互補方向的短片段，則3』短片段也將是反向互補方向的短片段。
[0030]「參考序列(reference sequence)"指可對利用所述短片段形成整個鹼基序列提供參考的鹼基序列。在鹼基序列分析中，通過將基因組測序儀所輸出的大量短片段參照參考序列進行映射而完成整個鹼基序列。在本發明中，所述參考序列既可以是鹼基序列分析時預先設定的序列(例如人類的整個鹼基序列等)，或者也可以將基因組測序儀中產生的鹼基序列使用為參考序列。
[0031]「鹼基(base)」為構成參考序列及短片段的最小單位。如前所述，構成DNA的鹼基可由A、C、G、T等四個字母表示的鹼基構成，將這些分別稱為鹼基。換言之，對於DNA而言，可用四種鹼基表示，短片段也是如此。
[0032]「片段序列(fragment sequence)"(或者簡稱為「片段」)指一種序列，該序列成為為了短片段的映射而比較短片段與參考序列時的單位。從理論上講，為了將短片段映射於參考序列，需要把整個短片段從參考序列的最前端部分開始依次比較並計算短片段的映射位置。然而，由於這種方法在映射一個短片段時消耗過多的時間並要求過高的計算能力，因此實際上要先把短片段的一部分所構成的片，即片段映射於參考序列而找出整個短片段的映射候選位置，然後將整個短片段映射於對應候選位置(Global Alignment)。
[0033]「種子序列(seed sequence)」指由短片段產生的片段中與參考序列相匹配的片段。即，在本發明的實施例中將由短片段產生的各片段分別與參考序列進行精確匹配(exact matching)，並進行用於將其中不與所述參考序列精確匹配的片段排除在外的篩選過程，且將所述篩選過程中精確匹配的片段作為種子序列，而將這些種子序列的集合稱為種子序列集合。此時，與所述參考序列相匹配的片段指與所述參考序列進行精確匹配(exact matching)時不一致的鹼基數為預先設定的允許值以下的片段。此時，如果所述允許值為O，則種子序列集合中只包含與所述參考序列精確匹配(S卩，沒有不一致的鹼基)的片段。
[0034]圖1為用於說明根據本發明一個實施例的鹼基序列比對方法100的圖。在本發明的實施例中，鹼基序列比對方法100指通過將基因組測序儀(genome sequencer)中輸出的一對短片段(雙末端短片段或配對短片段)與參考序列進行比較而確定對應短片段在所述參考序列中的映射(或比對)位置的一系列過程。在以下實施例中，將構成所述一對短片段的兩個短片段(5』短片段及3』短片段)分別稱為第一短片段及第二短片段。
[0035]首先，當從基因組測序儀(genome sequencer)接收到第一短片段及第二短片段(步驟102)時，分別針對輸入的兩個短片段的正向序列及反向互補序列而計算最小誤差估計值(MEB;Minimum Error Bound)(步驟104)。S卩，在本步驟中將分別計算包含第一短片段的正向序列、第一短片段的反向互補序列、第二短片段的正向序列、第二短片段的反向互補序列在內的四個序列的最小誤差估計值。此時，所述最小誤差估計值是指將所述各序列映射於參考序列時可能發生的誤差的最小值。
[0036]圖2為用於舉例表示所述步驟104中的MEB計算過程的圖。首先，如圖2的(a)所示，將初始MEB設定為0，並從對象序列的第一個鹼基開始向右側逐個移動而嘗試精確匹配。此時，如圖2 (b)所示，假定在對象序列的特定鹼基(在圖中以第二個T標記的地方)開始無法再實現精確匹配，則這種情況說明從序列的匹配起始位置到當前位置之間的區間中的某處出現了誤差。因此，在這種情況下將MEB值增加I (MEB=I)之後在下一個位置上重新開始精確匹配(在圖中標記為(c))。如果在以後又遇到判斷為無法精確匹配的情況，則是說明從重新開始精確匹配的位置到當前位置之間的區間某處又出現了誤差，因此又將MEB值增加I (MEB=2)之後在下一個位置上重新開始精確匹配(在圖中標記為(d))。通過這種過程，到達序列末端時的MEB值將成為對應序列的MEB值。
[0037]通過如上所述的過程，將分別計算出包括第一短片段的正向序列、第一短片段的反向互補序列、第二短片段的正向序列、第二短片段的反向互補序列在內的共計4個序列各自的MEB值。
[0038]然後將計算出的4個MEB值與預先設定的最大誤差允許值(maxError)進行比較(步驟106)。此時，如果計算出的4個MEB值均超過所述最大誤差允許值，則判定針對對應短片段的比對失敗。
[0039]與之相反，如果在所述步驟106中判斷的結果至少有部分序列的MEB為所述最大誤差允許值以下，則選擇計算出的MEB為最大誤差允許值以下的序列(步驟108)，並構造出所選擇序列各自的種子序列集合(步驟110)。然後將所述參考序列劃分為多個區間，並按所述多個區間分別計算所述所選擇序列的總映射值而生成映射直方圖(步驟112)，且利用所述映射直方圖而將所述一對短片段比對到所述參考序列(步驟114)。
[0040]以下詳細說明所述步驟110至步驟114的具體過程。[0041]由選擇的序列構成種子序列集合(步驟110)
[0042]該步驟為利用從所述步驟108中選擇的短片段序列生成一個以上的種子序列的步驟。首先，考慮所選擇序列的一部分或全部而生成多個片段。例如，可通過將所述序列的全部或特定區間分割為多個片或者將分割的片進行組合而生成片段。這種情況下，生成的片段可以連續地相連，然而並非一定要如此，也可以用序列內分離的片的組合構成片段。並且，生成的片段並非一定要具有相同的長度，在一個短片段內也可以生成多種長度的片段。總而言之，本發明中的由短片段序列生成片段的方法並不受到特別的局限，從短片段序列的一部分或全體中提取片段的各種算法均可不受限制地使用。
[0043]如果通過上述過程生成了與所選擇的各序列分別對應的片段，接著便通過從生成的片段中除去與參考序列不匹配的片段的篩選過程而構成種子序列集合。即，嘗試所生成的片段與所述參考序列之間的精確匹配(exact matching)，然後用不一致的鹼基數為預先設定的允許值以下的片段(種子序列)構成種子序列集合。此時，可通過適當考慮序列長度及從中提取的片段長度等而確定所述允許值。例如，序列長度較短的情況下(約為50bp以下)，優選為只對與所述參考序列精確匹配的片段予以考慮，此時所述允許值可以為O。而隨著序列長度變長，可將所述允許值提高為I或2等，從而防止映射的準確率過低。
[0044]牛成映射官方圖(步驟112)
[0045]如果通過上述過程構成了種子序列集合，接著便構成對應於各序列的映射直方圖(histogram)。在本發明中，映射直方圖為具有預定大小的整數陣列(integer array),整數陣列的值對應於將參考序列劃分為具有相同大小的多個區間時的各區間。例如，將參考序列劃分為具有65536(=216)bp大小的區間時,參考序列的O?65535bp的區間對應於映射直方圖h的第一個值h [O]，而65536?131071的區間對應於映射直方圖h的第二個值h [I]。可通過這種方式使參考序列的各劃分區間對應於映射直方圖。
[0046]而且，映射直方圖的各值h[i]中存儲有在對應的參考序列區間當中按各短片段序列分別提取的種子序列的總映射值。此時，所述映射值可以是對應參考序列區間中的所述種子序列的總映射長度。例如，假設從特定短片段序列中提取的種子序列中53-67種子序列(從所述短片段序列的第53?第67個鹼基提取的種子序列)及61-75種子序列被映射於映射直方圖的第一個區間，則在這種情況下對應區間的直方圖值將是23(=75-53+1)。
[0047]另外，所述映射值也可以是對應參考序列區間中的所述種子序列的總映射個數。在上面的例中，由於映射於映射直方圖的第一個區間的種子序列個數為2，因此對應區間的直方圖值將是2。並且，根據實施例的不同，也可以將多個區間各自的總映射長度及總映射個數作為所述映射值一併存儲。
[0048]一對短片段的比對(步驟114)
[0049]如果通過上述過程生成了第一短片段及第二短片段的序列各自的映射直方圖，則利用所生成的映射直方圖而將所述一對短片段比對到所述參考序列。
[0050]圖3為用於詳細說明根據本發明一個實施例的比對步驟(步驟114)的順序圖。
[0051]首先，判斷是否可以利用所述步驟106中選擇的短片段序列構成序列對(Sequence pair)(步驟 300)。
[0052]例如，在所述一對短片段為雙末端短片段的情況下，判斷MEB值為基準值(即最大誤差允許值)以下的序列是否可以構成如下配對中的至少一個。[0053](第一短片段的正向序列-第二短片段的反向互補序列)
[0054](第一短片段的反向互補序列-第二短片段的正向序列)
[0055]而如果所述一對短片段為配對短片段，則判斷MEB值為基準值(即最大誤差允許值)以下的序列是否可以構成如下配對中的至少一個。
[0056](第一短片段的正向序列-第二短片段的正向序列)
[0057](第一短片段的反向互補序列-第二短片段的反向互補序列)
[0058]如果在所述步驟300中判斷的結果為至少可以在上述配對中實現一種構成，則對構成序列對的兩個短片段序列的直方圖值進行比較，從而判斷是否存在兩個序列的直方圖值均為直方圖切值(Histogram Cut)以上的參考序列區間(步驟302)。
[0059]如果在所述步驟302中判斷的結果為存在兩個序列的直方圖值(映射值)均為直方圖切值(Histogram Cut, H)以上的參考序列區間的情況下，將對應區間選擇為映射對象區間(步驟304)，並在所選區間內對構成所述序列對的兩個短片段序列執行一次比對(步驟306、步驟308)。具體而言，在所述步驟306中，在所述映射對象區間內執行分別針對構成序列對的兩個短片段序列的全局比對(global alignment)，並將根據所述全局比對結果計算出的兩個短片段序列的比對位置對當中滿足預先設定的短片段間距範圍(插入大小，insert size)的比對位置對(有效配對,valid pair)選擇為所述第一短片段及所述第二短片段的比對位置。此時，所述有效配對要滿足如下三個條件。
[0060](I)兩個序列的比對方向要與初始輸入的一對短片段相同或對應。當輸入的一對短片段為雙末端短片段時各序列應具有反向互補關係。即，如果一個序列為正向序列，則另一個序列應該是反向互補序列。而且，當輸入的一對短片段為配對短片段時，兩個序列的比對方向應該相同。
[0061 ] (2)兩個序列中的至少一個應具有最大誤差允許值以下的誤差。
[0062](3)兩個序列的比對位置間距應該在預先設定的可映射範圍之內。此時，所述可映射範圍可用如下數學式I確定。
[0063][數學式I]
[0064]L「k.D ≤ L2 ≤ L^k.D
[0065](L1為構成序列對的第一個序列的映射位置，L2為第二個序列的映射位置，作為加權值的k具有大於O而小於1.8的值，D為預先設定的序列間的距離差(插入大小))
[0066]此時，之所以給所述插入大小賦予加權值k，是因為鹼基序列的特性決定了一些鹼基的插入或刪除可能導致序列間的距離發生變化，故使用加權值k來反映。
[0067]如果舉例說明搜尋所述有效配對的過程則如圖4所示。假定在圖示的映射對象區間內構成序列對的兩個序列中的第一序列映射於位置A和B，而第二序列映射於位置C。這種情況下將生成如下的兩個比對位置對。
[0068](A, C)
[0069](B, C)
[0070]假定所述A與C之間的插入大小(I1為1500bp,B與C之間的插入大小d2為650bp,基於所述數學式I的可映射範圍為_750bp~750bp。這種情況下，由於在兩個比對位置對當中滿足前述可映射範圍的是(B，C)，因此所述第一短片段和第二短片段的比對位置將是B和C。[0071]如上所述，將在所選區間內滿足前述範圍的比對位置對稱為有效配對(validpair)。即，在上述例中有效配對是(B，C)，如果找到該有效配對，則對應雙末端短片段的比對即成功。
[0072]而如果與之相反，在所述步驟304中在所選區間內進行一次比對的結果不存在有效配對，或者在所述步驟302中判斷的結果為不存在兩個序列的直方圖值均為H以上的區間，則將構成序列對的兩個序列中的某一個序列的直方圖值為H以上的區間選擇為映射對象區間(步驟310 )，並在所選映射對象區間當中執行二次比對(步驟312、步驟314 )。
[0073]對所述二次比對過程進行如下更為詳細的說明。首先，在兩個序列中選擇一個序列，並在所選序列的所述映射區間內計算比對位置。此時，所選的序列可以是兩個序列當中在對應映射對象區間內的直方圖值為H以上的序列。
[0074]然後，判斷剩餘序列在以計算出的所述比對位置為基準而設定的可映射範圍內是否得到映射(局部比對，local alignment)。即，判斷所述可映射範圍內是否存在滿足前述三個條件的有效配對。此時，所述可映射範圍如同前述數學式I所示。即，在該二次比對過程中是將直方圖值較大的序列作為一種錨點(anchor)來判斷剩餘序列在對應序列的周邊是否得到映射。
[0075]如果進行所述映射的結果存在有效配對，則對應的一對短片段的比對成功完成。而如果與之不同，即執行所述步驟312、步驟314的結果發現不存在有效配對，則所述短片段的比對失敗，在這種情況下將所述第一短片段及第二短片段分別在參考序列中進行全局比對，並輸出進行所述全局比對的結果當中比對分數(alignment score)最高的比對位置(步驟322)。此時，由於與各短片段的全局比對及比對分數的計算有關的內容為本發明所屬【技術領域】中的普通知識，因此省略其詳細說明。
[0076]另外，如果在所述步驟300中判斷的結果為無法構成兩個序列的MEB均為最大誤差允許值以下的序列對，接著便判斷兩個當中哪一序列的MEB是否為最大誤差允許值以下(步驟316)。此時，當在所述步驟316中判斷的結果某一序列的MEB為最大誤差允許值以下的情況下，計算MEB為最大誤差允許值以下的序列的相對所述參考序列的比對位置(步驟318, single end alignment)。然後,判斷在以計算出的所述比對位置為基準設定的可映射範圍內，是否存在剩餘序列滿足前述三個條件的有效配對(步驟320, local alignment)。此時，所述可映射範圍如同前述數學式I所示。即，在該二次比對過程中是將MEB為最大誤差允許值以下的序列作為一種錨點(anchor)來判斷剩餘序列在對應的序列周圍是否得到映射。
[0077]如果進行所述映射的結果存在有效配對，則對應的一對短片段的比對成功完成。而如果與之不同，即執行所述步驟318、320的結果不存在有效配對，則所述一對短片段的比對失敗，在這種情況下將所述第一短片段及第二短片段分別在參考序列中進行全局比對，並輸出進行所述全局比對的結果當中比對分數(alignment score)最高的比對位置(步驟322)。而且，在所述步驟316中判斷的結果所有序列的MEB值均超過最大誤差允許值的情形也進行與此相同的處理。
[0078]計算盲方圖切倌(Histogram Cut)
[0079]在上述實施例中，可通過如下方式計算直方圖切值。
[0080]首先，以映射於對應區間的種子序列的個數定義所述直方圖值(S卩，各區間內的映射值)的情況下，所述直方圖切值至少應該是2。這是由於映射的基本單位為種子序列，因此只映射一個種子序列的區間發生短片段映射的可能性很低。即，以映射於各區間的種子序列的個數定義所述直方圖值的情況下，可通過恰當考慮短片段長度、種子序列長度等而從2以上的整數值中確定所述直方圖切值。
[0081]其次，以映射於對應區間的種子序列的長度定義所述直方圖值的情況下，以如下方式計算直方圖切值。在f表示片段大小、s表示為了生成片段而在短片段中移動的距離、L表示短片段長度、e表示短片段中允許的最大誤差個數、H表示直方圖切值的情況下，可用如下數學式求出短片段中不受誤差影響的區域長度T。
[0082]T=L - f.e~s
[0083]此時，由於L和e為執行本發明時預先確定的值，因此由f和s的值決定T。即，算法的性能變化取決於f和s的值如何變化。
[0084]首先，在確定H的值時考慮以下兩個條件。其中，必須條件為必須要滿足的條件，而附加條件只在可能的情況下予以考慮。
[0085]必須條件:由於映射的基本單位為片段，因此無論直方圖切值多小，至少要具有能夠包含重疊(overlap)的兩個以上片段的大小。例如圖2所示，在f=15、s=4的情況下，由於重疊的兩個片段的最小長度為15+4=19，因此H值應該至少為19。而且，由於要將所述H值設定為至少包含兩個片段，因此H值至少要比f+s更大或相等。如後所述，f值應至少為15，因此將s值假定為其最小值I的情況下，H值至少為16 (=15+1)。
[0086]附加條件:在理想情況下，通過設定H=T並搜尋映射了 T以上的序列的直方圖，便可以找到對應於給定誤差的所有映射。然而，如前所述，在參考序列本身包含許多重複的情況下，根據情況可能會遇到需要增大片段長度的情形。因此，考慮到這一點，在確定H值時使用比T略小的T-s可能有利`於映射率。如果假定H=T，則H=L-f.e-s，如果假定其中的e取最小值I (由於e=0的情況為與參考序列精確匹配的情形，因此此時將在前述步驟104中映射完畢)，則有H=L-f-s。該值將是直方圖值的最大值。如果假定L=75bp、f=15bp、s=l，H的最大值便成為75-15-1=59。
[0087]綜上,所述H值應該滿足如下範圍。
[0088]f+s ^ H ^ L - (f+s)
[0089]然後，在滿足以下兩個條件的值當中選擇較大值作為f值。必須條件仍然要必須滿足，而附加條件只在可能的情況下考慮。
[0090]必須條件:f應該取15以上，這是由於如果片段長度為14以下，則參考序列中的映射位置的個數將急劇增加。
[0091]如下的表1表示根據片段長度的在人類基因組中的片段平均出現頻率。
[0092][表 I]
[0093]
【權利要求】
1.一種鹼基序列比對系統，用於將包含第一序列及第二序列的一對鹼基序列比對到參考序列，包括: 種子序列生成單元，從所述第一序列及所述第二序列中分別生成一個以上的片段，並由此構成第一種子序列集合及第二種子序列集合； 映射值計算單元，將所述參考序列劃分為多個區間，並按所述多個區間分別計算包含於所述第一種子序列集合中的種子序列在對應區間內的映射值即第一映射值、以及包含於所述第二種子序列集合中的種子序列在對應區間內的映射值即第二映射值； 比對單元，從所述多個區間中選擇計算出的所述第一映射值及所述第二映射值均為基準值以上的第一區間，並在所述第一區間內搜尋所述第一序列及所述第二序列的映射位置。
2.如權利要求1所述的鹼基序列比對系統，其特徵在於，所述第一種子序列集合只包括從所述第一序列提取的一個以上的片段中與所述參考序列相匹配的片段，而所述第二種子序列集合只包括從所述第二序列提取的一個以上的片段中與所述參考序列相匹配的片段。
3.如權利要求2所述的鹼基序列比對系統，其特徵在於，與所述參考序列相匹配的片段為進行與所述參考序列之間的精確匹配的結果，不一致的鹼基數為設定個數以下的片段。
4.如權利要求1所述的鹼基序列比對系統，其特徵在於，所述映射值計算單元基於包含在所述第一種子序列集合中的種子序列在對應區間內的總映射長度而計算所述第一映射值，並基於包含在所述第二種子序列集合中的種子序列在對應區間內的總映射長度而計算所述第二映射值。
5.如權利要求1所述的鹼基序列比對系統，其特徵在於，所述映射值計算單元基於包含在所述第一種子序列集合中的種子序列在對應區間內的總映射個數而計算所述第一映射值，並基於包含在所述第二種`子序列集合中的種子序列在對應區間內的總映射個數而計算所述第二映射值。
6.如權利要求1所述的鹼基序列比對系統，其特徵在於，所述比對單元在所述第一區間內執行針對所述第一序列及所述第二序列的全局比對，並將進行所述全局比對的結果計算出的所述第一序列及所述第二序列的比對位置對當中滿足預先設定的序列之間的距離範圍的比對位置對選擇為所述第一序列及所述第二序列的比對位置。
7.如權利要求1所述的鹼基序列比對系統，其特徵在於，當無法選擇所述第一區間時，所述比對單元便選擇所述第一映射值與所述第二映射值中的某一映射值為基準值以上的第二區間，並在選擇的所述第二區間內搜尋所述第一序列及所述第二序列的映射位置。
8.如權利要求7所述的鹼基序列比對系統，其特徵在於，所述比對單元在所述第二區間內計算相對從所述第一序列與所述第二序列中選擇的序列的比對位置，並在以計算出的所述比對位置為基準而設定的可映射範圍內執行針對剩餘序列的全局比對。
9.如權利要求8所述的鹼基序列比對系統，其特徵在於，所述選擇的序列為所述第一序列與所述第二序列當中在所述第二區間內的映射值更大的序列。
10.如權利要求8所述的鹼基序列比對系統，其特徵在於，所述可映射範圍是以所述選擇的序列的映射位置為基準向所述參考序列的前後端各延伸kXD長度的區間，其中，k為加權值，D為預先設定的序列之間的距離。
11.如權利要求10所述的鹼基序列比對系統，其特徵在於，所述加權值k為1.8以下。
12.—種鹼基序列比對系統，用於將包含第一序列及第二序列的一對鹼基序列比對到參考序列，包括: 誤差估計單元，分別計算所述第一序列及所述第二序列的最小誤差估計值； 比對單元，從所述第一序列與所述第二序列中選擇計算出的所述最小誤差估計值較小的序列，並計算該序列的相對所述參考序列的比對位置，且在以計算出的所述比對位置為基準而設定的可映射範圍內執行針對剩餘序列的全局比對， 其中，所述誤差估計單元由從所述第一序列與所述第二序列中選擇的序列的第一個鹼基開始以一個鹼基為單位逐個移動而將所述選擇的序列與所述參考序列進行精確匹配，而在所述選擇的序列的特定位置上無法實現精確匹配的情況下，從對應位置的下一個鹼基開始以一個鹼基為單位逐個移動的同時重新執行精確匹配，並在到達所述選擇的序列的末尾鹼基時，將判斷為無法實現精確匹配的位置的個數設定為所述選擇的序列的最小誤差估計值。
13.—種鹼基序列比對方法，用於在鹼基序列比對系統中將包含第一序列及第二序列的一對鹼基序列比對到參考序列，包括如下步驟: 在種子序列生成單元中，從所述第一序列及所述第二序列中分別生成一個以上的片段，並由此構成第一種子序列集合及第二種子序列集合； 在映射值計算單元中，將所述參考序列劃分為多個區間，並按所述多個區間分別計算包含於所述第一種子序列集合中的種子序列在對應區間內的映射值即第一映射值、以及包含於所述第二種子序列集合中的種子序列在對應區間內的映射值即第二映射值； 在比對單元中，從所述多個區`間中選擇計算出的所述第一映射值及所述第二映射值均為基準值以上的第一區間，並在所述第一區間內搜尋所述第一序列及所述第二序列的映射位置。
14.如權利要求13所述的鹼基序列比對方法，其特徵在於，所述第一種子序列集合只包括從所述第一序列提取的一個以上的片段中與所述參考序列相匹配的片段，而所述第二種子序列集合只包括從所述第二序列提取的一個以上的片段中與所述參考序列相匹配的片段。
15.如權利要求14所述的鹼基序列比對方法，其特徵在於，與所述參考序列相匹配的片段為進行與所述參考序列之間的精確匹配的結果，不一致的鹼基數為設定個數以下的片段。
16.如權利要求13所述的鹼基序列比對方法，其特徵在於，在進行所述計算的步驟中，基於包含在所述第一種子序列集合中的種子序列在對應區間內的總映射長度而計算所述第一映射值，並基於包含在所述第二種子序列集合中的種子序列在對應區間內的總映射長度而計算所述第二映射值。
17.如權利要求13所述的鹼基序列比對方法，其特徵在於，在進行所述計算的步驟中，基於包含在所述第一種子序列集合中的種子序列在對應區間內的總映射個數而計算所述第一映射值，並基於包含在所述第二種子序列集合中的種子序列在對應區間內的總映射個數而計算所述第二映射值。
18.如權利要求13所述的鹼基序列比對方法，其特徵在於，搜尋所述映射位置的過程包括如下步驟: 在所述第一區間內執行針對所述第一序列及所述第二序列的全局比對； 將進行所述全局比對的結果計算出的所述第一序列及所述第二序列的比對位置對當中滿足預先設定的序列之間的距離範圍的比對位置對選擇為所述第一序列及所述第二序列的比對位置。
19.如權利要求13所述的鹼基序列比對方法，其特徵在於，搜尋所述映射位置的過程還包括如下步驟: 當無法選擇所述第一區間時，選擇所述第一映射值與所述第二映射值中的某一映射值為基準值以上的第二區間，並在選擇的所述第二區間內搜尋所述第一序列及所述第二序列的映射位置。
20.如權利要求19所述的鹼基序列比對方法，其特徵在於，在搜尋所述映射位置的過程中，在所述第二區間內計算相對從所述第一序列與所述第二序列中選擇的序列的比對位置，並在以計算出的所述比對位置為基準而設定的可映射範圍內執行針對剩餘序列的全局比對，且所述選擇的序列為所述第一序列與所述第二序列當中在所述第二區間內的映射值更大的序列。
【文檔編號】G06F19/22GK103793626SQ201310367008
【公開日】2014年5月14日 申請日期:2013年8月21日 優先權日:2012年10月29日
【發明者】樸旻胥, 樸商賢, 呂潤九 申請人:三星Sds株式會社, 延世大學校產學協力團