1
$\begingroup$

I want to find and compare the results from STAR and BWA- MEM mapped. I have 150bp paired end reads in sorted.bam files in each case and i want to find in each case uniquely mapped reads, number of reads mapped to multiple loci and number of unmappped reads.

Especially for uniquely mapped reads i tried different solutions from biostars and this Obtaining uniquely mapped reads from BWA mem alignment but i dont know what to trust,

Could you please suggest the best way to find mapping results from bam paired end files ?

$\endgroup$
3
  • 1
    $\begingroup$ Please clarify your specific problem or provide additional details to highlight exactly what you need. As it's currently written, it's hard to tell exactly what you're asking. $\endgroup$
    – Community Bot
    Commented Feb 14, 2022 at 14:22
  • $\begingroup$ It would be helpful to describe the genome and a bit of project background regarding what you have already tried. Biostars is a different site to Stackexchange. $\endgroup$
    – M__
    Commented Feb 14, 2022 at 14:28
  • $\begingroup$ I have paired end mapped reads from human in sorted bam files and i want to find the unique mapped , the multimapped and the unmapped reads from BWA and from STAR. I already have the .bam files from each mapper. i Tried samtools flagstat and the solution above but i am not sure if this is correct especially for unique reads $\endgroup$ Commented Feb 14, 2022 at 15:08

1 Answer 1

2
$\begingroup$

Bam files contain a MAPQ field that informs on whether an alignment is unique or not. With STAR, a MAPQ value of 255 means the read is mapped uniquely and a score lower than 255 means not unique.

You can use samtools on a bam file to count or output the uniquely or multimapped reads by filtering for the MAPQ as described here: https://groups.google.com/forum/#!topic/rna-star/ipHc9FrgaPU

For example, count uniquely mapped reads: samtools view -c -q 255 myBamFile.bam
Output uniquely mapped reads: samtools view -q 255 myBamFile.bam

There are differences between STAR and BWA MAPQ output though. see here for details: https://sequencing.qcfail.com/articles/mapq-values-are-really-useful-but-their-implementation-is-a-mess/

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.