Question: Using samtools to assess the number of reads that cover a variant position, as well as the ratio of ref/alt among those reads
gravatar for claudiadast
19 months ago by
claudiadast0 wrote:

I currently have variant data per sample and the corresponding RNA seq BAMs per sample. What I'm trying to do is tap into the corresponding bam files, determine the number of reads associated with the variant, and then deduce the %ref and %alt alleles in those reads.

So for instance, I have a variant of interest at position chr2:86042011, where the ref allele is CCAGA and the alt allele is C, and let's say there are 10 reads that cover that variant position - I want to know how many of those reads have the ref allele and how many have the alt allele.

I'm currently trying to use samtools mpileup -I [bed file] [bam]. So taking the above example again, I feed a BED file that looks like:

chr2      86042010              86042012

The output I then get from running that command is:

chr2    86042011        N       2       CC      sJ
chr2    86042012        N       2       CC      sJ

One of my confusions is how you know how many reads cover your variant. For instance, should my BED file intervals be greater? Am I doing this right?

ADD COMMENTlink modified 19 months ago • written 19 months ago by claudiadast0

You can use samtools depth to get the number of reads covering a position.

    samtools depth [options] [in1.sam|in1.bam|in1.cram [in2.sam|in2.bam|in2.cram] [...]]

    Computes the depth at each position or region.
ADD REPLYlink modified 19 months ago • written 19 months ago by genomax83k

Thanks! So going off of my example above, I did samtools depth -r chr2:86042011-86042011 myBAM.bam and got:

chr2    86042011        3

So how does that read depth number then get factored in to the samtools mpileup command? Does that mean the interval has to be 86042011- 86042014 or how do I know the corresponding position of the 3 reads covering the variant?

Also, why is there a discrepancy between what mpileup (2) gives me versus depth (3) in the number of reads?

ADD REPLYlink modified 19 months ago • written 19 months ago by claudiadast0

See: Difference between samtools DEPTH and MPILEUP

ADD REPLYlink written 19 months ago by genomax83k
Please log in to add an answer.


Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1460 users visited in the last hour