rRNA去除之后就开始进行数据比对了，这一步骤作者使用了三个比对软件：Tophat2，STAR，HISAT2，BWA。
相应代码抠出来：
HISAT2：
image-20220326230833394.png
继续扣出来运行~
一、下载人的参考基因组序列
去ensemble数据库：
# 使用axel多线程下载数据，速度可到10M/s
mkdir GRCh38
cd GRCh38
# 下载fa文件
axel -n 100 http://ftp.ensembl.org/pub/release-105/fasta/homo_sapiens/dna/Homo_sapiens.GRCh38.dna.primary_assembly.fa.gz
# 下载gtf文件
axel -n 100 http://ftp.ensembl.org/pub/release-105/gtf/homo_sapiens/Homo_sapiens.GRCh38.105.chr.gtf.gz
# 解压
gunzip Homo_sapiens.GRCh38.dna.primary_assembly.fa.gz Homo_sapiens.GRCh38.105.chr.gtf.gz
二、构建索引
# 激活小环境
conda activate rna
# Human数据
# HISAT2构建索引，Hisat2Index.sh内容为下
mkdir Hisat2Index
fasta=Homo_sapiens.GRCh38.dna.primary_assembly.fa
gtf=Homo_sapiens.GRCh38.105.chr.gtf
fasta_baseName=GRCh38
# 提取外显子位点
hisat2_extract_exons.py $gtf > Hisat2Index/${fasta_baseName}.exon
# 提取可变剪切信息
hisat2_extract_splice_sites.py $gtf > Hisat2Index/${fasta_baseName}.ss
# 构建索引
hisat2-build -p 12 -f $fasta --exon Hisat2Index/${fasta_baseName}.exon --ss Hisat2Index/${fasta_baseName}.ss Hisat2Index/${fasta_baseName}
# 建议上面内容写成sh脚本，后台运行，好费时间稍微久一点
nohup sh Hisat2Index.sh >Hisat2Index.sh.log &
索引内容如下：
├── GRCh38.0.rf
├── GRCh38.10.rf
├── GRCh38.11.rf
├── GRCh38.1.ht2
├── GRCh38.1.rf
├── GRCh38.2.ht2
├── GRCh38.2.rf
├── GRCh38.3.ht2
├── GRCh38.3.rf
├── GRCh38.4.ht2
├── GRCh38.4.rf
├── GRCh38.5.rf
├── GRCh38.6.rf
├── GRCh38.7.ht2
├── GRCh38.7.rf
├── GRCh38.8.ht2
├── GRCh38.8.rf
├── GRCh38.9.rf
├── GRCh38.exon
├── GRCh38.rf
└── GRCh38.ss
三、数据比对
- --dta：reports alignments tailored for transcript assemblers
- --summary-file：print alignment summary to this file
- -U：单端fq数据
# 激活小环境
conda activate rna
# 创建文件夹
mkdir -p alignment/hisat2
index_base=../GRCh38/Hisat2Index/GRCh38
outdir=alignment/hisat2
ls alignment/rRNA_dup/*gz |while read id
do
sample_name=${id##*/}
sample_name=${sample_name%%.*}
echo "hisat2 --summary-file ${outdir}/${sample_name}_hisat2_summary.txt -p 12 --dta -x $index_base -U $id | samtools view -@ 12 -hbS - >${outdir}/${sample_name}_hisat2.bam"
done >hisat2.sh
# 运行
nohup sh hisat2.sh>hisat2.sh.log&
hisat2.sh的内容：
image-20220327170306770.png
运行完之后目录下每个样本会生成一个比对统计结果文件*_hisat2_summary.txt与比对文件*_hisat2.bam
其中，*_hisat2_summary.txt有总比对率，唯一比对率等比对的重要指标。
网友评论