RNA-seq-Liver

质控

input_dir=/home/gongyuqi/project/RNA-seq/liver-mouse/clean
output_dir=/home/gongyuqi/project/RNA-seq/liver-mouse/clean/QC
#对每个样本进行质控
nohup fastqc -t 18 $input_dir/*fq.gz -o $output_dir &
#统计合并所有指控结果
dir=/home/gongyuqi/project/RNA-seq/liver-mouse/clean/QC
multiqc $dir/*fastqc.zip --ignore $dir/*.html

Hisat2 进行比对

#因为不想让本机太辛苦，要让它歇息一会儿,所以将KO和WT分批次跑
ls *1.clean.fq.gz >fq1.txt
ls *2.clean.fq.gz >fq2.txt
paste fq1.txt fq2.txt >sample.txt 
#将下列脚本存放在hisat2-alignment.sh文件中
ref=/home/gongyuqi/ref/mm10/index/hisat2/mm10
cat sample.txt|while read id
do
arr=(${id})
fq1=${arr[0]}
fq2=${arr[1]}
sample=${fq1%%_*}.sort.bam
hisat2 -p 8 -x $ref -1 $fq1 -2 $fq2 | samtools sort -@ 8 -o /home/gongyuqi/project/RNA-seq/liver-mouse/align/$sample -
done
#运行比对程序
chmod +x hisat2-alignment.sh
nohup bash hisat2-alignment.sh &

统计比对情况

#建立索引
ls *.bam|xargs -i samtools index {}
#生成状态文件，统计比对情况
ls *.bam|xargs -i samtools flagstat -@ 8 {} | grep "mapped (" > mapping_ratio.txt
ls *.bam > bam.txt
paste bam.txt mapping_ratio.txt|gawk 'BEGIN{FS=" "}{print $1,$6}'|sed 's/(/: /g' > mapping_results.txt

比对后质量控制

qualimap bamqc & qualimap multi-bamqc

qualimap --java-mem-size=16G multi-bamqc -d multi-bamqc.txt -outdir qualimap_multi-bamqc -outformat PDF:HTML

qualimap rnaseq & multiqc

cat qualimap_rnaseq.sh
gtf=/home/gongyuqi/ref/mm10/GTF/gencode.vM24.annotation.gtf
ls *.bam | while read id; do qualimap --java-mem-size=16G rnaseq -pe -bam $id -gtf $gtf -outdir ${id%%.*}_qualimap_rnaseq -outformat PDF:HTML   -oc ${id%%.*}; done
nohup ./qualimap_rnaseq.sh &
#对qualimap rnaseq生成的文件进行multiqc
multiqc *_rnaseq

featureCounts生成表达矩阵

gtf=/home/gongyuqi/ref/mm10/GTF/gencode.vM24.annotation.gtf
nohup featureCounts -T 8 -p -t exon -g gene_id -a $gtf -o ../matrix/RNA-matrix.txt *.bam 1>counts.id.log 2>&1 &
#生成counts的质控文件
multiqc RNA-matrix.txt.summary