A Practical Guide to Statistical Analysis of Metagenomic Count Data: From Normalization to Differential Abundance Testing
By Dr. Zubair Khalid, DVM, MS, PhD ·

Key Takeaways
- Metagenomic count data is inherently compositional, sparse, and overdispersed, necessitating specialized statistical approaches beyond standard methods designed for other data types. Raw counts do not directly reflect biological abundance due to the fixed total number of reads per sample, where an increase in one feature necessitates a decrease in others.
- Normalization methods like TMM, RLE, and CSS are crucial for adjusting for technical variation and library size differences, but their suitability depends on data characteristics; TMM and RLE assume a stable majority of features, while CSS is designed for sparse data with many zeros. Total count normalization is insufficient as it can distort the relative abundances of rare taxa.
- Differential abundance testing requires models that account for overdispersion and sparsity, such as the negative binomial distribution used by DESeq2 and edgeR, or bias-corrected methods like ANCOM-BC that explicitly address compositional effects. Zero-inflated mixed models are recommended for data with excess zeros.
- Reproducibility in metagenomic analysis hinges on meticulous record-keeping, including a detailed decision log documenting every analytical choice (e.g., filtering thresholds, normalization method, model parameters) and its rationale, alongside a pre-registered analysis plan specifying the entire workflow before data exploration.
- Common failure patterns include ignoring the compositional nature of data, applying RNA-seq methods without adjustment for microbiome-specific characteristics, overlooking batch effects, inadequate multiple testing correction, and misinterpreting relative abundance as absolute abundance, all of which can lead to spurious findings.
- Clinical applications of metagenomics, such as pathogen detection, must account for host DNA and compositional structure, and quantitative measures like mapped read numbers may have prognostic value, as demonstrated in studies of pulmonary tuberculosis and severe acute pancreatitis.
Metagenomic count data analysis requires a deliberate workflow that moves from raw taxonomic or functional abundance tables through normalization to differential abundance testing. The central problem researchers face is that count tables from metagenomic sequencing are compositional, sparse, and overdispersed, so standard statistical methods designed for other data types will produce misleading results. This guide provides a concrete decision framework for choosing normalization methods and differential abundance tools, with practical R code examples and troubleshooting criteria for common failures.
At a Glance
The table below summarizes the main normalization methods and differential abundance tools discussed in this guide, along with their primary use cases and limitations.
| Method or Tool | Primary Use Case | Key Assumption | Common Limitation |
|---|---|---|---|
| TMM (Trimmed Mean of M-values) | Comparing relative abundance between conditions when most taxa do not change | Most features are not differentially abundant | Can be biased when a large proportion of features change |
| RLE (Relative Log Expression) | Scaling samples by median fold change relative to a reference | Majority of features are stable across samples | Sensitive to extreme outlier features |
| CSS (Cumulative Sum Scaling) | Normalizing sparse microbiome data with many zeros | A quantile of low-abundance features represents the scaling factor | May not perform well with very low sequencing depth |
| DESeq2 | Differential abundance testing with count-based models | Negative binomial distribution with shrinkage estimation | Requires integer counts and may be conservative with sparse data |
| edgeR | Differential expression and abundance analysis | Negative binomial distribution with empirical Bayes moderation | Designed for RNA-seq and may need adjustment for microbiome sparsity |
| ANCOM-BC | Differential abundance testing with bias correction | Accounts for compositional structure and sampling fraction | Computationally intensive with many features |
Understanding Metagenomic Count Data Structure
Metagenomic count tables record the number of sequencing reads assigned to each taxonomic or functional feature in each sample. These tables come from profiling pipelines that process raw sequencing data through quality filtering, taxonomic assignment, and abundance estimation. The structure of these tables creates specific statistical challenges that shape every downstream decision.
Compositional Nature of Count Data
A metagenomic count table is compositional because each sample has a fixed total number of reads, so the abundance of any feature is relative to all other features in that sample. If one taxon increases in abundance, the relative counts of all other taxa must decrease even when their absolute abundance is unchanged. This constraint means that raw counts do not directly reflect biological abundance and that standard correlation or distance measures can produce spurious results.
The compositional structure also affects how normalization methods work. Methods that scale each sample by a single factor assume that the total count reflects sequencing effort instead of biological signal. When the total count is dominated by a few highly abundant taxa, scaling by total count can distort the apparent abundance of rare taxa.
Sparsity and Zero Inflation
Metagenomic count tables typically contain a high proportion of zeros. Many taxa are absent from most samples, and rare taxa may be missed entirely due to limited sequencing depth. This sparsity creates two problems for statistical analysis. First, zero counts do not distinguish between biological absence and technical failure to detect a present taxon. Second, the distribution of counts is zero-inflated, meaning that standard count distributions such as the Poisson or negative binomial may not fit the data well.
Zero-inflated mixed models have been developed specifically for microbiome and metagenomics data analysis to address this issue. These models combine a count distribution with a separate process that generates excess zeros, allowing researchers to model both the presence or absence of a taxon and its abundance when present. The NBZIMM package implements negative binomial and zero-inflated mixed models for this purpose [<a href="#ref-1">1</a>].
Overdispersion Relative to Poisson Expectations
Count data from metagenomic sequencing show more variability than expected under a simple Poisson model, where the variance equals the mean. This overdispersion arises from biological variability between samples, technical variability in library preparation, and the compositional structure of the data. Differential abundance testing methods must account for this extra variability to avoid false positives.
Negative binomial models are the standard approach for overdispersed count data because they include a dispersion parameter that allows the variance to exceed the mean. Both DESeq2 and edgeR use negative binomial models with different strategies for estimating dispersion, and both require careful attention to the sparsity and compositional structure of metagenomic data.
Data Input Requirements and Quality Checks
Before any statistical analysis, the count table must be assembled and checked for quality issues that could invalidate downstream results. The quality of the input data determines the reliability of every subsequent step, so this stage deserves deliberate attention.
Required Input Files
The primary input is a count matrix with features in rows and samples in columns. Each cell contains the number of sequencing reads assigned to that feature in that sample. Feature identifiers may be taxonomic labels, operational taxonomic unit identifiers, or functional gene identifiers depending on the profiling pipeline used.
A metadata table is also required, with samples in rows and experimental variables in columns. This table must include the grouping variable for differential abundance testing, along with any covariates such as sequencing batch, extraction batch, or clinical variables that should be included in the model.
The count matrix and metadata table must use matching sample identifiers. Mismatched identifiers are a common source of errors that are difficult to detect after analysis has begun.
Sequencing Depth Assessment
Sequencing depth varies between samples, and this variation must be accounted for in normalization. The total number of reads per sample should be examined before analysis to identify samples with very low depth that may not represent the underlying community adequately.
Samples with extremely low sequencing depth may need to be removed or analyzed separately. There is no universal threshold for minimum depth, but samples that fall far below the distribution of other samples should be flagged for review. The decision to exclude a sample should be documented in the analysis record.
Library Size and Composition Checks
The distribution of library sizes, meaning the total counts per sample, should be visualized to identify outliers. Samples with unusually high or low total counts may indicate technical problems such as contamination, failed amplification, or sequencing errors.
The composition of each sample should also be examined. Samples dominated by a single taxon may indicate contamination or a technical artifact instead of a true biological signal. Comparing the relative abundance of known control features across samples can help identify technical variation.
Batch Effect Detection
Batch effects arise when samples are processed in different groups, such as different sequencing runs, different extraction batches, or different operators. These technical differences can create systematic variation that mimics or masks biological signals.
Batch effects should be assessed before differential abundance testing. Principal component analysis or clustering of samples by batch can reveal whether technical factors explain a large proportion of the variation. If batch effects are present, they should be included as covariates in the statistical model or addressed through batch correction methods.
Normalization Methods for Metagenomic Count Data
Normalization adjusts count data to make samples comparable by removing technical variation while preserving biological signal. The choice of normalization method has a substantial impact on downstream differential abundance results, so the decision should be based on the data structure and the research question.
Why Total Count Normalization Is Insufficient
The simplest normalization approach divides each count by the total number of reads in that sample, producing relative abundances. This method assumes that the total count is purely technical and that all samples have the same biological composition.
Total count normalization fails when a few highly abundant taxa dominate the library. If one taxon increases dramatically in one sample, the relative abundances of all other taxa decrease even when their absolute counts are unchanged. This distortion can create false differential abundance signals and mask true signals.
TMM Normalization
Trimmed Mean of M-values normalization estimates a scaling factor for each sample based on the weighted mean of log fold changes between that sample and a reference sample, after trimming the most extreme values. The method assumes that most features are not differentially abundant, so the trimmed mean of fold changes represents the technical scaling factor.
TMM is implemented in the edgeR package and is widely used for metagenomic count data. The method performs well when the assumption of a stable majority of features holds, but it can be biased when a large proportion of features change between conditions.
RLE Normalization
Relative Log Expression normalization scales each sample by the median of the ratios of each feature count to the geometric mean of that feature across samples. The method assumes that most features are stable across samples and that the median ratio represents the technical scaling factor.
RLE is implemented in the DESeq2 package and is the default normalization for that tool. The method is robust to outliers when the majority of features are stable, but it can be sensitive to features with very high or very low counts.
CSS Normalization
Cumulative Sum Scaling normalizes samples by the sum of counts up to a quantile that is determined from the data. The method is designed for sparse microbiome data where many features have zero counts and the distribution of counts is heavily skewed.
CSS is implemented in the metagenomeSeq package and is particularly useful for data with many rare features. The method selects a quantile that captures the majority of the signal while excluding the highly abundant features that dominate the total count.
Choosing a Normalization Method
The choice of normalization method should be guided by the data structure and the downstream analysis. For data where most features are expected to be stable, TMM or RLE are reasonable choices. For data with extreme sparsity and skewed abundance distributions, CSS may be more appropriate.
The normalization method should be reported in the methods section of any publication, along with the rationale for the choice. Sensitivity analyses comparing results across normalization methods can help determine whether conclusions are robust to the normalization choice.
Differential Abundance Testing Methods
Differential abundance testing identifies features whose abundance differs between conditions. The choice of method depends on the data structure, the experimental design, and the assumptions that can be reasonably made about the data.
DESeq2 for Count-Based Differential Abundance
DESeq2 models count data with a negative binomial distribution and uses empirical Bayes shrinkage to estimate dispersion and fold changes. The method is designed for RNA-seq data but is widely applied to metagenomic count data.
DESeq2 requires integer counts and works best when the number of samples per group is reasonably balanced. The method estimates size factors for normalization and then fits a generalized linear model for each feature. The shrinkage of dispersion estimates helps stabilize results for features with low counts.
For metagenomic data, DESeq2 can be conservative because the method assumes that most features are not differentially abundant and shrinks fold changes toward zero. This assumption may not hold for microbiome data where many taxa change between conditions.
edgeR for Empirical Bayes Moderation
edgeR also uses a negative binomial model with empirical Bayes moderation of dispersion. The method is implemented in the edgeR package and provides a range of statistical tests for differential abundance.
edgeR requires a normalization step before testing, typically using TMM. The method estimates dispersion for each feature and moderates extreme dispersion estimates toward a common value. This moderation helps improve power for features with low counts.
For metagenomic data, edgeR can be sensitive to the choice of normalization method and to the presence of features with very low counts. The method assumes that the negative binomial model fits the data, which may not hold for zero-inflated microbiome data.
ANCOM-BC for Bias Correction
ANCOM-BC (Analysis of Composition of Microbiomes with Bias Correction) accounts for the compositional structure of metagenomic data by estimating and correcting for sampling fractions. The method uses a linear regression framework with bias correction to identify differentially abundant features.
ANCOM-BC does not require a separate normalization step because it estimates the sampling fraction from the data. The method is designed specifically for microbiome data and addresses the compositional constraint that affects other methods.
The main limitation of ANCOM-BC is computational intensity when the number of features is large. The method also requires careful specification of the model formula and may be sensitive to the choice of covariates.
Zero-Inflated Mixed Models
Zero-inflated mixed models address the excess zeros in metagenomic count data by modeling the presence or absence of a feature separately from its abundance when present. The NBZIMM package implements these models for microbiome and metagenomics data analysis [<a href="#ref-1">1</a>].
These models are particularly useful for longitudinal or clustered data where repeated measurements from the same subject create correlation. The mixed model framework allows random effects for subjects or other grouping factors.
The main limitation of zero-inflated mixed models is complexity. These models require careful specification and may be difficult to fit when the number of features is large. The interpretation of results also requires attention to the distinction between presence or absence effects and abundance effects.
Class Prediction and Feature Selection Approaches
Beyond differential abundance testing, some analyses aim to predict class membership from metagenomic count data or to select features that discriminate between classes. Linear optimization approaches have been developed for class prediction and feature selection with metagenomic count data [<a href="#ref-2">2</a>].
These methods identify a subset of features that together predict the class label, instead of testing each feature individually. This approach can be useful for biomarker discovery or for building diagnostic classifiers.
The main limitation of class prediction approaches is the risk of overfitting, particularly when the number of features greatly exceeds the number of samples. Cross-validation and external validation are essential to assess the reliability of any classifier.
Practical Workflow for Differential Abundance Analysis
The workflow below provides a step-by-step approach to differential abundance analysis of metagenomic count data. Each step includes concrete decisions and checks that should be documented in the analysis record.
Step 1: Load and Validate Input Data
Load the count matrix and metadata table into R. Verify that sample identifiers match between the two files and that the count matrix contains integer values. Check for missing values and decide how to handle features with very low counts across all samples.
counts <- read.csv("count_matrix.csv", row.names = 1)
metadata <- read.csv("metadata.csv", row.names = 1)
stopifnot(all(colnames(counts) %in% rownames(metadata)))
Step 2: Filter Low-Count Features
Filter features with very low counts across most samples. A common approach is to retain features with at least a minimum count in a minimum number of samples. The filtering threshold should be documented and justified.
keep <- rowSums(counts >= 10) >= 2
counts_filtered <- counts[keep, ]
Step 3: Assess Sequencing Depth and Sample Quality
Calculate total counts per sample and visualize the distribution. Flag samples with unusually low depth for review. Examine the composition of each sample to identify potential contamination or technical artifacts.
Step 4: Choose and Apply Normalization
Select a normalization method based on the data structure and research question. Apply the normalization and verify that the scaling factors are reasonable. Document the normalization method and any sensitivity analyses.
Step 5: Fit the Statistical Model
Specify the model formula with the grouping variable and any covariates. Fit the differential abundance model using the chosen method. Check that the model converges and that dispersion estimates are reasonable.
Step 6: Extract and Interpret Results
Extract the results table with fold changes, p-values, and adjusted p-values. Apply multiple testing correction to control the false discovery rate. Interpret the results in the context of the biological question and the limitations of the data.
Step 7: Validate and Report
Validate the results through sensitivity analyses, such as comparing results across normalization methods or testing different filtering thresholds. Report the full analysis workflow, including software versions, parameters, and quality checks.
Records and Measurements for Reproducible Analysis
Reproducibility requires careful record keeping at every stage of the analysis. The records should be sufficient for another researcher to repeat the analysis and obtain the same results.
Analysis Log
Maintain an analysis log that records the date, software versions, parameters, and decisions for each analysis step. The log should include the rationale for each decision, such as the choice of normalization method or filtering threshold.
Version Control
Use version control for analysis scripts and documentation. The Carpentries lessons provide foundational training in version control with Git, which is essential for tracking changes to analysis code and collaborating with other researchers [<a href="#ref-3">3</a>].
Containerized Workflows
Containerized workflows, such as those provided by nf-core, ensure that the analysis environment is reproducible across different computing systems. These workflows package the software, dependencies, and configuration in a standardized format [<a href="#ref-4">4</a>].
Data Management
Store raw sequencing data, count tables, metadata, and analysis outputs in organized directories with clear naming conventions. The NCBI provides data resources for depositing and accessing sequencing data, which supports data sharing and reproducibility [<a href="#ref-5">5</a>].
Common Failure Patterns and Troubleshooting
Several recurring problems affect metagenomic count data analysis. Recognizing these patterns early can prevent wasted effort and invalid conclusions.
Failure Pattern 1: Ignoring Compositional Structure
Researchers who treat raw counts as absolute abundances and apply standard statistical methods without normalization will obtain results that reflect library size differences instead of biological differences. The solution is to use appropriate normalization and to interpret results as relative abundance changes.
Failure Pattern 2: Applying RNA-seq Methods Without Adjustment
DESeq2 and edgeR are designed for RNA-seq data and make assumptions that may not hold for metagenomic data. The sparsity and compositional structure of microbiome data require careful evaluation of whether these methods are appropriate. Sensitivity analyses comparing multiple methods can help identify robust results.
Failure Pattern 3: Overlooking Batch Effects
Batch effects can create spurious differential abundance signals when samples from different conditions are processed in different batches. The solution is to include batch as a covariate in the model or to use batch correction methods before analysis.
Failure Pattern 4: Inadequate Multiple Testing Correction
Testing thousands of features simultaneously requires multiple testing correction to control the false discovery rate. Failure to apply correction will produce many false positives. The Benjamini-Hochberg procedure is the standard approach for controlling the false discovery rate.
Failure Pattern 5: Ignoring Zero Inflation
Standard negative binomial models may not fit zero-inflated data well, leading to biased estimates and incorrect p-values. Zero-inflated mixed models or other methods that account for excess zeros should be considered when sparsity is extreme [<a href="#ref-1">1</a>].
Failure Pattern 6: Overinterpreting Results From Small Samples
Metagenomic studies often have small sample sizes, which limits statistical power and increases the risk of false positives and false negatives. Results from small studies should be interpreted cautiously and validated in independent cohorts.
Limitations and Interpretation Boundaries
Statistical analysis of metagenomic count data has inherent limitations that should be acknowledged in any report or publication.
Relative Abundance Versus Absolute Abundance
Normalization methods produce relative abundance estimates, not absolute abundance. A feature can appear to decrease in relative abundance when its absolute abundance is unchanged or even increased, if other features increase more. This limitation affects the interpretation of differential abundance results.
Detection Limits and Sequencing Depth
Features present at very low abundance may not be detected, particularly in samples with low sequencing depth. The absence of a feature from the count table does not prove its biological absence. This limitation is particularly relevant for clinical applications where pathogen detection is the goal.
Taxonomic Resolution
The resolution of taxonomic assignment depends on the profiling method and the reference database. Features assigned at different taxonomic levels may not be directly comparable. The choice of taxonomic level for analysis should be documented and justified.
Reference Database Dependence
Taxonomic and functional assignment depends on the reference database used. Different databases may produce different assignments for the same sequencing data. The database version should be reported, and results should be interpreted in the context of the database coverage.
Clinical and Diagnostic Context
Metagenomic next-generation sequencing is increasingly used in clinical settings for pathogen detection and diagnosis. The statistical considerations discussed in this guide have direct implications for clinical interpretation.
Pathogen Detection in Clinical Samples
Metagenomic sequencing can identify pathogens in clinical samples, including cases where conventional culture is negative. A study of pulmonary tuberculosis diagnosis in patients with acute exacerbations of chronic obstructive pulmonary disease found that metagenomic sequencing demonstrated higher sensitivity than TB-DNA detection and Xpert MTB/RIF assay, while also identifying bacterial or viral co-infections in a proportion of tuberculosis cases [<a href="#ref-6">6</a>].
The statistical analysis of clinical metagenomic data must account for the compositional structure of the sample and the presence of host DNA. The stringently mapped read number showed a positive correlation with intensive care unit admission rate and in-hospital mortality in this study, suggesting that quantitative measures from metagenomic data may have prognostic value [<a href="#ref-6">6</a>].
Time-Dependent Microbiology in Disease Progression
The microbiological characteristics of clinical samples can evolve over time. A prospective observational study of peripancreatic drainage fluid in severe acute pancreatitis found that microbiological positivity was low within 14 days of disease onset but increased in patients undergoing drainage beyond 14 days. Metagenomic sequencing identified a broader spectrum of pathogens than conventional culture, particularly polymicrobial, anaerobic, and fungal organisms [<a href="#ref-7">7</a>].
These findings illustrate the importance of considering disease timing in the interpretation of metagenomic results. Statistical analyses that ignore temporal dynamics may miss clinically relevant patterns.
Professional Escalation Criteria
When metagenomic analysis is used for clinical decision-making, certain findings warrant escalation to specialized personnel. These include detection of pathogens with public health implications, unusual resistance patterns, or discrepancies between metagenomic results and clinical presentation.
The decision to escalate should be guided by institutional protocols and clinical judgment. Statistical results should be interpreted in the context of the patient's clinical status and other diagnostic information.
Training and Reproducibility Resources
Developing the skills for metagenomic count data analysis requires training in bioinformatics, statistics, and reproducible research practices.
Bioinformatics Training Pathways
The EMBL-EBI Training program provides learning pathways for bioinformatics data resources and practical analysis education [<a href="#ref-8">8</a>]. These resources cover the fundamentals of working with biological sequence data and the tools used for metagenomic analysis.
Workflow Training
The Galaxy Training Network offers accessible workflow training and analysis tutorials for a range of bioinformatics applications [<a href="#ref-9">9</a>]. These tutorials provide hands-on experience with metagenomic analysis workflows and reproducible analysis practices.
Foundational Computing Skills
The Carpentries lessons provide foundational training in computing, data, shell, Git, and programming [<a href="#ref-3">3</a>]. These skills are essential for managing the large data files and complex analysis workflows involved in metagenomic research.
Community Pipeline Standards
The nf-core documentation describes community pipeline standards for reproducible workflow usage and configuration [<a href="#ref-4">4</a>]. These standards support the development of reliable, maintainable analysis pipelines.
Bioconductor Resources
The Bioconductor project provides official package documentation, workflow guidance, and installation instructions for reproducible genomic analysis [<a href="#ref-10">10</a>]. Many of the tools discussed in this guide, including DESeq2 and edgeR, are available through Bioconductor.
Building a Decision Log and Pre-Registered Analysis Plan for Metagenomic Count Data
The most common cause of irreproducible metagenomic results is not the choice of a specific normalization method or differential abundance tool, but the absence of a documented decision trail that explains why those choices were made. Researchers often select methods based on convenience, habit, or what a colleague used, then struggle to reconstruct the rationale when reviewers ask questions or when results fail to replicate. A decision log and pre-registered analysis plan address this problem directly by forcing explicit documentation of every analytical choice before results are examined, reducing the risk of unconscious bias and providing a clear record for troubleshooting when results are unexpected.
Why a Decision Log Matters for Metagenomic Analysis
Metagenomic count data analysis involves dozens of discrete decisions, each of which can alter the final list of differentially abundant features. The filtering threshold, the normalization method, the dispersion estimation approach, the multiple testing correction procedure, and the handling of covariates all interact in ways that are difficult to predict without systematic testing. When these decisions are made ad hoc during analysis, the researcher may unconsciously choose options that produce favorable results, a problem known as researcher degrees of freedom.
A decision log addresses this problem by creating a timestamped record of each analytical choice and its justification. The log should be started before any analysis begins and updated whenever a decision is made. This practice is particularly important for metagenomic data because the compositional structure and sparsity create many legitimate analytical options, and the results can vary substantially depending on which options are selected.
The decision log also serves a practical troubleshooting function. When a downstream analysis produces unexpected results, the log allows the researcher to trace the problem back to a specific decision point. For example, if a differential abundance test returns an unusually large number of significant features, the log can reveal whether the filtering threshold was too lenient or whether the normalization method was inappropriate for the data structure.
Components of a Pre-Registered Analysis Plan
A pre-registered analysis plan extends the decision log by specifying the complete analytical workflow before any data exploration begins. The plan should be written as a standalone document that another researcher could follow to reproduce the analysis without additional guidance.
The first component is a clear statement of the primary research question and the specific hypothesis being tested. This statement should define the comparison groups, the primary outcome measure, and the features of interest. For metagenomic data, the hypothesis might specify a particular taxonomic level, such as genus or species, and a minimum effect size that would be considered biologically meaningful.
The second component is a complete specification of the data processing steps. This includes the filtering criteria for low-count features, the normalization method and its parameters, and the differential abundance testing approach. Each step should include the exact thresholds or parameters that will be used, along with the rationale for those choices.
The third component is a plan for sensitivity analyses. These are additional analyses that test whether the main results are robust to reasonable changes in the analytical approach. For example, the plan might specify that results will be compared across two normalization methods or across different filtering thresholds. The sensitivity analysis plan should be specified in advance to prevent selective reporting of favorable results.
The fourth component is a set of criteria for interpreting results. This includes the significance threshold, the multiple testing correction method, and the minimum effect size for biological relevance. The plan should also specify how outliers and technical failures will be handled.
Creating the Decision Log Template
A practical decision log template for metagenomic analysis should include the following fields for each decision point: the date, the decision made, the rationale, the alternatives considered, and the expected impact on results. The template should be simple enough to maintain consistently but detailed enough to be useful for troubleshooting.
The decision log should be maintained as a plain text file or a spreadsheet that is stored alongside the analysis scripts. Version control systems such as Git provide a natural framework for tracking changes to the decision log over time. The Carpentries lessons provide foundational training in version control with Git, which is essential for tracking changes to analysis code and documentation [<a href="#ref-3">3</a>].
A practical template for each decision log entry might look like this:
| Date | Decision Point | Decision Made | Rationale | Alternatives Considered | Expected Impact |
|---|---|---|---|---|---|
| 2026-01-15 | Filtering threshold | Retain features with count >= 10 in at least 2 samples | Reduces noise from rare features while preserving moderately abundant taxa | Count >= 5 in 1 sample, count >= 20 in 3 samples | Fewer features tested, reduced multiple testing burden |
| 2026-01-15 | Normalization method | TMM | Most features expected to be stable between conditions | RLE, CSS | Scaling factors will adjust for library size differences |
| 2026-01-16 | Differential abundance test | DESeq2 with default settings | Negative binomial model appropriate for overdispersed counts | edgeR, ANCOM-BC | Shrinkage may reduce power for rare features |
The decision log should be updated whenever a decision is made, not retrospectively at the end of the analysis. This practice ensures that the rationale is recorded while the reasoning is still fresh and prevents the common problem of post hoc justification.
Pre-Registration in Practice
Pre-registration of the analysis plan can be done through formal platforms or through internal documentation. Formal pre-registration provides a timestamped public record of the analysis plan, which can be cited in publications. Internal pre-registration, such as a dated document stored in the project repository, serves a similar function for internal quality control.
The level of detail in the pre-registration should match the complexity of the analysis. A simple two-group comparison with a single normalization method requires less detail than a longitudinal study with multiple covariates and repeated measurements. The pre-registration should be specific enough to prevent major analytical deviations but flexible enough to accommodate necessary adjustments when data quality issues are discovered.
When deviations from the pre-registered plan are necessary, they should be documented in the decision log with a clear explanation. Common reasons for deviation include unexpected batch effects, extreme sparsity that makes the planned normalization method inappropriate, or the discovery of sample mislabeling. Each deviation should be recorded with the date, the reason, and the impact on the analysis.
Integrating the Decision Log with the Analysis Workflow
The decision log should be integrated into the analysis workflow at specific checkpoints. The first checkpoint occurs after data validation and before any filtering or normalization. At this point, the researcher should record the initial data quality assessment, including the number of samples, the number of features, the distribution of library sizes, and any samples flagged for review.
The second checkpoint occurs after filtering and before normalization. The researcher should record the filtering criteria, the number of features retained, and the proportion of features removed. This information is essential for interpreting the final results because the filtering step determines which features are eligible for differential abundance testing.
The third checkpoint occurs after normalization and before differential abundance testing. The researcher should record the normalization method, the range of scaling factors, and any samples with extreme scaling factors that might indicate technical problems.
The fourth checkpoint occurs after differential abundance testing and before interpretation. The researcher should record the number of significant features, the distribution of p-values, and any features with extreme fold changes that warrant manual review.
Common Failure Patterns in Decision Documentation
Several recurring problems undermine the usefulness of decision logs and pre-registered analysis plans. Recognizing these patterns can help researchers avoid them.
The first failure pattern is documenting decisions after results are known. This practice defeats the purpose of pre-registration because the researcher may unconsciously rationalize choices that produce favorable results. The solution is to enforce a strict rule that the decision log is updated before any results are examined.
The second failure pattern is recording decisions without rationale. A log that states what was done but not why is of limited value for troubleshooting or for convincing reviewers that the analysis was conducted appropriately. Each entry should include a substantive justification based on the data structure or the research question.
The third failure pattern is failing to update the log when the analysis plan changes. Analyses rarely proceed exactly as planned, and deviations are common. When the plan changes, the log should be updated immediately with the new decision and the reason for the change.
The fourth failure pattern is treating the decision log as a bureaucratic requirement instead of a scientific tool. The log is most valuable when it is used actively during troubleshooting and when it informs the interpretation of results. A log that is written once and never consulted provides little benefit.
Using the Decision Log for Troubleshooting
When a differential abundance analysis produces unexpected results, the decision log provides a systematic framework for identifying the cause. The troubleshooting process should proceed through the log in reverse chronological order, examining each decision point for potential problems.
The first step is to check whether the filtering criteria were appropriate for the data structure. If the filtering threshold was too lenient, the analysis may include many features with very low counts that produce unstable estimates. If the threshold was too strict, the analysis may exclude biologically relevant rare features.
The second step is to examine the normalization scaling factors. Extreme scaling factors can indicate that the normalization method was inappropriate for the data. For example, if one sample has a scaling factor that is an order of magnitude different from the others, the normalization method may be sensitive to a few highly abundant features in that sample.
The third step is to review the model specification. The choice of covariates, the handling of batch effects, and the dispersion estimation approach can all influence the results. The decision log should indicate whether batch effects were assessed and how they were handled.
The fourth step is to compare results across sensitivity analyses. If the pre-registered plan included sensitivity analyses, the results of these analyses can reveal whether the main findings are robust to analytical choices. If the results change substantially across sensitivity analyses, the conclusions should be interpreted with caution.
Records and Measurements for Decision Documentation
The decision log should be stored with the analysis scripts and data files in a structured project directory. A recommended structure includes separate directories for raw data, processed data, scripts, results, and documentation. The decision log and pre-registered analysis plan should be stored in the documentation directory with clear naming conventions that include the date.
The analysis scripts should be written so that each decision point is clearly marked with a comment that references the corresponding decision log entry. This practice creates a direct link between the code and the documentation, making it easier to verify that the analysis was conducted as planned.
Containerized workflows, such as those provided by nf-core, ensure that the analysis environment is reproducible across different computing systems [<a href="#ref-4">4</a>]. The decision log should record the container version and any configuration parameters that affect the analysis.
The NCBI provides data resources for depositing and accessing sequencing data, which supports data sharing and reproducibility [<a href="#ref-5">5</a>]. The decision log should reference the accession numbers for the raw data so that reviewers can verify the analysis against the original sequencing data.
Professional Escalation Criteria for Decision Documentation
Certain situations warrant escalation to a supervisor, statistician, or bioinformatics specialist. These include cases where the decision log reveals systematic problems with the analysis, where the pre-registered plan cannot be followed due to data quality issues, or where the results are highly sensitive to analytical choices.
If the decision log reveals that multiple normalization methods produce substantially different results, the analysis should be escalated to a statistician with experience in compositional data analysis. The statistician can help determine whether the data structure violates the assumptions of the chosen methods and whether alternative approaches are warranted.
If the pre-registered plan cannot be followed because of unexpected data quality issues, such as extreme batch effects or sample contamination, the analysis should be escalated to the study team to determine whether the data can be salvaged or whether additional data collection is needed.
If the results are highly sensitive to the filtering threshold or the normalization method, the analysis should be escalated to the study team to determine whether the biological conclusions are robust enough to report. In some cases, the sensitivity of the results to analytical choices indicates that the biological signal is weak and that additional data are needed.
Training Resources for Decision Documentation
Developing the skills for systematic decision documentation requires training in reproducible research practices and data management. The EMBL-EBI Training program provides learning pathways for bioinformatics data resources and practical analysis education [<a href="#ref-8">8</a>]. These resources cover the fundamentals of working with biological sequence data and the tools used for metagenomic analysis.
The Galaxy Training Network offers accessible workflow training and analysis tutorials for a range of bioinformatics applications [<a href="#ref-9">9</a>]. These tutorials provide hands-on experience with metagenomic analysis workflows and reproducible analysis practices, including the use of workflow histories that serve as a form of decision documentation.
The Carpentries lessons provide foundational training in computing, data, shell, Git, and programming [<a href="#ref-3">3</a>]. These skills are essential for managing the large data files and complex analysis workflows involved in metagenomic research, and for maintaining the version control systems that support decision documentation.
The Bioconductor project provides official package documentation, workflow guidance, and installation instructions for reproducible genomic analysis [<a href="#ref-10">10</a>]. Many of the tools discussed in this guide, including DESeq2 and edgeR, are available through Bioconductor, and their documentation includes guidance on parameter choices that can inform the decision log.
Frequently Asked Questions
What is the difference between normalization and differential abundance testing?
Normalization adjusts the count data to make samples comparable by removing technical variation such as differences in sequencing depth. Differential abundance testing uses the normalized data to identify features whose abundance differs between conditions. Normalization is a prerequisite for meaningful differential abundance testing because raw counts reflect both biological signal and technical variation.
Why cannot I use raw counts directly for differential abundance testing?
Raw counts reflect both biological abundance and sequencing effort. Samples with more total reads will have higher counts for all features, regardless of biological differences. This technical variation must be removed through normalization before comparing features between conditions. Additionally, the compositional structure of metagenomic data means that raw counts do not directly reflect absolute abundance.
How do I choose between TMM, RLE, and CSS normalization?
The choice depends on the data structure and the research question. TMM and RLE assume that most features are stable across samples and are appropriate when this assumption holds. CSS is designed for sparse data with many zeros and may perform better when the abundance distribution is heavily skewed. Sensitivity analyses comparing results across normalization methods can help determine whether conclusions are robust.
When should I use ANCOM-BC instead of DESeq2 or edgeR?
ANCOM-BC accounts for the compositional structure of metagenomic data by estimating and correcting for sampling fractions. This method is designed specifically for microbiome data and may be more appropriate than DESeq2 or edgeR when compositional effects are a concern. However, ANCOM-BC is computationally intensive and may not be practical for very large feature sets.
How do I handle zero counts in metagenomic data?
Zero counts can represent biological absence or technical failure to detect a present taxon. Filtering features with very low counts across most samples is a common first step. For features retained in the analysis, zero-inflated mixed models can model the excess zeros separately from the count distribution [<a href="#ref-1">1</a>]. The choice of approach depends on the research question and the proportion of zeros in the data.
What is the role of multiple testing correction in differential abundance analysis?
Testing thousands of features simultaneously creates a high probability of false positives if no correction is applied. Multiple testing correction controls the expected proportion of false positives among the features declared significant. The Benjamini-Hochberg procedure is the standard approach for controlling the false discovery rate in metagenomic studies.
How do I account for batch effects in metagenomic analysis?
Batch effects should be assessed before differential abundance testing by examining whether samples cluster by batch in principal component analysis or other ordination methods. If batch effects are present, include batch as a covariate in the statistical model. Alternatively, batch correction methods can be applied before analysis, but these methods have their own assumptions and limitations.
What should I report in the methods section of my publication?
Report the software versions, normalization method, filtering thresholds, statistical model, and multiple testing correction procedure. Include the rationale for each decision and any sensitivity analyses performed. This information is essential for reproducibility and for readers to assess the validity of the results.
Related Bioinformatics Guides
- Genomic Data Analysis Tools: A Comparative Guide for Researchers
- Metagenomics Functional Profiling: Tools and Databases for Pathway Analysis
- Proteomics Data Analysis in R: A Practical Workflow for Differential Expression and Visualization
- Metabolomics Data Analysis in R: A Practical Workflow
- Microbiome Data Analysis in R: A Practical Guide for Compositional Data
Related Clinical & Scientific Guides
- A Practical Guide to Detecting Antimicrobial Resistance Genes in Shotgun Metagenomic Data
- Computational Immunology: Modeling the Immune System
- How to Set Hard Filters for Germline Variant Calling: A Practical Guide to GATK Best Practices
References and Further Reading
[1] [NBZIMM: negative binomial and zero-inflated mixed models, with application to microbiome/metagenomics data analysis](https://doi.org/10.1186/s12859-020-03803-z). BMC Bioinformatics, 2020. [2] [Class Prediction and Feature Selection with Linear Optimization for Metagenomic Count Data](https://doi.org/10.1371/journal.pone.0053253). Plos One, 2013. [3] [The Carpentries Lessons](https://carpentries.org/lessons). The Carpentries. [4] [nf-core Documentation](https://nf-co.re/docs). nf-core. [5] [NCBI Data Resources](https://www.ncbi.nlm.nih.gov/). National Center for Biotechnology Information. [6] [Metagenomic Next-Generation Sequencing for Pulmonary Tuberculosis Diagnosis and Infection Risk Factor Analysis in AECOPD Patients: A Single-Center Retrospective Study.](https://doi.org/10.3390/jcm15124507). 2026. [7] [Time-dependent microbiology of peripancreatic drainage fluid in severe acute pancreatitis: a prospective real-world observational study using metagenomic sequencing and culture.](https://doi.org/10.3389/fmed.2026.1795250). 2026. [8] [EMBL-EBI Training](https://www.ebi.ac.uk/training). European Bioinformatics Institute. [9] [Galaxy Training Network](https://training.galaxyproject.org/). Galaxy Project. [10] [Bioconductor](https://bioconductor.org/). Bioconductor Project.This article is educational and does not replace validated analysis plans, institutional policy, clinical interpretation, or specialist review.