# Reads count in metagenomics

Background: I am developing a pipeline for metagenomic studies of human gut microbiote. In particular, I am mapping the reads data originating from shotgun whole genome sequencing onto a gene catalogue (similar to IGC) and counting reads mapped on each feature/gene. Till now I have been using NGLess, which is a rather convenient tool. However, we have some questions about how exactly it counts reads and it might be not the best option for scaling up our pipeline.

Question:
I am looking for alternative tools for gene count. A brief search has shown that htseq-count and featureCount are two frequently used tools. However they mainly appear in RNAseq context, and I am not certain whether they could be used in my case. I would also appreciate comments on the possibility of using samtools mpileup.

## 1 Answer

I think that there is a good argument for using technologies adapted from RNA-seq for metagenomic abundance quantification. The tool that I am most familiar with (that I use for this purpose) is kallisto, which uses pseudoalignment instead and is therefore less compute intensive.

For some background on the conceptual linkages between metagenomic abundance and transcriptomic abundance you can see the kallisto metagenome abundance paper or the explanatory blog post.

I think that it is preferable to use a tool that is specifically adapted for metagenomics, as there are always some numerical issues that might have been accounted for by the authors. I think that for this purpose it is maybe best not to naively use tools like mpileup, though I have no specific objection to it. However, if you look around you can see various papers using e.g. Cufflinks and Cuffdiff for effectively the same purpose, and Cufflinks is definitely intended for RNA-seq.

Depending on your specific application, you might want to use a different tool like mOTUS2, if you are interested in estimating abundance of specific genomes or lineages in a metagenome mixture. This method is probably more accurate as it is based on averaging across the whole genome, but it is not very useful if what you are interested in is more e.g. functional profiling of the metagenome.

• Thank you for answering. We did try motus2, but it is restricted to a custom-made gene database and the marker genes grounded in the associeted species classification. What I am looking for is a tool for creating a gene count table, which we then could use to develop our own methods of species identification. Jan 15 '21 at 12:24
• @Vadim In that case, I think that a tool like kallisto (or even cufflinks) might be most helpful. These will be providing counts as TPMs (though they also provide raw counts), of course, so if there is some statistical issue with that kind of compositional count then it might be a problem. However, metagenomics is intrinsically compositional so I'm not sure other methods will be much better. If you are ultimately going to do deconvolution of some kind, it might be helpful to look at what metabat2 is doing (more for ideas than for usefulness as a tool, as it's solving a different problem). Jan 15 '21 at 20:46