Hi, I am on qiime2-2025.7, and I would like a little more clarification on how the rarefaction process works with low-frequency samples in boots. I am planning to run the following command:
qiime boots core-metrics \
--i-table table-dada2-filtered.qza \
--i-phylogeny rooted-tree.qza \
--m-metadata-file meta.tsv \
--p-sampling-depth 5000 \
--p-n 100 \
--p-no-replacement \
--p-alpha-average-method median \
--p-beta-average-method medoid \
--output-dir boots-core-metrics \
However, my smallest sample has 5000 reads, so I'm not entirely sure what will happen to it after the first iteration of 5000 subsampling if --p-no-replacement is used. Is the sample filtered out from all subsequent iterations? I would imagine that could do some weird things to my average diversity metrics. I'm hoping to keep this sample involved in my analysis because it's part of a rarer treatment group.
Should I just use bootstrapping with --p-replacement instead?
Thank you in advance!
Hi @Emma_Gardner,
my smallest sample has 5000 reads, so I'm not entirely sure what will happen to it after the first iteration of 5000 subsampling if --p-no-replacement is used. Is the sample filtered out from all subsequent iterations?
Yes, it would be filtered from the feature table in each iteration, and would not be presented in the outputs (alpha diversity, beta diversity). It won't do anything weird to your average diversity metrics - it just won't be represented in them. (Is there some specific weirdness you're worried about - if so, let me know in case I'm misunderstanding the question.)
Should I just use bootstrapping with --p-replacement instead?
We don't have a lot of experience with using this for this purpose yet. I added it to boots for the purpose of experimentation, and in general we see very high correlation across the two approaches (with/without replacement). But, as the sampling depth gets lower, we need to start thinking of the low sampling depth as indicative of a failure of the sequencing of that sample, rather than simply lower coverage. That, in my mind, could lead to weirdness in the average diversity metrics. For example, think of an extreme case where we observed one sequence for a sample: sampling with replacement would always give the same feature table and alpha diversity, and probably make that sample an extreme outlier in beta diversity. So if you're going to use the with replacement option, I would do that after filtering samples with low sequence counts from your feature table. It might also be worth running the same analysis without replacement, and confirming that you don't see any major differences.
Do others have thoughts/experiences on with/without replacement rarefaction?
I am on qiime2-2025.7
If possible, I'd upgrade to 2026.7. We've made some big improvements to boots since 2025. The results won't change in ways that will impact your interpretation, so this isn't essential, but a few things have changed:
- we now default to
medoid for beta diversity averaging - we like this approach better, but prior to 2026.1 the implementation we were using was too memory intensive to use on real data
- and the random seed for rarefaction is now stored in data provenance, so the analysis itself will be more fully reproducible.
Hi @gregcaporaso! I appreciate you for taking the time to write out such a thorough reply.
I have a lot of follow-up questions, but I'll limit myself to two(ish 
) for now.
I'm just interested in how the code/math works. It probably doesn't matter in the grand scheme of things, but I don't want to report a significance value "based on 50 samples" if only 25 samples were actually included in the analysis based on the command parameters I used.
If I do subsampling at p-5,000 without replacement for 100 iterations, I imagine the samples included at each iteration would look something like this:
| Sample |
# of Reads |
Included in iteration 1? |
Included in iteration 2? |
Included in iteration 3? |
.... Included in iteration 100? |
| A |
5,000 |
yes |
no |
no |
no |
| B |
10,000 |
yes |
yes |
no |
no |
| C |
15,000 |
yes |
yes |
yes |
no |
| ...M |
499,000 |
yes |
yes |
yes |
no |
| ...Z |
3,000,000 |
yes |
yes |
yes |
yes |
From this result, how is the median/medioid found? Does qiime/boots completely exclude samples A through M if they don't make it through all 100 iterations, or is the median/medioid compiled from the results of each rarefaction?
I would love to, but my IT department requires me to go through a rather superfluous process to update my distribution on the HPC, so I'd rather not. I'd also need to update conda since I've seen some people were having issues with compatibility, so that would mean another collaboration with IT.
If I do go through the trouble of getting permission for an updated deployment, will I need to go back and re-do my whole pipeline (importing, denoising, etc.)? If I were to not re-do the pipeline, and I publish my methods with "QIIME2 2025.7 was used for some things" and "QIIME 2026.10 was used for other things", would you find that weird as a reviewer?
Again, I appreciate your time. I'll continue looking at literature/source code to see if I can find more answers.
Ah - I understand now.
The with/without replacement happens in the context of a single rarefying iteration, not across the iterations. So, as long as the total frequency (i.e., number of reads) in the input feature table for sample A is greater than or equal to the value provided for the sampling-depth parameter, it will be represented in every iteration's rarefied feature table (and the corresponding diversity calculations). In each iteration, what will differ is whether the sampling-depth reads that are selected to represent that sample are sampled with or without replacement from the input feature table.
Does that clarify?
Makes sense - in that case, you should be good to go with the version that you're using.