-
Notifications
You must be signed in to change notification settings - Fork 2
How to...
lomereiter edited this page Dec 17, 2012
·
2 revisions
This page presents how common operations on BAM file are done with BioD.
For more information, see module/function/method documentation in the source code.
See also examples.
- Get all reads in BAM file
import bio.bam.reader;
...
auto bam = new BamReader(filename);
foreach (read; bam.reads)
{
...
}- Get only reads overlapping a region (index must be available)
import bio.bam.reader;
...
auto bam = new BamReader(filename);
foreach (read; bam["chr17"][100000 .. 200000])
{
...
}- Make a pileup and iterate columns
import bio.bam.reader;
import bio.bam.pileup;
...
auto bam = new BamReader(filename);
auto reads = bam["chr11"][500000 .. 600000];
foreach (column; makePileup(reads))
{
...
}- Recover reference bases from read sequence, CIGAR, and MD tag (they start from the first mapped base of the read)
import bio.bam.md.reconstruct;
import std.conv; // optional, to convert the result to a string
...
auto reference_bases = to!string(dna(read));- Recover reference bases for a whole bunch of reads (they must be sorted by coordinate!)
import bio.bam.md.reconstruct;
import std.conv;
...
auto reference_bases = to!string(dna(bam["chr12"][850000 .. 1850000]));-
Get tag value
- Preferred way. Throws exception in case of error, does type conversions automatically.
import std.conv; auto value = to!ushort(read["ZF"]); // will check if the value is in 0..65535 range
- Faster but unsafe way. Exact type of tag must be known.
auto v = read["ZF"]; auto value = *cast(int*)(&v);