REMAG

Recovery of Eukaryotic Metagenome-Assembled Genomes using contrastive learning. A specialized metagenomic binning tool designed for recovering high-quality eukaryotic genomes from mixed prokaryotic-eukaryotic samples.

Index

• Installation
• Quick Start
• Usage
• Common Options
• How It Works
• Output
• Requirements
• Acknowledgments
• License
• Citation

Installation

Conda (recommended)

conda create -n remag -c bioconda -c conda-forge remag
conda activate remag

PyPI

Install miniprot separately first:

conda create -n remag python=3.9
conda activate remag
conda install -c bioconda miniprot
pip install remag

From source

conda create -n remag python=3.9
conda activate remag

git clone https://github.com/danielzmbp/remag.git
cd remag
conda install -c bioconda miniprot
pip install .

Docker

docker pull danielzmbp/remag:latest

Optional plotting dependencies

pip install "remag[plotting]"

GPU acceleration

REMAG automatically uses NVIDIA CUDA or Apple Silicon Metal (mps) when available, otherwise CPU. No extra flag is required.

NVIDIA GPUs (Linux)

Install CUDA-enabled PyTorch before REMAG. This example uses CUDA 12.8; use the PyTorch selector to choose a build matching your GPU and driver.

conda create -n remag-gpu -c conda-forge -c bioconda python=3.11 pip miniprot
conda activate remag-gpu
python -m pip install torch --index-url https://download.pytorch.org/whl/cu128
python -m pip install remag
python -c "import torch; print('CUDA:', torch.cuda.is_available())"

Apple Silicon (macOS)

In your native Apple Silicon REMAG environment:

python -m pip install torch remag
python -c "import torch; print('MPS:', torch.backends.mps.is_available())"

The relevant GPU check should report True.

Quick Start

Conda

remag contigs.fasta -c alignments.bam

Docker

docker run --rm -v "$(pwd):/data" danielzmbp/remag:latest \
  /data/contigs.fasta -c /data/alignments.bam -o /data/output

Singularity

singularity build remag.sif docker://danielzmbp/remag:latest
singularity run --bind "$(pwd):/data" remag.sif \
  /data/contigs.fasta -c /data/alignments.bam -o /data/output

Usage

Command line interface

After installation, you can use REMAG via the command line:

# Basic usage
remag contigs.fasta -c alignments.bam

# With explicit output directory
remag contigs.fasta -c alignments.bam -o output_directory

# Multiple samples
remag contigs.fasta -c sample1.bam -c sample2.bam

# Multiple samples using shell-expanded globs
remag contigs.fasta -c samples/*.bam

# Using precomputed coverage tables (one TSV per sample)
remag contigs.fasta -c sample1.tsv -c sample2.tsv

# Using interval coverage files (one COV/bedGraph per sample)
remag contigs.fasta -c sample1.bam.cov.gz -c sample2.bam.cov.gz

# Only run eukaryotic filtering (skip binning)
remag contigs.fasta --filter-only

# Keep intermediate files
remag contigs.fasta -c alignments.bam -k

# Recompute results in an existing output directory
remag contigs.fasta -c alignments.bam -o output_directory --force

Quote each coverage path containing spaces, for example -c "sample 1.bam" -c "sample 2.bam".

During feature generation, REMAG rejects repeated identifiers among retained contigs, using the first whitespace-separated token in each FASTA header.

Coverage is optional. To use sequence composition alone:

remag contigs.fasta -o output_directory

Reusing or replacing results

REMAG reuses saved results when their recorded model, settings and inputs match the current run. Verified embeddings can be reused without retraining, even without saved model weights. Clustering checks its settings and inputs and restarts rescue from saved pre-rescue assignments. Missing provenance or changed inputs/settings require a new output directory or --force. When reusing results, REMAG prints:

Existing REMAG results found. Reusing available outputs. Use --force to recompute.

Use --force to remove recognized REMAG outputs and rerun the requested workflow. This clears cached features, model weights, embeddings, gene mappings, previous bins, filtering outputs, temporary miniprot directories, and logs. Other files in the directory are left in place. If a supplied input would be removed, REMAG stops before deleting anything. --force also refuses a symlinked bins directory or temporary directories containing symbolic links; use a different output directory.

Use --force or a different output directory when changing inputs or analysis settings. --force --filter-only clears previous results and runs filtering only. Saved umap_coordinates.csv and umap_plot.pdf are also cleared; regenerate plots after the run.

Python module mode

python -m remag contigs.fasta -c alignments.bam

Precomputed Coverage Formats

Precomputed coverage files are supported as an alternative to BAM/CRAM. Use one file per sample. REMAG auto-detects the precomputed coverage layout from the file contents.

Contig-level TSV/TXT

Use this format when you have one average coverage value per contig.

• Column 1: contig ID
• Last column: coverage value for that contig
• No header row

Example:

contig_1	12.4
contig_2	3.8
contig_3	0.0

Contig-level TSV/TXT input provides contig-level coverage only. REMAG cannot infer fragment-specific coverage for augmented fragments from this layout, so every fragment from the same contig gets the same coverage value.

Interval COV/bedGraph

Use this format when you have per-interval coverage, such as output with run-length encoded coverage blocks.

• Column 1: contig ID
• Column 2: interval start, 0-based inclusive
• Column 3: interval end, 0-based exclusive
• Column 4: coverage value for that interval
• No header row

Example:

ctg0	0	256	1
ctg0	256	859	2
ctg0	859	861	1
ctg1	0	73	6

Supported extensions include .cov, .cov.gz, .bedgraph, .bedgraph.gz, .bg, and .bg.gz. Missing intervals are treated as zero coverage. Interval coverage is used to compute fragment-specific mean and standard deviation coverage, so it supports REMAG's random augmentations.

Do not mix BAM/CRAM inputs with precomputed coverage inputs in the same run.

BAM/CRAM coverage uses --cores readers and processes one alignment file at a time. More readers can improve speed but increase memory use, especially for CRAM. Queued work is bounded and can be reduced further with --coverage-batch-size. Memory or worker failures stop the run instead of producing zero coverage. If memory is exhausted, reduce --cores or allocate more RAM.

Common Options

• -c, --coverage: one or more BAM, CRAM, contig-level TSV/TXT, or interval COV/bedGraph coverage inputs
• -o, --output: output directory; defaults to remag_output next to the input FASTA
• -k, --keep-intermediate: retain features, model weights, encoder embeddings, and other optional intermediate files; embeddings.csv is saved without this flag
• --force: remove existing REMAG outputs and recompute results
• --filter-only: stop after eukaryotic filtering and write filtered FASTA output
• --save-filtered-contigs: also write classified contigs rejected by the eukaryotic filter
• --skip-bacterial-filter: disable the HyenaDNA filter
• --graph-min-contig-length: minimum length entering the Leiden graph; defaults to 1,000 bp with one coverage file, otherwise 3,000 bp
• --min-bin-size: minimum bin size written to FASTA after rescue; defaults to 500,000 bp

Use remag -h for a quick reference and remag --help for the full CLI, including training, clustering, filtering, and rescue options.

REMAG defaults to a minimum contig and training-fragment length of 1,000 bp for every run. With exactly one coverage file, the initial graph also includes contigs from 1,000 bp. With multiple coverage files or no coverage, the graph starts at 3,000 bp, while admitted shorter contigs remain available for rescue. Use --min-contig-length and --graph-min-contig-length to override these settings. A default graph cutoff is raised to the admission cutoff if needed. Effective settings are logged and saved in params.json; reusing outputs requires compatible admission and graph cutoffs. Use --force or a new output directory when changing them.

For example, the default settings with two coverage files are:

remag contigs.fasta -c sample1.bam -c sample2.bam

All admitted contigs share the same embedding space. The graph cutoff must be at least the admission cutoff; it does not change filtering or training-fragment length. After fragmented bins are merged, unassigned graph contigs and admitted shorter contigs are recruited together in embedding row order. --skip-rescue skips both merging and recruitment.

The k-NN graph uses 15 neighbors, and HyenaDNA filtering is enabled. Use --leiden-k-neighbors and --skip-bacterial-filter to adjust these settings. Multiple coverage files lower the default base learning rate from 0.005 to 0.0005; an explicitly supplied learning rate is preserved.

How It Works

REMAG recovers eukaryotic bins with a multi-stage pipeline:

1. Eukaryotic filtering: By default, REMAG filters contigs with the integrated HyenaDNA classifier. This step can be disabled with --skip-bacterial-filter.
2. Feature extraction: REMAG combines 4-mer composition with optional multi-sample coverage data. Contigs are augmented into fragments for training when their lengths permit; contigs longer than 50 kb receive augmentations from each half.
3. Contrastive learning: A Siamese network trained with Barlow Twins learns embeddings that place fragments from the same contig close together.
4. Core gene annotation: miniprot maps eukaryotic single-copy core genes to support clustering and quality checks. Multiple non-overlapping matches to a marker family within a contig are reported as separate copies; clustering and rescue still count that family once per contig.
5. Greedy clustering and rescue: REMAG applies greedy Leiden clustering using cosine-similarity edge weights across multiple resolutions, merges fragmented bins, then recruits all eligible unassigned contigs in one pass. Bins below the minimum size remain available as seeds until the final size filter.

With coverage, the network uses composition widths 136 → 128 → 64 and coverage widths p → 64 → 32, where p is the number of coverage features. Both encoders use BatchNorm, LeakyReLU and dropout 0.05. Coverage hidden widths are fixed for every positive p. Fusion uses bidirectional four-head cross-attention, three fusion scales, nonlinear compression, adaptive dropout and a final nonlinear MLP. Each modality uses a fixed 0.5 direct-versus-attended gate and a learned residual weight.

flowchart LR
    K["Composition: 136 → 128 → 64"] --> F["Cross-attention + fixed 0.5 gates<br/>Learned residual + multiscale nonlinear fusion"]
    C["Coverage: p → 64 → 32"] --> F
    F --> E["D-dimensional embedding<br/>default D = 256"]
    E --> G["Graph and clustering"]
    E --> P["Training projection: D → 128 → 512"]
    P --> L["Barlow Twins loss"]

The projection head is used for training; graph construction uses the pre-projection D-dimensional embedding. A dedicated generator, seeded with --random-seed, is shared by the DataLoader and its sampler to make batch order independent of model initialization.

Without coverage, the composition encoder uses 136 → 256 → 128, followed by a projection to D and a D → 512 → 512 training head.

Rescue uses length-weighted centroids and a cosine similarity threshold of 0.70. Centroids stay fixed during merging and are recomputed once before recruitment. Each candidate tries its nearest bin only. With the default marker limits, a marker-bearing candidate must leave at most 5% duplicated marker families and increase duplication by strictly less than 5 percentage points. Counts update after each accepted contig; marker-free candidates have no marker veto. A whole-bin merge has one additional allowance: at similarity of at least 0.95, a target already above the duplication ceiling may accept a merge that does not worsen its duplication. This allowance does not apply to individual-contig recruitment. --rescue-max-total-duplication controls the ceiling (capped at 10% for whole-bin merges), and --rescue-max-duplication-increase controls the strict increase limit.

Whole-bin merges also require at least one original clustering-graph edge between their current contigs, including contigs absorbed by earlier merges. If the nearest bin has no graph connection, the merge is rejected without trying another target. This safeguard also applies to the non-worsening marker allowance; individual-contig recruitment keeps its existing rule.

Output

All runs retain params.json with the effective settings, including admission and graph cutoffs and the rescue settings. Matching reruns preserve the original record.

Trained runs also save training_provenance.json with model dimensions, training/filter settings and input/output identities used to validate reuse. siamese_model.pt contains the PyTorch state dictionary; provenance is stored separately.

Binning outputs

• bins/: FASTA files for bins meeting --min-bin-size
• bins.csv: Final assignments for contigs in saved bins; excludes noise and bins below the minimum size
• pre_rescue.csv: Initial assignments for all admitted contigs; shorter contigs outside the graph are initially noise
• rescue_assignments.csv: Assignments before the final size filter, including noise and small bins, with an exported flag
• clustering_provenance.json: Clustering/rescue settings and input fingerprints used to check cache compatibility
• knn_graph_edges.csv: Original k-NN graph edge list used for clustering and the merge safeguard
• knn_graph_contigs.csv: Ordered graph contig IDs; edge indices refer to these rows, which may be a subset of embeddings.csv
• knn_graph_stats.json: Graph construction statistics and fingerprints for cache validation
• embeddings.csv: Embeddings for original contigs, including those that do not enter a saved bin
• fragments.pkl: Fragment sequences and coordinates used by the pipeline; currently written even without -k
• remag.log: Detailed log file
• gene_contig_mappings.json: Cached gene mappings and positions; older caches without positions are regenerated when annotations are loaded
• core_gene_duplication_results.json: Core gene duplication analysis for the final saved bins, after rescue and minimum-size filtering. Includes within_contig_duplications, identifying marker families with multiple locations on the same contig.

Binning stops with an error if required gene annotation fails. A successful search with no accepted matches still reports zero detected genes.

The original graph is retained even without --keep-intermediate and reused with compatible pre-rescue assignments. Older rescue versions or missing/changed graph evidence require a new output directory or --force; existing outputs are preserved until recomputation is requested.

Additional outputs with -k / --keep-intermediate

• siamese_model.pt: Trained Siamese neural network model
• kmer_embeddings.csv: K-mer encoder embeddings (before fusion)
• coverage_embeddings.csv: Coverage encoder embeddings, when coverage features are present (before fusion)
• features.csv: Extracted k-mer and coverage features
• temp_gene_mapping/: Miniprot files used to generate gene mappings
• temp_miniprot/: Per-bin miniprot files, when the fallback duplication-check path runs

Filtering outputs

• *_hyenadna_classification.tsv: Predictions for contigs meeting --min-contig-length
• *_eukaryotic_filtered.fasta: Contigs retained by the HyenaDNA filter
• *_non_eukaryotic.fasta: Rejected contigs, when --save-filtered-contigs is used and at least one contig is rejected

--filter-only stops after this stage and does not produce binning outputs. Use it with filtering enabled; --skip-bacterial-filter disables that stage. Contigs below the minimum length are excluded from classification and the rejected-contig FASTA.

Filtering uses the July 2026 single 4 kb HyenaDNA model for every eligible contig, including those shorter than 4 kb. Contigs up to 4,096 bp use one window. Longer contigs use up to 16 evenly spread 4,096 bp windows, including both endpoints. The window count is min(16, 1 + ceil((length - 4096) / 2048)); positions use Python's rounding. REMAG retains the entire original contig when the arithmetic mean of window probabilities is at least 0.30. There is no random resampling or prediction early stopping.

The tokenizer uses left padding and 4,096 tokens including CLS/SEP, so each window contributes at most 4,094 biological bases. Inference uppercases input and preserves unknown characters as UNK; retained and rejected FASTA records preserve the original headers, sequence content and order. The classification table retains its columns, reports probabilities at full precision, and always reports resampled=False; confidence is the native 0.5 window-agreement diagnostic, not a filtering rule.

Filtering uses fixed batches of 64 eligible contigs and up to 256 windows per inference call. The legacy --hyenadna-batch-size option is accepted but other values are ignored with a warning. Automatic inference uses CUDA when available, otherwise CPU. All probabilities come directly from the selected device, including those near the cutoff. MPS is not used because equivalence has not been established.

New runs record the policy, checkpoint SHA256, cutoff and execution settings in params.json. Reusing an older filtered FASTA does not rerun prediction; REMAG warns when its policy differs or is unknown and preserves the original provenance. Use --force or a new output directory to apply the current policy.

If no contigs pass the classifier, REMAG saves an empty filtered FASTA and finishes successfully without training or producing bins. Classifier initialization failures still fall back to the input; consult remag.log for details.

Visualization

The plotting helper currently runs from a REMAG source checkout. From the repository root, install the plotting dependencies and generate UMAP plots from embeddings.csv and bins.csv:

pip install -e ".[plotting]"
python scripts/plot_features.py --features output_directory/embeddings.csv --clusters output_directory/bins.csv --output output_directory

• umap_coordinates.csv: UMAP projections for visualization
• umap_plot.pdf: UMAP visualization plot with cluster assignments

Requirements

• Python 3.9+
• miniprot on PATH for core gene analysis during binning; --filter-only does not require it. The Dockerfile includes miniprot; install it separately for PyPI and source installations.
• Plotting extras are optional: pip install "remag[plotting]"

The package includes the pre-trained July 2026 HyenaDNA classifier (4,096-token context; SHA256 39d05fcecfdb0849fcf1822ffb0c816ee119e0436a19e0284e33e42f280f3689) and the protein reference database used for core gene mapping.

Acknowledgments

The integrated HyenaDNA classifier uses a pre-trained genomic foundation model:

• Repository: HazyResearch/hyena-dna
• Paper: Nguyen E, Poli M, Faizi M, et al. HyenaDNA: Long-Range Genomic Sequence Modeling at Single Nucleotide Resolution. NeurIPS 2023.

License

MIT License - see LICENSE file for details.

Citation

If you use REMAG in your research, please cite:

@article {G{\'o}mez-P{\'e}rez2026.03.05.709928,
	author = {G{\'o}mez-P{\'e}rez, Daniel and Raguideau, S{\'e}bastien and Warring, Sally and James, Robert and Hildebrand, Falk and Quince, Christopher},
	title = {REMAG: recovery of eukaryotic genomes from metagenomic data using contrastive learning},
	elocation-id = {2026.03.05.709928},
	year = {2026},
	doi = {10.64898/2026.03.05.709928},
	publisher = {Cold Spring Harbor Laboratory},
	URL = {https://www.biorxiv.org/content/early/2026/03/08/2026.03.05.709928},
	eprint = {https://www.biorxiv.org/content/early/2026/03/08/2026.03.05.709928.full.pdf},
	journal = {bioRxiv}
}