學(xué)習(xí)資料:B站視頻:生物技能樹 第八季轉(zhuǎn)錄組測序數(shù)據(jù)分析
一淋叶、轉(zhuǎn)錄組比對常用軟件
轉(zhuǎn)錄組比對常用:hisat2、subjunct(subread)伪阶、star煞檩。
基因組測序比對常用bowtie2处嫌、bwa。
要對參考基因組構(gòu)建索引斟湃,不同軟件熏迹,索引也有差異。ChIP-seq 流程里用的bowtie2生成sam凝赛,之后要用samtools轉(zhuǎn)成bam注暗。而subjunc直接就生成bam,但之后還需要用samtools進(jìn)行sort和index墓猎,才能用IGV來查看捆昏,或者別的分析。
二毙沾、subread
subread是一款能快速進(jìn)行序列比對和轉(zhuǎn)錄組定量的軟件屡立。其中subjunct命令用來比對,featureCounts進(jìn)行定量搀军。
1膨俐,下載fa
https://hgdownload.soe.ucsc.edu/goldenPath/hg38/bigZips/hg38.fa.gz
在電腦上下載后上傳到服務(wù)器。
先解壓基因組文件到我的個(gè)人目錄里
gunzip -c ~/reference/hg38.fa.gz > ~/reference/hg38.fa
2罩句,建立索引
然后使用subread軟件的subread-buildindex完成索引的建立
subread-buildindex -o ~/reference/subread_index/hg38 ~/reference/hg38.fa
報(bào)錯(cuò)說需要6G內(nèi)存焚刺,查看服務(wù)器上默認(rèn),用的c-4-1表示四核门烂,每核1G內(nèi)存的CPU計(jì)算節(jié)點(diǎn)乳愉。那么就是有4G的意思啦?明顯不夠啊屯远。
查詢有哪些可用CPU計(jì)算節(jié)點(diǎn)
sinfo
提交作業(yè)時(shí)蔓姚,用-p選定合適的節(jié)點(diǎn) -p c-16-1是不是可以啊。
sbatch -p c-16-1 subread_index.slurm
提交之后運(yùn)行了不短時(shí)間慨丐。
運(yùn)行后的記錄顯示20分鐘坡脐。
生成的index有5個(gè)文件
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環(huán)境
module add SAMtools/1.15.1-GCC-11.2.0 #加載samtools
project=~/rnaseq #項(xiàng)目文件夾
#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
提交的時(shí)候要選擇核數(shù)高的才能運(yùn)行房揭。
sbatch -p c-16-1 rnaseq_align.slurm
#16核备闲,每核1G內(nèi)存.
運(yùn)行時(shí)間大約每個(gè)樣本生成bam要10-20分鐘,sort要3-5分鐘左右捅暴。