I am analysing shotgun metagenomic data from approximately 9,300 stool samples using MaAsLin3. We have raw paired-end read pairs before host and quality filtering, as well as MGS count library size, defined as the sum of counts across the complete, unfiltered MGS count table after host and quality filtering, mapping and MGS quantification.
An analysis template supplied to us used MGS count library size as unscaled_abundance. However, our understanding from the MaAsLin3 documentation is that this argument requires an independent estimate of total microbial load, such as qPCR, spike-in data or cell counts, and that sequencing-derived library size should not be used. Is this interpretation correct?
I found the previous forum recommendation to use the untransformed number of post-QC, host-filtered reads entering microbial mapping as the read-depth covariate. Unfortunately, this number is not available in our data. In this situation, would you recommend using raw read pairs or MGS count library size as the covariate?
The distinction matters in our data because MGS count library size differs between the phenotype groups, whereas raw sequencing depth does not. Any advice on which measure would be most appropriate in this setting would be greatly appreciated.
First, your understanding of unscaled_abundance is correct in that there must be an independent estimate of microbial load to use this parameter correctly. The sequencing-based library size is chosen by the researcher and is independent of the actual cell counts, so it does not provide any useful information for determining the absolute abundances of the taxa in the samples. It’s probably worth mentioning this to whomever supplied the analysis template so that other people won’t be confused. (If you want to link this thread or point them to the forum for questions, feel free!)
Instead, like you said, the right thing to do with read depth is include it as an untransformed numeric covariate.
I don’t think I understand the difference between the raw read pairs and the MGS count size in your case if you’re saying you don’t have the post-QC, host-filtered reads. I thought at the top you defined MGS count library size as “the sum of counts across the complete, unfiltered MGS count table after host and quality filtering,” which should just be the post-QC host-filtered reads we recommend. Maybe you’re saying you only have one number across the whole dataset rather than per sample? If you could clarify which numbers you have, I could make a recommendation.
Thank you so much for your reply. That is very helpful, and I’m sorry that my description of the MGS count library size was unclear.
Both measures are available per sample:
Raw read pairs: The number of paired-end read pairs per sample before any filtering and removal of host DNA.
MGS count library size: The row sum across all 6,809 MGS in the count table for each sample. This measure is obtained after host and quality filtering – but also after mapping to the HMR05 gene catalogue and MGS quantification based don signature genes.
As far as I understand, the best measure would be the number of reads remaining after quality and host filtering, but before mapping to the HMR05 gene catalogue. Unfortunately, we don´t have this measurement. So, given the two available measures, would you recommend using raw read pairs (measured before quality and host filtering) or summed MGS counts (measured after mapping and MGS quantification) as the read-depth covariate?
For additional context, the two measures are only moderately correlated (Spearman’s ρ = 0.50), and summed MGS counts differ between the phenotype groups, while raw read pairs do not.
Ah, that makes sense. I think raw read pairs would be the better option since the read filtering and host DNA is probably relatively similar across the samples, so the raw read pairs are probably well correlated with the post-QC, host-depleted reads you’d want. By contrast, I think the MGS count library size might differ a lot sample-to-sample based on the mapping rate of what’s in the sample rather than the number of reads used to map in the first place.