BWA-mem and sambamba read group line error

Question

this question has been asked [and answered] on Stack Overflow

This is a two-part question:

1. help interpreting an error;
2. help with coding.

---

1. I'm trying to run bwa-mem and sambamba to align raw reads to a reference genome and to sort by position. These are the commands I'm using:

bwa mem  \
  -K 100000000 -v 3 -t 6 -Y \
  -R '\@S200031047L1C001R002\S*[1-2]' \
  /path/to/reference/GCF_009858895.2_ASM985889v3_genomic.fna \
  /path/to/raw-fastq/'S\S[^_]*_L01_[0-9]+-[0-9]+'_1.fq.gz \
  /path/to/raw-fastq/'S\S[^_]*_L01_[0-9]+-[0-9]+'_2.fq.gz | \
/path/to/genomics/sambamba-0.8.2 view -S -f bam \
  /dev/stdin | \
/path/to/genomics/sambamba-0.8.2 sort  \
  /dev/stdin \
  --out host_removal/${SAMPLE}/${SAMPLE}.hybrid.sorted.bam

This is the error message I'm getting: [E::bwa_set_rg] the read group line is not started with @RG.

My sequences were generated with an MGI sequencer and the readgroups are identified like this @S200031047L1C001R0020000243/1, i.e., they don't begin with an @RG. How can I specify to sambamba that my readgroups start with @S and not @RG?

---

2. The commands written above are a published pipeline I'm modifying for my own research. However, among several changes, I'm not confident on how to define sample id as such stated in the last line of the code: --out host_removal/${SAMPLE}/${SAMPLE}.hybrid.sorted.bam (I'm referring to ${SAMPLE}). Any insights?

---

Thank you very much!

Answer 1

Read groups must start with RG, that is a hardcoded requirement of the SAM format, nothing you can do about that.

See https://samtools.github.io/hts-specs/SAMv1.pdf

Answer 2

My sequences were generated with an MGI sequencer and the readgroups are identified like this @S200031047L1C001R0020000243/1

Those are MGI read names, not read groups. (And the @ is the FASTQ syntax that introduces a read name in a FASTQ file, not part of the read name or qname itself.)

Read groups more or less correspond to the sample that was sequenced. If a sample was sequenced across multiple sequencing machine runs and the resulting data merged into one file, there might be several read groups corresponding to the individual runs all pertaining to the same sample — though with increasing data volumes from individual sequencing runs this is perhaps less common than it used to be.

You get to make up the ids yourself, usually corresponding to what the sample is known as in your workflows. So your -R option should be something like

-R '@RG\tID:S123\tSM:S123\tDS:Sample 123\tPL:DNBSEQ'

(BWA's -R syntax is unfortunately not very easy to use. See addreplacerg for a similar -r option that is rather more user-friendly.)

The name you use for ${SAMPLE} would usually be the same as the value of the SM: field within that -R option — the identifier that is used in your LIMS and workflows to identify that particular sample.