前期文章
一個超簡單的轉錄組項目全過程--iMac+RNA-Seq（一）
一個超簡單的轉錄組項目全過程--iMac+RNA-Seq（二）QC
3 Alignment
有很多軟件都可以比對到參考基因組，hisat2，bowtie2，STAR等等（手上這臺iMac不夠內存跑STAR）
首先學一下hisat2的用法
hisat2 --help ## 主要要學會使用--help等參數調用幫助文檔（軟件說明書）
Usage:
hisat2 [options]* -x <ht2-idx> {-1 <m1> -2 <m2> |
-U <r> | --sra-acc <SRA accession number>} [-S <sam>]
-p/--threads <int> number of alignment threads to launch (1)
用法一，一步一個腳印
## 用“指鹿為馬“的方法，將會重復出現的路徑等縮略
wkd=/Users/Desktop/project
## 比對
cd $wkd/clean
ls *gz|cut -d"_" -f 1,2 |sort -u |\
while read id;
do
ls -lh ${id}R1_val_1.fq.gz ${id}R2_val_2.fq.gz
hisat2 -p 4 -x \
/Users/Desktop/project/reference/hg38/hisat2 \
-1 ${id}1_val_1.fq.gz \
-2 ${id}2_val_2.fq.gz -S ${id}.hisat.sam
done
## 轉為二進制文件sam to bam
ls *.bam| while read id;
do
(samtools sort -@ 4 -o $(basename ${id}".sam").bam ${id});
done
rm *.sam # 移除sam文件
用法二：一步到位
vim align.sh ## 創(chuàng)建一個腳本,寫入以下內容
#!/bin/bash
set -e
set -u
set -o pipefail
# source activate rna
# set PATH
HISAT2=/Users/Desktop/project/reference/hg38/hisat2
pwd
ls *gz | cut -d"_" -f 1,2 | sort -u |\
while read id
do
echo "processing ${id}_R1_val_1.fq.gz ${id}_R2_val_2.fq.gz"
hisat2 -p 4 -x \
$HISAT2/hg38.ht2 \
-1 ${id}_R1_val_1.fq.gz \
-2 ${id}_R2_val_2.fq.gz |\
samtools sort -@ 4 > ${id}.hisat.bam
done
# source deactivate
運行腳本
nohup bash align.sh &
查看結果
請關注比對報告中的aligned concordantly exactly 1 time這一項內容，至少要80 %以上。
小彩蛋：洲更老師的錦囊 jobs and disown的使用
當我忘記nohup腳本的時候，但我又要離開實驗室了，這時候使用disown幫助把跑到一半的任務丟到后臺，驚呆了！神仙走位！
小結：本章所需知識背景
（1）什么是比對，比對前后數據的變化；
（2）fasta，fastq文件的規(guī)則以及包含的信息；
（3）比對結果的查看，比對率的由來。
以上知識點均可以在網上查到，建議在進行著一部分項目的時候，一定要補充上最基本的知識點。
問題集錦
收到了幾個小問題，不過我?guī)缀醵紵o法解答，但我可以找到答案~
-
我的怎么QC怎么回事？跑著跑著跑著就停了呀。
QC報錯？
問了代碼才知道問題可能出在input上，過濾結束后的樣本名是val.fq.gz ，出現這樣結果的原因，應該是過濾出錯了。
for id in {xx1,xx2,xxx3};
do
echo "Processin sample ${id}"
hisat2 -p 5 -x /data1/reference/index/hisat2/hg19/hg19 \
-1 /data1/projects/clean/${id}_R1_trimmed.fq.gz \
-2 /data1projects/clean/${id}_R2_trimmed.fq.gz \
-S /data1/projects/align1/${id}.hisat.sam
done
-
為什么比對率這么低？
image.png
解答：
image.png -
比對率的計算好像不對?
怎么加起來都不等于96.7%呀？
這樣的問題，一般我都解決不了，我會直接放棄，然后尋求幫助。
解答：
具體請看簡書：2018-12-04 Hisat2 map 結果與 samtools flagstat 結果不一致
感謝洲更
