I'm trying to investigate whether there are errors within my BAM. After looking at the BAM header to see whether the read groups exists (using samtools, i.e.

samtools view -H file1.bam


Here, it appears that the header includes @RG tags. However, I suspect that the BAM contains reads that either lack the @RG tag or perhaps these read groups are not well-defined, i.e. it looks like some reads in the BAM may be missing an RG tag.

How do I explicitly check whether all reads have an RG tag, and whether these are well-defined?

According to the man page, running samtools stats --split RG <file1.bam> should produce summary statistics separated by read group. If it doesn't produce a list of ungrouped reads, counts/statistics can be compared to running without the --split RG argument.

Here's some example output from a BAM file with combined read groups:

samtools stats --split RG both.bam | grep '^SN' SN raw total sequences: 2392 SN filtered sequences: 0 SN sequences: 2392 SN is sorted: 1 SN 1st fragments: 2392 SN last fragments: 0 SN reads mapped: 2341 SN reads mapped and paired: 0 # paired-end technology bit set + both mates mapped SN reads unmapped: 51 SN reads properly paired: 0 # proper-pair bit set SN reads paired: 0 # paired-end technology bit set SN reads duplicated: 0 # PCR or optical duplicate bit set SN reads MQ0: 0 # mapped and MQ=0 SN reads QC failed: 0 SN non-primary alignments: 0 SN total length: 2292231 # ignores clipping SN bases mapped: 2255361 # ignores clipping SN bases mapped (cigar): 1538560 # more accurate SN bases trimmed: 0 SN bases duplicated: 0 SN mismatches: 254720 # from NM fields SN error rate: 1.655574e-01 # mismatches / bases mapped (cigar) SN average length: 958 SN maximum length: 9078 SN average quality: 20.0 SN insert size average: 0.0 SN insert size standard deviation: 0.0 SN inward oriented pairs: 0 SN outward oriented pairs: 0 SN pairs with other orientation: 0 SN pairs on different chromosomes: 0  The --split command has also produced additional .bamstat files, which can be looked at for the summary statistics for individual read groups. In this case I'm only interested in the "raw total sequences" line:  grep 'raw total sequences' *.bamstat
both.bam_minusOnly.bamstat:SN   raw total sequences:    1901
both.bam_plusOnly.bamstat:SN    raw total sequences:    491


There are 1901 sequences from the "minusOnly" group, and 491 sequences from the "plusOnly" group, and I know from the previous statistics that there are 2392 total sequences. 1901 + 491 = 2392, so I know there are no ungrouped reads in this file.

• I'm not sure how one extrapolates from this whether each of the reads have a RG tag. Could you extrapolate a bit how I would deduce from these statistics? Oct 1 '17 at 22:29
• I've added additional information to demonstrate how this can be done.
– gringer
Oct 1 '17 at 23:39

using samjdk and invoking the function getReadGroup()

getReadGroup() returns The SAMReadGroupRecord from the SAMFileHeader for this SAMRecord, or null if 1) this record has no RG tag, or 2) the header doesn't contain the read group with the given ID.or 3) this record has no SAMFileHeader

java -jar dist/samjdk.jar -e 'return record.getReadGroup()==null;'  input.bam


the output will be a sam/bam file with missing/bad RG .

A complete solution to validate a SAM/BAM file would be picard's ValidateSamFile. This will find many other issues with BAM files.