坐标转换困境
一些公开发表的论文中及很多数据库中经常提到变异,一般变异的表现形式有三种:1)基因组坐标:2)cDNA 坐标;3)蛋白氨基酸坐标。举个例子TP53上的某个变异的基因组坐标是g.chr17:74026C>A,cDNA坐标是c.1001G>T,蛋白氨基酸坐标是p.G334V。在数据分析的过程中经常会遇到这三种坐标相关转换的情况,例如你从文献或者某个数据库中收集到了几百个肿瘤靶向药的用药位点,而你在你样本中检测到了很多变异,想知道你的样本中包含多少收集到的已知的用药位点。但通常文献或者数据库会以第二种或者第三种形式表示变异,而我们自己检测的变异通常会以vcf格式存储,这样就无法直接匹配。当然可以对vcf格式的变异进行ANNOVAR注释,然后对cDNA或者蛋白氨基酸坐标形式的变异进行比较,但尝试过的人都表示特别痛苦:需要考虑的规则太多!尝试两次,还是放弃了:一是匹配规则不通用;二是总担心有没有考虑到过的情况。所以急需一个能完成这种坐标转换的工具。15年发表在NATURE METHODS上的题为:TransVar: a multilevel variant annotator for precision genomics的文章中推出了一款名为TransVar的软件成了解决不同层面变异坐标转换的神器。下面小编就介绍一下这款软件(Linux版),没有Linux基础的也不用担心,后续会写一篇基于Web版TransVar进行注释(坐标转换)的文章。
TransVar软件简介
Transvar 是一款多种方向的突变/坐标转换工具,它支持基因组坐标、cDNA 坐标以及蛋白氨基酸坐标之间的转换。
如上图所示,该软件的功能可细分为下面3种:
1)正向注释:对于基因组坐标的变异进行mRNA(cDNA)和蛋白注释,这款工具会提供所有的可能结果;
2)反向注释:将mRNA(cDNA)坐标和蛋白坐标的变异转换成所有可能基因组坐标形式的变异;
3)等价注释:对于某一给定的蛋白坐标的变异,搜索所有可能的与其为相同基因组坐标,但在不同转录本上的蛋白坐标变异。
软件下载和安装:
软件下载地址:
1,旧版(最近没有在更新):https://bitbucket.org/wanding/transvar/src/master/
2,新版(一直在更新):https://github.com/zwdzwd/transvar
按照方法如下:
sudo pip install transvar ## 全局安装,需要root权限
或者:
pip install --user transvar ##用户安装,没有root权限的用此方法
软件更新:
pip install -U transvar
这款软件在安装后要自己配置数据库操作起来也比较简单:
# set up databases
transvar config --download_anno --refversion hg19 #默认的hg19的 dbSNP 数据库是2016年的,部分数据库如dbSNP新版数据库收录内容有很大变化(主要是数量的提升),所以建议自行重新下载
# in case you don't have a reference
transvar config --download_ref --refversion hg19
# in case you do have a reference to link
transvar config -k reference -v [path_to_hg19.fa] --refversion hg19
需要注意的是直接使用Transvar的命令下载数据库容易因网络问题出错,导致下载的数据库是不完整的(不报错的,是个深坑!)可以到http://transvar.info/transvar_user/annotations/直接下载后进行配置。
软件的使用
这款软件即可以单点注释,也可以批量处理,下面分别介绍一下:
单点注释用 -i传入待注释位点,包括3种:
# 基因组正向注释
transvar ganno --ccds -i 'chr3:g.178936091G>A'
# cDNA反向注释
transvar canno --ccds -i 'PIK3CA:c.1633G>A'
# 氨基酸反向注释
transvar panno -i 'PIK3CA:p.E545K' --ensembl
# 其中--ccds、--ensembl为使用不同的数据库,如网页版,可以同时多选,\
# 如 --ccds --ensembl --refseq --ucsc 来进行多选
批量注释:
/*/software/anaconda3/bin/transvar canno -l mutiation.canno.list -m 1 -o 2 --refseq --longestcoding --gseq
###
canno:指cDNA反向注释,备选包括panno( 蛋白氨基酸反向注释)和ganno(基因组正向注释)
-l:输入文件,变异与canno、panno、ganno对应。格式示例如下:
![image.png](https://upload-images.jianshu.io/upload_images/22041438-ba466242c2050f60.png?imageMogr2/auto-orient/strip%7CimageView2/2/w/1240)
-m:-l指定的输入文件可以有多列,通过-m指定哪列是待注释列,不加-m参数默认是第一列
-o:同时可以通过-o来指定-l中的那一列作为输出文件的首列,不加-o,默认是第一列
--refseq:使用哪个数据库的转录本进行注释,还有其他数据库可选如 ensembl/gencode/ucsc/ccds/aceview等。
--longestcoding: 有多个转录本时,仅选择最长的转录本。如果不加这个参数会把涉及到的所有转录本都输出出来,这时候你就要自己制定标准进行筛选了
--gseq :在输出文件中增加类似VCF格式的变异信息,包括染色体,起始位置,终止位置,参考基因组序列,突变后的序列。
软件官方教程
官网:https://transvar.readthedocs.io/en/latest/
这里有对软件详细的介绍,这里就不赘述了,想深入研究的可以去官网看看。
说在最后的
transvar 在转换时总会有很多损失,个人经验损失主要来自于两部分:
1,输出结果中没有该变异,直接被丢掉了;
2,输出结果中有该变异,但在你选择的数据库中没有这个转录本,提示“no_valid_transcript_found”。
为了尽量提高成功转换的比例可以做如下尝试:
1,用所有能用的库去注释,不过还是建议以一个库的结果为准,把其它库包含但该库不包含的变异加上;
2,对于longestcoding没有成功去掉该参数后再尝试,然后自行选一个靠谱的转录本,如果不知道该怎么选就随机选一个;
3,如果你拿到的变异信息有对应的转录本,选取与所提供的转录本一致的数据库,分析时不加--longestcoding,然后根据转录本信息对转换结果进行匹配,这种是准确性最高的。
原创文字,如果觉得对你有帮助留下你的赞哦~