DNA序列的k-merindex問(wèn)題-基于Hash算法快速檢索

1、2015山東科技大學(xué)數(shù)學(xué)建模競(jìng)賽承 諾 書我們仔細(xì)閱讀了山東科技大學(xué)數(shù)學(xué)建模競(jìng)賽的競(jìng)賽規(guī)則.我們完全明白，在競(jìng)賽開始后參賽隊(duì)員不能以任何方式（包括電話、電子郵件、網(wǎng)上咨詢等）與隊(duì)外的任何人研究、討論與賽題有關(guān)的問(wèn)題。我們知道，抄襲別人的成果是違反競(jìng)賽規(guī)則的, 如果引用別人的成果或其他公開的資料（包括網(wǎng)上查到的資料），必須按照規(guī)定的參考文獻(xiàn)的表述方式在正文引用處和參考文獻(xiàn)中明確列出。我們鄭重承諾，嚴(yán)格遵守競(jìng)賽規(guī)則，以保證競(jìng)賽的公正、公平性。如有違反競(jìng)賽規(guī)則的行為，我們將受到嚴(yán)肅處理。我們參賽選擇的題號(hào)是（從A/B/C中選擇一項(xiàng)填寫）： 我們的參賽序號(hào)為： 所屬學(xué)院（請(qǐng)?zhí)顚懲暾娜?參賽隊(duì)員

2、 (打印并簽名) ：1. 2. 3. 日期： 年 月 日基于Hash表在大量DNA序列中快速索引查找k-mer 的算法摘要：為了解決在大量DNA數(shù)據(jù)中快速查找到k-mer所在位置，我們基于Hash算法思想建立適合此題的快速索引方法（桶式定址法），利用四進(jìn)制轉(zhuǎn)十進(jìn)制的方式（）取得關(guān)鍵碼值建立哈希表進(jìn)行索引查詢，8G內(nèi)存限制下在codeblocks集成開發(fā)環(huán)境中，借助C語(yǔ)言進(jìn)行編寫使k支持114。針對(duì)問(wèn)題將依次進(jìn)行下列敘述：對(duì)建立索引的算法進(jìn)行敘述和沖突分析；對(duì)建立索引算法的計(jì)算復(fù)雜度和空間復(fù)雜度進(jìn)行分析；對(duì)索引查詢進(jìn)行敘述及性能分析；分析整套算法程序在不同k值下內(nèi)存占用及極限分析?？偨Y(jié)分析整套索

3、引算法性能，對(duì)算法進(jìn)行缺陷分析及改進(jìn)方案。關(guān)鍵詞：索引算法、Hash表、k-mer快速索引、數(shù)據(jù)結(jié)構(gòu)一、問(wèn)題重述1.1背景 給定一個(gè)DNA序列，這個(gè)系列只含有4個(gè)字母ATCG，如 S =“CTGTACTGTAT”。給定一個(gè)整數(shù)值k，從S的第一個(gè)位置開始，取一連續(xù)k個(gè)字母的短串，稱之為k-mer（如k= 5，則此短串為CTGTA）， 然后從S的第二個(gè)位置， 取另一k-mer（如k= 5，則此短串為TGTAC），這樣直至S的末端，就得一個(gè)集合，包含全部k-mer 。 如對(duì)序列S來(lái)說(shuō)，所有5-mer為CTGTA，TGTAC，GTACT，TACTG，ACTGT，TGTAT通常這些k-mer需一種數(shù)據(jù)索

4、引方法，可被后面的操作快速訪問(wèn)。例如，對(duì)5-mer來(lái)說(shuō)，當(dāng)查詢CTGTA，通過(guò)這種數(shù)據(jù)索引方法，可返回其在DNA序列S中的位置為1，6。1.2問(wèn)題現(xiàn)在以文件形式給定 100萬(wàn)個(gè) DNA序列，序列編號(hào)為1-，每個(gè)基因序列長(zhǎng)度為100 。（1）要求對(duì)給定k， 給出并實(shí)現(xiàn)一種數(shù)據(jù)索引方法，可返回任意一個(gè)k-mer所在的DNA序列編號(hào)和相應(yīng)序列中出現(xiàn)的位置。每次建立索引，只需支持一個(gè)k值即可，不需要支持全部k值。（2）要求索引一旦建立，查詢速度盡量快，所用內(nèi)存盡量小。（3）給出建立索引所用的計(jì)算復(fù)雜度，和空間復(fù)雜度分析。 （4）給出使用索引查詢的計(jì)算復(fù)雜度，和空間復(fù)雜度分析。（5）假設(shè)內(nèi)存限制為8G，

5、分析所設(shè)計(jì)索引方法所能支持的最大k值和相應(yīng)數(shù)據(jù)查詢效率。（6）按重要性由高到低排列，將依據(jù)以下幾點(diǎn)，來(lái)評(píng)價(jià)索引方法性能 索引查詢速度 索引內(nèi)存使用 8G內(nèi)存下，所能支持的k值范圍 建立索引時(shí)間二、問(wèn)題分析在生物技術(shù)快速發(fā)展的今天，人類分析人類自身編碼的需求也越來(lái)越高，人們利用計(jì)算機(jī)來(lái)處理分析大量DNA序列，k-mer快速查找也是處理大量數(shù)據(jù)的問(wèn)題，所以必須依賴數(shù)據(jù)結(jié)構(gòu)原理，建立模型構(gòu)造算法，從而利用有限的資源解決復(fù)雜問(wèn)題。針對(duì)問(wèn)題一：按照給定k值，將所有數(shù)據(jù)按題目要求分組，求出每組數(shù)據(jù)的關(guān)鍵碼值，并將關(guān)鍵碼值與此組k-mer所在位置建立對(duì)應(yīng)關(guān)系并存儲(chǔ)到表中，從而建立哈希表。針對(duì)問(wèn)題二：將要查找

6、的k-mer序列求出其關(guān)鍵碼值，直接輸出其關(guān)鍵碼值在表中對(duì)應(yīng)位置信息，大大加快了索引查詢速度。針對(duì)問(wèn)題三：詳見四-（二）-2，3。針對(duì)問(wèn)題四：詳見四-（三）-2，3。針對(duì)問(wèn)題五：根據(jù)k值大小動(dòng)態(tài)分配內(nèi)存建立哈希表，最終實(shí)現(xiàn)k支持114的范圍，因?yàn)橹苯雨P(guān)鍵碼值尋址輸出，所以索引查詢速度非?？臁ａ槍?duì)問(wèn)題六：按照重要性首先考慮索引查詢速度，其次動(dòng)態(tài)內(nèi)存分配盡量減少索引對(duì)內(nèi)存的消耗，在8G內(nèi)存限制下，使k值支持114，最后優(yōu)化添加計(jì)數(shù)器記錄已經(jīng)存在地址的k-mer個(gè)數(shù)，倘若達(dá)到所有k-mer種類數(shù)，則停止建立索引，索引成功建立。三、符號(hào)說(shuō)明符號(hào)符號(hào)說(shuō)明H(x)關(guān)鍵碼值生成函數(shù)，其中x代表一個(gè)k-mer

7、Ci代表A，T，C，G中的任意一個(gè)Ii與Ci相對(duì)應(yīng)的四進(jìn)制數(shù)k一個(gè)k-mer的長(zhǎng)度M內(nèi)存空間占用（單位：GB）四、算法設(shè)計(jì)思路及性能分析（代碼見附錄一）(一) 哈希表設(shè)計(jì)：1、k-mer關(guān)鍵碼值生成函數(shù)H(x)由于DNA序列由4個(gè)字母排列而成，所以每個(gè)k-mer都是一個(gè)四進(jìn)制數(shù)，H(x)函數(shù)根據(jù)這個(gè)特征將四進(jìn)制數(shù)轉(zhuǎn)為十進(jìn)制數(shù)作為哈希表關(guān)鍵碼值。 x= “CTGTA”如上圖為每個(gè)字母代表的四進(jìn)制數(shù)字，例如一個(gè)5-mer “CTGTA” 可以表示為四進(jìn)制數(shù)21310，其十進(jìn)制表示為628，628即為5-mer “CTGTA”在哈希表中的關(guān)鍵碼值。關(guān)鍵碼值計(jì)算的一般公式為: （1）2、哈希表結(jié)構(gòu)將k

8、-mer通過(guò)公式（1）轉(zhuǎn)換為十進(jìn)制的關(guān)鍵碼值存入哈希表中的關(guān)鍵碼值一列，并將關(guān)鍵碼值與此k-mer所在位置建立對(duì)應(yīng)關(guān)系，從而便于索引尋址。這里采用的方法是：桶定址法。桶：一片足夠大的存儲(chǔ)空間。桶定址：為表中的每個(gè)地址關(guān)聯(lián)一個(gè)桶。如果桶已經(jīng)滿了，可以使用開放定址法來(lái)處理。如圖。3、沖突處理方法開放地址法有一個(gè)公式：Hi=(H(key)+di) MOD m i=1,2,.,k(k=m-1）其中，m為哈希表的表長(zhǎng)。di 是產(chǎn)生沖突的時(shí)候的增量序列。如果di值可能為1,2,3,.m-1，稱線性探測(cè)再散列。如果di取1，則每次沖突之后，向后移動(dòng)1個(gè)位置.如果di取值可能為1,-1,4,-4,9,-9,1

9、6,-16,.k*k,-k*k(k=m/2）稱二次探測(cè)再散列。(二) 建立索引的算法、計(jì)算復(fù)雜度及空間復(fù)雜度分析1、 算法分析hash表初始化時(shí)，根據(jù)用戶輸入的k值， 計(jì)算出存儲(chǔ)哈希表需要的空間，利用內(nèi)存動(dòng)態(tài)分配函數(shù)動(dòng)態(tài)分配內(nèi)存，k-mer所在行數(shù)為100 0000以內(nèi)所以我們采用int型（占用4字節(jié)）存儲(chǔ)，而其在行的第幾個(gè)位置數(shù)在100以內(nèi)，我們則采用char（占用1字節(jié)）型存儲(chǔ)，充分減少空間損耗。哈希表初始化代碼：求關(guān)鍵碼值如上敘述采用四進(jìn)制轉(zhuǎn)十進(jìn)制數(shù)作為關(guān)鍵碼值存入哈希表用于位置索引。求關(guān)鍵碼值代碼：建立哈希表過(guò)程中先將傳入的長(zhǎng)度為100的字符串分成每k個(gè)字符一段的串，每分好一個(gè)就求一個(gè)

10、關(guān)鍵碼值存儲(chǔ)哈希表，為了提高建立索引的速度，每個(gè)關(guān)鍵碼值只對(duì)應(yīng)一個(gè)位置信息，也就是說(shuō)索引查詢時(shí)只返回k-mer的一個(gè)位置信息。倘若每個(gè)關(guān)鍵碼值都有其對(duì)應(yīng)的位置信息了，就退出索引建立，索引建立完成。建立哈希表代碼：2、計(jì)算復(fù)雜度分析不考慮8G空間限制，k在150的范圍內(nèi)每條DNA序列分段的循環(huán)不斷增加，計(jì)算復(fù)雜度也隨之增加，而50100范圍內(nèi)循環(huán)則又開始減少。根據(jù)計(jì)算復(fù)雜度定義1，建立索引的計(jì)算復(fù)雜度為O(-k2+100k)3、空間復(fù)雜度分析邏輯上哈希表空間占用計(jì)算公式：M=(54k)/(102410241024)（每條索引信息邏輯空間占用為5Byte）k值建立索引后空間占用（GB）10. 20

11、. 30. 40. 50. 60. 70. 80. 90. 100. 110. 120. 130. 141. 155. 1620. 由以上圖表可以清晰的看到，索引建立的內(nèi)存消耗隨著k的增長(zhǎng)冪指數(shù)增長(zhǎng)，與k-mer排列組合種類數(shù)（4k）相對(duì)應(yīng)，所以減少每條索引數(shù)據(jù)的存儲(chǔ)空間，可大大增加k的取值范圍，k-mer所在行數(shù)為100 0000以內(nèi)所以我們采用int型（占用4字節(jié)）存儲(chǔ)，而其在行的第幾個(gè)位置數(shù)在100以內(nèi)，我們則采用char（占用1字節(jié)）型存儲(chǔ)，盡可能的縮小了每條索引數(shù)據(jù)的存儲(chǔ)空間，最終使k支持114。但由上表可以看出，邏輯上k=15也能支持，但設(shè)備資源有限，未能實(shí)際測(cè)試。k=14時(shí)，索引

12、建立實(shí)際內(nèi)存消耗為KB=1.GB，與邏輯分析基本一致。如下圖：(三) 索引查詢的算法、計(jì)算復(fù)雜度及空間復(fù)雜度分析1、算法分析索引查詢需要將要查找的k-mer的關(guān)鍵碼值求出來(lái)，然后直接從哈希表中輸出對(duì)應(yīng)關(guān)鍵碼值的位置信息，求關(guān)鍵碼值的函數(shù)同建立索引過(guò)程的函數(shù)，不再重復(fù)給出。索引查詢代碼如下：2、計(jì)算復(fù)雜度分析索引查詢的過(guò)程中只有求關(guān)鍵碼值時(shí)需要進(jìn)行計(jì)算，其計(jì)算復(fù)雜度為：O(k)3、空間復(fù)雜度分析索引查詢過(guò)程用到只一個(gè)存關(guān)鍵碼值的中間變量，根據(jù)空間復(fù)雜度定義2（臨時(shí)占用存儲(chǔ)空間大小的量度）可知，索引查詢過(guò)程中空間復(fù)雜度為：O(1)(四) 整套索引算法程序在不同k值下內(nèi)存占用及極限分析k值k-mer

13、種類數(shù)建立索引耗時(shí)（秒）索引查詢耗時(shí)（秒）1400216003640042560.0010510240.0010640960.00807163840.03108655360.17091.0680108.2901133.55801236.08301339.31601443.9470測(cè)量上表耗時(shí)數(shù)據(jù)設(shè)備參數(shù)：CPU：Intel 酷睿i3 3110M 2.4GHz 雙核心/四線程內(nèi)存（RAM）：DDR3 1600MHz 8G硬盤：500G 5400轉(zhuǎn)比對(duì)以上兩折線圖不難發(fā)現(xiàn)建立索引耗時(shí)與k-mer種類數(shù)目變化趨勢(shì)并不是一致的，原因有兩點(diǎn)，第一點(diǎn)建立索引用事在9秒前幾乎為零因?yàn)榇怂惴ㄔ诖_認(rèn)所有種類k-

14、mer都有對(duì)應(yīng)位置后索引就建立完成不管是否將所有數(shù)據(jù)都遍歷過(guò)，9秒前k-mer種類較少，所以沒有檢索所有數(shù)據(jù)就已經(jīng)完成了索引建立；第二點(diǎn)建立索引耗時(shí)并沒有跟著k-mer種類指數(shù)上升而是放緩，原因是k-mer種類太多，將所有數(shù)據(jù)都遍歷結(jié)束后，有的k-mer也沒有找到與其匹配的位置，所以耗時(shí)基本就是遍歷所有數(shù)據(jù)的時(shí)間。(五) 缺陷分析及改進(jìn)方案1、 缺陷分析建立索引過(guò)程中會(huì)有一個(gè)k-mer與多個(gè)位置對(duì)應(yīng)，由于hash表存儲(chǔ)限制發(fā)生了沖突，只能返回k-mer的一個(gè)位置信息，不能將k-mer的所有位置信息都返回。2、 改進(jìn)方案 針對(duì)這個(gè)缺陷，可以采取hash算法沖突處理方法開放定址法，或采用hash的改進(jìn)算法tri-樹（字典樹）算法，在這里我們采用沖突處理的方案對(duì)此算法進(jìn)行改進(jìn)，在內(nèi)存有限制的情況下，放棄一點(diǎn)k值的范圍，騰出部分空間來(lái)存儲(chǔ)更多k-mer位置信息，從而彌補(bǔ)之前的缺陷，代碼見附錄二。五、參考文獻(xiàn)1創(chuàng)建者：flyingmyself時(shí)間復(fù)雜度/link?url=4RJO-qHmJs867JEK