pastForward: a Snakemake pipeline for ancient and historical DNA with eukaryote-wide taxonomic screening and tracking of copy-number variation

Technical Summary of pastForward

Problem Statement
Ancient (aDNA) and historical (hDNA) sequencing offers direct windows into past populations and evolutionary dynamics, yet bioinformatic processing remains a significant bottleneck. Standard genomics pipelines are ill-suited for aDNA/hDNA due to the unique characteristics of degraded DNA: short fragment lengths (35–100 bp), post-mortem damage (C→T and G→A substitutions), and high rates of contamination. Existing specialized pipelines (e.g., PALEOMIX, nf-core/eager, mapache) often lack integrated taxonomic screening, fail to distinguish nuclear mitochondrial DNA segments (NUMTs) from authentic mitochondrial signals, do not support copy-number variation (CNV) tracking across time series, or rely on workflow managers requiring extensive prior knowledge. Furthermore, taxonomic classifiers often demand excessive memory or fail to classify short fragments effectively.

Methodology
pastForward is a fully automated, reproducible Snakemake pipeline designed to process raw reads into damage-rescaled BAM files and interactive quality reports. It is organized into four modular components, each configurable via a central configuration file:

1. Read Module: Handles raw FASTQ data (single- and paired-end). It performs adapter trimming and poly-X tail removal using fastp, merges overlapping pairs, and filters reads based on length and quality. Crucially, it performs taxonomic screening using two complementary tools:

   • ECMSD (Efficient Comprehensive Mitochondrial Sequence Detector): A novel tool that aligns reads against a curated NCBI RefSeq mitochondrial database. It distinguishes authentic eukaryotic signals from NUMTs by filtering for breadth coverage (default >25%), thereby avoiding false positives from nuclear insertions.
   • Centrifuge: Used for detecting bacteria, archaea, and viruses. It was selected over Kraken2 for its lower memory footprint and ability to classify shorter fragments (down to 22 bp) via exact matching rather than k-mer counting.

2. Reference Module: Maps reads to one or more reference sequences (genomes or marker genes like COI) using BWA-MEM2 (with BWA-aln or minimap2 as alternatives). It includes a parallelized, optimized version of DeDup for PCR duplicate removal and mapDamage2 for damage profiling and base-quality rescaling. The module calculates breadth and depth coverage statistics and estimates library complexity using Preseq.

3. REVEAL Module (Read-based Estimation and visualization of Element Abundance and Loci): A novel component for quantifying copy-number variation of user-defined genomic features (e.g., transposable elements, genes). It maps reads to target sequences and normalizes depth coverage against single-copy genes (SCGs), which are automatically identified via BUSCO if not provided. This normalization allows for comparable abundance estimates across samples with varying sequencing depths and endogenous DNA content.

4. Summary Module: Aggregates outputs from the Read and Reference modules into interactive MultiQC HTML reports. These reports visualize read counts, taxonomic composition, coverage statistics, library complexity, and damage patterns, facilitating rapid quality assessment.

Key Contributions

• Integrated Workflow: pastForward is the first pipeline to combine user-friendly, configurable aDNA/hDNA processing with automated taxonomic screening (including NUMT filtering) and targeted CNV tracking in a single Snakemake workflow.
• ECMSD: Introduces a memory-efficient method for eukaryotic taxonomic screening that specifically addresses the NUMT artifact problem by leveraging mitochondrial breadth coverage, filling a gap left by general-purpose metagenomic classifiers that exclude non-human metazoans.
• REVEAL: Provides a method to longitudinally track the abundance of genetic features (such as transposable elements or gene duplications) across ancient, historical, and modern datasets by normalizing read depth against SCGs.
• Accessibility: The pipeline manages all dependencies via Conda, automatically configures necessary databases (Centrifuge and ECMSD indices), and runs on standard desktop hardware, making it accessible to researchers without extensive bioinformatics infrastructure.

Results
Two case studies validated the pipeline:

1. Canine Domestication (aDNA): Applied to a genomic time series of dogs (Neolithic to modern) and wolves. pastForward successfully identified Canis lupus as the dominant taxon, filtering out NUMTs. The REVEAL module confirmed the known increase in AMY2B (starch-digesting enzyme) copy number during domestication, showing a shift from a single copy in wolves/Neolithic dogs to approximately 10 copies in modern breeds.
2. Drosophila melanogaster (hDNA): Applied to museum specimens spanning ~200 years. The pipeline detected low-level human contamination in 19th-century samples and successfully tracked the invasion of the retrotransposon opus. REVEAL visualized the absence of opus in 1800s samples and its continuous presence from 1933 onward, recapitulating previously published invasion timelines.

Significance
The authors position pastForward as a tool that shifts the focus from technical data processing to biological inquiry. By integrating taxonomic screening, damage assessment, and CNV tracking into a single, reproducible workflow, it enables researchers to efficiently process large-scale historical and ancient genomics datasets. The pipeline is particularly significant for "museomics," allowing natural history collections to be analyzed in-house to track evolutionary changes, such as gene amplification and transposable element invasions, directly over time. The authors emphasize that the tool is designed to be flexible and configurable, accommodating diverse sample types from bone to pinned insects, thereby facilitating the analysis of vast, previously underutilized specimen collections.