保姆级有参考RNA-seq分析基本流程1（表达矩阵）

我用的基本软件trim_galore数据过滤，fastqc质量控制，hisat2做mapping，Samtools做格式转换，htseq得到表达矩阵，python合并矩阵后，调用R做富集分析，

最后两步python和R做富集分析下一篇详细讲解，一步步到表达矩阵具体代码如下

首先启动miniconda3环境

$ source ~/anaconda3/bin/activate

解压缩.gz文件(日常处理可以不用解压缩的)
gzip -d *.gz


质控
fastqc *.gz


trim_galore数据过滤
#trim_galore -q 20 --phred33 --stringency 3 --length 20 -e 0.1 \
            --paired *_1.fq.gz *_2.fq.gz  \
            --gzip -o /Users/yqWang/Local/wangyiqi/ZL1_VS_ZL3/ZL1_2/ZL1_2clean_data


进入python2
source activate python2
然后启动hisat2
hisat2 -h


hisat2 -t -x /Users/yqWang/Local/wangyiqi/reference/index/hg38/genome -1 ZL1_2_1_val_1.fq.gz -2 ZL1_2_2_val_2.fq.gz -S ZL1_2.sam 

for i in `seq 1 1`
do
    samtools view -S ZL1_2.sam -b > ZL1_2.bam
    samtools sort ZL1_1.bam -o ZL1_1_sorted.bam
    samtools index ZL1_1_sorted.bam
done


对bam进行质控
bam_stat.py  -i ZL1_1_sorted.bam

下载hg38的基因覆盖率
wget https://sourceforge.net/projects/rseqc/files/BED/Human_Homo_sapiens/hg38_RefSeq.bed.gz/download

查看基因组覆盖率
read_distribution.py -i /Users/yqWang/Local/wangyiqi/ZL1_VS_ZL3/ZL1_1/ZL1_1clean_data/ZL1_1_sorted.bam -r /Users/yqWang/Local/wangyiqi/reference/index/hg38_RefSeq.bed

Reads计数，得到表达矩阵
htseq-count -r name -f bam /Users/yqWang/Local/wangyiqi/ZL1_VS_ZL3/ZL1_1/ZL1_1clean_data/ZL1_1_sorted.bam /Users/yqWang/Local/wangyiqi/reference/index/gencode.v32.gtf > /Users/yqWang/Local/wangyiqi/ZL1_VS_ZL3/ZL1_1/ZL1_1clean_data/ZL1_1.count

htseq-count -r pos -f bam /Users/yqWang/Local/wangyiqi/ZL1_VS_ZL3/ZL1_1/clean_data/ZL1_1_sorted.bam /Users/yqWang/Local/wangyiqi/reference/index/gencode.v32.gtf 

继续python3
source deactivate