最近在学习SV相关的资料,简而言之就是碰壁又碰壁,焦头烂额,不得不承认科研这条路如我想象中一样充满艰辛。所以这段时间争取能把发现的一些问题总结一下。
我在使用smartie-sv软件call基因组之间的SV位点过程中发现一个很有意思的问题:拿得到的BED格式结果文件用ANNOVAR软件进行注释,所有的Deletion都注释不出来,这是怎么回事呢?
A01,2555047,2555065,AATAATAAGAAGAAGAAG,-,.,.,.,.,.
A01,2558257,2558302,TAATAATAATAATAATAATAATAATAATAATAATAATAATAATAA,-,.,.,.,.,.
A01,2591196,2591252,TATATATATATATATATATATATATATATATATATATATATATATATATATATATA,-,.,.,.,.,.
A01,2603248,2603281,ATATATATATATATATAATATATATATATATAT,-,.,.,.,.,.
A01,2627079,2627115,TATATATATATATATATATATATATATATATATATA,-,.,.,.,.,.
...
结果如上,原本后面的.
处应该是有信息的,之前我在公司使用ANNOVAR时也没有出现过这种情况。
原来ANNOVAR软件的输入文件是avinput格式,也就是前五列分别为Chr|Start|End|Ref|Alt的文本,我当时是按照BED文件中的坐标直接进行转换,而refGene.txt则是用GFF3格式文件转换,结果一查资料……
参考网站:https://www.biostars.org/p/84686/
原来基因组坐标有两种不同表示方法,分别为0-based与1-based,同一个SV用这两种方法写出来是不一样的。
按照参考资料中给出的坐标转换伪代码,他们之间的关系应该是这样:
#0-based->1-based
if (type=SNV){start=start+1; end=end;}
if (type=DEL){start=start+1; end=end;}
if (type=INS){start=start; end=end+1;}
#1-based->0-based
if (type=SNV){start=start-1; end=end;}
if (type=DEL){start=start-1; end=end;}
if (type=INS){start=start; end=end-1;}
我恍然大悟,原来是因为GFF3文件与BED文件中基因坐标表示方法不一致导致我的Deletion全军覆没,差点翻了大车。
以后在对各种格式的生信文件互相转换时,一定要注意基因组坐标的表示方法:GFF、SAM、VCF等都是1-based,BED、BAM等都是0-based。
PS:最后提醒,如果要转换成ANNOVAR的avinput格式,若Insertion的ref为-,请务必保证start与end坐标一致,不要按照上述教程中end+1,否则也会导致ANNOVAR无法识别!