I'm struggling with finding a proper statistical test for my problem. My data is the following: imagine a set of 5 independent samples, and for each sample I can count 10 species populating it. Each sample will have different size.

This is the situation in absolute counts:

absolute counts

And this is the situation shown as percentages:


The question is: for each species (color) which samples significantly differ from which? E.g. I would like a test to tell me "The amount of Species 4 (purple) is significantly different between IX1 and VIII2, IX1 and VIII5, IX1 and VIII4, VIII1 and VIII2, VIII1 and VIII5..." ...etc, for all possible pairwise comparisons. Similarly to what a z-test would tell me for a pairwise comparison of proportions.

Now, given the difference in size between samples, my guess is that it would be fair to work with percentages. But maybe I could also randomly subsample (rarefy) all the samples to a size of ~800k (the size of the smallest sample, VIII4), and work with absolute counts...?

Thanks for your suggestions.

Edit: the samples are 5 "soups" of 10 microbes that are put in a DNA sequencer. For each sample/soup, the sequencer spits out hundreds of thousands of DNA reads, and the total amount of these reads can vary between the 5 samples for multiple reasons, it's never perfectly even. For this study's purpose, I have to classify these reads and assign each one of them to one of the 10 species, represented by colors (easy). Finally, the question is: given these assignments - represented as counts or percentages in the plots - is there a way to prove that "species purple is significantly different across samples; in particular between sample xxx and yyy"? Something like ANOVA and its post-hoc tests do.

  • 1
    $\begingroup$ You generally don't want to throw away data and it's generally best to work as close as possible to the original data (counts), but it's hard to know how best to proceed unless you can say why the samples have different sizes. The more details you can provide about the nature of the study and the sampling, the better. Please provide that information by editing the question, as comments are easy to overlook and can get deleted. $\endgroup$
    – EdM
    Commented Jan 10, 2022 at 20:00
  • $\begingroup$ Thanks for your comment. Edited the post, hopefully it's clearer. $\endgroup$
    – theDude83
    Commented Jan 10, 2022 at 20:51
  • $\begingroup$ Ps. I'm not excluding the possibility that I cannot draw any conclusion from this kind of data, and I just have to observe and describe the difference in counts as the main finding, that's it. $\endgroup$
    – theDude83
    Commented Jan 10, 2022 at 21:01

1 Answer 1


Before you jump into analyzing the data, do think more about what might be leading to the different numbers of mapped reads among the samples. Biological and biochemical differences that could bias your results come to mind.

For an example of the biological: what if some other microbial species are present, competing against the 10 you're evaluating? Then your unmapped reads might have been mapped against the other species, had you known to look.

For an example of the biochemical: GC-rich sequences pose difficulties for many sequencing methods, potentially leading to un-mappable reads. What if different analysis batches (for different samples) had different sensitivity to GC content, and your species differ in such GC-rich sequences? Then apparent differences in species prevalence among samples might represent measurement difficulties rather than true prevalence differences.

If you've ruled out such explanations adequately, then you could consider two types of generalized linear model: a multinomial regression or a count-based (Poisson or negative binomial) regression. Either analysis would be based on the number of counts for each species and sample.

A multinomial regression would label the 10 species as outcome classes. You would choose one as a reference species and then effectively examine a set of binomial comparisons of each other species against that reference species, with the samples included as levels of a single multi-level predictor variable. The regression coefficients would represent the log-odds, within a sample type, of a count being mapped to a corresponding non-reference species rather than to the reference species. (Prevalence of the reference species among samples is included in intercepts.)

A count-based regression would use counts as the outcome, and include the species, the sample, and sample:species interactions as predictors. The overall number of counts in a sample are accounted for by an offset term, the log of the number of total counts (as count modeling is done with what's called a log link between the predictors and the outcome). The regression coefficients are related to the log of the number of observed counts in the different combinations of sample and species, with the total number of counts for each sample taken into account.

In both types of analyses, the coefficient covariance matrix and the formula for the variance of a weighted sum of variables or Wald tests can be used to evaluate whatever post-hoc comparisons among species/samples you have in mind. Standard statistical software should be able to do such comparisons for you. This UCLA web page has links to how to perform multinomial regression or count-based (Poisson or negative-binomial) modeling with 5 different statistical software packages.

With the very large numbers of counts involved, it seems likely that many or all of the differences among species/samples will be "statistically significant" even after you take the multiple comparisons into account. With Poisson-distributed counts, the variance equals the mean. So with even as few as 10,000 counts the standard deviation due to sampling is only 100 counts, or 1%. I suspect that your major task will be trying to evaluate the biological significance of the results.

Finally, this seems to be a fairly simple example of metagenomic analysis. Although I don't have experience with that, I suspect that there are already well developed tools for this type of analysis, tools that can take things like batch effects and GC content into account. Check the Bioconductor website.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.