一、下载数据
GSE111229 :https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE111229
只下载7个SRR进行学习;
数据下载的话可以参考前面写的一篇笔记:
二、安装conda
参考前面的笔记:
安装软件:
conda install -y sra-tools multiqc fastp hisat2
三、SRA转fastq文件
ls * |while read id; do (nohup fastq-dump --gzip --split-3 -O ./ ${id} &); done
四、质控
ls *gz |xargs fastqc
五、过滤
- 新建几个文件夹
mkdir 01clean 02mapping 03count
- 批量fastp
ls *.fastq.gz | while read id; do (nohup fastp -i ${id} -o ../01clean/$(basename ${id} ".gz").clean);done
六、比对
下载参考基因组
hisat2官网有建好索引的基因组下载:
http://daehwankimlab.github.io/hisat2/download/#m-musculus
批量比对
cd 02mapping/
ls ../01clean/*clean |while read id;do (nohup hisat2 -p 2 -x /mnt/e/bioinfo/mmreference/mm10/genome -U $id -S ${id%%.*}.hisat.sam &);done
%%.* 语法参考:
https://blog.csdn.net/qq_30130417/article/details/80911989
file=/dir1/dir2/dir3/my.file.txt
{file#*/}: 删掉第一个 / 及其左边的字符串:dir1/dir2/dir3/my.file.txt
{file##*/}: 删掉最后一个 / 及其左边的字符串:my.file.txt
{file#*.}: 删掉第一个 . 及其左边的字符串:file.txt
{file##*.}: 删掉最后一个 . 及其左边的字符串:txt
{file%/*}: 删掉最后一个 / 及其右边的字符串:/dir1/dir2/dir3
{file%%/*}: 删掉第一个 / 及其右边的字符串:(空值)
{file%.*}: 删掉最后一个 . 及其右边的字符串:/dir1/dir2/dir3/my.file
{file%%.*}: 删掉第一个 . 及其右边的字符串:/dir1/dir2/dir3/my
由于是在自己电脑上跑数据,没用服务器,结果死机了,还是一个个数据单独比对吧,这里只演示一个:
nohup hisat2 -p 2 -x /mnt/e/bioinfo/mmreference/mm10/genome -U ../01clean/SRR6790711.fastq.clean -S SRR6790711.hisat.sam &
sam格式 转bam格式
ls *.sam | while read id;do (samtools sort -O bam -@ 2 -o $(basename ${id} ".sam").bam ${id});done
为bam文件构建index
ls *.bam |xargs -i samtools index {}
七、count
下载featureCounts
wget https://sourceforge.net/projects/subread/files/subread-2.0.2/subread-2.0.2-Linux-x86_64.tar.gz
下载gtf文件
genecode: https://www.gencodegenes.org/mouse/
全路径运行featureCounts
/mnt/d/bioinfo/biosoft/featureCounts/subread-2.0.2-Linux-x86_64/bin/featureCounts -T 5 -t exon -g gene_id -a /mnt/e/bioinfo/mmreference/gencode.vM30.annotation.gtf -o ../03count/all.id.txt *.bam 1>counts.id.log 2>&1 &
参考:
B站生信技能树: