modkit BigWig Generation Failure Using "modkit bedmethyl tobigwig" Despite Correct Sorting

After creating the bedmethyl file from modbam using the following command:

#######################################
## 3. MODIFIED-BASE BAM TO BEDMETHYL ##
#######################################

## CPG CONTEXT ##
# Extract 5mC methylation at CpG sites: ignore 5hmC calls (--ignore h option)
modkit pileup iPSC_MNCdLS1.sorted.bam iPSC_MNCdLS1_5mC_cpgcontext.bedMethyl --cpg --ignore h -t 30 --ref /home/concology/Dorado/reference_genome/hg38/GRCh38.primary_assembly.genome.fa

I want to try to get a big-wig using the "modkit bedmethyl tobigwig" function ([https://nanoporetech.github.io/modkit/intro_bedmethyl_merge.html#convert-bedmethyl-to-bigwig])):

#1. Convert bedMethyl file for 5mC ALL CONTEXT to BigWig:
sort -k1,1 -k2,2n iPSC_MNCdLS1_5mC_allcontext.bedMethyl > iPSC_MNCdLS1_5mC_allcontext.sorted.bedMethyl
modkit bedmethyl tobigwig -t 30 --negative-strand-values --sizes chrom.sizes --mod-codes m iPSC_MNCdLS1_5mC_allcontext.sorted.bedMethyl iPSC_MNCdLS1_5mC_allcontext.bigWig

The problem is that, despite sorting exactly as the command told me, I get this error:

> Error! Input bedGraph not sorted by chromosome. Sort with sort -k1,1 -k2,2n.

How can I solve the problem and then get the BigWig?

Thanks a lot in advance

Mar 12 '25 23:03 lucagu1965

Hello @lucagu1965.

That shouldn't happen, could you share the input file with me? If you send me an email at art.rand [at] nanoporetech.com, we can work out how to share it (assuming you're willing). You shouldn't need to sort the output bedMethyl (it should be sorted when it's written).

Could you try streaming the output of pileup through bedmethyl tobigwig?

modkit pileup \
  iPSC_MNCdLS1.sorted.bam \
  - \
  --cpg --ignore h -t 30 \
  --ref /home/concology/Dorado/reference_genome/hg38/GRCh38.primary_assembly.genome.fa \
  tee >(modkit bedmethyl tobigwig --negative-strand-values --sizes chrom.sizes --mod-codes m - iPSC_MNCdLS1_5mC_allcontext.bigWig) \
  > iPSC_MNCdLS1_5mC_cpgcontext.bedMethyl

It shouldn't really matter, but if this works it might help me debug it. Thanks.

Mar 13 '25 22:03 ArtRand

Hello ArtRand,

Many thanks for your help.

Sure, I'm sending the code for the pipeline I used below, and here is a link where you can download one of the bedMethyl files:

https://www.transfernow.net/dl/20250316ULiWl21i

####################################################
#### WGS-5mC (Nanopore) - iPSCs & G12 (MNCdLS1) ####
####################################################

## DOWNLOAD REFERENCE GENOME ##
cd /home/concology/Dorado/reference_genome/hg38
wget https://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_47/GRCh38.primary_assembly.genome.fa.gz
gunzip GRCh38.primary_assembly.genome.fa.gz


## DORADO INSTALLATION ##
cd /home/concology/Dorado
wget https://cdn.oxfordnanoportal.com/software/analysis/dorado-0.9.1-linux-x64.tar.gz
tar -xzvf dorado-0.9.1-linux-x64.tar.gz
echo 'export PATH="$PATH:/home/concology/Dorado/dorado-0.9.1-linux-x64/bin"' >> ~/.bashrc # Insertion into the path
source ~/.bashrc # Applying changes
echo $PATH # Check
tail -n 10 ~/.bashrc # Check


## MODKIT INSTALLATION ##
cd /home/concology/modkit
curl https://sh.rustup.rs -sSf | sh  # Cargo Installation
source "$HOME/.cargo/env"            # Insert Cargo into Path
cargo install --git https://github.com/nanoporetech/modkit.git  # InstallModKit Installation
which modkit                # See where modkit is located
modkit --version           # See the modkit version


# Generate FASTA index to create .fai file [DO ONCE]
conda activate samtools
samtools faidx /home/concology/Dorado/reference_genome/hg38/GRCh38.primary_assembly.genome.fa

# Creating chrom.sizes file from FASTA index [DO ONCE]
cut -f1,2 /home/concology/Dorado/reference_genome/hg38/GRCh38.primary_assembly.genome.fa.fai > chrom.sizes



##############################################
#### 1. BASE-CALLING & ALIGNMENT (DORADO) ####
##############################################

## BASECALLING EXECUTION (DORADO) - SUPER-ACCURATE - 5mC and 5hmC ALL CONTEXTS ##

# dorado basecaller --device cuda:all --emit-fastq sup --modified-bases 5mC_5hmC /mnt/big_storage/1_2024_10_30/iPSC_MNCdLS1/20241030_1215_P2S-01869-A_PBA38952_c66b0889/pod5 > iPSC_MNCdLS1.fastq
# dorado basecaller --device cuda:all --emit-fastq sup --modified-bases 5mC_5hmC --batchsize 16 /mnt/big_storage/1_2024_10_30/iPSC_MNCdLS1/20241030_1215_P2S-01869-A_PBA38952_c66b0889/pod5 > iPSC_MNCdLS1.fastq

dorado basecaller --device cuda:all --emit-sam sup \
  --modified-bases 5mC_5hmC \
  --reference /home/concology/Dorado/reference_genome/hg38/GRCh38.primary_assembly.genome.fa \
  /mnt/big_storage/1_2024_10_30/iPSC_MNCdLS1/20241030_1215_P2S-01869-A_PBA38952_c66b0889/pod5 \
  > iPSC_MNCdLS1.sam



##################################
#### 2. CONVERTING SAM TO BAM ####
##################################

conda activate samtools
cd /mnt/big_storage/1_2024_10_30/iPSC_MNCdLS1/20241030_1215_P2S-01869-A_PBA38952_c66b0889/Output

# 1. Converting SAM to BAM
samtools view -bS iPSC_MNCdLS1.sam > iPSC_MNCdLS1.bam

# 2. BAM sorting using multi-core
samtools sort -@ 30 iPSC_MNCdLS1.bam -o iPSC_MNCdLS1.sorted.bam

# 3. Indexing of the ordered BAM
samtools index iPSC_MNCdLS1.sorted.bam

# Execute in the same line, consecutively
# samtools view -bS iPSC_MNCdLS1.sam > iPSC_MNCdLS1.bam && samtools sort -@ 30 iPSC_MNCdLS1.bam -o iPSC_MNCdLS1.sorted.bam && samtools index iPSC_MNCdLS1.sorted.bam

# Alignment Statistics Report
samtools flagstat iPSC_MNCdLS1.sorted.bam



###########################################
#### 3. MODIFIED-BASE BAM TO BEDMETHYL ####
###########################################

## CPG CONTEXT ##

# Extract 5mC methylation at CpG sites: ignore 5hmC calls (--ignore h option)
modkit pileup iPSC_MNCdLS1.sorted.bam iPSC_MNCdLS1_5mC_cpgcontext.bedMethyl --cpg --ignore h -t 30 --ref /home/concology/Dorado/reference_genome/hg38/GRCh38.primary_assembly.genome.fa

# Extract 5hmC methylation at CpG sites: ignore 5mC calls (--ignore m option)
modkit pileup iPSC_MNCdLS1.sorted.bam iPSC_MNCdLS1_5hmC_cpgcontext.bedMethyl --cpg --ignore m -t 30 --ref /home/concology/Dorado/reference_genome/hg38/GRCh38.primary_assembly.genome.fa

## ALL CONTEXT ##

# Extract 5mC methylation in all contexts: ignore 5hmC calls
modkit pileup iPSC_MNCdLS1.sorted.bam iPSC_MNCdLS1_5mC_allcontext.bedMethyl --ignore h -t 30 --ref /home/concology/Dorado/reference_genome/hg38/GRCh38.primary_assembly.genome.fa

# Extract 5hmC methylation in all contexts: ignore 5mC calls
modkit pileup iPSC_MNCdLS1.sorted.bam iPSC_MNCdLS1_5hmC_allcontext.bedMethyl --ignore m -t 30 --ref /home/concology/Dorado/reference_genome/hg38/GRCh38.primary_assembly.genome.fa



#############################################
#### 4. BIG-WIG CREATION FROM BED-METHYL ####
#############################################

# Sorting of BedMethyl Files
sort -k1,1 -k2,2n iPSC_MNCdLS1_5mC_allcontext.bedMethyl > iPSC_MNCdLS1_5mC_allcontext.sorted.bedMethyl
sort -k1,1 -k2,2n iPSC_MNCdLS1_5hmC_allcontext.bedMethyl > iPSC_MNCdLS1_5hmC_allcontext.sorted.bedMethyl
sort -k1,1 -k2,2n iPSC_MNCdLS1_5mC_cpgcontext.bedMethyl > iPSC_MNCdLS1_5mC_cpgcontext.sorted.bedMethyl
sort -k1,1 -k2,2n iPSC_MNCdLS1_5hmC_cpgcontext.bedMethyl > iPSC_MNCdLS1_5hmC_cpgcontext.sorted.bedMethyl
ls


# 1. Convert bedMethyl file for 5mC ALL CONTEXT to BigWig
modkit bedmethyl tobigwig -t 30 --negative-strand-values --sizes chrom.sizes --mod-codes m iPSC_MNCdLS1_5mC_allcontext.sorted.bedMethyl iPSC_MNCdLS1_5mC_allcontext.bigWig 

# 2. Convert bedMethyl file for 5hmC ALL CONTEXT to BigWig
modkit bedmethyl tobigwig -t 30 --negative-strand-values --sizes chrom.sizes --mod-codes h iPSC_MNCdLS1_5hmC_allcontext.sorted.bedMethyl iPSC_MNCdLS1_5hmC_allcontext.bigWig

# 3. Convert bedMethyl file for 5mC CPG CONTEXT to BigWig
modkit bedmethyl tobigwig -t 30 --negative-strand-values --sizes chrom.sizes --mod-codes m iPSC_MNCdLS1_5mC_cpgcontext.sorted.bedMethyl iPSC_MNCdLS1_5mC_cpgcontext.bigWig

# 4. Convert bedMethyl file for 5hmC CPG CONTEXT to BigWig
modkit bedmethyl tobigwig -t 30 --negative-strand-values --sizes chrom.sizes --mod-codes h iPSC_MNCdLS1_5hmC_cpgcontext.sorted.bedMethyl iPSC_MNCdLS1_5hmC_cpgcontext.bigWig

Mar 16 '25 00:03 lucagu1965

Hello @lucagu1965,

Sorry for the delay, but I think I've tracked down the problem and made a fix. Could you try the build attached here?

Thanks.

Aug 15 '25 01:08 ArtRand

Hello @lucagu1965,

Sorry for the delay, but I think I've tracked down the problem and made a fix. You should not need this step:

sort -k1,1 -k2,2n iPSC_MNCdLS1_5mC_allcontext.bedMethyl

anymore. Could you try the build attached here?

Thanks.

modkit_dev442e6ba_u16_x86_64.tar.gz

Aug 15 '25 01:08 ArtRand

This change has been included in the latest release candidate. If you're still encountering a problem, please re-open this issue.

Sep 18 '25 16:09 ArtRand