New module: modkit/entropy

modules/nf-core/modkit/entropy/environment.yml

@@ -0,0 +1,7 @@
+---
+# yaml-language-server: $schema=https://raw.githubusercontent.com/nf-core/modules/master/modules/environment-schema.json
+channels:
+  - conda-forge
+  - bioconda
+dependencies:
+  - "bioconda::ont-modkit=0.6.1"

modules/nf-core/modkit/entropy/main.nf

@@ -0,0 +1,61 @@
+process MODKIT_ENTROPY {
+    tag "$meta.id"
+    label 'process_high'
+
+    conda "${moduleDir}/environment.yml"
+    container "${ workflow.containerEngine in ['singularity', 'apptainer'] && !task.ext.singularity_pull_docker_container ?
+        'https://depot.galaxyproject.org/singularity/ont-modkit:0.6.1--hcdda2d0_0':
+        'quay.io/biocontainers/ont-modkit:0.6.1--hcdda2d0_0' }"
+
+    input:
+    // stageAs '?/*' prevents filename collisions when multiple BAMs from the
+    // same sample (e.g. technical replicates) are passed to a single run.
+    tuple val(meta),  path(bams, stageAs: "in/?/*"), path(bais, stageAs: "in/?/*")
+    tuple val(meta2), path(fasta), path(fai)
+    tuple val(meta3), path(regions)
+
+    output:
+    tuple val(meta), path("*.bed")                                  , emit: bed           , optional: true
+    tuple val(meta), path("entropy_regions/*.bed")                  , emit: regions_bed   , optional: true
+    tuple val(meta), path("entropy_regions/*.bedgraph")             , emit: bedgraph      , optional: true
+    tuple val(meta), path("*.log")                                  , emit: log           , optional: true
+    tuple val("${task.process}"), val('modkit'), eval("modkit --version | sed 's/modkit //'"), emit: versions_modkit, topic: versions
+
+    when:
+    task.ext.when == null || task.ext.when
+
+    script:
+    def args     = task.ext.args ?: ''
+    def prefix   = task.ext.prefix ?: "${meta.id}"
+    def bam_args = bams instanceof List ? bams.collect { "--in-bam ${it}" }.join(' ') : "--in-bam ${bams}"
+    // modkit entropy's --out-bed expects a FILE without --regions, and a DIRECTORY with --regions
+    def out_arg  = regions ? "--regions ${regions} --out-bed entropy_regions --prefix ${prefix}" : "--out-bed ${prefix}.bed"
+    def mkdir    = regions ? "mkdir -p entropy_regions" : ""
+    """
+    ${mkdir}
+
+    modkit \\
+        entropy \\
+        $args \\
+        --threads ${task.cpus} \\
+        --ref ${fasta} \\
+        ${out_arg} \\
+        ${bam_args}
+    """
+
+    stub:
+    def prefix = task.ext.prefix ?: "${meta.id}"
+    // With --regions, modkit writes a directory containing <prefix>_regions.bed
+    // and <prefix>_windows.bedgraph instead of the single BED written otherwise.
+    if (regions) {
+        """
+        mkdir -p entropy_regions
+        touch entropy_regions/${prefix}_regions.bed
+        touch entropy_regions/${prefix}_windows.bedgraph
+        """
+    } else {
+        """
+        touch ${prefix}.bed
+        """
+    }
+}

modules/nf-core/modkit/entropy/meta.yml

@@ -0,0 +1,159 @@
+name: modkit_entropy
+description: |
+  Calculate methylation entropy over genomic windows from one or more mod-BAMs.
+  Entropy is a per-window measure of the diversity of methylation patterns
+  across reads covering the window and is complementary to mean methylation.
+  When `--regions` is supplied (via `ext.args` or the `regions` input), modkit
+  writes `<prefix>_regions.bed` and `<prefix>_windows.bedgraph` into a
+  directory; otherwise a single genome-wide BED file is produced.
+keywords:
+  - modkit
+  - methylation
+  - entropy
+  - nanopore
+  - ont
+  - modbam
+tools:
+  - "modkit":
+      description: A bioinformatics tool for working with modified bases in Oxford Nanopore
+        sequencing data.
+      homepage: https://github.com/nanoporetech/modkit
+      documentation: https://nanoporetech.github.io/modkit/
+      tool_dev_url: https://github.com/nanoporetech/modkit
+      licence:
+        - "Oxford Nanopore Technologies PLC. Public License Version 1.0"
+      identifier: ""
+input:
+  - - meta:
+        type: map
+        description: |
+          Groovy Map containing sample information
+          e.g. `[ id:'sample1' ]`. The output inherits this meta.
+    - bams:
+        type: list
+        description: |
+          One or more modBAM files to compute entropy from. Passing more than
+          one BAM aggregates counts across them (multi-sample entropy).
+        pattern: "*.{bam,cram}"
+        ontologies:
+          - edam: http://edamontology.org/format_2572
+    - bais:
+        type: list
+        description: |
+          BAM indices (`.bai` or `.csi`) matching each input BAM, one-to-one.
+        pattern: "*.{bai,csi}"
+        ontologies: []
+  - - meta2:
+        type: map
+        description: |
+          Groovy Map containing reference information
+          e.g. `[ id:'mm10' ]`.
+    - fasta:
+        type: file
+        description: Reference FASTA the BAM was aligned to.
+        pattern: "*.{fa,fasta,fna}"
+        ontologies:
+          - edam: http://edamontology.org/format_1929
+    - fai:
+        type: file
+        description: Samtools FASTA index for `fasta`.
+        pattern: "*.fai"
+        ontologies:
+          - edam: http://edamontology.org/format_3475
+  - - meta3:
+        type: map
+        description: |
+          Groovy Map containing region information
+          e.g. `[ id:'promoters' ]`. May be `[[], []]` to skip.
+    - regions:
+        type: file
+        description: |
+          Optional BED file of regions over which to compute per-region
+          descriptive statistics. When provided, modkit writes
+          `<prefix>_regions.bed` and `<prefix>_windows.bedgraph` into a
+          directory instead of a single genome-wide BED.
+        pattern: "*.{bed,bed.gz}"
+        ontologies:
+          - edam: http://edamontology.org/format_3003
+output:
+  bed:
+    - - meta:
+          type: map
+          description: |
+            Groovy Map containing sample information
+            e.g. `[ id:'sample1' ]`.
+      - "*.bed":
+          type: file
+          description: |
+            Genome-wide entropy BED file, produced when `regions` is not
+            provided. One row per window.
+          pattern: "*.bed"
+          ontologies:
+            - edam: http://edamontology.org/format_3003
+  regions_bed:
+    - - meta:
+          type: map
+          description: |
+            Groovy Map containing sample information
+            e.g. `[ id:'sample1' ]`.
+      - "entropy_regions/*.bed":
+          type: file
+          description: |
+            Per-region descriptive statistics, produced when `regions` is
+            provided (`entropy_regions/<prefix>_regions.bed`). One row per
+            region with mean/median/min/max entropy and window counts.
+          pattern: "*.bed"
+          ontologies:
+            - edam: http://edamontology.org/format_3003
+  bedgraph:
+    - - meta:
+          type: map
+          description: |
+            Groovy Map containing sample information
+            e.g. `[ id:'sample1' ]`.
+      - "entropy_regions/*.bedgraph":
+          type: file
+          description: |
+            Per-window entropy, produced when `regions` is provided
+            (`entropy_regions/<prefix>_windows.bedgraph`).
+          pattern: "*.bedgraph"
+          ontologies:
+            - edam: http://edamontology.org/format_3583
+  log:
+    - - meta:
+          type: map
+          description: |
+            Groovy Map containing sample information
+            e.g. `[ id:'sample1' ]`.
+      - "*.log":
+          type: file
+          description: |
+            Optional modkit debug log (only emitted when `--log-filepath
+            <name>.log` is passed via `ext.args`).
+          pattern: "*.log"
+          ontologies: []
+  versions_modkit:
+    - - ${task.process}:
+          type: string
+          description: The name of the process
+      - modkit:
+          type: string
+          description: The name of the tool
+      - modkit --version | sed 's/modkit //':
+          type: eval
+          description: The expression to obtain the version of the tool
+topics:
+  versions:
+    - - ${task.process}:
+          type: string
+          description: The name of the process
+      - modkit:
+          type: string
+          description: The name of the tool
+      - modkit --version | sed 's/modkit //':
+          type: eval
+          description: The expression to obtain the version of the tool
+authors:
+  - "@sahuno"
+maintainers:
+  - "@sahuno"

modules/nf-core/modkit/entropy/tests/main.nf.test

@@ -0,0 +1,167 @@
+nextflow_process {
+
+    name "Test Process MODKIT_ENTROPY"
+    script "../main.nf"
+    process "MODKIT_ENTROPY"
+    config "./nextflow.config"
+
+    tag "modules"
+    tag "modules_nfcore"
+    tag "modkit"
+    tag "modkit/entropy"
+
+    // modkit entropy's floating-point columns are not byte-deterministic across
+    // CPU counts / architectures, so the output files are asserted on structure
+    // (header line, column count) rather than md5, and only the versions topic
+    // is snapshotted.
+
+    test("homo sapiens - nanopore modbam - cpg") {
+
+        when {
+            process {
+                """
+                input[0] = [
+                    [ id: 'test' ],
+                    [ file(params.modules_testdata_base_path + 'genomics/homo_sapiens/nanopore/bam/test.sorted.phased.bam', checkIfExists: true) ],
+                    [ file(params.modules_testdata_base_path + 'genomics/homo_sapiens/nanopore/bam/test.sorted.phased.bam.bai', checkIfExists: true) ]
+                ]
+                input[1] = [
+                    [ id: 'genome' ],
+                    file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.fasta', checkIfExists: true),
+                    file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.fasta.fai', checkIfExists: true)
+                ]
+                input[2] = [[],[]]
+                """
+            }
+        }
+
+        then {
+            def bed_lines = path(process.out.bed[0][1]).readLines()
+            assertAll(
+                { assert process.success },
+                { assert process.out.bed },
+                // no --regions: a single genome-wide BED, no entropy_regions/ dir
+                { assert !process.out.regions_bed },
+                { assert !process.out.bedgraph },
+                { assert bed_lines[0] == '#chrom\tstart\tend\tentropy\tstrand\tnum_reads' },
+                { assert bed_lines.size() > 1 },
+                // every data row carries all six columns and a parseable entropy
+                { assert bed_lines[1..-1].every { it.split('\t').size() == 6 } },
+                { assert bed_lines[1..-1].every { it.split('\t')[3].isNumber() } },
+                { assert snapshot([versions: process.out.versions_modkit]).match() }
+            )
+        }
+    }
+
+    test("homo sapiens - nanopore modbam - cpg - regions") {
+
+        when {
+            process {
+                """
+                input[0] = [
+                    [ id: 'test' ],
+                    [ file(params.modules_testdata_base_path + 'genomics/homo_sapiens/nanopore/bam/test.sorted.phased.bam', checkIfExists: true) ],
+                    [ file(params.modules_testdata_base_path + 'genomics/homo_sapiens/nanopore/bam/test.sorted.phased.bam.bai', checkIfExists: true) ]
+                ]
+                input[1] = [
+                    [ id: 'genome' ],
+                    file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.fasta', checkIfExists: true),
+                    file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.fasta.fai', checkIfExists: true)
+                ]
+                input[2] = [
+                    [ id: 'regions' ],
+                    file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.bed', checkIfExists: true)
+                ]
+                """
+            }
+        }
+
+        then {
+            def regions_lines   = path(process.out.regions_bed[0][1]).readLines()
+            def bedgraph_lines  = path(process.out.bedgraph[0][1]).readLines()
+            assertAll(
+                { assert process.success },
+                // with --regions modkit writes a directory instead of a single BED
+                { assert !process.out.bed },
+                { assert process.out.regions_bed },
+                { assert process.out.bedgraph },
+                { assert file(process.out.regions_bed[0][1]).name == 'test_regions.bed' },
+                { assert file(process.out.bedgraph[0][1]).name == 'test_windows.bedgraph' },
+                { assert regions_lines[0].startsWith('chrom\tstart\tend\tregion_name\tmean_entropy') },
+                { assert regions_lines.size() > 1 },
+                // the region name column echoes the interval from the input BED
+                { assert regions_lines[1].split('\t')[3] == 'chr22:0-40001' },
+                { assert bedgraph_lines[0] == '#chrom\tstart\tend\tentropy\tstrand\tnum_reads' },
+                { assert bedgraph_lines.size() > 1 },
+                { assert snapshot([versions: process.out.versions_modkit]).match() }
+            )
+        }
+    }
+
+    test("homo sapiens - nanopore modbam - cpg - stub") {
+
+        options "-stub"
+
+        when {
+            process {
+                """
+                input[0] = [
+                    [ id: 'test' ],
+                    [ file(params.modules_testdata_base_path + 'genomics/homo_sapiens/nanopore/bam/test.sorted.phased.bam', checkIfExists: true) ],
+                    [ file(params.modules_testdata_base_path + 'genomics/homo_sapiens/nanopore/bam/test.sorted.phased.bam.bai', checkIfExists: true) ]
+                ]
+                input[1] = [
+                    [ id: 'genome' ],
+                    file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.fasta', checkIfExists: true),
+                    file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.fasta.fai', checkIfExists: true)
+                ]
+                input[2] = [[],[]]
+                """
+            }
+        }
+
+        then {
+            assertAll(
+                { assert process.success },
+                { assert snapshot(process.out).match() }
+            )
+        }
+    }
+
+    test("homo sapiens - nanopore modbam - cpg - regions - stub") {
+
+        options "-stub"
+
+        when {
+            process {
+                """
+                input[0] = [
+                    [ id: 'test' ],
+                    [ file(params.modules_testdata_base_path + 'genomics/homo_sapiens/nanopore/bam/test.sorted.phased.bam', checkIfExists: true) ],
+                    [ file(params.modules_testdata_base_path + 'genomics/homo_sapiens/nanopore/bam/test.sorted.phased.bam.bai', checkIfExists: true) ]
+                ]
+                input[1] = [
+                    [ id: 'genome' ],
+                    file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.fasta', checkIfExists: true),
+                    file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.fasta.fai', checkIfExists: true)
+                ]
+                input[2] = [
+                    [ id: 'regions' ],
+                    file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.bed', checkIfExists: true)
+                ]
+                """
+            }
+        }
+
+        then {
+            assertAll(
+                { assert process.success },
+                // the stub mirrors the --regions layout, not the single-BED one
+                { assert !process.out.bed },
+                { assert process.out.regions_bed },
+                { assert process.out.bedgraph },
+                { assert snapshot(process.out).match() }
+            )
+        }
+    }
+}