Rasmus Kirkegaard 22 December, 2023
With the release of R10.4.1 I wanted to check the quality of the raw reads and the consensus accuracy to ensure that we had the right approach for continuing our nanopore only assemblies from metagenomic samples (Mantas et al. 2022). So we sequenced the Zymo HMW mock DNA to evaluate the quality of the raw reads but more importantly the assembled genomes. With the introduction of dorado as a basecaller that should replace guppy it was also interesting to see how well that performed.
I think it is work taking a moment to think about how amazing it is that we now have several independent methods (sequencing tech+assembler) that can reproduce a bacterial genome of e.g. P aeruginosa (6832199 bp) with error rates of less than 0.01 per 100 kbp \<2 total errors.
Ending 2023 with nanopore data improving from their single strand data (simplex) barely hitting their Q20 goal to going well beyond it and when reading both stands (duplex) even racing beyond Q30 is quite amazing. Lots of innovation on run conditions and basecalling models. The promise of Q30 single strand reads in 2024 does not seem too far fetched.
The data has been added to the NCBI-SRA PRJNA934154. I have managed to upload the fastq and the fast5 files. The fastq and fast5 files should be available through the “cloud delivery service”.
Fastq data (fast,hac & Sup): SRR23563655
Fast5 data: SRR23437037
Fastq data: - fast: SRR24893246 - hac: SRR24893245 - sup: SRR24893244
Pod5 data: The data has been added to the ENA PRJEB64570. (NCBI did not accept pod5 and told me they do not want to do that)
The data has been added to the ENA PRJEB65462
Fastq data: - sup: ERR11901474 - sup duplex: ERR11901475
Pod5 data: - ERR11924124
Phred scores for perfect matching reads are calculated as recommended by Armin Topfer which takes length into account.
With the launch of 5 khz sampling rate around London Calling 2023 ONT was hoping that the GPU power needed for basecalling could be decreased dramatically as this should allow for the use of the faster HAC model to replace the compute intensive SUP model at 4khz. To test this we here compare the consensus assemblies with the 4 khz SUP model and the new 5khz HAC and SUP models.
The indel rate seems to be higher with 5khz HAC than both SUP regardless of sample rate. Interestingly the 5khz sample rate with SUP performs much better than 4 khz SUP for some organisms but not S enterica.
While HAC is on par with SUP for some organisms it is never the best option for mismatches.
The 4.3 model includes some modifications found in bacteria and it seems to pay off. The model gives a clear improvement for the consensus accuracy across many of the genomes. Indel rates without any polishing ranging from below 0.01 to 1 per 100 kbp translates to 1-100 errors for a 10 Mbp genome. Similar levels were achieved regarding mismatches. The Bacillus showed notably poorer consensus accuracy than the others which I assume could be due to some real differences or errors in the Pacbio HiFi based references.When we are getting in the 0.01-0.1 errors per 100 kbp I assume that it gets difficult to tell whether they are actual errors without hand curating reference genomes with multiple independent technologies as per Ryan Wicks “Perfect bacterial genome tutorial”.
Here is a brief description of the tools used. For the exact commands check out the Snakefile in this repository (Snakemake v. 7.18.2).
DNA sample was the Zymo Mock HMW standard. The DNA was prepared for sequencing using the nanopore ligation sequencing kit (SQK-LSK114) and sequenced on a R10.4.1 nanopore promethion flowcell (FLO-PRO114M) with the “400 bp/s” mode (4khz sampling). The DNA was prepared for sequencing using the nanopore ligation sequencing kit (SQK-LSK114) and sequenced on a R10.4.1 nanopore MinION flowcell (FLO-MIN114) with the “400 bp/s” mode (5khz sampling). The DNA was prepared for sequencing using the nanopore ligation sequencing kit (SQK-LSK114) and sequenced on a R10.4.1 nanopore PromethION flowcell (FLO-PRO114HD) with the “400 bp/s” mode (5khz sampling).
The reads were basecalled using dorado (v. 0.1.1) with fast, hac and sup accuracy mode using the 4.0.0 and 4.1.0 models.
The reads were basecalled using dorado (v. 0.3.0) with fast, hac and sup accuracy mode using the 4.2.0 models.
The reads were basecalled using dorado (v. 0.3.4) with fast, hac and sup accuracy mode using the 4.2.0 models. The reads were basecalled using dorado (v. 0.5.0) with fast, hac and sup accuracy mode using the 4.3.0 models.
Reads were mapped to the updated zymo reference genomes (hopefully goes public soon) using minimap2 (v. 2.24), and QC information was obtained using NanoPlot (v. 1.41.0).
The reads were subsampled using seqtk (v. 1.3) and assembled using flye (v. 2.9.1). The metagenome assemblies were then polished using medaka (v. 1.11.3).
The assembled contigs were compared to the reference contigs using QUAST (v. 5.2.0) and fastANI (v. 1.33).