We're working on a nasopharyngeal 16S rRNA study (V3-V4) comparing a pre-pandemic cohort with a post-pandemic cohort in healthy individuals and individuals with invasive pneumococcal disease. Our problem is a methodological/design issue we'd like a second opinion on.
All pre-pandemic samples were sequenced on MiSeq, and all post-pandemic samples on NextSeq 2000. Sequencing platform is therefore completely confounded with study period. We're aware this means no batch-correction method can statistically separate the platform effect from the pandemic effect, since they're perfectly collinear.
Here are some study details:
16S V3-V4, 2x300 bp on both platforms (so amplicon coverage should be comparable).
92 pre-pandemic samples across 4 MiSeq runs and 107 post-pandemic samples in 1 NextSeq 2000 run.
Groups (healthy/disease) are distributed across the 4 pre-pandemic runs.
What we're already planning is to run DADA2 with per-run error learning (the 4 MiSeq runs + the NextSeq run separately), using the same truncation length across runs so the ASVs can be merged. We're also considering re-sequencing a subset of pre-pandemic samples on the NextSeq 2000 as bridge samples, to directly measure the platform effect.
Our questions are:
Given the full confounding, is a pre/post comparison interpretable at all? How have others framed this limitation?
For the bridge samples, how many pairs would you consider a defensible minimum to characterize the platform effect?
Are there additional harmonization steps you'd recommend when merging MiSeq and NextSeq 2000 amplicon data?
What's the current recommended way to handle NextSeq binned quality scores within a QIIME2 workflow, especially when mixing binned (NextSeq) and non-binned (MiSeq) runs in the same study?
Any published examples of MiSeq vs NextSeq 2000 comparisons for 16S data (or any Illumina cross-platform 16S comparison) that we could reference?
We haven't found much on this, so this could be a really interesting one to discuss. Looking forward to your thoughts. Thanks so much in advance!
Hi @ggc ,
As sequencing run is confounded with study period, I don't think that you can meaningfully compare the study periods. As seqeuncing platform is also a confounder this will only make matters worse... basically, you will not be able to determine whether differentially abundant features are linked to batch effects, platform, study period...
you could mitigate the issue by using (i) negative controls/internal standards for contaminant detection, e.g., with decontam (in the q2-quality-control plugin), or (ii) positive controls/mock communities/replicates in each run to detect batch effects. The latter would be the best option to see if there is any difference at all between runs/platforms where the same identical samples were sequenced.
perfect, but it's not enough to see the differences between platforms, you will need these sequenced on the same run as the others if you specifically want to detect batch effects as opposed to platform effects.
Platform effects are certainly possible (e.g., as NextSeq uses the newer patterned flowcell and is more prone to, e.g., index hopping errors) but in general my experience is that any such effects are quite small... mock community samples sequenced on MiSeq and NextSeq seem quite comparable. But batch effects are a well recognized issue, and as you are comparing samples from different cohorts collected at totally different times this issue will be compounded because you will have possible batch effects in collection, DNA extraction, amplification, library prep, sequencing... as well as cohort differences that are unrelated to disease. So I suggest carefully checking for such effects wherever possible. Personally, I would just resequence everything on a single run to eliminate at least one covariate.
One clarification regarding the cohorts themselves: post-pandemic individuals were recruited using the same criteria as the pre-pandemic ones, and samples were processed following identical protocols (sample collection, DNA extraction...). So sequencing platform and batch effects are really the only major variables that differ structurally between the two cohorts.
Regarding batch effect monitoring, each of our runs already includes negative extraction and amplification controls, so decontam is already part of our planned pipeline. I hadn’t mentioned this previously, apologies.
Unfortunately, re-sequencing everything in a single run isn't feasible for us. We don't have enough remaining pre-pandemic sample material to redo it all, and funding wouldn't cover a full re-sequencing effort either. So we're trying to build the most defensible partial solution rather than eliminate the confound entirely.
With this additional information, if you have any further recommendations, we would be very grateful.
Dada2 reports base-transition-stats.qza and there's a new dada2 plot-base-transitions action. This shows you how different the error profiles are between the two runs, allowing you are argue that the error profiles is similar even though the sequencing chemistry is different. Hopefully!
There are other sources of batch effect beyond Illumina error profile, but it's something you can run right now on the data you have!
Somehow, these are never included on the runs where we need them most.
