从 vcf 文件准备 ANNOVAR 数据库

annovar可以说是最常用的变异注释软件了，除了基于基因位置进行注释，还有丰富的第三方数据库支持，clinvar, cosmic等等，但是annovar提供下载数据库版本较老，需要自行下载第三方 vcf 进行转换。

以 clinvar 为例说说 vcf 格式转换成 annovar 可识别表格格式，annovar 仅需要 vcf 中部分信息，为了脚本通用这里输出vcf所有信息到 tab 格式，详细处理脚本如下，同时实现脚本存放 github。

import gzip


def Parse_Variant_Annovar_Style(line):
	CHROM, POS, ID, REF, ALT, QUAL, FILTER, INFO = line.split('\t')[:8]
    # 一般公共数据下载的 vcf 文件仅包含8列
    POS = int(POS)
    dat = {}
    dat['CHROM'] = CHROM
    dat['ID'] = ID
    if len(REF) == len(ALT) == 1 and REF != ALT:    # SNV
        dat['START'], dat['END'] = POS, POS
        dat['REF'], dat['ALT'] = REF, ALT
    elif len(REF) == 1 and len(ALT) != 1 and ALT.startswith(REF):    # Insert
        dat['START'], dat['END'] = POS, POS
        dat['REF'], dat['ALT'] = '-', ALT[1:]
    elif len(REF) != 1 and len(ALT) == 1 and REF.startswith(ALT):    # Delete
        dat['START'], dat['END'] = POS+1, POS+len(REF)-1
        dat['REF'], dat['ALT'] = REF[1:], '-'
    else:    # complex variant
        dat['START'], dat['END'] = POS, POS+len(REF)-1
        dat['REF'], dat['ALT'] = REF, ALT
    dat.update(INFO2dict(INFO))
    return dat


def Parse_Variant(line):
	CHROM, POS, ID, REF, ALT, QUAL, FILTER, INFO = line.split('\t')[:8]
    # 一般公共数据下载的 vcf 文件仅包含8列
    POS = int(POS)
    dat = {}
    dat['CHROM'] = CHROM
    dat['ID'] = ID
    dat['START'], dat['END'] = POS, POS+len(REF)-1
    dat['REF'], dat['ALT'] = REF, ALT
    dat.update(INFO2dict(INFO))
    return dat


def INFO2dict(info):
    tmp = [i.strip() for i in info.split(';') if i.strip()]
    infor = {}
    for i in tmp:
        if '=' in i:
            # index = i.index('=')
            # infor[i[:index]] = i[index+1:]
            i = i.split('=')
            infor[i[0]] = i[1]
        else:
            infor[i] = 'true'
    return infor


fopen = gzip.open if clinvar.endswith('gz') else open
# 这种文件打开适用于 python2, python3 中存在 bits 和 string 类型的差别
with fopen(input_vcf) as f:
    tags = []
    for i in f:
        if not i.strip():
            continue

        if i.startswith('##INFO='):
             tags.append(i.split(',')[0].split('=')[-1])
        elif i.startswith('#CHROM'):
             print('#CHROM\tSTART\tEND\tID\tREF\tALT\t{}\n'.format('\t'.join(tags)))
        elif not i.startswith('#'):
             dat = Parse_Variant_Annovar_Style(i)
             # dat = Parse_Varitn(i)
             for x, label in enumerate(['CHROM', 'START', 'END', 'ID', 'REF', 'ALT']+tags):
                if x != 0:
                    print('\t')
                print(dat.get(label, '--'))    # INFO 列中某些 tag 并不存在于所有行
             print('\n')