学习资料:B站视频:生物技能树 第八季转录组测序数据分析
一、转录组比对常用软件
转录组比对常用:hisat2、subjunct(subread)、star。
基因组测序比对常用bowtie2、bwa。
要对参考基因组构建索引,不同软件,索引也有差异。ChIP-seq 流程里用的bowtie2生成sam,之后要用samtools转成bam。而subjunc直接就生成bam,但之后还需要用samtools进行sort和index,才能用IGV来查看,或者别的分析。
二、subread
subread是一款能快速进行序列比对和转录组定量的软件。其中subjunct命令用来比对,featureCounts进行定量。
1,下载fa
https://hgdownload.soe.ucsc.edu/goldenPath/hg38/bigZips/hg38.fa.gz
在电脑上下载后上传到服务器。
先解压基因组文件到我的个人目录里
gunzip -c ~/reference/hg38.fa.gz > ~/reference/hg38.fa
2,建立索引
然后使用subread软件的subread-buildindex完成索引的建立
subread-buildindex -o ~/reference/subread_index/hg38 ~/reference/hg38.fa
报错说需要6G内存,查看服务器上默认,用的c-4-1表示四核,每核1G内存的CPU计算节点。那么就是有4G的意思啦?明显不够啊。
查询有哪些可用CPU计算节点
sinfo
提交作业时,用-p选定合适的节点 -p c-16-1是不是可以啊。
sbatch -p c-16-1 subread_index.slurm
提交之后运行了不短时间。
运行后的记录显示20分钟。
生成的index有5个文件
3,比对
建立rnaseq_align.slurm
如下
#!/bin/bash
#SBATCH --output=rnaseq_align.out
#SBATCH --error=rnaseq_align.err
#SBATCH --mail-type=end
#SBATCH --mail-user=zmeraner@126.com
module add Anaconda3/2020.02
source activate
conda activate py2.7 #激活2.7环境
module add SAMtools/1.15.1-GCC-11.2.0 #加载samtools
project=~/rnaseq #项目文件夹
#step4:align & sort &index
mkdir $project/align
mkdir $project/bam
output_dir=$project/align
bam_dir=$project/bam
ls $project/clean/*.gz | while read id
do
file=$(basename $id) #文件名
sample=${file%%_trimmed*} #样本名
echo " start alignment for $sample" `date`
subjunc -i /home/cloudam/reference/subread_index/hg38 -r $id -o $output_dir/$sample.bam
echo " end alignment for $sample" `date`
echo " sort &index for $sample" `date`
samtools sort $output_dir/$sample.bam -o $bam_dir/$sample.sorted.bam
samtools index $bam_dir/$sample.sorted.bam
rm $output_dir/$sample.bam `date`
samtools flagstat $bam_dir/$sample.sorted.bam >$bam_dir/$sample.flagstat.txt
echo "end sort index for $sample"
done
提交的时候要选择核数高的才能运行。
sbatch -p c-16-1 rnaseq_align.slurm
#16核,每核1G内存.
运行时间大约每个样本生成bam要10-20分钟,sort要3-5分钟左右。