NGS基礎(chǔ)：測(cè)序原始數(shù)據(jù)批量下載
生物或醫(yī)學(xué)中涉及高通量測(cè)序的論文，一般會(huì)將原始測(cè)序數(shù)據(jù)上傳到公開的數(shù)據(jù)庫(kù)，上傳方式見測(cè)序文章數(shù)據(jù)上傳找哪里；并在文章末尾標(biāo)明數(shù)據(jù)存儲(chǔ)位置和登錄號(hào),如 The data from this study was deposited in NCBI Sequence Read Archive under accession SRA: SRP114962.。
NCBI的SRA (Sequence Read Archive) 數(shù)據(jù)庫(kù)(http://www.ncbi.nlm.nih.gov/sra/) 是最常用的存儲(chǔ)測(cè)序數(shù)據(jù)的數(shù)據(jù)庫(kù)。目前SRA數(shù)據(jù)的組織方式分為下面4個(gè)層次：
Studies—研究課題；
Experiments—實(shí)驗(yàn)設(shè)計(jì)；
Runs—測(cè)序結(jié)果集；
Samples—樣品信息。
進(jìn)入SRA官網(wǎng)：https://www.ncbi.nlm.nih.gov/sra, Search框中輸入SRA編號(hào)（SRP114962），獲得如下圖的界面：
點(diǎn)擊第一個(gè)樣品即可查看其詳細(xì)信息。
當(dāng)樣品比較多時(shí)，可以點(diǎn)擊Send results to Run selector（圖中畫圈的位置）進(jìn)入篩選頁(yè)面。
從圖中可發(fā)現(xiàn)，測(cè)序平臺(tái)是Illumina HiSeq 4000，5748個(gè)Runs，每個(gè)Run的名字、樣本名、測(cè)序類型（全基因組/外顯子組等）、tissue、treatment等。
在如此多的Runs中，假設(shè)我們想獲取其中兩個(gè)病人的化療前和化療后的外顯子組測(cè)序數(shù)據(jù)，觀察其化療前后究竟有哪些基因突變以及突變的頻率怎么樣。數(shù)據(jù)來自于文章 腫瘤化療無效是對(duì)預(yù)先存在的突變的選擇還是誘發(fā)新突變，Cell給你答案。
5748個(gè)Runs，有116Page，怎么找呢？
在Facets下拉框中先勾選Assay Type，等待頁(yè)面相應(yīng)后勾選wxs，即全外顯子組數(shù)據(jù)，等待頁(yè)面相應(yīng)。
在Facets下拉框中勾選Sample name，等待頁(yè)面相應(yīng)后勾選ktn102及ktn134兩個(gè)病人的分別四個(gè)樣本（四種treatment：pre、2cycleschemo、operative和blood），如圖。等待頁(yè)面相應(yīng)。獲得Run編號(hào)（藍(lán)色框）：SRR5908363、SRR5908362…
然后使用NCBI提供的工具SRAToolkit下載。
SRA toolkit https://trace.ncbi.nlm.nih.gov/Traces/sra/sra.cgi?view=software, 根據(jù)服務(wù)器操作系統(tǒng)類型下載對(duì)應(yīng)的二進(jìn)制編碼包，下載解壓放到環(huán)境變量即可使用。
使用NCBI提供的SRA-toolkit中的工具fastq-dump直接下載SRR文件，并轉(zhuǎn)換為FASTQ格式，--split-3參數(shù)表示如果是雙端測(cè)序就自動(dòng)拆分，如果是單端不受影響。--gzip轉(zhuǎn)換fastq為壓縮文件，節(jié)省空間。
下載的數(shù)據(jù)集一般比較大，放入后臺(tái)不中斷下載 (nohup cmd &)。
nohup fastq-dump -v --split-3 --gzip SRR5908360 &
nohup fastq-dump -v --split-3 --gzip SRR5908361 &
nohup fastq-dump -v --split-3 --gzip SRR5908362 &
nohup fastq-dump -v --split-3 --gzip SRR5908363 &
nohup fastq-dump -v --split-3 --gzip SRR5906250 &
nohup fastq-dump -v --split-3 --gzip SRR5906251 &
nohup fastq-dump -v --split-3 --gzip SRR5906252 &
nohup fastq-dump -v --split-3 --gzip SRR5906253 &注意：如果數(shù)據(jù)量很大可能需要下載1-2天。數(shù)據(jù)下載完會(huì)在~/ncbi下面存在緩存的sra文件，記得定時(shí)清空。
按照上述步驟下載完畢后可看到很多個(gè)fastq.gz格式測(cè)序文件。
數(shù)據(jù)比較多時(shí)，一個(gè)個(gè)手動(dòng)寫也比較麻煩？怎么處理呢？
下載上面的metadata后自己生成一個(gè)批量下載并重命名的腳本。下下來的metadat文件通常名字是SraRunTable.txt，列很多，且是CSV文件，分隔符是逗號(hào)，某一個(gè)字段里面還有逗號(hào)的存在。Linux命令直接不好處理。
先寫一個(gè)單行R腳本提取Run和Sample Name列，并另存到文件SraRunTable.tsv。
Rscript -e 'write.table(read.table("SraRunTable.txt",sep=",", header=T, row.names=NULL)[,c("Run","Sample.Name")],"SraRunTable.tsv",sep="\t", quote=F,col.names=T, row.names=F)'然后用awk就可以讀取第一列批量下載并結(jié)合第二列批量重命名就可以了。
awk 'FNR>1{system("fastq-dump -v --split-3 --gzip "$1"; rename "$1" "$2" "$1"*");}' SraRunTable.tsv
等同于手動(dòng)輸入了如下命令
fastq-dump -v --split-3 --gzip SRR12603383; rename SRR12603383 PFER12d1 SRR12603383*
fastq-dump -v --split-3 --gzip SRR12603384; rename SRR12603384 PFER12d3 SRR12603384*
fastq-dump -v --split-3 --gzip SRR12603385; rename SRR12603385 PFER12d2 SRR12603385*
fastq-dump -v --split-3 --gzip SRR12603386; rename SRR12603386 PFER12d1 SRR12603386*
fastq-dump -v --split-3 --gzip SRR12603387; rename SRR12603387 PFER9d3 SRR12603387*
fastq-dump -v --split-3 --gzip SRR12603388; rename SRR12603388 PFER9d2 SRR12603388*
fastq-dump -v --split-3 --gzip SRR12603389; rename SRR12603389 PFER9d1 SRR12603389*
fastq-dump -v --split-3 --gzip SRR12603390; rename SRR12603390 PFER9d3 SRR12603390*這個(gè)為生信學(xué)習(xí)和生信作圖打造的開源R教程真香?。?！
這個(gè)為生信學(xué)習(xí)打造的開源Linux教程真香?。?！
往期精品(點(diǎn)擊圖片直達(dá)文字對(duì)應(yīng)教程)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
后臺(tái)回復(fù)“生信寶典福利第一波”或點(diǎn)擊閱讀原文獲取教程合集
（請(qǐng)備注姓名-學(xué)校/企業(yè)-職務(wù)等）