HISAT2+StringTie+Ballgown安装及使用流程
HISAT2+StringTie+Ballgown安裝及使用流程
2015年Nature Methods上面發表了一款快速比對工具hisat,作為接替tophat和bowtie的比對工具,它具有更快的比對速度和更高的比對率,最近把這個流程走完一遍,感覺優勢還是很明顯的。?
一、HISAT2:?
1、下載安裝:?
hisat2下載地址:ftp://ftp.ccb.jhu.edu/pub/infphilo/hisat2/downloads/hisat2-2.1.0-Linux_x86_64.zip?
hisat2官方手冊:http://ccb.jhu.edu/software/hisat2/manual.shtml?
下載完成后解壓縮:?
unzip hisat2-2.0.5-Linux_x86_64.zip?
進入hisat2-2.0.5文件夾:?
這里面的綠色文件都是可執行文件,所以只需要把目錄添加到環境變量中即可:?
vim進入編輯bashrc文件,在文本中輸入紅色方框內的內容,保存退出,然后source ~/.bashrc 即可?
此時我們就可以直接調用hisat2命令了。?
2、建立索引:?
如同tophat一樣,比對之前需要利用bowtie建立index,hisat2同樣需要建立索引:?
首先提取gtf文件中的剪切位點和外顯子位置:?
extract_splice_sites.py gencode.vM4.annotation.gtf >gencode.vM4.annotation.for.hisat2.ss?
extract_exons.py gencode.vM4.annotation.gtf >gencode.vM4.annotation.for.hisat2.exon?
建立索引:?
hisat2-build -p 30 --ss gencode.vM4.annotation.for.hisat2.ss --exon gencode.vM4.annotation.for.hisat2.exon GRCm38.p3.genome.fa mouseGencodeIndex?
##如果電腦內存<200G,那么可以不用--ss/--exon參數,但是在比對的時候需要加?
--known-splicesite-infile參數。3、比對:?
我的數據是雙段的無鏈特異性數據,此處需要把sam文件轉化為bam文件,所以需要提前安裝samtools:?
??????? hisat2 --known-splicesite-infile gencode.vM4.annotation.for.hisat2.ss --dta -t -p 24 -x mouseGencodeIndex -1 samp_1.fq.gz -2 samp_2.fq.gz -S accepted_hits.sam &> alignment_summary.txt?
?????? samtools view -bS accepted_hits.sam > accepted_hits.bam?
?????? samtools sort accepted_hits.bam -o accepted_hits_sorted.bam?
?????? rm accepted_hits.bam?
?????? rm accepted_hits.sam?
二、StringTie:?
比對完生成了sam文件,我們利用samtools將它轉化為了排好序的bam文件,下一步就需要量化和確定表達值了,這里用到的StringTie相比之前的cufflinks來說功能強大了好多。?
1、下載安裝:?
stringtie下載地址:http://ccb.jhu.edu/software/stringtie/dl/stringtie-1.3.3b.Linux_x86_64.tar.gz?
stringtie官方手冊:http://ccb.jhu.edu/software/stringtie/index.shtml?t=manual?
直接下載解壓就可以用了,它是可執行文件,也可以按上述方法將路徑添加到環境變量中方便調用。?
2、運行:?
stringtie accepted_hits_sorted.bam -o outRes.gtf -p 28 -G gencode.vM4.annotation.gtf -A gene_abund.tab -B -e?
運行后每個樣本文件夾下結果如下:?
這里我生成了結果gtf文件outRes.gtf和ballgown需要的.ctab文件,還有基因的表達量文件gene_abund.tab,該文件包括基因的表達量FPKM以及TPM等。當然如果你想要轉錄本的表達量,直接打開t_data.ctab這個文件,這里面有轉錄本的FPKM值。?
當然如果我們想利用DESeq2或者edgeR等計算差異表達,那我們就需要得到原始counts值矩陣來作為輸入,此時我們需要利用StringTie自帶的腳本prepDE.py來計算counts值,它可以同時對多個樣本做:?
prepDE.py -i stringtieRes/ -g countsRes/gene_count_matrix.csv -t countsRes/transcript_count_matrix.csv?
stringtieRes/文件夾下面是我所有的樣本的文件夾。?
*這里我們能得到所有樣本的count matrix,但是只能拿到每個樣本對應的FPKM值,又有什么方法能得到合并在一起的FPKM matrix呢?這就需要借助ballgown了。?
三、Ballgown:?
1、安裝:?
首先你需要下載安裝R,我的是3.4.0版本。?
source("https://bioconductor.org/biocLite.R") biocLite("ballgown")?
這里可能提示你安裝XML包的時候會出現錯誤提示:Cannot find xml2-config?
這就需要你在自己電腦上安裝相應的模塊了,我的是centos7,于是安裝相應的模塊:
yum install libxml2-devel?
順利安裝上ballgown包。?
2、使用:?
讀取所有樣本到ballgown對象中:de>?
bg = ballgown(dataDir=de>de>de>YSde>, samplePattern='YT1', meas='all');?
#其中de>de>YS是我的所有樣本的父目錄,每個樣本文件夾名字都包含YT1。?
#計算轉錄本和基因的FPKM值?
de>de>transcript_fpkm = texpr(bg, 'FPKM')?
row.names(de>de>de>transcript_fpkmde>) = transcriptNames(bg)?
write.csv(de>de>transcript_fpkm,"de>de>de>transcript_fpkm_matrix.csvde>")
de>de>gene_expression = gexpr(bg)?
de>de>write.csv(de>de>de>gene_expressionde>,"de>de>de>de>gene_fpkmde>_matrix.csvde>")?
任務完成。?
3、差異表達分析:?
ballgown可以做case/control兩兩比較的差異表達,也可以做多組比較的差異表達(此時不能計算Fold Change值),?
當然也可以做時間序列的差異。?
de>de>de>pData(bg) = data.frame(id=sampleNames(bg), group=rep(c(1,0), each=10))?
#這里是條件矩陣,每行是一個樣本,第二列是條件,如果是case/control那么就是0/1.?
de>de>de>de>stat_results = stattest(bg, feature='transcript', meas='FPKM', getFC=TRUE, covariate='group')?
#注意getFC在多組比較時候不能用,feature參數可以對基因'gene'或者轉錄本'transcript'或者外顯子'exon'做?
差異表達分析。?
de>de>de>de>de>Data(bg) = data.frame(pData(bg), time=rep(1:10, 2)) #dummy time covariate timecourse_results = stattest(bg, feature='transcript', meas='FPKM', covariate='time', timecourse=TRUE)de>?
de>?
但是我個人不太推薦使用ballgown,喜歡使用DESeq2和edgeR來計算。?
de>
轉載于:https://www.cnblogs.com/wangprince2017/p/9937395.html
總結
以上是生活随笔為你收集整理的HISAT2+StringTie+Ballgown安装及使用流程的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 帝豪gs质量怎么样 探究帝豪gs产品的质
- 下一篇: 福特探险者昆仑巅峰版这车怎么样?