04 May 2026

Short-read Methylation Pipeline Using Nextflow

Giang Nguyen

We evaluated both BSBolt and Rastair for methylation calling. BSBolt provides accurate and fast bisulfite sequencing alignments and methylation calls, outperforming Bismark, BSSeeker2, BISCUIT, and BWA-Meth based on alignment accuracy and methylation calling accuracy. Rastair achieves F1 scores exceeding 0.99 for datasets above 30x depth while processing a 30x depth file in under 30 minutes with 32 CPU cores.

References:

nf-core/methylation is widely regarded as a leading pipeline for methylation calling for short-read sequencing, backed by a large and active community. However, it does come with several limitations:

Reputation bias: Widely adopted tools are often retained as defaults, even when they may no longer offer the best performance. This can lead to unnecessary resource consumption and suboptimal efficiency in certain steps.
Legacy features: Some integrated tools have not been actively maintained for years, resulting in performance gaps. Despite this, they are often preserved for backward compatibility, increasing codebase complexity and the maintenance burden when introducing new features.
Customization constraints: The pipeline is designed for general-purpose and research-oriented use. In industrial production settings, workflows typically need to be streamlined—retaining only essential steps and achieving near-linear scalability with resources. In such cases, a more minimal and purpose-built approach is often more effective.

Therefore, we developed the nf-short-read-methylation specifically to support large-scale methylation analysis from bisulfite sequencing data (BSBolt) and TAPs sequencing data (Rastair).

Architecture

Key Features

Two Methylation Pathways: BSBolt (traditional bisulfite) and Rastair (TAPS-based)
Multiple Input Formats: FASTQ (full pipeline), BAM, and CRAM (skip alignment, auto-convert)
Quality Control: Fastp
M-bias Analysis: Automatic M-bias calculation and trimming optimization
Flexible Deduplication: Samtools MarkDup (BSBolt), GATK MarkDuplicates (Rastair)
Cross-sample Aggregation: Methylation matrix generation for comparative analysis
CRAM Compression: Built-in CRAM→BAM conversion for efficient re-analysis

Methylation Callers

BSBolt

BiSulfite Bolt is a bisulfite sequencing analysis platform that provides accurate and fast bisulfite sequencing alignments and methylation calls. It outperforms Bismark, BSSeeker2, BISCUIT, and BWA-Meth based on alignment accuracy and methylation calling accuracy.

Key advantages:

Fast, accurate bisulfite-aware alignment
Supports both directional and undirectional libraries
CGmap and BedGraph output formats
Cross-sample matrix aggregation for comparative analysis

Rastair

Rastair is an integrated software toolkit for simultaneous SNP detection and methylation calling from mC→T sequencing data (such as TAPS+ and Illumina's 5-Base chemistries). It combines machine-learning-based variant detection with genotype-aware methylation estimation.

Key advantages:

F1 scores exceeding 0.99 for datasets above 30x depth
Processes a 30x depth file in under 30 minutes with 32 CPU cores
Integrated variant and methylation calling
Standard-compliant outputs in VCF, BAM, and BED formats
Identifies positions where variants disrupt or create CpG sites

Benchmark

BSBolt Performance

According to the BiSulfite Bolt publication, BSBolt outperforms existing bisulfite alignment tools:

Tool	Alignment Accuracy	Methylation Calling
BSBolt	Highest	Most accurate
Bismark	Lower	Good
BSSeeker2	Lower	Good
BISCUIT	Moderate	Moderate
BWA-Meth	Lower	Lower

Rastair Performance

According to the Rastair publication on NA12878 benchmark datasets:

F1 Score: >0.99 for datasets above 30x depth
Processing Time: <30 minutes for 30x depth file with 32 CPU cores
GPU Acceleration: ~2x faster with GPU available
Additional Detection: Reports ~500,000 additional positions where SNPs create "de-novo" CpGs

Quick Start

Configuration for Primary Use Case

Default configuration uses:

BSBolt as methylation caller (set taps: true for Rastair)
GRCh38 reference genome
FASTP quality filtering
Automatic deduplication
Per-cytosine methylation reports

1
pixi run nextflow run main.nf -profile docker -resume

Prepare a Samplesheet

Create a CSV samplesheet with your input. The pipeline supports three input modes:

- FASTQ Input (Full Pipeline)

1
sample,lane,fastq_1,fastq_2
2
sample1,L001,/path/to/sample1_R1.fastq.gz,/path/to/sample1_R2.fastq.gz
3
sample2,L001,/path/to/sample2_R1.fastq.gz,/path/to/sample2_R2.fastq.gz

- BAM Input (Skip Alignment)

1
sample,lane,bam,bai
2
sample1,L001,/path/to/sample1.bam,/path/to/sample1.bam.bai

- CRAM Input (Skip Alignment + Auto-Convert)

1
sample,lane,cram,crai
2
sample1,L001,/path/to/sample1.cram,/path/to/sample1.cram.crai

CRAM Benefits:

Compressed input: CRAM files are ~4x smaller than BAM
Faster pipeline: Skip alignment step when re-running methylation calling
Automatic conversion: CRAM→BAM conversion integrated into pipeline

Samplesheet Columns:

sample: Sample identifier
lane: Sequencing lane (optional, defaults to L001)
FASTQ mode: fastq_1, fastq_2 (gzipped FASTQ files)
BAM mode: bam, bai (aligned BAM + index)
CRAM mode: cram, crai (compressed alignment + index)

Run the Pipeline

It will automatically detect the input format (FASTQ, BAM, or CRAM) to run the appropriate steps:

1
nextflow run main.nf \
2
  --input samplesheet.csv \
3
  --profile docker \
4
  -resume

Advanced Options

1
# Run with Rastair (TAPS-based methylation)
2
nextflow run main.nf \
3
  --input samplesheet.csv \
4
  --taps true \
5
  --profile docker \
6
  -resume
7

8
# Run with custom trim parameters (Rastair)
9
nextflow run main.nf \
10
  --input samplesheet.csv \
11
  --taps true \
12
  --trim_OT 10 \
13
  --trim_OB 10 \
14
  --profile docker \
15
  -resume
16

17
# BSBolt with pre-built index
18
nextflow run main.nf \
19
  --input samplesheet.csv \
20
  --bsbolt_index /path/to/bsbolt_index \
21
  --profile docker \
22
  -resume
23

24
# Custom reference genome
25
nextflow run main.nf \
26
  --input samplesheet.csv \
27
  --reference /path/to/reference.fa \
28
  --profile docker \
29
  -resume

For test mode with sample data:

1
nextflow run main.nf -profile docker,test -resume

View Results

Output files will be generated in the results/ directory. File structure depends on pipeline mode:

-BSBolt Results:

results/alignment/*.bam - Aligned BAM files
results/bsbolt/methylation_calls/*.cgmap.gz - CGmap format methylation calls
results/bsbolt/methylation_calls/*.bedGraph.gz - BedGraph format for visualization
results/bsbolt/aggregate_matrix/*_matrix.txt - Cross-sample methylation matrix
results/deduplicated/*.bam - Deduplicated BAM files

-Rastair Results:

results/rastair/mbias/ - M-bias calculation results and plots
results/rastair/call/*.txt - Methylation calls
results/rastair/methylkit/*.methylkit.txt.gz - MethylKit format for R analysis

-Common output files: