柏程思

摘 ? 要：當(dāng)前生物信息學(xué)過濾測(cè)序接頭序列的軟件不能涵蓋所有Illumina測(cè)序平臺(tái)的接頭序列。這樣造成了分析NGS數(shù)據(jù)平臺(tái)的局限性。本文通過R語言編程利用爬蟲技術(shù)對(duì)Illumina質(zhì)控文件的分析，挖掘所有不能被過濾軟件識(shí)別的接頭序列。

關(guān)鍵詞：生物信息學(xué) ?R語言編程 ?爬蟲 ?Illumina測(cè)序

中圖分類號(hào)：Q811.4 ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?文獻(xiàn)標(biāo)識(shí)碼：A ? ? ? ? ? ? ? ? ? ? ? 文章編號(hào)：1674-098X（2019）02（b）-0136-02

當(dāng)前生物信息測(cè)序領(lǐng)域中，Illumina公司屬于二代測(cè)序的壟斷公司。其開發(fā)的Illumina Hiseq、Illumina Miseq、Illumina GAII等平臺(tái)已經(jīng)是流行于全世界。絕大多數(shù)分子生物學(xué)、基因組學(xué)和細(xì)胞生物學(xué)實(shí)驗(yàn)室都在使用Illumina平臺(tái)測(cè)序。

在二代測(cè)序的分析流程中，拿到的數(shù)據(jù)是FASTQ數(shù)據(jù)，需要先對(duì)數(shù)據(jù)進(jìn)行質(zhì)量控制。質(zhì)量控制通常是使用FastQC軟件對(duì)FASTQ數(shù)據(jù)進(jìn)行分析，判斷測(cè)序數(shù)據(jù)是否具有高質(zhì)量。如果質(zhì)量低，則不支持后續(xù)生物信息學(xué)分析，需要過濾。通常情況下，由于測(cè)序儀機(jī)器的誤差，從測(cè)序儀下機(jī)的數(shù)據(jù)都或多或少有低質(zhì)量的序列，這些低質(zhì)量需要有的是堿基質(zhì)量低，有的是測(cè)序接頭未去掉（盡管Illumina大多數(shù)測(cè)序平臺(tái)的測(cè)序儀在2013年以后能保證數(shù)據(jù)下機(jī)自動(dòng)去接頭，但是部分測(cè)序平臺(tái)依然不能自動(dòng)去接頭）。需要過濾，而過濾使用的軟件一般為Trimmomatic軟件。但是Trimmomatic軟件過濾使用的文件是自身adapter文件夾中自帶的Truseq文件過濾測(cè)序的接頭，而這些接頭文件只包含了Hiseq、Miseq和GAII測(cè)序平臺(tái)的接頭文件，沒有包括全部的接頭文件。未去接頭的序列在質(zhì)控結(jié)果中可以將接頭序列以O(shè)verrepresent形式表示出來。所以，如果我們測(cè)序時(shí)選擇的測(cè)序平臺(tái)不能自動(dòng)下機(jī)去接頭，我們需要手動(dòng)自己去接頭。爬蟲是利用計(jì)算機(jī)技術(shù)對(duì)網(wǎng)絡(luò)數(shù)據(jù)的挖掘，因?yàn)榛ヂ?lián)網(wǎng)數(shù)據(jù)基本都是儲(chǔ)存在網(wǎng)絡(luò)服務(wù)器中，網(wǎng)絡(luò)服務(wù)器末端端口是用戶。用戶可以通過網(wǎng)頁訪問網(wǎng)絡(luò)服務(wù)器，網(wǎng)頁是由HTML語言搭建的可視化端口。HTML是HyperText Markup Language（超文本標(biāo)記語言）的縮寫，這個(gè)語言使用<標(biāo)簽>內(nèi)容基本格式進(jìn)行網(wǎng)頁編輯[1]。例如

This is how to use HTML

The way to use HTML

You can learn it

將上面這個(gè)代碼復(fù)制到一個(gè)新建TXT文本中，并將后綴命名為.html，雙擊該文件打開會(huì)出現(xiàn)以下內(nèi)容（見圖1）。

1 ?分析方法

FASTQC軟件輸出的質(zhì)控結(jié)果就有HTML本地文件，如果出現(xiàn)了接頭序列就會(huì)在Overrepresent中出現(xiàn)，Overrepresent有其對(duì)應(yīng)的HTML標(biāo)簽。R語言[2-3]可以通過追溯內(nèi)容所在的標(biāo)簽追溯到內(nèi)容，這個(gè)追溯內(nèi)容的路徑稱之為xpath，最后通過正則表達(dá)式篩選我們要的內(nèi)容即可。首先在Linux系統(tǒng)上存放測(cè)序數(shù)據(jù)的路徑（該路徑中只能含有測(cè)序數(shù)據(jù)文件）下使用命令fastqc `ls $pwd`，然后下載輸出的HTML文件。先用網(wǎng)頁查看是否有接頭序列，再用Notepad++打開文件找接頭序列所在的xpath（見圖2）。

編寫如下R語言代碼

library（rvest）

library（stringi）

setwd（"D：/test/fastQC"）

myQCfile<-dir（"D：/test/fastQC"）

truseqindex<-NULL

for （i in 1：length（myQCfile）） {

qc<-read_html（myQCfile[i]）

a<-qc%>%

html_nodes（xpath = "http://tr/td"）%>%

html_text（）

b<-NULL

for （j in 1：length（a）） {

if（grepl（a[j]，pattern = "（TruSeq|Primer）"））{b<-c（b，a[j-3]，a[j]）} }

truseqindex<-c（truseqindex，b）}

2 ?結(jié)語

我們開發(fā)的挖掘當(dāng)前過濾軟件無法過濾的接頭腳本更具有實(shí)用性，可以適用于所有NGS數(shù)據(jù)分析過濾腳本。使分析結(jié)果更具有可靠性。

參考文獻(xiàn)

[1] 鄧子云.爬蟲系統(tǒng)中標(biāo)簽刪除功能的設(shè)計(jì)及優(yōu)化[J].計(jì)算機(jī)系統(tǒng)應(yīng)用，2019，28（1）：176-181.

[2] 許慶煒.B語言—生物信息學(xué)可視化流程語言[J].計(jì)算機(jī)與數(shù)字工程，2009，37（5）：90-93.

[3] 吳棟楊. 構(gòu)建基于R語言的生物信息學(xué)研究平臺(tái)[A].第十次中國(guó)生物物理學(xué)術(shù)大會(huì)論文摘要集[C].中國(guó)生物物理學(xué)會(huì)，2006：1.