Why are my reads filtered out (98%) when quality seems okay?

User Support
,
1

I am using qiime2 2026.7 on a WSL system in Ubuntu and working on paired-end 18S seq data (Illumina Miseq 2x300) from a mixed seawater community. My data is already imported and demultiplexed, but at the denoising step (using dada2), around 98% of reads are filtered out for some reason. How can I figure out the reason for this? In the quality plot, the quality of sequences seem normal to me, and I have already removed primers. I also tried increasing the error rate allowed to 10 with no difference seen.

Here is my command and output (I previously tried the same with no specification of p-max-ee-f and r):

qiime dada2 denoise-paired --i-demultiplexed-seqs demux_hvalfjordur.qza --p-trim-left-f 20 --p-trunc-len-f 300 --p-trim-left-r 18 --p-trunc-len-r 300 --o-representative-sequences hvalfjordur_asv-seqs.qza --o-table hvalfjordur_asv-table.qza --o-denoising-stats hvalfjordur_denoising-stats.qza --o-base-transition-stats hvalfjordur_base-transition-stats.qza --verbose --p-max-ee-f 10 --p-max-ee-r 10
Running external command line application(s). This may print messages to stdout and/or stderr.
The command(s) being run are below. These commands cannot be manually re-run as they will depend on temporary files that no longer exist.

Command: run_dada.R --input_directory /tmp/tmpxoechzqg/forward --input_directory_reverse /tmp/tmpxoechzqg/reverse --output_path /tmp/tmpxoechzqg/output.tsv.biom --output_track /tmp/tmpxoechzqg/track.tsv --output_err_track /tmp/tmpxoechzqg/err_track.tsv --filtered_directory /tmp/tmpxoechzqg/filt_f --filtered_directory_reverse /tmp/tmpxoechzqg/filt_r --truncation_length 300 --truncation_length_reverse 300 --trim_left 20 --trim_left_reverse 18 --max_expected_errors 10 --max_expected_errors_reverse 10 --truncation_quality_score 2 --min_overlap 12 --max_merge_mismatch 0 --trim_overhang False --pooling_method independent --chimera_method consensus --min_parental_fold 1.0 --allow_one_off False --num_threads 1 --learn_min_reads 1000000 --retain_unmerged False

R version 4.5.2 (2025-10-31)
Loading required package: Rcpp
DADA2: 1.38.0 / Rcpp: 1.1.2 / RcppParallel: 5.1.11.2
2) Filtering The filter removed all reads: /tmp/tmpxoechzqg/filt_f/sample86_L001_R1_001.fastq.gz and /tmp/tmpxoechzqg/filt_r/sample86_L001_R2_001.fastq.gz not written.
Some input samples had no reads pass the filter.
...............................................................................................x..............
3) Learning Error Rates
29801520 total bases in 106434 reads from 109 samples will be used for learning the error rates.
30014388 total bases in 106434 reads from 109 samples will be used for learning the error rates.
3) Denoise samples .............................................................................................................
.............................................................................................................
5) Remove chimeras (method = consensus)
6) Report read numbers through the pipeline
7) Write output
Saved FeatureTable[Frequency] to: hvalfjordur_asv-table.qza
Saved FeatureData[Sequence] to: hvalfjordur_asv-seqs.qza
Saved SampleData[DADA2Stats] to: hvalfjordur_denoising-stats.qza
Saved DADA2BaseTransitionStats to: hvalfjordur_base-transition-stats.qza

I am also attaching my qzv files for the quality of reads and denoising stats.

demux_hvalfjordur.qzv (333.5 KB)

hvalfjordur_denoising-stats.qzv (312.5 KB)

I considered using just the forward reads instead, but would really like to understand why this is happening before I pivot.

Hope you can help!

Best,

Line

2

Hi! Apparently, there are two likely reasons why your reads are getting filtered out. The first being from quality filtering, and the second from chimera removal. Did you check the stats to see if there were a lot of chimeras detected? If not, it's due to the quality filtering.

Having a lot of overlap will cause the merging to fail just like small reads that get over-truncated. Soilrotifer recommends in this discussion to adjust your truncation or adjust the parameter --p-min-fold-parent-over-abundance to 8 or 16.

3

Yes, I did already check the stats. As you can see below, the main reason for removal of reads is the filtering step (98% removed), not merging or chimera removal (some reads are also being removed because of this, but not as many). So it is the quality filtering.

image
image1795×723 59.3 KB

It is not entirely clear to me which criteria the filtering step applies and which of them my reads are failing, especially since it seems the overall quality of my reads seem fine in terms of Phred-scores.

4

Hi @lineroager,

Some thoughts...

The quality plots you've share look typical of many. The provenance also shows no explicit primer removal steps were ran, i.e. cutadapt. Only Importing and DADA2. Thus, it appears that you are using DADA2's trim functions to remove the primers rather than cutadapt, which is valid.

Anyway, I think the main issue is setting the --p-trunc-len-f 300 --p-trunc-len-r 300. Given the quality plot, this will not likely work. For the forward read, you might want to truncate at ~280 bases, and the reverse read you might want to truncate at ~230 bases.

What 18S rRNA gene primers are you using? What length are you expecting?

-Mike

5

Hi Mike,

Yes, I am removing the primers by trimming in DADA2.

I tried now truncating at 280 and 230 like you suggested and this increased the number of reads kept after filtering to 25-35%. However, I am still losing a lot, and now also a bunch to merging and chimera removal steps.

image

I am using V4 primers that should theoretically give an amplicon of 415 bp, but I have seen larger bands on my gel when running the PCR products.

Any ideas on what else I could try to retain more data? Or would it be better to work with forward reads only at this point?

6

Thanks @lineroager,

We'll that is quite an improvement!

For the chimera removal I'd definitely do:

Okay this is great, as you have lots of room to truncate and still merge reads.

We'll try and calculate some values:

Basic read length:
We do not really need to take into account the removal of the primers to calculate overlap, but I am doing this just to get an overall view of the actual read lengths we're working with. Otherwise just calculate using the 300 bases for each read throughout the calculations.
FW read - FW primer = 280 bases
RV read - RV primer = 282 bases

Calculate read overlap for merging
We'll subtract a 20 base overlap, then subtract the estimated length of your expected amplicon size:
280 + 282 - 20 - 415 = 147 bases of overlap.

In a nutshell 2x300 might be overkill for this... However, it gives you a lot of breathing room to liberally truncate the reads. This is important, as you are trying to help the denoiser, by removing the low-quality portions at the end of your reads. This is also, why there are still many failing the quality control filter early on.

If I've done my calculation correctly, then you should not trim more than 147 bases (total) from both reads combined (say 47 off the forward read and 100 off the reverse read). Given this much potential overlap, you're more likely to encounter more mismatches between the two reads, causing merge fails.

I'd think I'd start here:
--p-trunc-len-f 241 --p-trunc-len-r 199.

^ this should give you 241 + 199 - 20 = ~420 base amplicon after merging. assuming you are using the DADA2 parameters. If using cutadapt to trim primers, you'd have to remake the quality plot post cutadapt then you can pick your truncation parameters from there

Then add or subtract to the truncation parameters as needed to optimize. It might take a couple of iterations.

I am sure someone else will point out any errors in my calculations, or provide an easier way to go about it. :laughing: If you search the forum for "overlap calculation" you'll find several better explanations.

-Mike

7

Thanks so much for these recommendations! I was able to increase the retained reads after filtering to ~50%, merged reads ~40%, and non-chimeral reads ~35-40%. I did this through --p-trunc-len-f 260 --p-trunc-len-r 200 --p-min-fold-parent-over-abundance 16

While this is not amazing, it is much better than the starting point, and I can proceed with this. Thanks again!