Hi
What I want to do is calculate the minor allele frequency of a SNP/
small indel from bam files
I am not sure what is the PREFERRED way to do this
I have been running samtools and just do pileup in those regions of
interest... and read the pileup file with Rsamtools...
With scanBam on a single position I get back all the reads that span
that location. So I've tried using the start position and the cigar and
the DNAStringSet in the seq slot to get the alleles at the position I
want (works ok)
Is there a function written that already does that? seems like something
of general utility , the only complication is the cigar
Would also be nice to have a existing function with writes the output
of scanBam back into sam format (assuming necessary attributes are
available).
ALSO perhaps I am being daft:
but in scanBamFlag
isValidVendorRead: A logical(1) indicating whether invalid (FALSE),
valid (TRUE), or any (NA) read should be returned. A 'valid'
read is one flagged by the vendor as passing quality control
criteria.
so isValidVendorRead=TRUE would return VALID vendor reads?
at the moment it is returning the INVALID vendor reads
isValidVendorRead=FALSE returns the valid ones. (Solid reads)
Thanks
Paul
> sessionInfo()
R version 2.13.0 Under development (unstable) (2010-11-03 r53517)
Platform: x86_64-unknown-linux-gnu (64-bit)
locale:
[1] LC_CTYPE=en_AU.UTF-8 LC_NUMERIC=C
LC_TIME=en_AU.UTF-8 LC_COLLATE=en_AU.UTF-8
LC_MONETARY=C
[6] LC_MESSAGES=en_AU.UTF-8 LC_PAPER=en_AU.UTF-8 LC_NAME=C
LC_ADDRESS=C LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_AU.UTF-8 LC_IDENTIFICATION=C
attached base packages:
[1] stats graphics grDevices utils datasets methods
base
other attached packages:
[1] ShortRead_1.9.2 lattice_0.19-13 Biobase_2.11.2
Rsamtools_1.3.3 Biostrings_2.19.0 GenomicFeatures_1.3.5
[7] GenomicRanges_1.3.2 IRanges_1.9.6
loaded via a namespace (and not attached):
[1] biomaRt_2.7.0 BSgenome_1.19.0 DBI_0.2-5
grid_2.13.0 hwriter_1.2 RCurl_1.4-4
RSQLite_0.9-2
[8] rtracklayer_1.11.3 tcltk_2.13.0 tools_2.13.0
XML_3.2-0
>
[[alternative HTML version deleted]]
_______________________________________________
Bioc-sig-sequencing mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/bioc-sig-sequencing
