[Help] DADA2 error rates could not be estimated for iSeq data — how to fix in R?

Other Bioinformatics Tools
1

Hello all,

I’m running DADA2 on Illumina iSeq paired-end data outside of QIIME 2 (in R), because the plugin version keeps failing at the error-learning step.

Here’s what happens in R:

errF <- learnErrors(filtFs, multithread = TRUE)
105659550 total bases in 704397 reads from 18 samples will be used for learning the error rates.
Error rates could not be estimated (this is usually because of very few reads).
Error in getErrors(err, enforce = TRUE) : Error matrix is NULL.

Even though there are plenty of reads, the error model cannot be estimated and errF is never created. I understand from GitHub issue #791 that this problem is caused by binned quality scores on iSeq / NovaSeq platforms (e.g., Q = 2, 11, 25, 37), which break the loess smoothing used in learnErrors().

:puzzle_piece: What I’ve already tried

  1. Running learnErrors() with different settings:

    errF <- learnErrors(filtFs, nbases = 1e7, randomize = TRUE, multithread = TRUE)

    → still returns the same error.

  2. Only using forward reads (filtFs) — sometimes works, sometimes not.

  3. I’ve read Ben Callahan’s comment suggesting a workaround:

    “The workaround right now is to enforce monotonicity in the fitted error model.
    The error model (getErrors(errF)) is just a matrix, so you can assign the value at Q = 40 to all lower entries…”

However, since learnErrors() fails completely, there is no errF object to modify.

:red_question_mark:My Question

How can I manually or programmatically create or repair the error model in R so that DADA2 can proceed, following the suggestion in issue #791?

In other words:

  • How can I modify or initialize the errF matrix before or after learnErrors() to enforce monotonicity?

  • Is there an R-based workaround (e.g., getErrors() or manually setting err_out) that will allow the pipeline to continue even if learnErrors() fails due to binned Q-scores?

Any working example of how to implement this fix in R would be deeply appreciated.

:dna: Context

  • Platform: Illumina iSeq 100 (2×150 bp)

  • Target: 16S rRNA V3–V4

  • DADA2 version: 1.30.0 (R 4.4.0 / Bioconductor 3.20)

  • Filtering summary: ~7 × 10⁵ reads, 18 samples

:folded_hands: Goal

I just want to make the DADA2 workflow run entirely within R and produce valid error models for iSeq data — or any reproducible script that implements the monotonic “error-rate fix” suggested in issue #791.

Thank you so much in advance!
— Shuoyun Jiang

2

This issue has already been opened on GitHub:

Looks like Ben is helping you on GitHub, so I'll unqueue this thread for now.

(There is no need to post twice. The same researchers and developers watch the GitHub repos and the Qiime2 forums, so your message will be seen either way.)

1 Like