I have a fastq file from Nanopore sequencing data that contains reads from both the UK and South Africa variants of SARS-CoV-2.

The variants are identified by three key mutations in the receptor-binding domain of Spike.

Aim Firstly, I want to determine the number of reads then belongs to UK and South Africa variant. Then, I will create a bam file from the fastq file using the SARS-CoV-2 reference genome.


  1. How do I determine the reads belonging to each variants?
  2. Are there any tools available to classify the reads based on the pre-defined mutation criteria?
  • $\begingroup$ I don't think that's possible in general. You can obviously map the reads to the references and split but that's if you trust the references and in that case why do you even need the reads? Just thinking out loud. $\endgroup$
    – juanjo75es
    Commented Mar 2, 2022 at 0:39

1 Answer 1


Based on a quick request from someone else, I wrote some code a couple of years ago for assigning read groups based on variants at particular locations; this was set up for HIV variant splitting, so I expect it should also work on SARS-CoV-2 nanopore reads as well:


It processes a SAM stream [i.e. in text mode] using Perl, so might be a bit slow, but should be fine for a few hundred thousand sequences. You'd use it as follows:

samtools view -h mapped_reads.bam | \
  samVarSplitter.pl [-ref <refName>] [-pos <int>] | \
  samtools sort > grouped_reads.bam

Then the reads can be split into different files based on their group using samtools split:

Usage: samtools split [-u <unaccounted.bam>] [-h <unaccounted_header.sam>]
                      [-f <format_string>] [-v] <merged.bam>
  • 1
    $\begingroup$ Fantastic, good response. $\endgroup$
    – M__
    Commented Feb 25, 2022 at 21:34

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.