最好用的融合基因查找工具终于正式发表了

就是STAR-fusion啦，它可以直接基于STAR比对好的bam文件来做分析，而大多数其它融合基因查找工具，需要从fastq文件开始，不太方便。之前我在生信技能树公众号介绍过它，那个时候发表该工具的文章是：STAR-Fusion: Fast and Accurate Fusion Transcript Detection from RNA-Seq  在biorxiv预印本：

image-20191029093943451

两年前我在生信技能树就写过它的教程：使用STAR-fusion来对转录组数据找融合基因 也强调大家自行读STAR-Fusion详细说明：https://github.com/STAR-Fusion/STAR-Fusion/wiki 今天（2019-10-29  ）翻它的主页发现上面写着已经正式发表了：

image-20191029094043567

还设计了一个好看（简洁）的logo，我留意到作者少了很多：

image-20191029094535954

而且才意识到，居然是 Aviv Regev 实验室出品的工具！如果生信技能树的粉丝们有在 Aviv Regev 实验室做科研的（博士，博士后均可），希望可以联系我一下，认识一下。

再次回顾STAR-fusion对RNA-seq测序数据找融合基因的流程

生物信息学鉴定融合转录本的方法一般有两种：

  • ①将RNA-seq数据与Reference genome做alignment，鉴别可能发生重排的基因；

  • ②先直接将reads装配成更长的转录本序列，再鉴别与重排序列一致的融合转录本。

而我们选择的Broad Institute的 Brian J. Haas 和冷泉港实验室（CSHL）的 Alex Dobin 等人开发的工具STAR-Fusion，其工作原理分为三步：

  • ①先将reads通过STAR比对到参考基因组，筛选出split和discordant reads作为候选的融合基因序列；

  • ②将候选融合基因序列与参考基因序列进行比对，根据overlaps预测出融合基因；

  • ③对预测结果做过滤，去除假阳性结果。

首先对RNA-seq测序数据进行star的two-pass比对：

批量运行，其中star_index 需要指定自己的star索引路径文件夹，而bin_star需要指定自己的star软件的可执行程序路径哦。

核心代码是：

start=$(date +%s.%N)
echo star start `date`
$bin_star --runThreadN  4  --genomeDir $star_index  \
--twopassMode Basic --outReadsUnmapped None --chimSegmentMin 12  \
--alignIntronMax 100000 --chimSegmentReadGapMax parameter 3  \
--alignSJstitchMismatchNmax 5 -1 5 5  \
--readFilesCommand zcat --readFilesIn $fq1 $fq2 --outFileNamePrefix  ${sample}_
mv ${sample}_Aligned.out.sam $sample.sam
$bin_samtools sort -o $sample.bam  $sample.sam
$bin_samtools index $sample.bam
$bin_samtools flagstat $sample.bam  > $sample.flagstat
touch  ok.star.$sample.status
rm  $sample.sam
echo star  end  `date`
dur=$(echo "$(date +%s.%N) - $start" | bc)
printf "Execution time for star : %.6f seconds" $dur

实际上就是一行命令在运行比对过程，但是呢，参数太多了，调起来很麻烦，通常如果不理解的话就不建议修改参数。 有一个链接需要注意：https://github.com/STAR-Fusion/STAR-Fusion/issues/104

值得注意的是，star这个比对软件，也需要选择比较新的版本哦：Please be sure to use the STAR aligner with min version: 2.7.0f 毕竟我上次写教程介绍它的时候是2017年。软件需要最新版，下载的数据库也需要最新版！

image-20191029104053046

比对本来是为了得到bam文件，但是我们这个流程是为了 得到'Chimeric.out.junction'这个文件，供STAR-fusion使用。最新版软件，关于得到'Chimeric.out.junction'这个文件，有一个参数：chimOutJunctionFormat 需要注意。

然后使用conda安装STAR-fusion

因为这个工具大量依赖perl模块，一般人不会这个语言，弄起来也麻烦，所以推荐conda一键式安装

一步法安装代码如下：

conda install -c bioconda star-fusion

可以看到下面那么多的perl模块都会被自动安装：

image-20191029094324826

一般来说，安装成功后，可以使用perl代码检测一下：

perl -e 'use Set::IntervalTree'
# 值得注意的是
STAR-Fusion 
#软件命令是 STAR-Fusion

如果没有明显的报错，说明你安装成功啦。如果你不想使用conda一步到位的安装，那么你会面临perl模块的巨坑：

  A typical perl module installation may involve:
   perl -MCPAN -e shell
   install DB_File
   install URI::Escape
   install Set::IntervalTree
   install Carp::Assert
   install JSON::XS
   install PerlIO::gzip

所以自己衡量一下时间和精力哦，不要逞强！

通常是软件版本问题

比如star版本就很奇怪，在GitHub是最新的：

# Get latest STAR source from releases
wget https://github.com/alexdobin/STAR/archive/2.7.3a.tar.gz
tar -xzf 2.7.3a.tar.gz
cd STAR-2.7.3a

# Alternatively, get STAR source using git
git clone https://github.com/alexdobin/STAR.git

但是conda安装的是2.7.0d , 就有可能报错：

qiEXITING because of FATAL ERROR: Genome version: 2.7.1a is INCOMPATIBLE with running STAR version: 2.7.0f
SOLUTION: please re-generate genome from scratch with running version of STAR, or with version: 2.7.0d
's
Oct 29 20:10:37 ...... FATAL ERROR, exiting

其实conda可以安装指定版本的软件哦！

而且上面的那个报错其实也不是字面上的意思。

需要下载指定参考基因组的数据库文件供STAR-fusion使用

这个需要考验大家网速了，其wiki说明：https://github.com/STAR-Fusion/STAR-Fusion/wiki 写的很清楚。需要注意的是保证参考基因组的一致性哦！

如果你有自己的star软件索引，就下载2.9G的，如果没有，就下载29G的，如果你是人类，下载37和38均可，取决于你对参考基因组的熟悉程度！https://data.broadinstitute.org/Trinity/CTAT_RESOURCE_LIB/  如果你不熟悉，最好是选择最新版，不要选择我下面截图的那些！！！

反正我下载的是：GRCh38_gencode_v31_CTAT_lib_Oct012019.plug-n-play.tar.gz ，足足31个G！

image-20191029112135492

而且需要下载 ; be sure to use a more modern version of the companion CTAT_GENOME_LIB  如果你报错的话，通常就是这些问题。

比如旧版的如下（文件夹里面是star的index文件）：

image-20191029105012523

新版数据库如下：

image-20191029142024193

运行STAR-fusion

同样也是批量运行：

lib='/home/yb77613/biosoft/starFusion/db/GRCh38_gencode_v31_CTAT_lib_Oct012019.plug-n-play/ctat_genome_lib_build_dir/'

# Please be sure to use the STAR aligner with min version: 2.7.0f 
ls ../alignment/*.junction  |while read id
do
file=$(basename $id )
sample=${file%%.*}
~/biosoft/starFusion/STAR-Fusion/STAR-Fusion --genome_lib_dir $lib -J $id --output_dir $sample
done

或者单独运行：

lib='/home/yb77613/biosoft/starFusion/db/GRCh38_gencode_v31_CTAT_lib_Oct012019.plug-n-play/ctat_genome_lib_build_dir/'
STAR-Fusion --genome_lib_dir $lib -J Lib_FUSCCTNBC001_Chimeric.out.junction  --output_dir s1

实际上，是非常多命令的组合，只不过你自己只需要运行一句话命令即可，比如：

CMD: mkdir -p /home/yb77613/data/public/tnbc/RNA-seq/star/s1
CMD: mkdir -p /home/yb77613/data/public/tnbc/RNA-seq/star/s1/_starF_checkpoints
CMD: mkdir -p /home/yb77613/data/public/tnbc/RNA-seq/star/s1/star-fusion.preliminary
-sample contains 41047586
* Running CMD: /home/yb77613/miniconda3/envs/rna/lib/STAR-Fusion/util/STAR-Fusion.map_chimeric_reads_to_genes  --genome_lib_dir /home/yb77613/biosoft/starFusion/db/GRCh38_gencode_v31_CTAT_lib_Oct012019.plug-n-play/ctat_genome_lib_build_dir/  -J /home/yb77613/data/public/tnbc/RNA-seq/star/Lib_FUSCCTNBC001_Chimeric.out.junction  > /home/yb77613/data/public/tnbc/RNA-seq/star/s1/star-fusion.preliminary/star-fusion.junction_breakpts_to_genes.txt

具体大家可以去查看log日志哈。

输出文件的解读

其实说明书也写的很清楚了，主要就是看JunctionReads列，看看有多少条序列支持这样的融合，然后也会列出融合基因的断点出，左右两个基因的衔接位置的碱基等等。

#FusionName    JunctionReadCount   SpanningFragCount   SpliceType  LeftGene    LeftBreakpoint  RightGene   RightBreakpoint LargeAnchorSupport  FFPM    LeftBreakDinuc  LeftBreakEntropy    RightBreakDinuc RightBreakEntropy   annots
IGH-@--MALAT1    4   16  INCL_NON_REF_SPLICE IGH-@^IGH-.g@   chr14:105708080:+   MALAT1^ENSG00000251562.8    chr11:65499045:+    YES_LDAS    0.4385  TG  1.7232  GT  1.8892  ["TCGA_StarF2019","INTERCHROMOSOMAL[chr14--chr11]"]
IGH-@-ext--MALAT1    4   16  INCL_NON_REF_SPLICE IGH-@-ext^IGH-.g@-ext   chr14:105708080:+   MALAT1^ENSG00000251562.8    chr11:65499045:+    YES_LDAS    0.4385  TG  1.7232  GT  1.8892  ["INTERCHROMOSOMAL[chr14--chr11]"]
IGH-@--KRT19    5   9   INCL_NON_REF_SPLICE IGH-@^IGH-.g@   chr14:105742042:+   KRT19^ENSG00000171345.13    chr17:41528296:-    YES_LDAS    0.3069  GC  1.6729  CG  1.9899  ["TCGA_StarF2019","INTERCHROMOSOMAL[chr14--chr17]"]
IGH-@-ext--KRT19    5   9   INCL_NON_REF_SPLICE IGH-@-ext^IGH-.g@-ext   chr14:105742042:+   KRT19^ENSG00000171345.13    chr17:41528296:-    YES_LDAS    0.3069  GC  1.6729  CG  1.9899  ["INTERCHROMOSOMAL[chr14--chr17]"]
AFF4--MAPK8    2   10  ONLY_REF_SPLICE AFF4^ENSG00000072364.13 chr5:132937067:-    MAPK8^ENSG00000107643.16    chr10:48401612:+    YES_LDAS    0.2632  GT  1.9219  AG  1.7232  ["INTERCHROMOSOMAL[chr5--chr10]"]
AFF4--MAPK8    1   10  ONLY_REF_SPLICE AFF4^ENSG00000072364.13 chr5:132937067:-    MAPK8^ENSG00000107643.16    chr10:48325897:+    YES_LDAS    0.2412  GT  1.9219  AG  1.9086  ["INTERCHROMOSOMAL[chr5--chr10]"]
RRM2B--AC016074.2    2   4   ONLY_REF_SPLICE RRM2B^ENSG00000048392.11    chr8:102224046:-    AC016074.2^ENSG00000286122.1    chr8:125650233:+    YES_LDAS    0.1316  GT  1.8323  AG  1.7465  ["INTRACHROMOSOMAL[chr8:23.41Mb]"]

然后新版的STAR-fusion增加了一个过滤的功能，主要是过滤掉那些线粒体基因融合好HLA相关基因融合，说明书一直提到了：Red Herrings: Fusion pairs that may not be relevant to cancer, and potential false positives. 主要是因为这个软件开发就是为癌症研究服务的。

如果你确实不需要这个默认的过滤功能，就需要修改参数啦。

STAR-fusion还有两个好哥们

If you plan on using the included FusionInspector for 'inspect' or 'validate' modes, please install the FusionInspector dependencies. 主要是可视化你的融合事件，通常是基于IGV。还有个是 FusionAnnotator 主要是基于癌症数据库，来为你的融合事件推荐注释。

后面我们根据粉丝的学习情况来看看是否需要加餐继续介绍！

赞 (0)

相关推荐