參考http://www.reibang.com/p/681e02e7f9af
http://www.reibang.com/p/93f96e7538da
任務：
- 比對軟件很多表悬，首先大家去收集一下溪厘，因為我們是帶大家入門紊册，請統(tǒng)一用hisat2，并且搞懂它的用法。
- 直接去hisat2的主頁下載index文件即可杰赛，然后把fastq格式的reads比對上去得到sam文件又活。
- 接著用samtools把它轉為bam文件，并且排序(注意N和P兩種排序區(qū)別)索引好搀突，載入IGV刀闷，再截圖幾個基因看看！
- 順便對bam文件進行簡單QC仰迁，參考直播我的基因組系列甸昏。
1. 比對軟件
- HISAT2：http://ccb.jhu.edu/software/hisat2/index.shtml
參考資料：http://blog.biochen.com/archives/337 - STAR：https://codeload.github.com/alexdobin/STAR/zip/master
參考資料：http://www.bio-info-trainee.com/727.html - TopHat：http://ccb.jhu.edu/software/tophat/index.shtml
參考資料：http://blog.sina.com.cn/s/blog_8808cae20101amqp.html - RapMap：https://github.com/COMBINE-lab/RapMap
參考文獻:https://academic.oup.com/bioinformatics/article/32/12/i192/2288985/RapMap-a-rapid-sensitive-and-accurate-tool-for - CIDANE：http://ccb.jhu.edu/software/cidane/
參考文獻：https://genomebiology.biomedcentral.com/articles/10.1186/s13059-015-0865-0 - CLASS2 ：https://sourceforge.net/projects/splicebox/files/?source=navbar
參考文獻：https://academic.oup.com/nar/article/44/10/e98/2516329/CLASS2-accurate-and-efficient-splice-variant
2. HISAT2的使用
為什么要用index？官網(wǎng)有描述：為了用整個index代表整個基因組徐许，HISAT2 用小的index覆蓋了整個基因組施蜜，每個index覆蓋了56 Kbp的范圍，覆蓋整個人類基因組需要55,000 indexes雌隅，這些index結合其他策略可以快速準確的比對序列翻默。
#index 下載
nohup wget ftp://ftp.ccb.jhu.edu/pub/infphilo/hisat2/data/hg19.tar.gz &
tar -zxvf *.tar.gz #解壓
# 刪除壓縮包
rm -rf *.tar.gz
- hisat2 -h查看幫助文件：
Usage:
hisat [options]* -x <bt2-idx> {-1 <m1> -2 <m2> | -U <r> | --sra-acc <SRA accession number>} [-S <sam>]
<bt2-idx> Index filename prefix (minus trailing .X.ht2).
<m1> Files with #1 mates, paired with files in <m2>.
Could be gzip'ed (extension: .gz) or bzip2'ed (extension: .bz2).
<m2> Files with #2 mates, paired with files in <m1>.
Could be gzip'ed (extension: .gz) or bzip2'ed (extension: .bz2).
<r> Files with unpaired reads.
Could be gzip'ed (extension: .gz) or bzip2'ed (extension: .bz2).
<SRA accession number> Comma-separated list of SRA accession numbers, e.g. --sra-acc SRR353653,SRR353654.
<sam> File for SAM output (default: stdout)
<m1>, <m2>, <r> can be comma-separated lists (no whitespace) and can be specified many times. E.g. '-U file1.fq,file2.fq -U file3.fq'.
- 參數(shù)說明：
-x 指定index文件名
-1 <m1> 雙端測序第一個文件
-2 <m2> 雙端測序第二個文件
-U 單端測序文件
--sra-acc 登錄號
-S 輸出文件為sam格式
-q 輸入文件為FASTQ .fq/.fastq格式
-比對到參考基因組,RNA-Seq數(shù)據(jù)是從 SRR3589957 ~ SRR3589962，6個樣本的PE數(shù)據(jù)恰起，也就是有6次循環(huán)修械，可以寫腳本，也可以直接在命令行里運行
for i in `seq 56 62`
do
hisat2 -t -x /opt/NfsDir/UserDir/qin/qin/Data/annotation/hg19/genome -1 /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/SRR35899${i}.sra_1.fastq.gz -2 /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/SRR35899${i}.sra_2.fastq.gz -S /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR35899${i}.sam &
done
運行時間如下：
Paste_Image.png
- 很耗CPU凹炫巍肯污！用的服務器！
Paste_Image.png
- 結果：
Paste_Image.png
SAM（sequence Alignment/mapping)數(shù)據(jù)格式是目前高通量測序中存放比對數(shù)據(jù)的標準格式吨枉，當然他可以用于存放未比對的數(shù)據(jù)仇箱。目前處理SAM格式的工具主要是SAMTools。samtools功能眾多东羹，在本次作業(yè)中剂桥，我們主要學會將sam文件轉換為bam文件，并對bam文件進行sorted（其中有兩種排序方式N和P）属提，最后建立索引权逗。
Program: samtools (Tools for alignments in the SAM format)
Version: 1.3.1-58-gcbee45e (using htslib 1.3.2-228-g0c32631)
Usage: samtools <command> [options]
Commands:
-- Indexing
dict create a sequence dictionary file
faidx index/extract FASTA
index index alignment
-- Editing
calmd recalculate MD/NM tags and '=' bases
fixmate fix mate information
reheader replace BAM header
rmdup remove PCR duplicates #移除PCR產(chǎn)生的重復序列
targetcut cut fosmid regions (for fosmid pool only)
addreplacerg adds or replaces RG tags
-- File operations
collate shuffle and group alignments by name
cat concatenate BAMs
merge merge sorted alignments
mpileup multi-way pileup
sort sort alignment file
split splits a file by read group
quickcheck quickly check if SAM/BAM/CRAM file appears intact
fastq converts a BAM to a FASTQ #格式轉換
fasta converts a BAM to a FASTA
-- Statistics
bedcov read depth per BED region #bed文件的測序深度
depth compute the depth
flagstat simple stats
idxstats BAM index stats
phase phase heterozygotes
stats generate stats (former bamcheck)
-- Viewing
flags explain BAM flags
tview text alignment viewer
view SAM<->BAM<->CRAM conversion
depad convert padded BAM to unpadded BAM
主要功能：
view: BAM-SAM/SAM-BAM 轉換和提取部分比對
sort: 比對排序
merge: 聚合多個排序比對
index: 索引排序比對
faidx: 建立FASTA索引，提取部分序列
tview: 文本格式查看序列
pileup: 產(chǎn)生基于位置的結果和 consensus/indel calling
# 首先將比對后的sam文件轉換成bam文件
# 利用的是samtools的view選項冤议，參數(shù)-S 輸入sam文件斟薇；參數(shù)-b 指定輸出的文件為bam；最后重定向寫入bam文件
$ for ((i=56;i<=62;i++));do samtools view -S /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR35899${i}.sam -b > /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR35899${i}.bam;done
# 將所有的bam文件進行排序
$ nohup for (( i=58;i<=62;i++ )); do samtools sort /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR35899${i}.bam -o /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR35899${i}.sorted.bam;done
# 將所有的排序文件建立索引恕酸，索引文件.bai后綴
$ for ((i=56;i<=62;i++));do samtools index /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR35899${i}.sorted.bam;done
合在一塊跑：
for i in `seq 56 58`
do
samtools view -S SRR35899${i}.sam -b > SRR35899${i}.bam
samtools sort SRR35899${i}.bam -o SRR35899${i}_sorted.bam
samtools index SRR35899${i}_sorted.bam
done
Paste_Image.png
比對質控(QC)
常用工具有
- Picard https://broadinstitute.github.io/picard/
- RSeQC http://rseqc.sourceforge.net/
- Qualimap http://qualimap.bioinfo.cipf.es/
java -jar /opt/NfsDir/BioDir/picard-tools-1.119/picard.jar CollectMultipleMetrics \
I=/opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR3589956.sorted.bam \
O=/opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/multiple_metrics \
R=/opt/NfsDir/PublicDir/reference/ucsc.hg19.fasta
統(tǒng)計bam文件：
/opt/NfsDir/BioDir/RSeQC/RSeQC-2.6.4/scripts/bam_stat.py -i /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR3589956.sorted.bam
提示出錯：
Paste_Image.png
檢查發(fā)現(xiàn)是路徑和名稱寫錯了,修改后就可以了
對bam文件進行統(tǒng)計分析
Paste_Image.png
#下載hg19_RefSeq.bed文件
https://sourceforge.net/projects/rseqc/files/BED/Human_Homo_sapiens/hg19_RefSeq.bed.gz/download
#查看基因組覆蓋率
/opt/NfsDir/BioDir/RSeQC/RSeQC-2.6.4/scripts/read_distribution.py -i /opt/NfsDir/UserDir/qin/qin/Data/RNAseq/align/SRR3589956.sorted.bam -r /opt/NfsDir/UserDir/qin/qin/Data/annotation/hg19/hg19_RefSeq.bed
Paste_Image.png
