zaeleus / noodles

Bioinformatics I/O libraries in Rust
MIT License
482 stars 52 forks source link

How to extract pileup information from a bam file? #205

Closed huangnengCSU closed 11 months ago

huangnengCSU commented 11 months ago

Hi, I am attempting to replace rust-htslib with noodles in my code. But I am unsure of how to extract the pileup information from a bam file. In rust-htslib, the method is as follows:

use rust_htslib::{bam, bam::Read};

let mut bam = bam::Reader::from_path(&"test/test.bam").unwrap();

// pileup over all covered sites
for p in bam.pileup() {
    let pileup = p.unwrap();
    println!("{}:{} depth {}", pileup.tid(), pileup.pos(), pileup.depth());

    for alignment in pileup.alignments() {
        if !alignment.is_del() && !alignment.is_refskip() {
            println!("Base {}", alignment.record().seq()[alignment.qpos().unwrap()]);
        }
        // mark indel start
        match alignment.indel() {
            bam::pileup::Indel::Ins(len) => println!("Insertion of length {} between this and next position.", len),
            bam::pileup::Indel::Del(len) => println!("Deletion of length {} between this and next position.", len),
            bam::pileup::Indel::None => ()
        }
    }
}

best, Neng

zaeleus commented 11 months ago

While noodles does have a simple pileup iterator to calculate depth at each position (see noodles-bam/examples/bam_depth.rs), it does not provide a more sophisticated engine like htslib.

You may be interested in tools like perbase.

huangnengCSU commented 11 months ago

Thank you!