【问题标题】:How to extract genotype information for each sample as a string from a VCF file using htslib?如何使用 htslib 从 VCF 文件中提取每个样本的基因型信息作为字符串?
【发布时间】:2021-05-10 23:11:04
【问题描述】:

我正在使用 htslib 提取 C++ 中的 VCF 文件中包含的所有信息。

目前,感谢VCF specification 和文件vcf.h 中的文档,我已经成功提取了标题(Meta-Information Lines)中的所有元数据信息,并且大部分信息包含在文件的主体(数据行)。

但是,我不知道如何提取基因型信息(样本列)。

我正在使用来自 1000G 项目的示例文件。这是文件两行的示例,它显示了格式字段和两个样本(该文件每行有超过 1000 个样本,我想提取所有样本的数据):

FORMAT      HG00096                         HG00097
GT:DS:GL    0|0:0.050:-0.48,-0.48,-0.48     0|0:0.050:-0.24,-0.40,-1.49
GT:DS:GL    0|0:0.000:-0.10,-0.69,-4.70     0|0:0.000:-0.05,-0.94,-5.00

我知道这是一项繁重的任务,需要一些计算时间。我已经提取了每列的名称(HG00096、HG00077...),但我不知道如何将每个样本的信息提取为完整的字符串(例如,“0|0:0.050:-0.48,- 0.48,-0.48"),作为键值对的集合(数组、映射、向量...)(例如,[("GT", "0|0"), ("DS", "0.050") , ("GL", "-0.48,-0.48,-0.48")) 或简单地作为值数组(例如,["0|0", "0.050", "-0.48,-0.48,-0.48" ]。我想对每个样本都这样做。

我一直在阅读vcf.h 文件中的文档,我认为函数 bcf_get_genotypes(hdr,line,dst,ndst) 可能适用于此,但我不确定如何使用它将值提取为字符串。另外,我认为这些信息可能存储在“bcf_fmt_t”的“p”指针中,但我不确定,它只包含一组 uint8_t 值的数组,我不知道是否为字符串(或char array) 可以按我想要的方式提取。

typedef struct bcf_fmt_t {
    int id;             
    int n, size, type;  
    uint8_t *p;        
    uint32_t p_len;
    uint32_t p_off:31, p_free:1;
} bcf_fmt_t;

有没有我想做的方法?

【问题讨论】:

    标签: c++ vcf-vcard vcf-variant-call-format vcftools htslib


    【解决方案1】:

    我终于明白了。有一些函数可以执行此操作,具体取决于格式 ID 标头中指定的类型:这些函数位于 htslib 中的 vcf.h 文件中:

    #define bcf_get_format_int32(hdr,line,tag,dst,ndst)  
    bcf_get_format_values(hdr,line,tag,(void**)(dst),ndst,BCF_HT_INT)
    #define bcf_get_format_float(hdr,line,tag,dst,ndst)  
    bcf_get_format_values(hdr,line,tag,(void**)(dst),ndst,BCF_HT_REAL)
    #define bcf_get_format_char(hdr,line,tag,dst,ndst)   
    bcf_get_format_values(hdr,line,tag,(void**)(dst),ndst,BCF_HT_STR)
    #define bcf_get_genotypes(hdr,line,dst,ndst)         
    bcf_get_format_values(hdr,line,"GT",(void**)(dst),ndst,BCF_HT_INT)
    
    HTSLIB_EXPORT
    int bcf_get_format_string(const bcf_hdr_t *hdr, bcf1_t *line, const char *tag, char ***dst, int *ndst);
    
    HTSLIB_EXPORT
    int bcf_get_format_values(const bcf_hdr_t *hdr, bcf1_t *line, const char *tag, void **dst, int *ndst, int type);
    

    【讨论】:

      猜你喜欢
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 2014-03-03
      • 1970-01-01
      • 2015-01-25
      • 1970-01-01
      相关资源
      最近更新 更多