下載
創(chuàng)建hisat目錄
mkdir hisat
進入
cd hisat
下載
wget http://ccb.jhu.edu/software/hisat/downloads/hisat-0.1.5-beta-source.zip
建議:去官網(wǎng)下2.0的更好旗唁,能多線程跑
解壓
unzip hisat-0.1.5-beta-source.zip
配置環(huán)境(網(wǎng)上教程)
echo 'PATH=$PATH:/home/wbq/hisat/hisat-0.1.5-beta/' >> ~/.pro
因為我是xshell交互式login shell外构,所以用的是profile块茁。
非交互式的用
echo 'PATH=$PATH:/home/wbq/hisat/hisat-0.1.5-beta/' >> ~/.bashrc
二者區(qū)別參考https://www.cnblogs.com/Chary/p/No0000177.html
比較好的配置環(huán)境操作
在~目錄中
vim .bash_profile
E #編輯改文件课幕。電腦無反應酒朵,沒關(guān)系,正常译暂。再輸入
PATH=$PATH:文件的絕對路徑
#輸入好了之后檢查無錯誤后
按鍵盤左上角Esc鍵 #電腦無反應抠忘,沒關(guān)系,正常秧秉。再輸入
:wq #保存并退出
source .bash_profile #使其生效
配置環(huán)境時一定要用絕對路徑
重啟后生效
檢查有沒有成功
echo $PATH
顯示多出一個hisat-0.1.5-beta 即為成功
運用
構(gòu)建索引文件
hisat-build /home/xxxx.fa /home/genome
genome是指給index的命名名字褐桌,而不是指文件夾。
比對
沒有注釋文件的比對方法
hisat2 -p 18 --dta -x ~/genome/HJX74 -1 /home/wbq/liuqing/rawdata/fq.gz/Rice_D908-T02_good_1.fq.gz -2 /home/wbq/liuqing/rawdata/fq.gz/Rice_D908-T02_good_2.fq.gz -S HJX74.sam
有2個輸入文件象迎,是因為雙端測序的結(jié)果有左端和右端。如果是單端測序呛踊,就只有一個輸入文件砾淌。
不知道每個option的意義的可以輸入
hisat2 -h #看說明書
21.60%的比對沒有一次對齊,61.81%的有一次對齊谭网,16.59%的超過一次對齊汪厨。
超過一次對齊的原因:存在重復序列,有基因不止一個拷貝愉择。
整體比對在最后一行劫乱,87.68%的序列能夠比對成功。
.sam文件全稱是 sequence alignment map format, 專門儲存高通量測序文件锥涕。保存了序列名字衷戈、長度、序列比對到染色體的什么位置(并校驗數(shù)據(jù)可靠性)层坠、program命令殖妇。
各列的意義:
序列,以什么狀態(tài)(沒)比對上(SAM Flags)破花,染色體谦趣,比對具體位置,比對質(zhì)量(the more the better)座每,比對差異情況(插入(i)和缺少(d)算miss match)
有注釋文件的比對方法
hisat2_extract_splice_site.py gene.gtf > genes.ss #把剪切位點提取出來
hisat2 -p 4 --known-splicesite-infile gene.ss -x genome(index文件) -1 SRRxxxx_1.fastq -2 SRRxxx_2.fastq -S SRRxxxx.sam
samtools
下載安裝
mkdir samtools
cd samtools
wget https://github.com/samtools/samtools/releases/download/1.10/samtools-1.10.tar.bz2
tar -jxvf samtools-1.10
./configure --prefix=/home/wbq/samtools/samtools-1.10
make
make install
vim .bash_profile
E
PATH=$PATH:/home/wbq/samtools/samtools-.1.10
source .bash_profile
echo $PATH
使用
就是把sam文件改成bam文件前鹅,壓縮成二進制文件
samtools view -bS HJX74.sam > /home/wbq/HJX74/HJX74.bam #查看bam文件內(nèi)容
samtools sort -@ 2 -o HJX74.sort.bam HJX74.bam #按比對位置排序+格式轉(zhuǎn)換
samtools index HJX74.bam #建立bam文件索引
samtools merge -@ 4 -h SRR1582649.bam merged.bam SRR1582646.bam SRR1582647.bam SRR1582648.bam SRR1582649.bam SRR1582650.bam SRR1582651.bam #把生成的bam文件合并為一個文件。因為每個文件的sam文件表頭都一樣峭梳,所以用-h指定某一個文件的表頭作為總文件的表頭舰绘。
以上命令不是全都需要打,samtools sort的命令必須打
StringTie
下載安裝
mkdir stringtie
cd stringtie
wget http://ccb.jhu.edu/software/stringtie/dl/stringtie-2.1.4.tar.gz
tar -zxvf stringtie-2.1.4.tar.gz
make
vim .bash_profile
i
PATH=$PATH:/home/wbq/stringtie/stringtie-2.1.4
source .bash_profile
echo $PATH
使用
有參考基因組文件
stringtie /home/wbq/HJX74/HJX74.sort.bam -e -G /home/wbq/liuqing/newW0/Oryza_HJX74_top_level.v2.2real.gff3 -o /home/wbq/HJX74/HJX74.gtf
結(jié)果
我們拼出來的有1947條轉(zhuǎn)錄本,但實際原本有5400+基因除盏,差異原因:
- 有些基因表達量比較低叉橱;
- 基因間距很小,有重疊者蠕,軟件會把他們拼接為一個轉(zhuǎn)錄本窃祝,但里面不止一個基因
- 鏈特異性
- ...
可視化結(jié)果
IGV軟件
表達量定量
三種定量指標
stringtie -p 4 -G genes.gtf -e -A SRR1582646.gene_exp \ -o SRR1582646.out SRR1582646.bam
FPKM和TPM都是樣本內(nèi)可比,樣本間不可比踱侣。
差異表達分析
有重復的分析
cuffdiff -p 4 genes.gtf SRR1582646.bam,SRR1582647.bam,SRR1582648.bam SRR1582649.bam,SRR1582650.bam,SRR1582651.bam
#重復用逗號隔開粪小,比較差異的組用空格隔開。
如果它報錯:BAM record error: found spliced alignment without XS attribute
沒關(guān)系抡句,是因為cufflinks版本太老探膊,與hisat2格式有點不合,等它報錯報完還是會繼續(xù)跑的待榔。但是缺點是檢查不到可變剪切逞壁。
結(jié)果:
fold_change大于1的就是顯著上調(diào),0~1的是不顯著變化锐锣,小于0的為顯著下調(diào)腌闯。