Question

locateVariants with sample names

0

Entering edit mode

by0 • 0

@by0-7984

Last seen 10.5 years ago

United Kingdom

I'm using the `VariantAnnotation` package to extract variants from a multi-sample vcf file like this:

library(VariantAnnotation)
library(TxDb.Scerevisiae.UCSC.sacCer3.sgdGene)

vcf    <- readVcf('~/Desktop/yeast.vcf', genome="sacCer3")
target <- rowData(vcf)
loc <- locateVariants(target, TxDb.Scerevisiae.UCSC.sacCer3.sgdGene, AllVariants())
names(loc) <- NULL
out <- as.data.frame(loc)

Is there anyway to have, for each variant, the sample or samples that have that variant in the output data?

Thank you

annotation variantannotation sample • 2.0k views

ADD COMMENT • link 10.5 years ago • updated 10.4 years ago by0 • 0

score 1 · Answer 1 · 2015-06-01

1

Entering edit mode

Michael Lawrence ★ 11k

@michael-lawrence-3846

Last seen 3.9 years ago

United States

The result of locateVariants has a "QUERYID" column that indexes back into the input object. So it's possible to join the input with the output. It would help to know more about your use case. One way to represent the information on a per-sample basis is to use a "long-form" table (VCF is "wide-form"), where each row corresponds to a position and sample. You could get there by coercing your vcf object into a VRanges, and pass it to locateVariants instead. You'll still need to join via the query ID, but it's easier.

ADD COMMENT • link 10.5 years ago Michael Lawrence ★ 11k

0

Entering edit mode

Thank you that's helpful. In my case I just have a vcf with called SNPs in multiple samples. What i'd like is to annotate the coding/noncoding variants, but for each variant also have which sample it exists in.

What would be the best way to get a sample to QUERYID map? I guess I could use something like

sample_geno = geno(vcf)$GQ
sample2queryid = apply(sample_geno, 2, function(z) which(!is.na(z)))

ADD REPLY • link 10.4 years ago by0 • 0

1

Entering edit mode

I was referring to the query ID map that locateVariants() returns. Since the return value of locateVariants() is not 1-to-1, one needs to somehow join the results with the input. locateVariants() makes this easier by returning a QUERYID column that indexes into the original object. Let's say you want a logical column in the locateVariants result for each sample, indicating whether the sample had a call there. You might do something like:

lv <- locateVariants(vcf)
matchedVCF <- vcf[lv$QUERYID,]
sample_geno <- !is.na(geno(matchedVCF)$GQ) 
cbind(lv, sample_geno)

ADD REPLY • link 10.4 years ago Michael Lawrence ★ 11k