RNA-Seq STAR 比对学习1

作者: 7f0a92cda77c | 来源:发表于2021-05-12 16:58 被阅读0次

RNA-seq 检测变异之 GATK 最佳实践流程
STAR的使用说明
RNA-Seq STAR 比对学习1
使用STAR-fusion进行融合基因的分析
比对软件STAR的使用
「RNA-seq分析软件」RNA-seq比对工具STAR学习笔记
STAR 使用手册
RNA-seq数据分析—STAR比对
RNA-seq：STAR 软件比对（SLURM递交系统）
RNAseq相关

最近使用STAR比对出现了状况，软件版本和对应的Index不相容，报错信息如下

#EXITING because of FATAL ERROR: Genome version: 2.7.0d is INCOMPATIBLE with running STAR version: 2.7.9a
#SOLUTION: please re-generate genome from scratch with running version of STAR, or with version: 2.7.4a

原因是之前构建好的index是使用以前版本的STAR软件构建好的；现在因为重新安装了STAR软件，出现了Genome version和STAR version不兼容的问题

第一步-直接下载了STAR软件

由于conda 环境出现问题，所以就直接下载了STAR软件

#https://github.com/alexdobin/STAR 提供版本是2.7.9a
wget https://github.com/alexdobin/STAR/archive/2.7.9a.tar.gz
tar -xzf 2.7.9a.tar.gz
cd /public/vip/biosoft/STAR-2.7.9a/source
make STAR #全部大写，命令是大写的，对其进行编译

第二步-构建这个版本对应的Index

mkdir  mm10_2_7_9a
/public/vip/biosoft/STAR-2.7.9a/source/STAR --runMode genomeGenerate \
    --genomeDir /data/server/reference/index/star/mm10_2_7_9a \ #将准备构建好的index存放在哪个文件夹，对应第一行构建的那个文件夹的路径
    --genomeFastaFiles /data/server/reference/genome/mm10/mm10.fa \ #参考基因组
    --sjdbGTFfile /data/server/reference/gtf/gencode/gencode.vM27.annotation.gtf  \ 
    --sjdbOverhang 149 --runThreadN 4
#starting to sort Suffix Array. This may take a long time

构建耗时比较长

最后生成的文件是这个：

sudo chmod 0755 * 
#把权限给改了下，可以执行
ls -hlt
# total 26G
# -rwxr-xr-x 1 root root 318K 5月  13 03:23 Log.out
# -rwxr-xr-x 1 root root 1.5G 5月  13 03:23 SAindex
# -rwxr-xr-x 1 root root  22G 5月  13 03:23 SA
# -rwxr-xr-x 1 root root 2.7G 5月  13 03:21 Genome
# -rwxr-xr-x 1 root root  842 5月  13 03:21 genomeParameters.txt
# -rwxr-xr-x 1 root root 8.0M 5月  13 03:12 sjdbInfo.txt
# -rwxr-xr-x 1 root root 7.1M 5月  13 03:12 sjdbList.out.tab
# -rwxr-xr-x 1 root root  493 5月  13 02:10 chrLength.txt
# -rwxr-xr-x 1 root root 1.4K 5月  13 02:10 chrNameLength.txt
# -rwxr-xr-x 1 root root  912 5月  13 02:10 chrName.txt
# -rwxr-xr-x 1 root root  720 5月  13 02:10 chrStart.txt
# -rwxr-xr-x 1 root root 8.7M 5月  13 02:10 sjdbList.fromGTF.out.tab
# -rwxr-xr-x 1 root root  12M 5月  13 02:10 exonInfo.tab
# -rwxr-xr-x 1 root root 9.5M 5月  13 02:10 transcriptInfo.tab
# -rwxr-xr-x 1 root root 2.3M 5月  13 02:10 geneInfo.tab
# -rwxr-xr-x 1 root root  29M 5月  13 02:10 exonGeTrInfo.tab

第三步-建立STAR比对

3.1在比对之前先创建一个文本文件，存放测序的文件名前缀，保证它可以跟后续的文件联系上

ls /public/vip/project/clean/*val_1.fq.gz|cut -d "/" -f 6|cut -d "_" -f 1 
#1-PC-LCJ9004
#2-CP-LCJ9001
#4-CP-LCJ9003
#对应的是它相关的数据
#1-PC-LCJ9004_combined_R1_val_1.fq.gz
#2-CP-LCJ9001_combined_R1_val_1.fq.gz
#4-CP-LCJ9003_combined_R1_val_1.fq.gz
ls /public/vip/project/clean/*val_1.fq.gz|cut -d "/" -f 6|cut -d "_" -f 1 > list.txt

3.2做个测试，将要比对的文件进行循环输出，并检查是否正确

cat list.txt |while read id; do echo -n "/public/vip/project/clean/${id}_combined_R1_val_1.fq.gz /public/vip/project/clean/${id}_combined_R2_val_2.fq.gz "; done

3.3正确后，再进行下面的操作，将要执行的命令进行了循环输出到文件夹

要改的就只有比对的文件的差异，其它参数都一样，所以可以这样操作

cat list.txt |while read id;
do echo -n "/public/vip/biosoft/STAR-2.7.9a/source/STAR --runThreadN 5 "; echo -n "--genomeDir /data/server/reference/index/star/mm10_2_7_9a "; 
echo -n "--outSAMtype BAM SortedByCoordinate --outReadsUnmapped Fastx "; 
echo -n "--quantMode GeneCounts --readFilesCommand zcat --twopassMode Basic "; 
echo -n "--outFilterMismatchNmax 999 --outFilterMismatchNoverReadLmax 0.04 "; 
echo -n "--outFilterType BySJout --outFilterMultimapNmax 20 "; 
echo -n "--alignSJoverhangMin 8 --alignSJDBoverhangMin 1 "; 
echo -n "--chimSegmentMin 20 --chimJunctionOverhangMin 20 --chimOutJunctionFormat 1 "; 
echo -n "--alignIntronMin 20 --alignIntronMax 1000000 --alignMatesGapMax 1000000 "; 
echo -n "--chimSegmentReadGapMax 0 --alignSJstitchMismatchNmax 0 -1 0 0 "; 
echo -n  "--readFilesIn /public/vip/project/clean/${id}_combined_R1_val_1.fq.gz /public/vip/project/clean/${id}_combined_R2_val_2.fq.gz --outFileNamePrefix ${id}"; done > star.sh

上述步骤均是参考了生信菜鸟团的帖子

less star.sh #循环生成的STAR比对命令
#/public/vip/biosoft/STAR-2.7.9a/source/STAR 
#--runThreadN 5 
#--genomeDir /data/server/reference/index/star/mm10_2_7_9a --outSAMtype BAM SortedByCoordinate 
#--outReadsUnmapped Fastx 
#--quantMode GeneCounts 
#--readFilesCommand zcat 
#--twopassMode Basic 
#--outFilterMismatchNmax 999 
#--outFilterMismatchNoverReadLmax 0.04 
#--outFilterType BySJout 
#--outFilterMultimapNmax 20 
#--alignSJoverhangMin 8 
#--alignSJDBoverhangMin 1 
#--chimSegmentMin 20 
#--chimJunctionOverhangMin 20 
#--chimOutJunctionFormat 1 
#--alignIntronMin 20 
#--alignIntronMax 1000000 
#--alignMatesGapMax 1000000 
#--chimSegmentReadGapMax 0 
#--alignSJstitchMismatchNmax 0 -1 0 0 
#--readFilesIn /public/vip/project/clean/1-PC-LCJ9004_combined_R1_val_1.fq.gz /public/vip/project/clean/1-PC-LCJ9004_combined_R2_val_2.fq.gz 
#--outFileNamePrefix 1-PC-LCJ9004

3.3正确后，再进行下面的操作

/public/vip/biosoft/STAR-2.7.9a/source/STAR --runThreadN 5 --genomeDir /data/server/reference/index/star/mm10_2_7_9a --outSAMtype BAM SortedByCoordinate --outReadsUnmapped Fastx --quantMode GeneCounts --readFilesCommand zcat --twopassMode Basic --outFilterMismatchNmax 999 --outFilterMismatchNoverReadLmax 0.04 --outFilterType BySJout --outFilterMultimapNmax 20 --alignSJoverhangMin 8 --alignSJDBoverhangMin 1 --chimSegmentMin 20 --chimJunctionOverhangMin 20 --chimOutJunctionFormat 1 --alignIntronMin 20 --alignIntronMax 1000000 --alignMatesGapMax 1000000 --chimSegmentReadGapMax 0 --alignSJstitchMismatchNmax 0 -1 0 0 --readFilesIn /public/vip/project/clean/1-PC-LCJ9004_combined_R1_val_1.fq.gz /public/vip/project/clean/1-PC-LCJ9004_combined_R2_val_2.fq.gz --outFileNamePrefix 1-PC-LCJ9004

如果这个运行成功，可以再批量后台运行这个

nohup bash star.sh &
top 
 #PID USER      PR  NI    VIRT    RES    SHR S  %CPU %MEM     TIME+ COMMAND                                                                                           
#27025 vip       20   0 28.500g 0.025t   5556 R 493.7 40.9   5:15.71 STAR                                                                                              
#27033 vip       20   0    4920    852    788 S   7.6  0.0   0:06.46 gzip                                                                                              
#27029 vip       20   0    4920    884    820 S   7.0  0.0   0:06.27 gzip

对于所选参数，再看下