陳梅麗，馬英克，李茹姣，鮑一明

1.中國科學(xué)院北京基因組研究所（國家生物信息中心），國家基因組科學(xué)數(shù)據(jù)中心和中國科學(xué)院基因組科學(xué)與信息重點(diǎn)實(shí)驗(yàn)室，北京 100101

2.中國科學(xué)院大學(xué)，未來技術(shù)學(xué)院，北京 100049

引言

隨著人類基因組測序計(jì)劃的完成，基因組學(xué)的影響力迅速擴(kuò)大，數(shù)以萬計(jì)的動(dòng)物、植物、微生物基因組被組裝[1-4]，游離DNA 應(yīng)用于無創(chuàng)產(chǎn)前檢測、通過基因檢測指導(dǎo)靶向藥物治療成為可能[5-6]、單細(xì)胞技術(shù)應(yīng)用于輔助生殖[7]、DNA 編輯技術(shù)廣泛應(yīng)用[8]、體外提供造血干細(xì)胞[9]、長壽基因被找到[10]、抗病蟲害和高產(chǎn)的農(nóng)業(yè)作物新品種被不斷培育[11-12]。這些基因組科學(xué)領(lǐng)域的巨大進(jìn)展，一方面來自于測序和實(shí)驗(yàn)技術(shù)的革新，同時(shí)也依賴于為適應(yīng)實(shí)驗(yàn)和測序技術(shù)進(jìn)步而不斷發(fā)展的分析手段和方法[4,13-14]。

基因組學(xué)測序數(shù)據(jù)增長迅猛，加快數(shù)據(jù)分析速度，提高數(shù)據(jù)處理效率，是對大數(shù)據(jù)整合分析工具和算法開發(fā)的迫切需求。如何用好多源海量基因組學(xué)數(shù)據(jù)，去除異質(zhì)性，并對其進(jìn)行整合分析和深度挖掘，是對數(shù)據(jù)分析工具和算法開發(fā)另一個(gè)層面的新要求?？茖W(xué)家們也在不斷地開發(fā)各種算法和工具來提高計(jì)算效率，比如第三代測序數(shù)據(jù)組裝算法wtdbg2，將拼接的分析速度提高5倍，少于數(shù)據(jù)產(chǎn)出時(shí)間[4]。而在此基礎(chǔ)上，人工智能等先進(jìn)技術(shù)也被廣泛應(yīng)用于基因組學(xué)大數(shù)據(jù)分析的工具開發(fā)[15]。本文將詳細(xì)地闡述隨著測序技術(shù)的發(fā)展，基因組、轉(zhuǎn)錄組、表觀組數(shù)據(jù)的分析算法和工具開發(fā)的現(xiàn)狀，以及大數(shù)據(jù)時(shí)代基因組學(xué)數(shù)據(jù)算法和工具開發(fā)在未來將面臨的問題和挑戰(zhàn)。

1 基因組測序數(shù)據(jù)分析

隨著測序技術(shù)的發(fā)展，各類基因組相關(guān)研究計(jì)劃接踵而來[16-17]，為生物多樣性、物種進(jìn)化、分子育種、臨床治療等研究提供了寶貴的數(shù)據(jù)資源?；蚪M數(shù)據(jù)分析主要包括基因組組裝、基因組注釋、基因組變異分析等。基因組序列和注釋信息包含了生物體的所有遺傳信息和功能信息，是多種組學(xué)研究的重要基礎(chǔ)數(shù)據(jù)[18]。通過基因組變異分析可以解析基因組變異對表型、物種進(jìn)化、疾病等的影響[19]。但是測序數(shù)據(jù)分析時(shí)間往往遠(yuǎn)大于測序數(shù)據(jù)產(chǎn)出時(shí)間，無法匹配數(shù)據(jù)爆發(fā)式增長的趨勢，是基因組數(shù)據(jù)分析面臨的最大挑戰(zhàn)之一。

1.1 基因組組裝

基因組組裝是將測序產(chǎn)生的讀段（read）片段經(jīng)過序列組裝成完整的基因組序列?；蚪M組裝方法主要有兩類：（1）基于參考基因組序列比對的有參組裝，常用于重測序和線粒體/葉綠體等保守細(xì)胞器基因組[20-21]；（2）從頭（de novo）組裝，目前主要有純二代測序組裝[22]、二代和三代測序混合組裝[23]、純?nèi)鷾y序組裝[24]等組裝策略。對于動(dòng)植物等大基因組可以利用遺傳圖譜、Bionano 光學(xué)圖譜、Hi-C、10X Genomics 等圖譜信息進(jìn)行整合輔助組裝，將基因組組裝提升到染色體級別[25-28]。動(dòng)植物基因組因存在基因組大、雜合度高、GC 含量高、重復(fù)序列多和多倍化水平較高等復(fù)雜因素，給組裝帶來了很大的挑戰(zhàn)，需要開發(fā)出兼顧效率和計(jì)算資源消耗，并且在重復(fù)區(qū)域獲得連續(xù)性和完整性都表現(xiàn)很好的高質(zhì)量基因組組裝結(jié)果的算法和工具。目前常用的軟件有SOAPdenovo2[22]、ALLPATHS-LG[29]、Canu[30]、Falcon[31]等。

為了解決大型基因組組裝效率、準(zhǔn)確性要求高和計(jì)算資源消耗大的問題，阮玨團(tuán)隊(duì)提出了三代數(shù)據(jù)組裝wtdbg2 算法，該算法遵循overlap-layoutconsensus 模式，以快速的讀段“全部對全部”比對實(shí)現(xiàn)方式和基于模糊布魯因圖這種新組裝圖理論——一種與稀疏布魯因圖和其變體A-Bruijn 圖有關(guān)的序列組裝的新數(shù)據(jù)結(jié)構(gòu)，來改進(jìn)現(xiàn)有的組裝程序并提高組裝效率[4]。在一臺計(jì)算機(jī)上，可在2 天內(nèi)完成4個(gè)～30X的人類基因組數(shù)據(jù)集的組裝，極大提高三代測序數(shù)據(jù)的分析效率。與此前提出的Flye 算法[32]相比，wtdbg2的分析速度提升了5倍，內(nèi)存使用僅為Flye的一半，組裝連續(xù)性和精度可與Canu、Falcon、Flye 等其他算法相媲美。wtdbg2 首次將測序數(shù)據(jù)分析時(shí)間降低到少于測序數(shù)據(jù)產(chǎn)出時(shí)間，是一種兼具高容錯(cuò)和高準(zhǔn)確的高效算法，并可擴(kuò)展到超大基因組，如32 Gb 大小的蠑螈基因組[4]。

當(dāng)前Hi-C 技術(shù)越來越多地應(yīng)用于輔助染色體水平二倍體基因組組裝[25]，Hi-C 技術(shù)是基于將線性距離遠(yuǎn)、空間結(jié)構(gòu)近的DNA 片段進(jìn)行交聯(lián)，并將交聯(lián)的DNA 片段富集后進(jìn)行高通量測序，對測序數(shù)據(jù)分析揭示染色質(zhì)的遠(yuǎn)程相互作用，明確基因組草圖中scaffold和染色體對應(yīng)關(guān)系、scaffold之間的連接順序和方向，從而將scaffold 掛載到染色體上。但是對于同源多倍體和近期加倍的異源多倍體來說，其同源染色體之間的Hi-C 交聯(lián)信號會將序列相似的等位片段連接在一起，導(dǎo)致同源染色體被錯(cuò)誤地連接到一起，形成大量嵌合的組裝，給組裝造成了較大困難。針對該問題明瑞光團(tuán)隊(duì)提出了ALLHiC 算法，該算法包括pruning、partition、rescue、optimization、building 5個(gè)步驟。ALLHiC 算法一方面通過修剪Hi-C 平行信號和弱信號進(jìn)行等位基因分型，減少了同源染色體間的嵌合連接；另一方面通過遺傳算法隨機(jī)優(yōu)化，極大地提高了contig 序列的排序和定向準(zhǔn)確性，成功解決了同源高倍體甘蔗、菠蘿和異源四倍體栽培花生等多倍體組裝難題[33-34]。

而針對于二倍體或者是多倍體，單體型基因組組裝是基因組組裝的最終目標(biāo)，目前已提出單體型區(qū)塊組裝策略，如FALCON-Unzip[31]、trio binning[35]等方法，嘗試將兩個(gè)或多個(gè)配子來源的染色體進(jìn)行分別組裝。但是當(dāng)前已有的單體型區(qū)塊組裝策略依然無法很好地解決單條染色體長度的單體型基因組組裝問題。針對該問題，Kronenberg 團(tuán)隊(duì)提出了一種單體型基因組組裝的新技術(shù)FALCON-Phase，可把雜合基因組的父本母本分型組裝，并可將其應(yīng)用于野外采集的樣本或缺乏譜系信息的生物。其主要原理和流程為：將基因組區(qū)域中顯示出高水平雜合性的區(qū)域鑒定為單體型區(qū)塊contig，基于鑒定結(jié)果對所有contig 拆分和打斷，通過將Hi-C 數(shù)據(jù)集成到圖形數(shù)據(jù)結(jié)構(gòu)中構(gòu)建標(biāo)準(zhǔn)化的互作矩陣，最后經(jīng)過Hi-C 輔助組裝掛載獲得單條染色體長度的單體型基因組組裝[36]。

優(yōu)化算法解決非橋式重復(fù)等重復(fù)序列導(dǎo)致的組裝連續(xù)性問題[32]，通過開發(fā)出對雜合子感知的一致性算法和分型組裝方法，兼具高容錯(cuò)、高準(zhǔn)確、高效、計(jì)算資源低耗的優(yōu)點(diǎn)，解決超大基因組組裝和大規(guī)?；蚪M組裝的計(jì)算困境，是未來基因組組裝工具開發(fā)需要繼續(xù)改進(jìn)優(yōu)化的方向。同時(shí)高效低耗的組裝算法也能解決真核生物泛基因組（pan-genome）構(gòu)建的困境[37]。在單體型基因組組裝方面，未來的發(fā)展可以優(yōu)化算法將單體型組裝擴(kuò)大到高雜合性樣品、更高倍性物種的應(yīng)用上，將分型算法整合到組裝中，并通過基于組裝圖的方法提高單體型區(qū)塊contig 放置的準(zhǔn)確率，提高算法性能[36]。

1.2 基因組注釋

基因組注釋是對組裝的基因組進(jìn)行基因結(jié)構(gòu)和功能注釋?；蝾A(yù)測方法依賴于利用兩種類型的基因信息：（1）內(nèi)容——局部位點(diǎn)，如剪接位點(diǎn)、起始密碼子、終止密碼子等，預(yù)測編碼和非編碼區(qū)域；（2）信號——蛋白質(zhì)功能位點(diǎn)，檢測功能性位點(diǎn)是否存在。原核生物的基因結(jié)構(gòu)簡單，各種局部位點(diǎn)特異性強(qiáng)、易識別，基因預(yù)測方法基本成熟，如GeneMarkS[38]、Glimmer[39]等。而真核生物的基因存在內(nèi)含子結(jié)構(gòu)，要經(jīng)過RNA 剪接步驟才能形成成熟的RNA，使得真核基因的預(yù)測更為復(fù)雜和困難，需要采用更加復(fù)雜的算法來預(yù)測真核生物的基因結(jié)構(gòu)[18]。如何準(zhǔn)確地識別可變剪切位點(diǎn)是真核基因預(yù)測的技術(shù)難點(diǎn)，目前有三種策略來預(yù)測真核基因結(jié)構(gòu)：（1）同源預(yù)測：有一些基因蛋白在相近物種間是高度保守的，可以使用已有的高質(zhì)量近緣物種注釋信息通過比對的方式確定外顯子邊界和剪切位點(diǎn)，該方法最可靠，但是對數(shù)據(jù)的依賴性很高，所選取物種間的親緣關(guān)系和進(jìn)化距離對結(jié)果影響較大[40]；（2）從頭預(yù)測：通過已有的概率模型（如隱馬爾科夫模型）來預(yù)測基因結(jié)構(gòu)，再預(yù)測剪切位點(diǎn)和非翻譯區(qū)，這類方法需要預(yù)先訓(xùn)練獲得預(yù)測參數(shù)，準(zhǔn)確性較低，尤其是對于研究較少的物種，這種情況建議利用近緣物種的同源基因數(shù)據(jù)以訓(xùn)練基因結(jié)構(gòu)預(yù)測參數(shù)[18]；（3）基于轉(zhuǎn)錄組預(yù)測：通過物種的RNASeq/EST 等轉(zhuǎn)錄組數(shù)據(jù)輔助注釋，能夠較為準(zhǔn)確地確定剪切位點(diǎn)和外顯子區(qū)域[41]，但是低表達(dá)基因和轉(zhuǎn)錄組組裝錯(cuò)誤也是該方法的主要瓶頸?？梢蚤_發(fā)出綜合上述三種策略的算法，兼具準(zhǔn)確性和特異性地識別可變剪切位點(diǎn)?；诳煽康幕蚪Y(jié)構(gòu)注釋，后續(xù)可以進(jìn)行基因功能注釋，目前普遍采用比對的方法[40,42]，但是序列相似與實(shí)際生物學(xué)功能相似無法完全等價(jià)，需考慮引入其它方法，如引入基因表達(dá)、調(diào)控網(wǎng)絡(luò)等信息，開發(fā)算法進(jìn)行信息整合進(jìn)一步完善基因功能注釋工作。

1.3 基因組變異分析

基因組變異是指與參考序列相比，基因組中發(fā)生的單堿基變異、DNA 序列片段插入、缺失、擴(kuò)增和復(fù)雜結(jié)構(gòu)變異等。目前基于測序方法進(jìn)行單核苷酸 多 態(tài) 性（Single Nucleotide Polymorphism,SNP）、短的插入缺失（Insertion-Deletion,InDel）等變異檢測的策略主要有：（1）基于比對結(jié)果直接檢測變異信息；（2）基于從頭組裝結(jié)果比對檢測變異信息。而在結(jié)構(gòu)變異（Structural Variation,SV）檢測方面目前主要有四種策略：（1）基于插入片段長度的異常分布（Read Pair,RP）；（2）基于讀段覆蓋度的異常分布（Read Depth,RD）；（3）基于未比對上的讀段進(jìn)行分割比對（Split Read,SR）；（4）長讀段或者從頭組裝[43]。不同策略可檢測到的變異類型不同，可將多個(gè)策略整合在一起，從而大大降低假陽性率[43]?；谄唇拥姆椒ㄊ钱?dāng)前獲得全基因組范圍內(nèi)各類型變異最好的方法，但是該方法的結(jié)果準(zhǔn)確性取決于高質(zhì)量組裝，有待高效準(zhǔn)確組裝方法的發(fā)展。GATK是當(dāng)前公認(rèn)的最主流的基因組變異分析軟件之一，可應(yīng)用于人和其它物種[44]。但是當(dāng)前基因組數(shù)據(jù)變異分析在計(jì)算處理能力的瓶頸限制了這些技術(shù)的廣泛使用，針對該問題已有用FPGA 進(jìn)行硬件加速、搭載GPU 加速的Parabricks 軟件等方案來加速GATK的運(yùn)行效率。測序讀長加長可以提高變異分析的準(zhǔn)確度，尤其是大的結(jié)構(gòu)變異分析。未來測序技術(shù)的發(fā)展，兼具測序讀長長和測序準(zhǔn)確性，也將極大的改變變異檢測方法的現(xiàn)狀[45]。

如何解讀這些變異信息，從海量的變異位點(diǎn)中篩選出真正具有生物學(xué)意義的基因位點(diǎn)是變異組學(xué)運(yùn)用的主要瓶頸之一。目前可以通過GWAS 分析獲得基因型和表型之間的關(guān)系[46]。在GWAS 分析中需要考慮適合樣本的關(guān)聯(lián)模型。樣本的表型分兩種：（1）數(shù)量性狀的表型，這種關(guān)聯(lián)分析主要采用線性回歸模型，包括一般線性模型和混合線性模型，對于受到多因素共同影響的復(fù)雜數(shù)量性狀常采用混合線性模型或者改進(jìn)后的模型進(jìn)行分析，仔細(xì)選擇模型或者整合互補(bǔ)模型可以對復(fù)雜性狀的遺傳基礎(chǔ)有更深入的了解[47-48]；（2）case-control 二元的表型，case和control 比例失衡時(shí)，可能會導(dǎo)致較高的假陽性率，針對該問題提出的SAIGE模型顯示出在樣品比例極其失衡時(shí)結(jié)果仍較為正常的特點(diǎn)[49]。當(dāng)前GWAS 分析對于低頻基因的挖掘能力有待提高，可以通過改進(jìn)開發(fā)新的模型和新的群體設(shè)計(jì)揭示表型變異的確定原因[47]。而將GWAS 關(guān)聯(lián)結(jié)合轉(zhuǎn)錄組、甲基化組等組學(xué)研究的數(shù)據(jù)分析，即eQTL、meQTL 等分析，則可以將表型與基因組變異，再到基因表達(dá)、甲基化等之間的關(guān)系聯(lián)系起來。整合多組學(xué)數(shù)據(jù)進(jìn)行聯(lián)合分析揭示復(fù)雜性狀的分子機(jī)制，這也是目前變異解析的主要挑戰(zhàn)[50]。在識別各類QTLs 位點(diǎn)的時(shí)候，需要進(jìn)行多重假設(shè)檢驗(yàn)和校正，常用的多重假設(shè)檢驗(yàn)校正方法有Bonferroni 校正、BH 校正和隨機(jī)打亂后進(jìn)行Storey 校正等[51]。當(dāng)前的多重假設(shè)檢驗(yàn)校正方法難以兼具準(zhǔn)確性和速度，保證準(zhǔn)確度需要耗費(fèi)大量的內(nèi)存和時(shí)間，尤其是對于trans-eQTLs 位點(diǎn)識別時(shí)用隨機(jī)打亂后進(jìn)行Storey 校正方法，通常的超級計(jì)算機(jī)集群是基本不可能完成的，亟需新的模型和算法優(yōu)化來攻克這道難關(guān)[52]。

2 轉(zhuǎn)錄組測序數(shù)據(jù)分析

轉(zhuǎn)錄組是細(xì)胞內(nèi)所有RNA 分子的總和。細(xì)胞中的DNA 記錄的遺傳信息通過轉(zhuǎn)錄表達(dá)，決定生物體的性狀。轉(zhuǎn)錄組中，mRNA是DNA和蛋白質(zhì)之間的中間產(chǎn)物，而非編碼RNA 則有多種多樣的功能。轉(zhuǎn)錄組技術(shù)可以捕捉到細(xì)胞內(nèi)所有RNA 分子瞬時(shí)的存在情況。轉(zhuǎn)錄組數(shù)據(jù)分析技術(shù)對轉(zhuǎn)錄組技術(shù)產(chǎn)生的數(shù)據(jù)進(jìn)行分析和挖掘，得到的結(jié)果和發(fā)現(xiàn)對基礎(chǔ)研究、臨床診斷等都有很大的推動(dòng)。比如Hammond 等[53]對小鼠發(fā)育過程中、老年和大腦受傷后的76 000個(gè)小神經(jīng)膠質(zhì)細(xì)胞進(jìn)行了單細(xì)胞測序，發(fā)現(xiàn)年輕小鼠的小神經(jīng)膠質(zhì)細(xì)胞的異質(zhì)性最高，老年小鼠小神經(jīng)膠質(zhì)細(xì)胞中的炎性細(xì)胞最多，并鑒定出受傷后發(fā)生反應(yīng)的小神經(jīng)膠質(zhì)細(xì)胞亞型，為今后鑒定和操作健康和疾病狀態(tài)的各種亞型的小神經(jīng)膠質(zhì)細(xì)胞奠定了基礎(chǔ)。而在臨床診斷中，轉(zhuǎn)錄組數(shù)據(jù)分析結(jié)果可以得出超出生理范圍的基因表達(dá)情況、等位基因的特異表達(dá)和異常的選擇性剪切等[54]，對基因組測序數(shù)據(jù)是一個(gè)重要的補(bǔ)充。

2.1 RNA-Seq 數(shù)據(jù)分析

RNA-Seq 技術(shù)是目前科研中最常用的轉(zhuǎn)錄組技術(shù)[55]。由于RNA-Seq 讀長短，一般無法測出全長的轉(zhuǎn)錄本，為了檢測哪些轉(zhuǎn)錄本在樣品中表達(dá)，在所研究物種有高質(zhì)量的參考基因組的情況下，通常會先把測序短讀段比對到參考基因組，再組裝成全長轉(zhuǎn)錄本[56]。在研究的物種沒有參考基因組或參考基因組質(zhì)量較差的情況下，需要對轉(zhuǎn)錄本進(jìn)行從頭組裝[41]。如果物種的參考轉(zhuǎn)錄組注釋質(zhì)量高或使用從頭組裝的轉(zhuǎn)錄組，則可以將讀段比對到轉(zhuǎn)錄組上直接定量[57]。對基因表達(dá)的定量一般是通過直接計(jì)算比對到基因上的讀段數(shù)，然后對基因長度和測序深度做歸一化，得到基因表達(dá)的RPKM（Reads Per Kilobase per Million mapped reads），F(xiàn)PKM（Fragments Per Kilobase per Million mapped reads）或TPM（Transcripts Per Million mapped reads）值[56]。也有一些定量方法不需要將讀段比對到基因序列上，比如通過對讀段中k-mer 計(jì)數(shù)的方法獲得基因的表達(dá)量[58]。在對基因的表達(dá)數(shù)據(jù)做樣本間、樣本內(nèi)歸一化、去除批次效應(yīng)和考慮其它可能產(chǎn)生偏差的因素后，可以基于對表達(dá)量統(tǒng)計(jì)分布的推斷（一般假設(shè)基因的表達(dá)量服從負(fù)二項(xiàng)分布或泊松分布），通過統(tǒng)計(jì)檢驗(yàn)確定在分組間差異表達(dá)的基因[59]。對轉(zhuǎn)錄本水平上的差異表達(dá)分析有兩種策略：基于轉(zhuǎn)錄異構(gòu)體（isoform）的策略和基于外顯子的策略。基于轉(zhuǎn)錄異構(gòu)體的策略對轉(zhuǎn)錄異構(gòu)體的表達(dá)進(jìn)行定量，然后檢測差異表達(dá)[56]；基于外顯子的策略通過對比對到外顯子上和跨剪切位點(diǎn)的讀段進(jìn)行計(jì)數(shù)來檢測選擇性剪切的信號[60]。融合基因的鑒定類似于新轉(zhuǎn)錄異構(gòu)體的鑒定，但由于融合的兩個(gè)基因可能位于不同的染色體上，搜索空間更大，并且要考慮到讀段比對的錯(cuò)誤和生物學(xué)上的重要性對找到的融合基因進(jìn)行篩選[61]。對小RNA 測序數(shù)據(jù)的分析一般是把讀段直接比對到已知的小RNA 序列上進(jìn)行定量，或者把讀段比對到參考基因組，提取上下游序列進(jìn)行莖環(huán)結(jié)構(gòu)預(yù)測，預(yù)測miRNA 前體序列和可能的新miRNA[62]。

2.2 三代轉(zhuǎn)錄組測序數(shù)據(jù)分析

三代轉(zhuǎn)錄組測序數(shù)據(jù)的特點(diǎn)是讀長長、錯(cuò)誤率高和通量較低。因其讀段長，大多數(shù)情況下可以直接得到全長的轉(zhuǎn)錄本序列，不需要組裝。但因?yàn)槠溴e(cuò)誤率較高，需要在轉(zhuǎn)錄本鑒定前對讀段進(jìn)行錯(cuò)誤校正[63]。一般錯(cuò)誤校正的方法分為兩類：（1）基于混合讀段的錯(cuò)誤校正，即利用對同一個(gè)樣品的RNASeq 短讀段序列對易出錯(cuò)的長讀段進(jìn)行校正[64]；（2）基于讀段重疊區(qū)的錯(cuò)誤校正，即通過比對不同讀段間重疊的序列來校正這段被多次測到的序列[65]。最新的三代轉(zhuǎn)錄組的建庫和測序策略（如PacBio 公司的ISO-Seq 技術(shù)[66]和Oxford Nanopore公司的R2C2 方法[67]），充分利用目前測序讀段長度遠(yuǎn)長于轉(zhuǎn)錄本長度的特點(diǎn)，在一個(gè)長讀段中對同一轉(zhuǎn)錄本連續(xù)進(jìn)行多次測序，將每次測到的全長轉(zhuǎn)錄本（稱為subread）合在一起比對得到轉(zhuǎn)錄本的一致序列（consensus sequence），有效排除了測序錯(cuò)誤的影響。由于之前三代測序技術(shù)的通量較低，一般使用二代和三代混合測序數(shù)據(jù)進(jìn)行轉(zhuǎn)錄本定量。但目前主流三代測序技術(shù)的通量都有很大提高，如PacBio Sequel II 測序儀的一個(gè)SMRT Cell的通量已達(dá)到8 百萬個(gè)長讀段，達(dá)到了可以對轉(zhuǎn)錄本進(jìn)行定量的標(biāo)準(zhǔn)，此時(shí)只需要對相應(yīng)的長讀段進(jìn)行簡單計(jì)數(shù)就可以得到基因或轉(zhuǎn)錄本的TPM值。三代測序數(shù)據(jù)由于讀長長的優(yōu)勢，在轉(zhuǎn)錄異構(gòu)體鑒定[68]、融合基因轉(zhuǎn)錄本鑒定[69]、轉(zhuǎn)錄本單體型鑒定（transcript haplotype phasing）[70]、轉(zhuǎn)錄本重復(fù)序列鑒定等應(yīng)用上無論是準(zhǔn)確度還是分析方法的簡便性對比RNASeq 技術(shù)都有很大優(yōu)勢，加上三代的RNA 直測技術(shù)[71]的興起可以去除建庫時(shí)反轉(zhuǎn)錄步驟帶來的人為干擾，使轉(zhuǎn)錄本的鑒定和定量更準(zhǔn)確，三代轉(zhuǎn)錄組測序技術(shù)是未來轉(zhuǎn)錄組技術(shù)進(jìn)步的一個(gè)方向。

2.3 單細(xì)胞RNA-Seq 數(shù)據(jù)分析

單細(xì)胞RNA-Seq 技術(shù)由于把瞬時(shí)基因表達(dá)譜的刻畫提升到了單個(gè)細(xì)胞的精度，極大地推進(jìn)了人們對生命系統(tǒng)的認(rèn)識，促成了很多新的科學(xué)發(fā)現(xiàn)，因此是當(dāng)下最熱門和發(fā)展最快的轉(zhuǎn)錄組技術(shù)。單細(xì)胞RNA-Seq 產(chǎn)生的讀段除了轉(zhuǎn)錄本的片段序列外，一般會帶有一段細(xì)胞特異的（實(shí)際上是井或者液滴特異的）條形碼和一段轉(zhuǎn)錄本特異的唯一分子標(biāo)識（Unique Molecular Identifier，UMI），首先需要根據(jù)這兩段標(biāo)記序列給讀段分類，得到一個(gè)計(jì)數(shù)矩陣（count matrix），每一行代表一個(gè)細(xì)胞，每一列是一個(gè)轉(zhuǎn)錄本在各個(gè)細(xì)胞中表達(dá)的讀段數(shù)[13]。在做基因表達(dá)分析之前，需要綜合考慮每個(gè)條形碼下的讀段數(shù)（計(jì)數(shù)深度，count depth）、每個(gè)條形碼下表達(dá)的基因數(shù)和每個(gè)條形碼下的線粒體基因的讀段數(shù)三個(gè)特征來對測序數(shù)據(jù)進(jìn)行質(zhì)量控制，去除死細(xì)胞、細(xì)胞膜破裂的細(xì)胞和一個(gè)井/液滴中含有兩個(gè)或多個(gè)細(xì)胞的情況[72]。由于建庫和測序過程中化學(xué)反應(yīng)的復(fù)雜性，即使對于完全相同的細(xì)胞，都可能獲得不同的計(jì)數(shù)深度，因此在比較細(xì)胞間基因表達(dá)差異之前，需要在細(xì)胞水平上對轉(zhuǎn)錄本表達(dá)計(jì)數(shù)進(jìn)行歸一化。單細(xì)胞RNA-Seq 計(jì)數(shù)矩陣中由于零值較多，需要專門設(shè)計(jì)歸一化的方法。單細(xì)胞RNA-Seq 計(jì)數(shù)矩陣歸一化的方法分為線性和非線性兩種：線性方法對一個(gè)細(xì)胞內(nèi)的所有轉(zhuǎn)錄本使用全局統(tǒng)一的縮放因子[73]；而非線性方法則可以去掉其它的偏差，比如批次因素[74]。在歸一化之后，計(jì)數(shù)矩陣會進(jìn)行l(wèi)og(x+1)變換。

為了減少下游分析的復(fù)雜度和減少數(shù)據(jù)中的噪聲，避免“維數(shù)災(zāi)難”，需要對計(jì)數(shù)矩陣進(jìn)行特征選擇和降維。特征選擇步驟選出對于下游分析信息量最大的基因，通常會選擇表達(dá)量波動(dòng)最大的1 000-5 000個(gè)基因（Highly Variable Genes,HVGs），作為降維算法的輸入[75]。降維算法在盡量保持?jǐn)?shù)據(jù)結(jié)構(gòu)的基礎(chǔ)上，進(jìn)一步將高維空間的數(shù)據(jù)映射到低維空間。最常用的降維方法分為兩類：（1）線性降維方法，主要關(guān)注在映射后保持表達(dá)差異點(diǎn)的距離，如主成分分析法（Principal Component Analysis,PCA）和多維標(biāo)度法（Multidimensional Scaling,MDS）；（2）非線性降維方法，主要關(guān)注在映射后保持表達(dá)模式相似的點(diǎn)足夠接近、保持?jǐn)?shù)據(jù)的局部結(jié)構(gòu)，同時(shí)也兼顧保持全局結(jié)構(gòu)，比如曲線成分分析（Curvilinear Components Analysis,CCA）、t-SNE（t-Distributed Stochastic Neighbor Embedding）、UMAP（Uniform Manifold Approximation and Projection）、擴(kuò)散映射（Diffusion maps）等。在數(shù)據(jù)降維后，可以實(shí)現(xiàn)對單細(xì)胞RNA-Seq 數(shù)據(jù)的可視化。

對經(jīng)上述步驟處理好的表達(dá)數(shù)據(jù)進(jìn)行模式識別是單細(xì)胞RNA-Seq 數(shù)據(jù)分析中的核心步驟。單細(xì)胞RNA-Seq 常見的下游分析有聚類分析、聚類注釋、組成分析、軌跡推斷和差異表達(dá)分析等。聚類分析根據(jù)基因表達(dá)譜定義的細(xì)胞間的距離對細(xì)胞聚類。細(xì)胞間距離的度量方式包括歐幾里得距離、余弦相似度、基于相關(guān)性的距離度量和通過對每個(gè)數(shù)據(jù)集用高斯核學(xué)習(xí)得到的距離度量等，后三種方式由于具有標(biāo)度不變性，考慮值之間的相對差異，對于建庫大小、細(xì)胞大小不同的細(xì)胞之間的比較更健壯。傳統(tǒng)上，聚類分析方法主要有：（1）k-均值法，其中最常用的是RaceID和SIMLR，這兩種方法針對識別罕見細(xì)胞類型做了特殊處理；（2）層次聚類法，常用的有CIDR和PCAReduce，分別針對測序深度低導(dǎo)致的dropout 事件和罕見細(xì)胞類型做了優(yōu)化；（3）社區(qū)發(fā)現(xiàn)法（community detection），該類方法的思路是把密集連接的點(diǎn)，而非距離近的點(diǎn)歸為一類，這類方法被認(rèn)為針對大數(shù)據(jù)集有更快的速度和更好的聚類效果，其中基于k-近鄰圖[76]的Louvain算法[77]在被Seurat和SCANPY 工具包采用后，在單細(xì)胞RNA-Seq 數(shù)據(jù)聚類分析中得到了最廣泛的應(yīng)用。近期針對單細(xì)胞RNA-Seq 數(shù)據(jù)，出現(xiàn)了一些新的聚類分析思路，比如基于概率模型的CellAssign、BAMM-SC、基于深度神經(jīng)網(wǎng)絡(luò)的GOAE、GONN、ACTINN、scScope和scDeepCluster，這些方法在一些數(shù)據(jù)集上顯示出具有更高的聚類準(zhǔn)確度，未來有望得到更廣泛的應(yīng)用。對聚類的注釋需要依靠Mouse Brain Atlas和Human Cell Atlas[78]這樣的參考數(shù)據(jù)庫，或者先找出研究的一類細(xì)胞與所有其它細(xì)胞差異表達(dá)的基因（稱為標(biāo)記基因）與數(shù)據(jù)庫中細(xì)胞的標(biāo)記基因?qū)φ眨蛘哂眠@類細(xì)胞的整個(gè)基因表達(dá)譜與數(shù)據(jù)庫中的細(xì)胞表達(dá)譜對照，推斷這類細(xì)胞的類型。對聚類進(jìn)行注釋之后，可以對細(xì)胞組成進(jìn)行統(tǒng)計(jì)建模，用統(tǒng)計(jì)檢驗(yàn)對不同數(shù)據(jù)集間的細(xì)胞組成差異進(jìn)行分析。測序樣品中的細(xì)胞如果處于連續(xù)變化的過程中，比如樣品中同時(shí)存在處于發(fā)育、疾病進(jìn)展、對處理?xiàng)l件響應(yīng)等過程中不同階段的細(xì)胞時(shí)，旨在尋找離散的分組的聚類分析技術(shù)無法準(zhǔn)確刻畫樣品中細(xì)胞的組成結(jié)構(gòu)，而軌跡推斷或偽時(shí)序分析技術(shù)則可以對細(xì)胞變化的路徑進(jìn)行解析，把每類細(xì)胞在進(jìn)程中的先后順序表示為線性的、二叉的、復(fù)雜圖的、樹的或多叉的拓?fù)浣Y(jié)構(gòu)[79]。在基因?qū)用嫔?，可以對?xì)胞間差異表達(dá)的基因進(jìn)行分析。根據(jù)測試，為普通RNA-Seq 設(shè)計(jì)的差異表達(dá)分析方法加上對基因進(jìn)行加權(quán)的方法的效果要好于目前專門為單細(xì)胞RNA-Seq 設(shè)計(jì)的差異表達(dá)分析方法[80]。

轉(zhuǎn)錄組分析技術(shù)已成為科學(xué)研究、臨床應(yīng)用中不可缺少的關(guān)鍵技術(shù)環(huán)節(jié)。一方面，現(xiàn)有轉(zhuǎn)錄組技術(shù)的不斷完善和新興的單細(xì)胞三代測序技術(shù)、空間轉(zhuǎn)錄組測序技術(shù)等新轉(zhuǎn)錄組技術(shù)的不斷涌現(xiàn)和發(fā)展為轉(zhuǎn)錄組數(shù)據(jù)分析技術(shù)的發(fā)展提出了更高的要求；另一方面，人們對科學(xué)問題理解的逐步加深又要求轉(zhuǎn)錄組分析技術(shù)不斷創(chuàng)新，能更充分更深入地挖掘基因表達(dá)數(shù)據(jù)中蘊(yùn)含的現(xiàn)象和規(guī)律?？梢灶A(yù)見，轉(zhuǎn)錄組分析技術(shù)在未來將發(fā)揮更大的作用。

3 表觀組測序數(shù)據(jù)分析

隨著測序技術(shù)的發(fā)展，表觀組學(xué)相關(guān)研究已經(jīng)成為生命科學(xué)領(lǐng)域炙手可熱的研究方向。科學(xué)家們發(fā)現(xiàn)，表觀組學(xué)對許多重要的細(xì)胞過程的調(diào)控起著關(guān)鍵作用，在精準(zhǔn)醫(yī)學(xué)、生長發(fā)育、分子育種、公共安全等領(lǐng)域研究中占有不可忽視的地位。表觀組研究主要包括組蛋白修飾、DNA 甲基化、染色質(zhì)可及性，以及染色質(zhì)三維結(jié)構(gòu)等。表觀組測序?qū)嶒?yàn)技術(shù)的特殊性決定了用于數(shù)據(jù)分析的各種計(jì)算方法的特征和特異性。隨著表觀組研究技術(shù)發(fā)展的突飛猛進(jìn)，必將開發(fā)出更多新的算法，產(chǎn)生更多解決問題的新工具用于數(shù)據(jù)整合和深度挖掘，也必將使得更多的生命科學(xué)秘密被發(fā)現(xiàn)，更多的改變生命過程的表觀遺傳標(biāo)記被發(fā)現(xiàn)并應(yīng)用于提高國計(jì)民生。

3.1 組蛋白修飾

目前，組蛋白修飾研究主要通過免疫共沉淀技術(shù)與深度測序技術(shù)相結(jié)合的辦法（ChIP-Seq）獲取實(shí)驗(yàn)數(shù)據(jù)。ChIP-Seq 數(shù)據(jù)分析中的核心步驟是富集峰識別（peak calling）。首先是通過整合比對到特定基因組區(qū)域的讀段數(shù)來生成信號圖譜。然后，依靠滑動(dòng)窗口方法將讀段計(jì)數(shù)的離散分布平滑為連續(xù)的信號輪廓分布，估計(jì)超出預(yù)定義的峰的讀段數(shù)量[81]，進(jìn)一步考慮正鏈和負(fù)鏈中讀段計(jì)數(shù)的對應(yīng)關(guān)系[82-83]，以提高峰的分離度。還有其他一些工具使用更復(fù)雜的方法在序列窗口中整合信號，例如使用局部泊松模型來識別基因組位置的局部偏差[84]、依賴核密度估計(jì)[85]、使用貝葉斯分層t-混合模型用于平滑基因組信號圖譜中的讀段計(jì)數(shù)[86]。在與樣本的比較分析中，背景分布的選擇也是識別富集峰中必不可少的步驟，多使用非結(jié)合蛋白質(zhì)的對照抗體實(shí)驗(yàn)的ChIPSeq 數(shù)據(jù)。大多數(shù)峰識別算法通過估計(jì)其相應(yīng)的p值和q值來對更重要的峰進(jìn)行排序和選擇[87]。通常，不同的峰識別工具產(chǎn)生的結(jié)果會有所不同，因此，使用者需根據(jù)研究的靶向蛋白和實(shí)驗(yàn)的情況，比如陰性對照樣本的選擇、是否有生物學(xué)重復(fù)等來選擇峰識別的工具。

3.2 DNA 甲基化

通常情況下，DNA 甲基化是指胞嘧啶甲基化。檢測DNA 甲基化狀態(tài)的實(shí)驗(yàn)主要包括重亞硫酸鹽芯片、MeDIP-Seq、重亞硫酸鹽測序BS-Seq 等等。Illumina 公司的Infinium Methylation Assay 能夠?qū)θ蚪M范圍內(nèi)的單堿基分辨度的甲基化水平進(jìn)行定量測量。在處理Infinium 芯片測定的數(shù)據(jù)時(shí)，要進(jìn)行背景校正，以消除非特異性信號和重復(fù)樣品之間的差異[88]，完成歸一化[88]和進(jìn)行批次效應(yīng)校正[89]等。MeDIP-Seq是利用5-甲基胞嘧啶特異性抗體識別甲基化DNA的免疫沉淀技術(shù)結(jié)合高通量測序的方法。MeDIP-Seq 數(shù)據(jù)可以采用與ChIP-Seq 數(shù)據(jù)相同的生物信息學(xué)方法識別富集區(qū)域，也可以使用專門為甲基化富集數(shù)據(jù)量身定制的方法[89]。在重亞硫酸鹽測序數(shù)據(jù)中，甲基化的胞嘧啶不與重亞硫酸鹽發(fā)生化學(xué)反應(yīng)，而未甲基化的胞嘧啶則在重亞硫酸鹽處理和測序后變成胸腺嘧啶。對于BS-Seq 數(shù)據(jù)，將已經(jīng)發(fā)生堿基變化的讀段比對到基因組，是甲基化數(shù)據(jù)分析的第一個(gè)難點(diǎn)。目前，已開發(fā)出多種比對策略來解決這個(gè)問題。一類是利用已有通用比對軟件[90-91]，采用三堿基（即T，G和A 或者A，T和C）方式進(jìn)行序列比對[92-94]。由于三堿基比對降低了序列的特異性，在比對到大的參考基因組序列時(shí)，一方面增加了比對時(shí)間，另一方面減少了唯一比對的讀段數(shù)量，造成測序讀段利用率變低。一旦重亞硫酸鹽測序讀段與參考基因組唯一比對，就可以估計(jì)特定基因組區(qū)域的甲基化水平，定量Cs和Ts的頻率。另一種比對策略是在比對過程中，不斷調(diào)整標(biāo)記比對評分的矩陣，以適應(yīng)堿基不匹配的情況[94]。這類比對工具往往高估了高甲基化的區(qū)域，但比對效率較高。在為BS-Seq 數(shù)據(jù)選擇比對工具時(shí)，要考慮測序數(shù)據(jù)是來自于定向還是非定向建庫測序。在BSSeq 建庫時(shí)往往會加入非甲基化的Lambda DNA，用于比對后評估重亞硫酸鹽反應(yīng)的轉(zhuǎn)化效率。近年來，單細(xì)胞實(shí)驗(yàn)技術(shù)和測序技術(shù)迅猛發(fā)展，單細(xì)胞甲基化研究多采用重亞硫酸鹽建庫方式，測序數(shù)據(jù)的比對率一直都比較低。2019年一項(xiàng)研究發(fā)現(xiàn)[95]，在單細(xì)胞BS-Seq 文庫制備時(shí)，通過重亞硫酸鹽轉(zhuǎn)化后基因組近端序列與微同源區(qū)（MR）重組產(chǎn)生嵌合分子，而且MR 內(nèi)的DNA 甲基化高度可變，必須刪除這些區(qū)域以準(zhǔn)確估算DNA 甲基化水平。而常規(guī)比對工具采用end-to-end 比對，難以達(dá)到較高的比對率。于是科學(xué)家們開發(fā)出支持單細(xì)胞甲基化測序數(shù)據(jù)局部比對的比對工具，在執(zhí)行質(zhì)控和重亞硫酸鹽測序數(shù)據(jù)局部比對的同時(shí)，進(jìn)行嵌合分子測定和MR 去除，大大提高了單細(xì)胞重亞硫酸鹽測序數(shù)據(jù)的比對率和功能元件的回收率。

對于DNA 甲基化測序數(shù)據(jù)的分析，如何在全基因組水平上準(zhǔn)確識別DNA 甲基化差異區(qū)域（Differentially Methylated Region,DMR）仍然是一個(gè)挑戰(zhàn)。DMR的識別通常采用以下幾種方法，第一種使用二項(xiàng)分布、負(fù)二項(xiàng)分布或具有過度分散參數(shù)的離散分布對甲基化/未甲基化讀段的數(shù)量進(jìn)行建模[96-98]。一種是應(yīng)用平滑算子（smoothing operator），考慮相鄰CpG 位點(diǎn)之間甲基化模式的相關(guān)性[98]。Metilene[99]采用分割算法來檢測單個(gè)或者一組重復(fù)實(shí)驗(yàn)之間的DMR，不需要對數(shù)據(jù)生成機(jī)制進(jìn)行任何模型假設(shè)。雖然識別DMR的方法很多，但是往往存在參數(shù)的依賴性，選擇不同的參數(shù)，會得到不同的DMR 識別結(jié)果。為了解決參數(shù)依賴和缺乏通用性的問題，科學(xué)家們采用完全貝葉斯方法開發(fā)工具ABBA[14]自動(dòng)平滑甲基化圖譜，并可靠地識別DMR。ABBA 基于潛在的高斯模型（Latent Gaussian Model,LGM）和集成嵌套拉普拉斯近似（Integrated Nested Laplace Approximation,INLA），對WGBS 實(shí)驗(yàn)的隨機(jī)采樣過程（甲基化和非甲基化的讀段數(shù)量分布作為非高斯響應(yīng)變量）進(jìn)行建模，采用概率分布指定所有未知量。ABBA 通過對平滑的未觀測到的甲基化圖譜的后驗(yàn)概率的評估來計(jì)算每個(gè)CpG 位點(diǎn)的后驗(yàn)概率，并對兩組數(shù)據(jù)之間的全基因組后驗(yàn)甲基化概率進(jìn)行差異比較，進(jìn)而識別DMRs。ABBA 考慮了WGBS 數(shù)據(jù)的一些內(nèi)在特征，比如通過具有特定組內(nèi)差異的隨機(jī)效應(yīng)進(jìn)行建模以消除組內(nèi)實(shí)驗(yàn)重復(fù)之間的DNA 甲基化差異；通過潛在的高斯場方程反映模型的鄰域結(jié)構(gòu)，并自動(dòng)適應(yīng)數(shù)據(jù)的變化，進(jìn)而明確DNA 甲基化模式的相關(guān)性。ABBA為潛在的高斯場方程的參數(shù)分配先驗(yàn)分布，從而充分考慮了這些量的不確定性。所有這些功能使ABBA的模型能夠根據(jù)實(shí)際情況進(jìn)行調(diào)整，而無需任何用戶定義的參數(shù)。盡管此方法實(shí)現(xiàn)了DMR 識別的無參數(shù)，但是算法繁瑣復(fù)雜，并且對于大樣本間的差異甲基化識別仍有很大的困難和挑戰(zhàn)。

3.3 染色質(zhì)可及性

目前，研究染色質(zhì)可及性的實(shí)驗(yàn)手段主要是通過酶解或者物理化學(xué)手段對開放區(qū)域的DNA 進(jìn)行片段化處理，進(jìn)一步對DNA 片段進(jìn)行基因組定位，明確富集區(qū)域，并對富集區(qū)域進(jìn)行功能分析。目前研究染色質(zhì)可及性的方法主要包括DNase-Seq、FAIRE-Seq、MNase-Seq和ATAC-Seq。MNase-Seq是研究核小體定位的方法，通過對核小體保護(hù)區(qū)域的DNA 序列進(jìn)行定位，進(jìn)而反映出失去核小體保護(hù)的DNA的區(qū)域，其他三種方法則是通過對染色質(zhì)開放區(qū)域的直接定位，反映染色質(zhì)開放性。分析這幾種實(shí)驗(yàn)類型測序數(shù)據(jù)的思路與ChIP-Seq 測序數(shù)據(jù)的分析思路基本上是一致的，即識別富集區(qū)域，并對富集區(qū)域進(jìn)行功能分析。但是染色質(zhì)開放性峰信號通常是寬序列讀取峰，需要區(qū)分真實(shí)峰和假象峰，也要對峰進(jìn)行分類[100-102]?？茖W(xué)家們不斷開發(fā)和優(yōu)化算法，使染色質(zhì)可及性測序數(shù)據(jù)的分析越來越準(zhǔn)確。近年來，單細(xì)胞染色質(zhì)可及性研究越來越多[103]，單細(xì)胞染色質(zhì)可及性對于在單細(xì)胞水平了解細(xì)胞核內(nèi)染色質(zhì)動(dòng)態(tài)變化、揭示染色質(zhì)可及性在細(xì)胞間的異質(zhì)性、鑒定細(xì)胞特異性的染色質(zhì)活性區(qū)域等具有重要意義。

3.4 染色質(zhì)三維結(jié)構(gòu)

三維基因組結(jié)構(gòu)研究利用的染色質(zhì)構(gòu)象技術(shù)在經(jīng)歷了3C、4C和5C 技術(shù)后，又開發(fā)出了Hi-C 技術(shù)。通過Hi-C 數(shù)據(jù)分析揭示染色質(zhì)的遠(yuǎn)程相互作用，推導(dǎo)出基因組的三維空間結(jié)構(gòu)和可能的基因之間的調(diào)控關(guān)系。Hi-C 數(shù)據(jù)的上游分析包括：（1）比對基因組，Hi-C 數(shù)據(jù)需要用特定的策略來比對基因組[104-105]；（2）基因組的一定區(qū)間范圍（bin）的選擇，bin的大小直接決定了最終分析結(jié)果的分辨度。由此，也產(chǎn)生了一些確定最佳bin 大小的方法，使選擇的bin 盡可能的小[106]；（3）歸一化，如迭代校正和特征向量分解歸一化，順序分量歸一化等[104-105,107-108]，并生成數(shù)據(jù)矩陣。Hi-C 數(shù)據(jù)的下游分析是從Hi-C 數(shù)據(jù)矩陣中提取有生物學(xué)意義的結(jié)果，包括：（1）隔室（A/B Compartments）識別，隔室識別多通過主成分分析獲得，由第一主成分值的正負(fù)來表示A Compartment或者B Compartment，也有依靠一種更快且內(nèi)存效率更高的方法來定義隔室得分，以反映出任何給定的bin 進(jìn)入“A”隔室的可能性[105,109]；（2）拓?fù)湎嚓P(guān)結(jié)構(gòu)域（TAD）識別，即識別沿著接觸矩陣的對角線顯示為高度自關(guān)聯(lián)的連續(xù)區(qū)域[106,110-112]；（3）相互作用點(diǎn)識別，也就是識別距離較遠(yuǎn)的染色質(zhì)區(qū)域之間相互聯(lián)系的特定位點(diǎn)。相互作用的計(jì)算識別需要定義背景模型，以便于識別相互作用頻率高于預(yù)期的相互作用聯(lián)系[113-114]。另外，數(shù)據(jù)的可視化在Hi-C 數(shù)據(jù)的分析中顯得尤為重要，可視化可以直觀地幫助研究者描述TAD 模式，也可以揭示分析產(chǎn)生的基因組結(jié)構(gòu)特征性信息。因此，科學(xué)家們不斷地開發(fā)一些可視化工具用于Hi-C 數(shù)據(jù)的分析，比如

Hi-Browse[115]、HicPlotter[116]、3D-GNOME[117]、HiC-3Dviewer[118]和3D genome browser[119]等以熱圖的方式可視化染色質(zhì)相互作用矩陣數(shù)據(jù)，并大多采用web 形式實(shí)現(xiàn)人機(jī)交互。近兩年開發(fā)的，如Delta[120]通過使用4個(gè)可視化視圖模塊，即Virtual 4C 繪圖、線性基因組視圖、基因組環(huán)形視圖和物理視圖，實(shí)現(xiàn)對不同類型的特征數(shù)據(jù)進(jìn)行全面的可視化展示；GITAR[121]除了可實(shí)現(xiàn)Hi-C 數(shù)據(jù)的全面分析和可視化，還提供大量人類和小鼠的數(shù)據(jù)集，便于用戶進(jìn)行數(shù)據(jù)集之間的比較；HiCcompare[122]能夠可視化分析兩個(gè)Hi-C 數(shù)據(jù)集之間的差異；GenomeFlow[123]在實(shí)時(shí)可視化Hi-C 數(shù)據(jù)3D 基因組模型的同時(shí)，還允許用戶在3D 基因組模型上附上基因注釋、基因表達(dá)數(shù)據(jù)和基因組甲基化數(shù)據(jù)；SilkDB 3.0[124]專門用于家蠶數(shù)據(jù)的交互式可視化分析?？傮w而言，Hi-C 技術(shù)的迅速發(fā)展，帶來了染色質(zhì)三維結(jié)構(gòu)數(shù)據(jù)集的快速增加，這也加快了數(shù)據(jù)分析方法的迅速興起，然而，雖然可以使用大量的方法學(xué)解決方案來識別TAD，但尚未完全了解TAD的結(jié)構(gòu)和功能，尤其是TAD的內(nèi)部結(jié)構(gòu)仍然難以捉摸[125]。因此，全面了解染色質(zhì)三維結(jié)構(gòu)的功能和作用，面臨著生物學(xué)和技術(shù)上的雙重挑戰(zhàn)。

4 挑戰(zhàn)和展望

從本世紀(jì)初高通量組學(xué)測序技術(shù)開始得到廣泛應(yīng)用至今，隨著測序花費(fèi)的下降和測序通量的不斷提高，各國政府和企業(yè)資助的大型基因組、轉(zhuǎn)錄組、表觀組計(jì)劃相繼開展[16-17,78,126]，對基因組學(xué)數(shù)據(jù)分析技術(shù)處理海量數(shù)據(jù)的速度提出了極高的要求。為了應(yīng)對基因組學(xué)數(shù)據(jù)的飛速增長，一方面，算法技術(shù)不斷改進(jìn)，如數(shù)據(jù)分析中耗時(shí)步驟是將讀段比對到基因組的算法，從線性比對[127]，到對基因組建立散列索引[128]，到建立后綴數(shù)組/Burrows-Wheeler 變換索引，到建立更加復(fù)雜的FM 索引[129]，到更高效的層次化圖結(jié)構(gòu)的FM 索引（Hierarchical Graph FM index，HGFM）[130]，一系列改進(jìn)使得序列比對速度有了上千倍的提高；另一方面，專門針對基因組數(shù)據(jù)分析設(shè)計(jì)的并行計(jì)算和高性能計(jì)算方案也開始涌現(xiàn)，比如使用多核共享內(nèi)存的超級計(jì)算機(jī)、特殊的硬件設(shè)備（FPGA、GPU、TPU 等）、多節(jié)點(diǎn)的高性能計(jì)算機(jī)、云計(jì)算、容器部署等[45]。已有文章針對上述硬件開發(fā)了應(yīng)用于多種基因組學(xué)數(shù)據(jù)分析的并行算法，比如基于Apache Spark的GATK RNA-Seq流程、基于迭代隨機(jī)森林的基因表達(dá)網(wǎng)絡(luò)構(gòu)建、基于MapReduce的應(yīng)用于年齡預(yù)測的CpG 島篩選、基于Apache Spark的對BWA-MEM的加速等，都達(dá)到了很好的效果。盡管有上述這些進(jìn)展，現(xiàn)有的算法和硬件技術(shù)發(fā)展對處理數(shù)據(jù)效率的提升仍然遠(yuǎn)遠(yuǎn)趕不上組學(xué)數(shù)據(jù)增長的速度，很多分析流程的復(fù)雜性決定了計(jì)算過程必須是并行和串行混合的，并且有很難突破的限速步驟，并行計(jì)算和分布式計(jì)算本身也存在理論上的極限[131]。因此，如何開發(fā)新的方法，提高組學(xué)大數(shù)據(jù)的處理效率，仍是目前亟待解決的問題。另外，實(shí)踐“功能即服務(wù)”的理念——把復(fù)雜的數(shù)據(jù)分析流程分割成微服務(wù)，部署在商業(yè)云或?qū)I(yè)云上，可以增加每部分服務(wù)的可伸縮性，對于世界各地需要進(jìn)行組學(xué)數(shù)據(jù)分析的中小實(shí)驗(yàn)室，可以減去其購買高性能服務(wù)器等基礎(chǔ)設(shè)施建設(shè)和搭建軟件環(huán)境的壓力，提高其數(shù)據(jù)分析的效率，也是未來的發(fā)展趨勢之一。

海量組學(xué)數(shù)據(jù)分析除了數(shù)據(jù)量大帶來的計(jì)算問題，數(shù)據(jù)來源廣產(chǎn)生的異構(gòu)性、不完整、尺度大、質(zhì)量參差不齊，也給組學(xué)大數(shù)據(jù)分析帶來了極大的挑戰(zhàn)。目前已有很多去除組學(xué)大數(shù)據(jù)批次效應(yīng)的方法和模型[132-133]，但是統(tǒng)計(jì)建模過程從特定數(shù)據(jù)類型，轉(zhuǎn)變?yōu)槎嘟M學(xué)整合交叉研究的自動(dòng)化建模是未來發(fā)展的方向。針對異構(gòu)性和數(shù)據(jù)缺失問題，雖然已提出了多核學(xué)習(xí)（multiple kernel learning）算法、DARTS 計(jì)算框架、有參模型、無參模型等大數(shù)據(jù)機(jī)器學(xué)習(xí)、人工智能、統(tǒng)計(jì)建模技術(shù)，為多源異構(gòu)組學(xué)數(shù)據(jù)構(gòu)建預(yù)測模型，但是算法的性能、精度和可靠性仍有很大的提升空間[15,134-135]。測序技術(shù)的發(fā)展也將推動(dòng)組學(xué)數(shù)據(jù)分析方法的發(fā)展，但是如何整合分析多組學(xué)數(shù)據(jù)全面解讀生物系統(tǒng)也成了一項(xiàng)挑戰(zhàn)[13]。多組學(xué)數(shù)據(jù)整合分析將帶來維度災(zāi)難、擴(kuò)展性、歸一化等問題，有待優(yōu)化和開發(fā)一些深度學(xué)習(xí)方法促進(jìn)多維組學(xué)大數(shù)據(jù)整合分析[136-137]。此外，臨床上將多組學(xué)數(shù)據(jù)與環(huán)境暴露數(shù)據(jù)、健康醫(yī)療檔案信息、臨床影像信息聯(lián)合分析將能更好的服務(wù)于精準(zhǔn)醫(yī)學(xué)的發(fā)展，雖然已經(jīng)有了一些研究范式[17,135]，但是制定數(shù)據(jù)標(biāo)準(zhǔn)化體系、建設(shè)生物醫(yī)學(xué)大數(shù)據(jù)共享平臺、開發(fā)多維生物醫(yī)學(xué)數(shù)據(jù)分析的新算法、建立多維知識圖譜等都是亟待解決的問題[138]。我國于2019年成立國家基因組科學(xué)數(shù)據(jù)中心來存儲、質(zhì)控、管理、發(fā)布、共享這些多維、異構(gòu)基因組學(xué)數(shù)據(jù)[139]，期待該中心未來也能開發(fā)出能解決數(shù)據(jù)多維、異構(gòu)問題的算法和工具。

利益沖突聲明

所有作者聲明不存在利益沖突關(guān)系。