Thanks for providing additional data and plots, @llenzi,
I would remove the --p-trunc-q
parameter and see if that impacts results.
Also run dada2 with the --verbose
flag, so that the denoising summary is printed out (maybe you are already doing this but it's not in that command — this would be useful information for us to diagnose). Could you share those logs?
Your read qualities look pretty good, and those trim lengths seem reasonable, so it is a bit mysterious that dada2 is not working. For V4 (515f/806r primers, at least), 220 + 220 nt should be adequate.
V4 is short enough that just using the 250 bp forward reads should nearly cover it. That is another option to compare — just use the forward reads and discard the reverse.
It sounds like your reads are being dropped at the merge step (we would know for sure if you show us the log printed by dada2 with the --verbose
flag on), not at the denoising step. If reads were being dropped at denoising, it would help to trim further in case there are low-quality segments near the 3' end. If these are being dropped at the merge stage (the log would tell us for sure), then that is suggesting that the reads are not long enough! But that does not make sense for V4 with these trimming parameters.
One more idea: are primer sequences trimmed from these reads? This should always be done for dada2/deblur and in dada2 can be done with the --p-trim-left
parameters. For deblur, you can use q2-cutadapt to trim primers.
I agree, it sounds like joining reads and then using deblur is better for you. Note that you can join reads with qiime vsearch join-reads
, if you wanted to keep it all in QIIME2 (thus preserving that information in provenance). I do not know the technical differences between PEAR and vsearch for read joining, though.
deblur does not use Q scores. It uses a static error model to predict read errors. You will need to see the original publication for more details.
So it looks like RDP and QIIME2's sklearn classifier give similar results for all taxa except for Listeria. It is not entirely clear, because you show genus-level results for the first two plots, then species-level with RDP. But Listeria definitely gets better classification results with RDP, which is interesting because these are two very similar methods (both naive Bayes classifiers). You may want to adjust the confidence
parameter to see if this improves your results, since you have mock communities to validate your results. The default is a good average across many different sequencing lengths/conditions but since you have a mock community you can re-tune this for your sequencing runs. Notably, QIIME2's default confidence
== 0.7, whereas RDP's recommendation for short reads is 0.5, 0.8 for longer reads.
See this thread for a little more discussion — I think they are using the same mock community that you are.
Good luck!