教程及软件有错误是常态

有些时候，我们可能真的很努力了，但事情就是解决不了，这个时候仍然是自己的问题，努力能解决一些问题，但不是所有。

最近在用star找融合基因，对人类物种的数据来说，没啥问题，很轻松就运行成功了，合作方给我了一些小鼠的测序数据，居然出现了问题，当然，我下载数据库肯定不会出错的啦。

├── [ 26G]  GRCh37_gencode_v19_CTAT_lib_Nov012017.plug-n-play.tar.gz
├── [ 20G]  GRCh38_gencode_v26_CTAT_lib_Nov012017.plug-n-play.tar.gz
├── [ 24G]  Mouse_M15_CTAT_lib_Nov012017.plug-n-play.tar.gz

报错如下：

-parsing GTF file: /home/jianmingzeng/biosoft/starFusion/db/Mouse_M15_CTAT_lib_Nov012017/ctat_genome_lib_build_dir//ref_annot.gtf
Error, cannot get transcript_id from gene_id "ENSMUSG00000102693.1"; gene_type "TEC"; gene_name "4933401J01Rik"; level 2; havana_gene "OTTMUSG00000049935.1"; of line
chr1    HAVANA  gene    3073253 3074322 .   +   .   gene_id "ENSMUSG00000102693.1"; gene_type "TEC"; gene_name "4933401J01Rik"; level 2; havana_gene "OTTMUSG00000049935.1"; at /home/jianmingzeng/miniconda3/lib/STAR-Fusion/util/../lib/GTF_utils.pm line 111, <$fh> line 6.
    GTF_utils::GTF_to_gene_objs("/home/jianmingzeng/biosoft/starFusion/db/Mouse_M15_CTAT_lib_N"...) called at /home/jianmingzeng/miniconda3/lib/STAR-Fusion/util/../lib/GTF_utils.pm line 30
    GTF_utils::index_GTF_gene_objs_from_GTF("/home/jianmingzeng/biosoft/starFusion/db/Mouse_M15_CTAT_lib_N"..., HASH(0x1998ac8)) called at /home/jianmingzeng/miniconda3/lib/STAR-Fusion/util/../lib/GTF_utils.pm line 20
    GTF_utils::index_GTF_gene_objs("/home/jianmingzeng/biosoft/starFusion/db/Mouse_M15_CTAT_lib_N"..., HASH(0x1998ac8)) called at /home/jianmingzeng/miniconda3/lib/STAR-Fusion/util/STAR-Fusion.predict line 331
    main::parse_GTF_features("/home/jianmingzeng/biosoft/starFusion/db/Mouse_M15_CTAT_lib_N"..., HASH(0x23174a8), HASH(0x23174d8)) called at /home/jianmingzeng/miniconda3/lib/STAR-Fusion/util/STAR-Fusion.predict line 100

看起来很复杂，又是perl又是模块，还有gtf的问题，事实上只需要看第一句话，就是star-fusion自带的基因注释数据库文件里面的gtf问题，有一些基因是没有对应转录本的。

我去gencode数据库里面搜索了同样的ID，发现这个问题一模一样，很有可能是gencode那边的人并没有做好这个数据库的维护，而star-fusion的开发团队直接沿用了这个数据库文件。

但是为什么star-fusion的开发团队并没有去仔细检查呢？

难道他们开发了数据库缺没有实际上运行任何一个小鼠的例子吗？

当然了，解决方案非常简单：

grep gene_id ref_annot.gtf |grep transcript_id

就是挑选那些既有基因ID，又有转录本ID的那些咯，可以看到原来是1745955行，经过过滤只剩下1693314行。

当然了，软件也成功运行了。

├── [4.0K]  star-fusion.filter.intermediates_dir
│   ├── [ 47K]  star-fusion.pre_blast_filter
│   ├── [ 19K]  star-fusion.pre_blast_filter.abridged
│   ├── [105K]  star-fusion.pre_blast_filter.filt_info
│   ├── [ 65K]  star-fusion.pre_blast_filter.filt_info.abridged
│   ├── [ 48K]  star-fusion.pre_blast_filter.post_blast_n_promisc_filter
│   └── [ 21K]  star-fusion.pre_blast_filter.post_blast_n_promisc_filter.abridged
├── [ 42K]  star-fusion.fusion_candidates.final
├── [ 17K]  star-fusion.fusion_candidates.final.abridged
├── [ 87K]  star-fusion.fusion_candidates.preliminary
├── [ 131]  star-fusion.predict.intermediates_dir
│   ├── [ 22M]  star-fusion.junction_breakpts_to_genes.txt
│   ├── [125K]  star-fusion.junction_read_names
│   └── [1.6M]  star-fusion.spanning_frag_names
├── [   0]  star-fusion.STAR-Fusion.filter.ok
└── [   0]  star-fusion.STAR-Fusion.predict.ok