[BioC] bug report/repair in vcf2sm of GGtools

Vincent Carey stvjc at channing.harvard.edu
Mon Jun 13 12:44:58 CEST 2011


It was called to my attention that base calls denoted "." in VCF are
mapped to REF instead of missing.
This has been repaired in release and devel.  The solution is to map
the genotype to missing (raw 0) if
any base is denoted ".".

Owing to issues with dependencies and the build system, it may take
time for the devel version repair to
emerge.  For those working with source, the patch is

--- vcfutils.R  (revision 56128)
+++ vcfutils.R  (working copy)
@@ -45,9 +45,10 @@
  meta = vec[1:nmetacol]
  calls = vec[-c(1:nmetacol)]
  nalt = strsplit(calls, "")
- nums = lapply(nalt, "[", c(1,3))
- nalt = sapply(nums, function(x) 2-sum(x=="0"))
- if (any(is.na(nalt))) nalt[which(is.na(nalt))]=-1
+ nums = lapply(nalt, "[", c(1,3))  # extract the call components
+ hasmiss = which(sapply(nums, function(x) any(x == ".")))
+ nalt = sapply(nums, function(x) 2-sum(x=="0"))  # this is correct
only for diallelic locus; note in doc
+ if (length(hasmiss)>0) nalt[hasmiss] = -1
  nalt = nalt+1
  chr = meta[1]
  id = meta[3]

At this time I personally am not making extensive use of VCF so
observations of shortcomings of this
functionality are welcome.  At present the function makes no use of
call uncertainty information, and it is
hoped this will be remedied before too long.



More information about the Bioconductor mailing list