diff --git a/conf/igenomes.config b/conf/igenomes.config index 71ff9f59..db1f1391 100644 --- a/conf/igenomes.config +++ b/conf/igenomes.config @@ -24,6 +24,8 @@ params.genomes = [ vep_species : "homo_sapiens", savana_contigs : "https://raw.githubusercontent.com/cortes-ciriano-lab/savana/main/example/contigs.chr.hg38.txt", savana_g1000_vcf : "1000g_hg38", + dmr_cpg_islands_bed : "https://raw.githubusercontent.com/IntGenomicsLab/test-datasets/main/references/dmr/GRCh38.cpg_islands.bed.gz", + dmr_gencode_gene_bed : "https://raw.githubusercontent.com/IntGenomicsLab/test-datasets/main/references/dmr/GRCh38.gencode_gene.bed.gz", vep_alphamissense : "https://storage.googleapis.com/dm_alphamissense/AlphaMissense_hg38.tsv.gz", vep_alphamissense_tbi : "https://g-608c0c.273595.03c0.data.globus.org/VEP_plugins/AlphaMissense_hg38.tsv.gz.tbi", // A dated release rather than the rolling vcf_GRCh38/clinvar.vcf.gz, whose VCF and @@ -54,6 +56,8 @@ params.genomes = [ vep_species : "homo_sapiens_gca009914755v4", savana_contigs : "https://raw.githubusercontent.com/IntGenomicsLab/test-datasets/main/references/savana/contigs.chr.chm13.txt", savana_g1000_vcf : "1000g_t2t", + dmr_cpg_islands_bed : "https://raw.githubusercontent.com/IntGenomicsLab/test-datasets/main/references/dmr/CHM13.cpg_islands.bed.gz", + dmr_gencode_gene_bed : "https://raw.githubusercontent.com/IntGenomicsLab/test-datasets/main/references/dmr/CHM13.gencode_gene.bed.gz", vep_alphamissense_aa : "https://g-608c0c.273595.03c0.data.globus.org/VEP_plugins/alphamissense_protein_v2023_uniprot-2026_03.tsv.gz", vep_alphamissense_aa_tbi : "https://g-608c0c.273595.03c0.data.globus.org/VEP_plugins/alphamissense_protein_v2023_uniprot-2026_03.tsv.gz.tbi", // Pinned to a release rather than current_variation/, which moves at every Ensembl release diff --git a/conf/modules.config b/conf/modules.config index 8b20e4fd..7aafefb0 100644 --- a/conf/modules.config +++ b/conf/modules.config @@ -345,6 +345,40 @@ process { ] } + // + // SUBWORKFLOW: DMR -- intermediate steps (tabix-indexing, CpG-island restriction, the raw + // modkit dmr pair output) are not published; DMR_NEAREST_GENE's annotated bed is the final + // deliverable. + // + withName: '.*:DMR:HTSLIB_BGZIPTABIX' { + // meta.haplotype ('hp1'/'hp2') is stripped from meta again right after this step (both + // haplotypes are rejoined under the plain sample meta), so it has to be in the filename + // instead -- otherwise hp1 and hp2 both stage as .bed.gz into DMR_HAPLOTYPE_REGIONS. + ext.prefix = { "${meta.id}_${meta.haplotype}" } + ext.args2 = '-p bed' + publishDir = [ + enabled: false + ] + } + withName: '.*:DMR:DMR_HAPLOTYPE_REGIONS' { + publishDir = [ + enabled: false + ] + } + withName: '.*:DMR:MODKIT_DMR' { + ext.args = '--base C' + publishDir = [ + enabled: false + ] + } + withName: '.*:DMR:DMR_NEAREST_GENE' { + publishDir = [ + path: { "${params.outdir}/${meta.id}/methylation/${meta.type}/dmr" }, + mode: params.publish_dir_mode, + saveAs: { filename -> filename.equals('versions.yml') ? null : filename } + ] + } + withName: '.*:FIBERTOOLSRS_PREDICTM6A' { ext.args = { [ diff --git a/modules.json b/modules.json index ab8746ac..c50d0b08 100644 --- a/modules.json +++ b/modules.json @@ -2,220 +2,321 @@ "name": "IntGenomicsLab/lrsomatic", "homePage": "https://github.com/IntGenomicsLab/lrsomatic", "repos": { + "https://github.com/YannVRB/modules.git": { + "modules": { + "nf-core": { + "modkit/dmr": { + "branch": "add-modkit-dmr", + "git_sha": "fbbd19f488dcc8cb0318e92c895a431b791a237f", + "installed_by": [ + "modules" + ], + "patch": "modules/nf-core/modkit/dmr/modkit-dmr.diff" + } + } + } + }, "https://github.com/nf-core/modules.git": { "modules": { "nf-core": { "ascat": { "branch": "master", "git_sha": "e753770db613ce014b3c4bc94f6cba443427b726", - "installed_by": ["modules"], + "installed_by": [ + "modules" + ], "patch": "modules/nf-core/ascat/ascat.diff" }, "bcftools/annotate": { "branch": "master", "git_sha": "3d9c2f4beaa4f62b3f006928fd9095a496d1e5a8", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "bcftools/concat": { "branch": "master", "git_sha": "6383d8fe58f9498eecd5aa303e71a4a932d1e9f6", - "installed_by": ["modules", "vcf_gather_bcftools"] + "installed_by": [ + "modules", + "vcf_gather_bcftools" + ] }, "bcftools/isec": { "branch": "master", "git_sha": "3b2c3559699a7bca6a7c2b220695a072e030e17d", - "installed_by": ["modules"], + "installed_by": [ + "modules" + ], "patch": "modules/nf-core/bcftools/isec/bcftools-isec.diff" }, "bcftools/merge": { "branch": "master", "git_sha": "3d9c2f4beaa4f62b3f006928fd9095a496d1e5a8", - "installed_by": ["modules"], + "installed_by": [ + "modules" + ], "patch": "modules/nf-core/bcftools/merge/bcftools-merge.diff" }, "bcftools/norm": { "branch": "master", "git_sha": "6383d8fe58f9498eecd5aa303e71a4a932d1e9f6", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "bcftools/query": { "branch": "master", "git_sha": "6383d8fe58f9498eecd5aa303e71a4a932d1e9f6", - "installed_by": ["modules"], + "installed_by": [ + "modules" + ], "patch": "modules/nf-core/bcftools/query/bcftools-query.diff" }, "bcftools/sort": { "branch": "master", "git_sha": "6383d8fe58f9498eecd5aa303e71a4a932d1e9f6", - "installed_by": ["modules", "vcf_gather_bcftools"], + "installed_by": [ + "modules", + "vcf_gather_bcftools" + ], "patch": "modules/nf-core/bcftools/sort/bcftools-sort.diff" }, "bcftools/view": { "branch": "master", "git_sha": "feef37435aea56816adf4b3bde1fc76aac327a8d", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "deepvariant/callvariants": { "branch": "master", "git_sha": "f2b138ee1d91f67d31c187317d7e83e429bf0309", - "installed_by": ["deepvariant"], + "installed_by": [ + "deepvariant" + ], "patch": "modules/nf-core/deepvariant/callvariants/deepvariant-callvariants.diff" }, "deepvariant/makeexamples": { "branch": "master", "git_sha": "f2b138ee1d91f67d31c187317d7e83e429bf0309", - "installed_by": ["deepvariant"], + "installed_by": [ + "deepvariant" + ], "patch": "modules/nf-core/deepvariant/makeexamples/deepvariant-makeexamples.diff" }, "deepvariant/postprocessvariants": { "branch": "master", "git_sha": "f2b138ee1d91f67d31c187317d7e83e429bf0309", - "installed_by": ["deepvariant"], + "installed_by": [ + "deepvariant" + ], "patch": "modules/nf-core/deepvariant/postprocessvariants/deepvariant-postprocessvariants.diff" }, "ensemblvep/download": { "branch": "master", "git_sha": "90cdd21fd96ccbdb3bc90797ca69570d18391055", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "ensemblvep/vep": { "branch": "master", "git_sha": "890fdcff71928fc1470d3e669d4c430c8c770297", - "installed_by": ["modules"], + "installed_by": [ + "modules" + ], "patch": "modules/nf-core/ensemblvep/vep/ensemblvep-vep.diff" }, + "htslib/bgziptabix": { + "branch": "master", + "git_sha": "cbe6025183336010a9b716e463f897b0c0f70f8a", + "installed_by": [ + "modules" + ] + }, "longphase/haplotag": { "branch": "master", "git_sha": "b8d30a43f33aee3148b0e9e9f00587984a4ac195", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "longphase/phase": { "branch": "master", "git_sha": "b8d30a43f33aee3148b0e9e9f00587984a4ac195", - "installed_by": ["modules"], + "installed_by": [ + "modules" + ], "patch": "modules/nf-core/longphase/phase/longphase-phase.diff" }, "minimap2/align": { "branch": "master", "git_sha": "5c9f8d5b7671237c906abadc9ff732b301ca15ca", - "installed_by": ["modules"], + "installed_by": [ + "modules" + ], "patch": "modules/nf-core/minimap2/align/minimap2-align.diff" }, "minimap2/index": { "branch": "master", "git_sha": "14980f759266eec42dac401fcafeb83d6c957b41", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "modkit/pileup": { "branch": "master", "git_sha": "3d81317a30d1016b533982d6b84df07713ae520a", - "installed_by": ["modules"], + "installed_by": [ + "modules" + ], "patch": "modules/nf-core/modkit/pileup/modkit-pileup.diff" }, "mosdepth": { "branch": "master", "git_sha": "6832b69ef7f98c54876d6436360b6b945370c615", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "multiqc": { "branch": "master", "git_sha": "98403d15b0e50edae1f3fec5eae5e24982f1fade", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "nanoplot": { "branch": "master", "git_sha": "682f789f93070bd047868300dd018faf3d434e7c", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "pigz/uncompress": { "branch": "master", "git_sha": "f84336b7fa91a65aa61d215b8c109fbb8e4b4ac6", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "samtools/cat": { "branch": "master", "git_sha": "f9edc59be2fe25bb6fc73ca4dfc0d28246f2a2d6", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "samtools/faidx": { "branch": "master", "git_sha": "b2e78932ef01165fd85829513eaca29eff8e640a", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "samtools/flagstat": { "branch": "master", "git_sha": "1d2fbdcbca677bbe8da0f9d0d2bb7c02f2cab1c9", - "installed_by": ["bam_stats_samtools"] + "installed_by": [ + "bam_stats_samtools" + ] }, "samtools/idxstats": { "branch": "master", "git_sha": "1d2fbdcbca677bbe8da0f9d0d2bb7c02f2cab1c9", - "installed_by": ["bam_stats_samtools"] + "installed_by": [ + "bam_stats_samtools" + ] }, "samtools/index": { "branch": "master", "git_sha": "1d2fbdcbca677bbe8da0f9d0d2bb7c02f2cab1c9", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "samtools/merge": { "branch": "master", "git_sha": "6d46786420b4d7bc88eba026eb389c0c5535d120", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "samtools/stats": { "branch": "master", "git_sha": "fe93fde0845f907fc91ad7cc7d797930408824df", - "installed_by": ["bam_stats_samtools"], + "installed_by": [ + "bam_stats_samtools" + ], "patch": "modules/nf-core/samtools/stats/samtools-stats.diff" }, "savana/classify": { "branch": "master", "git_sha": "7e612460bcb5aaf9c6f146cfb6475064c2d062ee", - "installed_by": ["modules"], + "installed_by": [ + "modules" + ], "patch": "modules/nf-core/savana/classify/savana-classify.diff" }, "savana/cna": { "branch": "master", "git_sha": "cc9117a41906d96832594d98f23da7805a87eca2", - "installed_by": ["modules"], + "installed_by": [ + "modules" + ], "patch": "modules/nf-core/savana/cna/savana-cna.diff" }, "savana/run": { "branch": "master", "git_sha": "8226aa5748298be6c00661fe611e231c64f88e71", - "installed_by": ["modules"], + "installed_by": [ + "modules" + ], "patch": "modules/nf-core/savana/run/savana-run.diff" }, "savana/to": { "branch": "master", "git_sha": "ad25ca6c4d27f57740e739886000439d8dff2061", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "severus": { "branch": "master", "git_sha": "4dd9d8439a429c7ee566e0e2347f76ddeef27e66", - "installed_by": ["modules"], + "installed_by": [ + "modules" + ], "patch": "modules/nf-core/severus/severus.diff" }, "untar": { "branch": "master", "git_sha": "447f7bc0fa41dfc2400c8cad4c0291880dc060cf", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "unzip": { "branch": "master", "git_sha": "4dd9d8439a429c7ee566e0e2347f76ddeef27e66", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "wget": { "branch": "master", "git_sha": "41dfa3f7c0ffabb96a6a813fe321c6d1cc5b6e46", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] }, "whatshap/stats": { "branch": "master", "git_sha": "bfab71f4d68c1aaff09335a3433e7b2836918b2a", - "installed_by": ["modules"] + "installed_by": [ + "modules" + ] } } }, @@ -224,31 +325,41 @@ "bam_stats_samtools": { "branch": "master", "git_sha": "7ac6cbe7c17c2dad685da7f70496c8f48ea48687", - "installed_by": ["subworkflows"] + "installed_by": [ + "subworkflows" + ] }, "deepvariant": { "branch": "master", "git_sha": "f2b138ee1d91f67d31c187317d7e83e429bf0309", - "installed_by": ["subworkflows"], + "installed_by": [ + "subworkflows" + ], "patch": "subworkflows/nf-core/deepvariant/deepvariant.diff" }, "utils_nextflow_pipeline": { "branch": "master", "git_sha": "c2b22d85f30a706a3073387f30380704fcae013b", - "installed_by": ["subworkflows"] + "installed_by": [ + "subworkflows" + ] }, "utils_nfcore_pipeline": { "branch": "master", "git_sha": "a3fb7351b1fdb2b1de282b765816bbea190e86a8", - "installed_by": ["subworkflows"] + "installed_by": [ + "subworkflows" + ] }, "utils_nfschema_plugin": { "branch": "master", "git_sha": "a7b27fd25bfa8dcc07d299e88bd790585901a436", - "installed_by": ["subworkflows"] + "installed_by": [ + "subworkflows" + ] } } } } } -} +} \ No newline at end of file diff --git a/modules/local/dmr/haplotype_regions/environment.yml b/modules/local/dmr/haplotype_regions/environment.yml new file mode 100644 index 00000000..8b02e7a2 --- /dev/null +++ b/modules/local/dmr/haplotype_regions/environment.yml @@ -0,0 +1,8 @@ +--- +# yaml-language-server: $schema=https://raw.githubusercontent.com/nf-core/modules/master/modules/environment-schema.json +channels: + - conda-forge + - bioconda +dependencies: + - bioconda::bedtools=2.31.1 + - bioconda::htslib=1.22.1 diff --git a/modules/local/dmr/haplotype_regions/main.nf b/modules/local/dmr/haplotype_regions/main.nf new file mode 100644 index 00000000..c149081d --- /dev/null +++ b/modules/local/dmr/haplotype_regions/main.nf @@ -0,0 +1,73 @@ +process DMR_HAPLOTYPE_REGIONS { + tag "$meta.id" + label 'process_single' + + conda "${moduleDir}/environment.yml" + // The plain bedtools:2.31.1 biocontainer has no tabix binary at all -- needed since + // haplotype_regions/main.nf queries the bgzip+tabix-indexed bedMethyl inputs directly. + // Reusing nf-core/modules' pints/caller container here (real, already in production use, + // confirmed via its own environment.yml to bundle both bedtools and htslib) rather than + // hand-constructing an unverifiable Wave/mulled tag; it carries unrelated pybedtools/pypints + // baggage this module doesn't use, but that's preferable to guessing a container hash. + container "${workflow.containerEngine in ['singularity', 'apptainer'] && !task.ext.singularity_pull_docker_container + ? 'https://community-cr-prod.seqera.io/docker/registry/v2/blobs/sha256/f1/f1a9e30012e1b41baf9acd1ff94e01161138d8aa17f4e97aa32f2dc4effafcd1/data' + : 'community.wave.seqera.io/library/pybedtools_bedtools_htslib_pip_pypints:39699b96998ec5f6'}" + + input: + // hp1/hp2 bedMethyl are not read for methylation values here, only for which CpG positions + // they cover -- restricting cpg_islands_bed to islands covered in BOTH haplotypes is what + // keeps modkit dmr pair from comparing real data on one haplotype against the other + // haplotype's absence of data at the same island. + tuple val(meta), path(hp1_bedmethyl), path(hp1_tbi), path(hp2_bedmethyl), path(hp2_tbi) + path cpg_islands_bed + tuple val(meta2), path(fai) + + output: + tuple val(meta), path("*.regions.bed"), emit: regions_bed + path "versions.yml" , emit: versions + + when: + task.ext.when == null || task.ext.when + + script: + def prefix = task.ext.prefix ?: "${meta.id}" + // Minimum valid coverage (bedMethyl column 10, Nvalid_cov) for a position to count as + // "covered" -- 1 preserves the previous any-depth behaviour; override via ext.args to + // require more. + def min_valid_coverage = task.ext.args ?: 1 + """ + cut -f1,2 ${fai} > genome.txt + + # Query directly via the bgzip+tabix index for positions inside the CpG islands, instead + # of decompressing/sorting every bedMethyl row genome-wide -- hp*_bedmethyl are already + # coordinate-sorted and tabix-indexed upstream, so this only reads the (small) subset of + # rows that can possibly matter here. + tabix -R ${cpg_islands_bed} ${hp1_bedmethyl} \\ + | awk -v min_cov=${min_valid_coverage} 'BEGIN{OFS="\\t"} \$10>=min_cov {print \$1,\$2,\$3}' \\ + | sort -k1,1 -k2,2n -u > hp1.cpg.bed + tabix -R ${cpg_islands_bed} ${hp2_bedmethyl} \\ + | awk -v min_cov=${min_valid_coverage} 'BEGIN{OFS="\\t"} \$10>=min_cov {print \$1,\$2,\$3}' \\ + | sort -k1,1 -k2,2n -u > hp2.cpg.bed + + bedtools intersect -u -a ${cpg_islands_bed} -b hp1.cpg.bed | sort -k1,1 -k2,2n > islands.hp1.bed + bedtools intersect -u -a ${cpg_islands_bed} -b hp2.cpg.bed | sort -k1,1 -k2,2n > islands.hp2.bed + bedtools intersect -u -a islands.hp1.bed -b islands.hp2.bed | bedtools sort -g genome.txt -i - > ${prefix}.regions.bed + + cat <<-END_VERSIONS > versions.yml + "${task.process}": + bedtools: \$(bedtools --version | sed 's/bedtools v//g') + tabix: \$(tabix --version 2>&1 | head -n1 | sed 's/tabix (htslib) //') + END_VERSIONS + """ + + stub: + def prefix = task.ext.prefix ?: "${meta.id}" + """ + touch ${prefix}.regions.bed + + cat <<-END_VERSIONS > versions.yml + "${task.process}": + bedtools: 2.31.1 + END_VERSIONS + """ +} diff --git a/modules/local/dmr/haplotype_regions/meta.yml b/modules/local/dmr/haplotype_regions/meta.yml new file mode 100644 index 00000000..c3f80db6 --- /dev/null +++ b/modules/local/dmr/haplotype_regions/meta.yml @@ -0,0 +1,76 @@ +# yaml-language-server: $schema=https://raw.githubusercontent.com/nf-core/modules/master/modules/meta-schema.json +name: "dmr_haplotype_regions" +description: | + Restrict a CpG-islands BED to islands with real read coverage in both of a sample's two + haplotype bedMethyl files, for use as modkit dmr pair's --regions-bed. Without this, an island + covered in one haplotype but not the other would be compared against absence of data rather than + a real methylation difference. +keywords: + - methylation + - dmr + - haplotype + - cpg islands + - bedtools +tools: + - "bedtools": + description: "A powerful toolset for genome arithmetic" + homepage: "https://bedtools.readthedocs.io/" + documentation: "https://bedtools.readthedocs.io/" + licence: ["GPL-2.0-only"] + +input: + - - meta: + type: map + description: | + Groovy Map containing sample information, e.g. `[ id:'sample1' ]` + - hp1_bedmethyl: + type: file + description: Bgzip-compressed bedMethyl for haplotype 1 + pattern: "*.bed.gz" + - hp1_tbi: + type: file + description: Tabix index for hp1_bedmethyl + pattern: "*.bed.gz.tbi" + - hp2_bedmethyl: + type: file + description: Bgzip-compressed bedMethyl for haplotype 2 + pattern: "*.bed.gz" + - hp2_tbi: + type: file + description: Tabix index for hp2_bedmethyl + pattern: "*.bed.gz.tbi" + - - cpg_islands_bed: + type: file + description: Genome-wide CpG islands BED to restrict + pattern: "*.bed" + - - meta2: + type: map + description: | + Groovy Map containing reference information, e.g. `[ id:'GRCh38' ]` + - fai: + type: file + description: FASTA index for the reference genome (used only for its chrom/length columns) + pattern: "*.fai" + +output: + regions_bed: + - meta: + type: map + description: | + Groovy Map containing sample information, e.g. `[ id:'sample1' ]` + - "*.regions.bed": + type: file + description: | + cpg_islands_bed restricted to islands covered in both haplotypes, sorted by + genome/coordinate order, ready to pass to modkit dmr pair's --regions-bed + pattern: "*.regions.bed" + versions: + - "versions.yml": + type: file + description: File containing software versions + pattern: "versions.yml" + +authors: + - "@YannVRB" +maintainers: + - "@YannVRB" diff --git a/modules/local/dmr/nearest_gene/environment.yml b/modules/local/dmr/nearest_gene/environment.yml new file mode 100644 index 00000000..45c307b0 --- /dev/null +++ b/modules/local/dmr/nearest_gene/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::bedtools=2.31.1 diff --git a/modules/local/dmr/nearest_gene/main.nf b/modules/local/dmr/nearest_gene/main.nf new file mode 100644 index 00000000..f2f4350b --- /dev/null +++ b/modules/local/dmr/nearest_gene/main.nf @@ -0,0 +1,47 @@ +process DMR_NEAREST_GENE { + tag "$meta.id" + label 'process_single' + + conda "${moduleDir}/environment.yml" + container "${workflow.containerEngine in ['singularity', 'apptainer'] && !task.ext.singularity_pull_docker_container + ? 'https://depot.galaxyproject.org/singularity/bedtools:2.31.1--hf5e1c6e_0' + : 'quay.io/biocontainers/bedtools:2.31.1--hf5e1c6e_0'}" + + input: + tuple val(meta), path(dmr_bed) + path gencode_gene_bed + + output: + tuple val(meta), path("*.nearestGene.bed"), emit: nearest_gene_bed + path "versions.yml" , emit: versions + + when: + task.ext.when == null || task.ext.when + + script: + // -D b: distance is signed relative to the gene's (B's) strand, which is meaningful -- + // DMR regions (A) are unstranded, so -D a's "distance relative to A's strand" is not. + def args = task.ext.args ?: '-D b -t first' + def prefix = task.ext.prefix ?: "${meta.id}" + """ + sort -k1,1 -k2,2n ${dmr_bed} > dmr.sorted.bed + + bedtools closest -a dmr.sorted.bed -b ${gencode_gene_bed} ${args} > ${prefix}.nearestGene.bed + + cat <<-END_VERSIONS > versions.yml + "${task.process}": + bedtools: \$(bedtools --version | sed 's/bedtools v//g') + END_VERSIONS + """ + + stub: + def prefix = task.ext.prefix ?: "${meta.id}" + """ + touch ${prefix}.nearestGene.bed + + cat <<-END_VERSIONS > versions.yml + "${task.process}": + bedtools: 2.31.1 + END_VERSIONS + """ +} diff --git a/modules/local/dmr/nearest_gene/meta.yml b/modules/local/dmr/nearest_gene/meta.yml new file mode 100644 index 00000000..39853526 --- /dev/null +++ b/modules/local/dmr/nearest_gene/meta.yml @@ -0,0 +1,51 @@ +# yaml-language-server: $schema=https://raw.githubusercontent.com/nf-core/modules/master/modules/meta-schema.json +name: "dmr_nearest_gene" +description: | + Sort a DMR bed (modkit dmr pair's output) by position and annotate each region with its nearest + gene from a GENCODE-derived gene BED, using bedtools closest. +keywords: + - methylation + - dmr + - gene annotation + - bedtools +tools: + - "bedtools": + description: "A powerful toolset for genome arithmetic" + homepage: "https://bedtools.readthedocs.io/" + documentation: "https://bedtools.readthedocs.io/" + licence: ["GPL-2.0-only"] + +input: + - - meta: + type: map + description: | + Groovy Map containing sample information, e.g. `[ id:'sample1' ]` + - dmr_bed: + type: file + description: modkit dmr pair's output BED + pattern: "*.bed" + - - gencode_gene_bed: + type: file + description: Genome-wide, position-sorted gene-level BED (chrom, start, end, gene_name, score, strand) + pattern: "*.bed" + +output: + nearest_gene_bed: + - meta: + type: map + description: | + Groovy Map containing sample information, e.g. `[ id:'sample1' ]` + - "*.nearestGene.bed": + type: file + description: dmr_bed with each region's nearest gene and signed distance appended + pattern: "*.nearestGene.bed" + versions: + - "versions.yml": + type: file + description: File containing software versions + pattern: "versions.yml" + +authors: + - "@YannVRB" +maintainers: + - "@YannVRB" diff --git a/modules/nf-core/htslib/bgziptabix/environment.yml b/modules/nf-core/htslib/bgziptabix/environment.yml new file mode 100644 index 00000000..ec62c057 --- /dev/null +++ b/modules/nf-core/htslib/bgziptabix/environment.yml @@ -0,0 +1,8 @@ +--- +# yaml-language-server: $schema=https://raw.githubusercontent.com/nf-core/modules/master/modules/environment-schema.json +channels: + - conda-forge + - bioconda +dependencies: + - bioconda::htslib=1.24 + - conda-forge::xz=5.8.3 diff --git a/modules/nf-core/htslib/bgziptabix/main.nf b/modules/nf-core/htslib/bgziptabix/main.nf new file mode 100644 index 00000000..573c05fc --- /dev/null +++ b/modules/nf-core/htslib/bgziptabix/main.nf @@ -0,0 +1,88 @@ +process HTSLIB_BGZIPTABIX { + tag "${meta.id}" + label 'process_low' + + conda "${moduleDir}/environment.yml" + container "${workflow.containerEngine in ['singularity', 'apptainer'] && !task.ext.singularity_pull_docker_container + ? 'https://community-cr-prod.seqera.io/docker/registry/v2/blobs/sha256/86/863ca0dbbba30c8367fa4fbd3fa3a84393532fb7b300a5c5c2e70f0dfc475bbf/data' + : 'community.wave.seqera.io/library/htslib_xz:32f2772a564b3cd2'}" + + input: + tuple val(meta), path(infile), path(infile_tbi), path(regions) + val action + val make_index + val out_ext + + output: + tuple val(meta), path("${outfile}"), emit: output + tuple val(meta), path("${outfile}.{tbi,csi}"), emit: index, optional: true + // all htslib tools have the same version, we use bgzip + tuple val("${task.process}"), val('htslib'), eval("bgzip --version | sed '1! d; s/bgzip (htslib) //'"), topic: versions, emit: versions_htslib + tuple val("${task.process}"), val('xz'), eval("xz --version | sed '1! d; s/xz (XZ Utils) //'"), topic: versions, emit: versions_xz + + when: + task.ext.when == null || task.ext.when + + script: + def allowed_actions = ["compress", "decompress"] + if (action !in allowed_actions) { + error("htslib/bgziptabix: Invalid action: ${action}. Allowed actions are: ${allowed_actions.join(', ')}") + } + + if (action == "decompress" && make_index) { + log.warn("htslib/bgziptabix: Cannot create index when decompressing. Ignoring make_index option.") + } + + def args = task.ext.args ?: '' + def args2 = task.ext.args2 ?: '' + prefix = task.ext.prefix ?: "${meta.id}" + outfile = action == "compress" ? (out_ext ? "${prefix}.${out_ext}.gz" : "${prefix}.gz") : (out_ext ? "${prefix}.${out_ext}" : "${prefix}") + + def compress_cmd = action == "compress" ? "bgzip -c ${args} -@ ${task.cpus}" : "cat" + def bgzip_cmd = action == "compress" ? "[ '\$(basename ${infile})' != '\$(basename ${outfile})' ] && ln -s ${infile} ${outfile}" : "bgzip -c -d ${args} -@ ${task.cpus} ${infile} > ${outfile}" + + def regions_arg = regions ? "-R ${regions}" : "" + def tabix_cmd = (make_index && !infile_tbi) ? "tabix -@ ${task.cpus} ${regions_arg} ${args2} -f ${outfile}" : "" + def link_tabix_cmd = make_index && infile_tbi ? "ln -s ${infile_tbi} ${outfile}.${infile_tbi.extension}" : "" + def uncompressed_cmd = action == "compress" ? "${compress_cmd} ${infile} > ${outfile}" : (infile.getName() == outfile ? "" : "ln -s ${infile} ${outfile}") + """ + ${link_tabix_cmd} + + FILE_TYPE=\$(htsfile ${infile}) + + case "\$FILE_TYPE" in + *BGZF-compressed*) + ${bgzip_cmd} ;; + *gzip-compressed*) + [ "\$(basename ${infile})" == "\$(basename ${outfile})" ] && echo "Input and output names cannot be the same" && exit 1 + bgzip -d -c -@ ${task.cpus} ${infile} | ${compress_cmd} > ${outfile} ;; + *bzip2-compressed*) + bzcat ${infile} | ${compress_cmd} > ${outfile} ;; + *XZ-compressed*) + xzcat ${infile} | ${compress_cmd} > ${outfile} ;; + *) + ${uncompressed_cmd} ;; + esac + + ${tabix_cmd} + """ + + stub: + def args = task.ext.args ?: '' + def args2 = task.ext.args2 ?: '' + prefix = task.ext.prefix ?: "${meta.id}" + outfile = action == "compress" ? (out_ext ? "${prefix}.${out_ext}.gz" : "${prefix}.gz") : (out_ext ? "${prefix}.${out_ext}" : "${prefix}") + + def touch_cmd = action == "compress" ? "echo | bgzip -c" : "echo" + def index_fmt = args2.contains('-C') ? 'csi' : 'tbi' + def tabix_cmd = make_index ? "touch ${outfile}.${index_fmt}" : "" + def link_tabix_cmd = make_index && infile_tbi ? "ln -s ${infile_tbi} ${outfile}.${infile_tbi.extension}" : "" + """ + echo ${args} + + ${touch_cmd} > ${outfile} + + ${tabix_cmd} + ${link_tabix_cmd} + """ +} diff --git a/modules/nf-core/htslib/bgziptabix/meta.yml b/modules/nf-core/htslib/bgziptabix/meta.yml new file mode 100644 index 00000000..4cdefd0e --- /dev/null +++ b/modules/nf-core/htslib/bgziptabix/meta.yml @@ -0,0 +1,125 @@ +name: "htslib_bgziptabix" +description: "Multi-purpose module to compress, decompress and index files using bgzip + and tabix." +keywords: + - compress + - decompress + - index + - bgzip + - tabix + - gzip + - bzip + - xz +tools: + - "htslib": + description: "C library for high-throughput sequencing data formats." + homepage: "http://www.htslib.org/" + documentation: "http://www.htslib.org/doc/" + tool_dev_url: "https://github.com/samtools/htslib" + doi: "10.1093/gigascience/giab007" + licence: + - "MIT" + identifier: biotools:htslib +input: + - - meta: + type: map + description: | + Groovy Map containing sample information + e.g. [ id:'sample1' ] + - infile: + type: file + description: Input file to compress or decompress + pattern: "*" + ontologies: [] + - infile_tbi: + type: file + description: Optional tabix index for the input file. + pattern: "*.{tbi,csi}" + ontologies: + - edam: http://edamontology.org/format_3616 # tabix + - regions: + type: file + description: Optional file of regions to extract (BED or chr:start-end format). + Only used when creating an index for the output file. + pattern: "*.{bed,txt,tsv}" + ontologies: + - edam: http://edamontology.org/format_3475 # TSV + - edam: http://edamontology.org/format_3003 # BED + - action: + type: string + description: Action to perform, either `compress` or `decompress` + - make_index: + type: boolean + description: Whether to create a tabix index for the output file; only used + if `action` is `compress` + - out_ext: + type: string + description: Output file extension without `.gz` suffix (for example `vcf`) +output: + output: + - - meta: + type: map + description: | + Groovy Map containing sample information + e.g. [ id:'sample1' ] + - ${outfile}: + type: file + description: Compressed or decompressed output file + pattern: "*" + ontologies: [] + index: + - - meta: + type: map + description: | + Groovy Map containing sample information + e.g. [ id:'sample1' ] + - ${outfile}.{tbi,csi}: + type: file + description: Tabix index file for the compressed output file + pattern: "*.{tbi,csi}" + ontologies: + - edam: http://edamontology.org/format_3616 # tabix + versions_htslib: + - - ${task.process}: + type: string + description: The name of the process + - htslib: + type: string + description: The name of the tool + - bgzip --version | sed '1! d; s/bgzip (htslib) //': + type: eval + description: The expression to obtain the version of the tool + versions_xz: + - - ${task.process}: + type: string + description: The name of the process + - xz: + type: string + description: The name of the tool + - xz --version | sed '1! d; s/xz (XZ Utils) //': + type: eval + description: The expression to obtain the version of the tool +topics: + versions: + - - ${task.process}: + type: string + description: The name of the process + - htslib: + type: string + description: The name of the tool + - bgzip --version | sed '1! d; s/bgzip (htslib) //': + type: eval + description: The expression to obtain the version of the tool + - - ${task.process}: + type: string + description: The name of the process + - xz: + type: string + description: The name of the tool + - xz --version | sed '1! d; s/xz (XZ Utils) //': + type: eval + description: The expression to obtain the version of the tool +authors: + - "@itrujnara" +maintainers: + - "@itrujnara" diff --git a/modules/nf-core/htslib/bgziptabix/tests/main.nf.test b/modules/nf-core/htslib/bgziptabix/tests/main.nf.test new file mode 100644 index 00000000..a7346506 --- /dev/null +++ b/modules/nf-core/htslib/bgziptabix/tests/main.nf.test @@ -0,0 +1,435 @@ +nextflow_process { + + name "Test Process HTSLIB_BGZIPTABIX" + script "../main.nf" + process "HTSLIB_BGZIPTABIX" + + tag "modules" + tag "modules_nfcore" + tag "htslib" + tag "htslib/bgziptabix" + + test("sarscov2 - vcf - decompress") { + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/sarscov2/illumina/vcf/test.vcf', checkIfExists: true), + [], + [] + ] + input[1] = 'decompress' // action + input[2] = false // make_index + input[3] = 'vcf' // out_ext + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(sanitizeOutput(process.out)).match() }, + { assert process.out.output.get(0).get(1).endsWith('.vcf') }, + { assert process.out.index.size() == 0 } + ) + } + + } + + test("sarscov2 - vcf - compress - index") { + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/sarscov2/illumina/vcf/test.vcf', checkIfExists: true), + [], + [] + ] + input[1] = 'compress' // action + input[2] = true // make_index + input[3] = 'vcf' // out_ext + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(sanitizeOutput(process.out)).match() }, + { assert process.out.output.get(0).get(1).endsWith('.vcf.gz') }, + { assert process.out.index.get(0).get(1).endsWith('.vcf.gz.tbi') } + ) + } + + } + + test("sarscov2 - vcf + regions - compress - index") { + when { + process { + """ + input[0] = [ + [ id:'example' ], + file(params.modules_testdata_base_path + 'genomics/sarscov2/illumina/vcf/test.vcf.gz', checkIfExists: true), + file(params.modules_testdata_base_path + 'genomics/sarscov2/illumina/vcf/test.vcf.gz.tbi', checkIfExists: true), + file('https://raw.githubusercontent.com/luisas/test-datasets/refs/heads/add-bedgraph-subset-illumina/data/genomics/sarscov2/illumina/bed/test.bed', checkIfExists: true) + ] + input[1] = 'compress' // action + input[2] = true // make_index + input[3] = 'vcf' // out_ext + """ + } + } + + then { + assertAll ( + { assert process.success }, + { assert snapshot( + sanitizeOutput(process.out), + path(process.out.output[0][1]).vcf.getVariantsMD5(), + ).match() } + ) + } + } + + test("sarscov2 - bgzip - decompress") { + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/sarscov2/illumina/vcf/test.vcf.gz', checkIfExists: true), + [], + [] + ] + input[1] = 'decompress' // action + input[2] = false // make_index + input[3] = 'vcf' // out_ext + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(sanitizeOutput(process.out)).match() }, + { assert process.out.output.get(0).get(1).endsWith('.vcf') }, + { assert process.out.index.size() == 0 } + ) + } + } + + test("sarscov2 - bgzip - compress - no index") { + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/sarscov2/illumina/vcf/test.vcf.gz', checkIfExists: true), + [], + [] + ] + input[1] = 'compress' // action + input[2] = false // make_index + input[3] = 'vcf' // out_ext + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(sanitizeOutput(process.out)).match() }, + { assert process.out.output.get(0).get(1).endsWith('.vcf.gz') }, + { assert process.out.index.size() == 0 } + ) + } + + } + + test("sarscov2 - gzip - decompress") { + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/sarscov2/illumina/fastq/test_1.fastq.gz', checkIfExists: true), + [], + [] + ] + input[1] = 'decompress' // action + input[2] = false // make_index + input[3] = 'fastq' // out_ext + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(sanitizeOutput(process.out)).match() }, + { assert process.out.output.get(0).get(1).endsWith('.fastq') }, + { assert process.out.index.size() == 0 } + ) + } + } + + test("sarscov2 - gzip - (re)compress - no index") { + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/sarscov2/illumina/fastq/test_1.fastq.gz', checkIfExists: true), + [], + [] + ] + input[1] = 'compress' // action + input[2] = false // make_index + input[3] = 'fastq' // out_ext + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(sanitizeOutput(process.out)).match() }, + { assert process.out.output.get(0).get(1).endsWith('.fastq.gz') }, + { assert process.out.index.size() == 0 } + ) + } + } + + test("sarscov2 - gzip - name clash") { + + when { + process { + """ + input[0] = [ + [ id:'test_1' ], + file(params.modules_testdata_base_path + 'genomics/sarscov2/illumina/fastq/test_1.fastq.gz', checkIfExists: true), + [], + [] + ] + input[1] = 'compress' // action + input[2] = false // make_index + input[3] = 'fastq' // out_ext + """ + } + } + + then { + assert process.failed + assertAll( + { assert process.errorReport.contains("Input and output names cannot be the same") } + ) + } + } + + test("metagenome - bz2 - decompress") { + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/prokaryotes/metagenome/rgi/card-data.tar.bz2', checkIfExists: true), + [], + [] + ] + input[1] = 'decompress' // action + input[2] = false // make_index + input[3] = 'tar' // out_ext + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(sanitizeOutput(process.out)).match() }, + { assert process.out.output.get(0).get(1).endsWith('.tar') }, + { assert process.out.index.size() == 0 } + ) + } + } + + test("metagenome - bz2 - (re)compress - no index") { + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/prokaryotes/metagenome/rgi/card-data.tar.bz2', checkIfExists: true), + [], + [] + ] + input[1] = 'compress' // action + input[2] = false // make_index + input[3] = 'tar' // out_ext + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(sanitizeOutput(process.out)).match() }, + { assert process.out.output.get(0).get(1).endsWith('.tar.gz') }, + { assert process.out.index.size() == 0 } + ) + } + } + + test("metagenome - xz - decompress") { + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/prokaryotes/metagenome/taxonomy/misc/taxa_sqlite.xz', checkIfExists: true), + [], + [] + ] + input[1] = 'decompress' // action + input[2] = false // make_index + input[3] = '' // out_ext + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot( + process.out, + process.out.findAll { key, val -> key.startsWith('versions') } + ).match() }, + { assert process.out.output.get(0).get(1).endsWith('test') }, + { assert process.out.index.size() == 0 } + ) + } + } + + test("metagenome - xz - (re)compress - no index") { + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/prokaryotes/metagenome/taxonomy/misc/taxa_sqlite.xz', checkIfExists: true), + [], + [] + ] + input[1] = 'compress' // action + input[2] = false // make_index + input[3] = '' // out_ext + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(sanitizeOutput(process.out)).match() }, + { assert process.out.output.get(0).get(1).endsWith('.gz') }, + { assert process.out.index.size() == 0 } + ) + } + } + + test("sarscov2 - vcf - compress - index - stub") { + + options "-stub" + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/sarscov2/illumina/vcf/test.vcf', checkIfExists: true), + [], + [] + ] + input[1] = 'compress' // action + input[2] = true // make_index + input[3] = 'vcf' // out_ext + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(sanitizeOutput(process.out)).match() } + ) + } + + } + + test("sarscov2 - vcf - decompress - stub") { + + options "-stub" + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/sarscov2/illumina/vcf/test.vcf.gz', checkIfExists: true), + [], + [] + ] + input[1] = 'decompress' // action + input[2] = false // make_index + input[3] = 'vcf' // out_ext + """ + } + } + + then { + assert process.success + assertAll( + { assert snapshot(sanitizeOutput(process.out)).match() } + ) + } + + } + + test("illegal action") { + + when { + process { + """ + input[0] = [ + [ id:'test' ], + file(params.modules_testdata_base_path + 'genomics/sarscov2/illumina/vcf/test.vcf', checkIfExists: true), + [], + [] + ] + input[1] = 'invalid_action' // action + input[2] = true // make_index + input[3] = 'vcf' // out_ext + """ + } + } + + then { + assert process.failed + assert process.errorReport.contains("Invalid action: invalid_action. Allowed actions are: compress, decompress") + } + + } + +} diff --git a/modules/nf-core/htslib/bgziptabix/tests/main.nf.test.snap b/modules/nf-core/htslib/bgziptabix/tests/main.nf.test.snap new file mode 100644 index 00000000..4d8b3d25 --- /dev/null +++ b/modules/nf-core/htslib/bgziptabix/tests/main.nf.test.snap @@ -0,0 +1,527 @@ +{ + "sarscov2 - gzip - (re)compress - no index": { + "content": [ + { + "index": [ + + ], + "output": [ + [ + { + "id": "test" + }, + "test.fastq.gz:md5,4161df271f9bfcd25d5845a1e220dbec" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + } + ], + "timestamp": "2026-07-10T08:07:32.981182", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + }, + "metagenome - xz - (re)compress - no index": { + "content": [ + { + "index": [ + + ], + "output": [ + [ + { + "id": "test" + }, + "test.gz:md5,b8d852a2b1ee52ed64d83046dcdb9de2" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + } + ], + "timestamp": "2026-07-10T08:08:57.748138", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + }, + "metagenome - bz2 - decompress": { + "content": [ + { + "index": [ + + ], + "output": [ + [ + { + "id": "test" + }, + "test.tar:md5,39e9e71fd16cfd09ceca12cd46e6abce" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + } + ], + "timestamp": "2026-07-10T08:07:44.215383", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + }, + "sarscov2 - vcf - decompress - stub": { + "content": [ + { + "index": [ + + ], + "output": [ + [ + { + "id": "test" + }, + "test.vcf:md5,68b329da9893e34099c7d8ad5cb9c940" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + } + ], + "timestamp": "2026-07-10T08:10:11.941036", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + }, + "sarscov2 - gzip - decompress": { + "content": [ + { + "index": [ + + ], + "output": [ + [ + { + "id": "test" + }, + "test.fastq:md5,4161df271f9bfcd25d5845a1e220dbec" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + } + ], + "timestamp": "2026-07-10T08:07:28.301585", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + }, + "sarscov2 - vcf - compress - index - stub": { + "content": [ + { + "index": [ + [ + { + "id": "test" + }, + "test.vcf.gz.tbi:md5,d41d8cd98f00b204e9800998ecf8427e" + ] + ], + "output": [ + [ + { + "id": "test" + }, + "test.vcf.gz:md5,68b329da9893e34099c7d8ad5cb9c940" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + } + ], + "timestamp": "2026-07-10T08:09:24.961486", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + }, + "sarscov2 - vcf - decompress": { + "content": [ + { + "index": [ + + ], + "output": [ + [ + { + "id": "test" + }, + "test.vcf:md5,8e722884ffb75155212a3fc053918766" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + } + ], + "timestamp": "2026-07-10T08:07:05.41219", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + }, + "metagenome - bz2 - (re)compress - no index": { + "content": [ + { + "index": [ + + ], + "output": [ + [ + { + "id": "test" + }, + "test.tar.gz:md5,39e9e71fd16cfd09ceca12cd46e6abce" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + } + ], + "timestamp": "2026-07-10T08:07:51.772557", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + }, + "sarscov2 - vcf - compress - index": { + "content": [ + { + "index": [ + [ + { + "id": "test" + }, + "test.vcf.gz.tbi:md5,7f005943c935f2b55ba3f9d4802aa09f" + ] + ], + "output": [ + [ + { + "id": "test" + }, + "test.vcf.gz:md5,8e722884ffb75155212a3fc053918766" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + } + ], + "timestamp": "2026-07-10T08:07:09.68169", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + }, + "metagenome - xz - decompress": { + "content": [ + { + "0": [ + [ + { + "id": "test" + }, + "test:md5,b8d852a2b1ee52ed64d83046dcdb9de2" + ] + ], + "1": [ + + ], + "2": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "3": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ], + "index": [ + + ], + "output": [ + [ + { + "id": "test" + }, + "test:md5,b8d852a2b1ee52ed64d83046dcdb9de2" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + }, + { + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + } + ], + "timestamp": "2026-07-10T08:08:19.920765", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + }, + "sarscov2 - bgzip - compress - no index": { + "content": [ + { + "index": [ + + ], + "output": [ + [ + { + "id": "test" + }, + "test.vcf.gz:md5,8e722884ffb75155212a3fc053918766" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + } + ], + "timestamp": "2026-07-10T08:07:23.405306", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + }, + "sarscov2 - vcf + regions - compress - index": { + "content": [ + { + "index": [ + [ + { + "id": "example" + }, + "example.vcf.gz.tbi:md5,d22e5b84e4fcd18792179f72e6da702e" + ] + ], + "output": [ + [ + { + "id": "example" + }, + "example.vcf.gz:md5,8e722884ffb75155212a3fc053918766" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + }, + "bc7bf3ee9e8430e064c539eb81e59bf9" + ], + "timestamp": "2026-07-10T08:07:15.187648", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + }, + "sarscov2 - bgzip - decompress": { + "content": [ + { + "index": [ + + ], + "output": [ + [ + { + "id": "test" + }, + "test.vcf:md5,8e722884ffb75155212a3fc053918766" + ] + ], + "versions_htslib": [ + [ + "HTSLIB_BGZIPTABIX", + "htslib", + "1.24" + ] + ], + "versions_xz": [ + [ + "HTSLIB_BGZIPTABIX", + "xz", + "5.8.3" + ] + ] + } + ], + "timestamp": "2026-07-10T08:07:19.313945", + "meta": { + "nf-test": "0.9.5", + "nextflow": "26.04.6" + } + } +} \ No newline at end of file diff --git a/modules/nf-core/modkit/dmr/environment.yml b/modules/nf-core/modkit/dmr/environment.yml new file mode 100644 index 00000000..62b97863 --- /dev/null +++ b/modules/nf-core/modkit/dmr/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: + - ont-modkit=0.6.1 diff --git a/modules/nf-core/modkit/dmr/main.nf b/modules/nf-core/modkit/dmr/main.nf new file mode 100644 index 00000000..41ccec1c --- /dev/null +++ b/modules/nf-core/modkit/dmr/main.nf @@ -0,0 +1,60 @@ +process MODKIT_DMR { + tag "${meta.id}" + label 'process_medium' + + // Conda is not supported: matches MODKIT_PILEUP's own patch (the guard in `script:` stops + // conda/mamba runs) -- kept on the same patched build as pileup for a single consistent + // modkit binary across the whole DMR chain, rather than mixing stock 0.6.1 (dmr pair) with + // the patched 0.6.4 (pileup). Revert to the biocontainer once the patch is upstreamed/released. + container "${(workflow.containerEngine == 'singularity' || workflow.containerEngine == 'apptainer') && !task.ext.singularity_pull_docker_container + ? 'oras://ghcr.io/ljwharbers/modkit-sif:0.6.4-pacbiofix-6e0afa2' + : 'ghcr.io/ljwharbers/modkit:0.6.4-pacbiofix-6e0afa2'}" + + input: + tuple val(meta), path(bedmethyl_a), path(bedmethyl_a_tbi) + tuple val(meta2), path(bedmethyl_b), path(bedmethyl_b_tbi) + tuple val(meta3), path(regions_bed) + tuple val(meta4), path(fasta), path(fai) + + output: + tuple val(meta), path("*.bed"), emit: bed + tuple val(meta), path("*.log"), emit: log + 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: + // Exit if running this module with -profile conda / -profile mamba + if (workflow.profile.tokenize(',').intersect(['conda', 'mamba']).size() >= 1) { + error "MODKIT_DMR does not support Conda: the patched container is Docker/Singularity/Apptainer-only. Use one of those, or --skip_dmr." + } + def args = task.ext.args ?: '' + def prefix = task.ext.prefix ?: "${meta.id}" + def regions = regions_bed ? "-r ${regions_bed}" : '' + // fai is never referenced directly -- modkit dmr pair has no --fai flag and looks for + // .fai next to --ref itself. Declaring it as an input only stages it alongside + // fasta so modkit finds and reuses it instead of building one from scratch, which for a + // multi-GB genome is a real, repeated cost across every task invocation. Confirmed by + // direct testing: with no .fai present, modkit dmr pair creates one; with one already + // staged next to fasta, it's read (unmodified mtime) rather than rebuilt. + """ + modkit \\ + dmr pair \\ + -a ${bedmethyl_a} \\ + -b ${bedmethyl_b} \\ + ${regions} \\ + --ref ${fasta} \\ + -o ${prefix}.bed \\ + --log-filepath ${prefix}.log \\ + -t ${task.cpus} \\ + ${args} + """ + + stub: + def prefix = task.ext.prefix ?: "${meta.id}" + """ + touch ${prefix}.bed + touch ${prefix}.log + """ +} diff --git a/modules/nf-core/modkit/dmr/meta.yml b/modules/nf-core/modkit/dmr/meta.yml new file mode 100644 index 00000000..63c958da --- /dev/null +++ b/modules/nf-core/modkit/dmr/meta.yml @@ -0,0 +1,140 @@ +name: modkit_dmr +description: Compare regions between a pair of bedMethyl samples (e.g. two haplotypes, + or tumor and normal) and report differentially methylated regions +keywords: + - methylation + - ont + - long-read + - differential methylation +tools: + - "modkit": + description: A bioinformatics tool for working with modified bases in Oxford + Nanopore sequencing data + homepage: https://github.com/nanoporetech/modkit + documentation: https://github.com/nanoporetech/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 for the first (control) input + e.g. `[ id:'test' ]` + - bedmethyl_a: + type: file + description: Bgzipped bedMethyl file for the first (control) sample, as + produced by modkit pileup with --bgzf + pattern: "*.bed.gz" + ontologies: + - edam: "http://edamontology.org/format_3003" # BED + - edam: "http://edamontology.org/format_3989" # GZIP format + - bedmethyl_a_tbi: + type: file + description: Tabix index for bedmethyl_a. Must have the same basename with + a .tbi extension next to bedmethyl_a + pattern: "*.bed.gz.tbi" + ontologies: + - edam: "http://edamontology.org/format_3700" # Tabix index format + - - meta2: + type: map + description: | + Groovy Map containing sample information for the second (experimental) + input e.g. `[ id:'test' ]` + - bedmethyl_b: + type: file + description: Bgzipped bedMethyl file for the second (experimental) sample, + as produced by modkit pileup with --bgzf + pattern: "*.bed.gz" + ontologies: + - edam: "http://edamontology.org/format_3003" # BED + - edam: "http://edamontology.org/format_3989" # GZIP format + - bedmethyl_b_tbi: + type: file + description: Tabix index for bedmethyl_b. Must have the same basename with + a .tbi extension next to bedmethyl_b + pattern: "*.bed.gz.tbi" + ontologies: + - edam: "http://edamontology.org/format_3700" # Tabix index format + - - meta3: + type: map + description: | + Groovy Map containing regions information + e.g. `[ id:'regions' ]` + - regions_bed: + type: file + description: Optional BED file of regions over which to compare methylation + levels. When omitted, methylation levels are compared at each site instead + pattern: "*.bed" + ontologies: + - edam: "http://edamontology.org/format_3003" # BED + - - meta4: + type: map + description: | + Groovy Map containing reference information + e.g. `[ id:'hg38' ]` + - fasta: + type: file + description: Reference sequence in FASTA format, used during pileup/alignment + pattern: "*.{fa,fasta}" + ontologies: + - edam: "http://edamontology.org/format_1929" # FASTA + - fai: + type: file + description: | + Optional FASTA index for fasta. modkit dmr pair has no --fai flag and looks for + .fai next to --ref itself; passing one here only stages it alongside fasta + so modkit reuses it instead of building one from scratch on every invocation, which + for a multi-GB genome is a real, repeated cost. Pass `[]` to omit. + pattern: "*.fai" + ontologies: [] +output: + bed: + - - meta: + type: map + description: | + Groovy Map containing sample information + e.g. `[ id:'test' ]` + - "*.bed": + type: file + description: BED file with the score column indicating the magnitude of + the difference in methylation between the two samples + pattern: "*.bed" + ontologies: + - edam: "http://edamontology.org/format_3003" # BED + log: + - - meta: + type: map + description: | + Groovy Map containing sample information + e.g. `[ id:'test' ]` + - "*.log": + type: file + description: File for debug logs to be written to + 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: + - "@YannVRB" +maintainers: + - "@YannVRB" diff --git a/modules/nf-core/modkit/dmr/modkit-dmr.diff b/modules/nf-core/modkit/dmr/modkit-dmr.diff new file mode 100644 index 00000000..bd6a816a --- /dev/null +++ b/modules/nf-core/modkit/dmr/modkit-dmr.diff @@ -0,0 +1,36 @@ +Changes in component 'nf-core/modkit/dmr' +'modules/nf-core/modkit/dmr/environment.yml' is unchanged +'modules/nf-core/modkit/dmr/meta.yml' is unchanged +Changes in 'modkit/dmr/main.nf': +--- modules/nf-core/modkit/dmr/main.nf ++++ modules/nf-core/modkit/dmr/main.nf +@@ -2,10 +2,13 @@ + tag "${meta.id}" + label 'process_medium' + +- 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'}" ++ // Conda is not supported: matches MODKIT_PILEUP's own patch (the guard in `script:` stops ++ // conda/mamba runs) -- kept on the same patched build as pileup for a single consistent ++ // modkit binary across the whole DMR chain, rather than mixing stock 0.6.1 (dmr pair) with ++ // the patched 0.6.4 (pileup). Revert to the biocontainer once the patch is upstreamed/released. ++ container "${(workflow.containerEngine == 'singularity' || workflow.containerEngine == 'apptainer') && !task.ext.singularity_pull_docker_container ++ ? 'oras://ghcr.io/ljwharbers/modkit-sif:0.6.4-pacbiofix-6e0afa2' ++ : 'ghcr.io/ljwharbers/modkit:0.6.4-pacbiofix-6e0afa2'}" + + input: + tuple val(meta), path(bedmethyl_a), path(bedmethyl_a_tbi) +@@ -22,6 +25,10 @@ + task.ext.when == null || task.ext.when + + script: ++ // Exit if running this module with -profile conda / -profile mamba ++ if (workflow.profile.tokenize(',').intersect(['conda', 'mamba']).size() >= 1) { ++ error "MODKIT_DMR does not support Conda: the patched container is Docker/Singularity/Apptainer-only. Use one of those, or --skip_dmr." ++ } + def args = task.ext.args ?: '' + def prefix = task.ext.prefix ?: "${meta.id}" + def regions = regions_bed ? "-r ${regions_bed}" : '' +************************************************************ diff --git a/modules/nf-core/modkit/dmr/tests/main.nf.test b/modules/nf-core/modkit/dmr/tests/main.nf.test new file mode 100644 index 00000000..d7d4e6de --- /dev/null +++ b/modules/nf-core/modkit/dmr/tests/main.nf.test @@ -0,0 +1,194 @@ +nextflow_process { + + name "Test Process MODKIT_DMR" + script "../main.nf" + tag "modules" + tag "modules_nfcore" + tag "modkit" + tag "modkit/dmr" + tag "modkit/pileup" + tag "htslib/bgziptabix" + process "MODKIT_DMR" + config "./nextflow.config" + + setup { + run("MODKIT_PILEUP") { + script "../../pileup/main.nf" + 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: 'test_ref' ], + 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] = [[],[]] + """ + } + } + + run("HTSLIB_BGZIPTABIX") { + script "../../../htslib/bgziptabix/main.nf" + process { + """ + // MODKIT_PILEUP --phased emits one [meta, [combined, hp1, hp2]] tuple per sample + // ("bedgz_files" is a single Path rather than a List when only one file matches + // the module's output glob -- true for -stub, since modkit pileup's stub script + // never splits by haplotype -- so normalize via [x].flatten() first). Split into + // one [meta, bedgz] tuple per haplotype so each can be tabix-indexed (modkit + // pileup --bgzf already bgzips, but never indexes) and fed to MODKIT_DMR as its + // own bedmethyl_a/bedmethyl_b input, exactly like real --modkit_phased usage. + input[0] = MODKIT_PILEUP.out.bedgz + .flatMap { meta, bedgz_files -> + def files = [ bedgz_files ].flatten() + def phased = files.findAll { f -> f.name =~ /_hp[12]\\.bed\\.gz\$/ } + if (phased) { + phased.collect { f -> + def haplotype = (f.name =~ /_(hp[12])\\.bed\\.gz\$/)[0][1] + [ meta + [ id: "\${meta.id}_\${haplotype}" ], f, [], [] ] + } + } else { + // -stub: modkit pileup's stub always emits one generic, unsplit file + // (it never reads --phased in stub mode) -- reuse it for both + // haplotypes, since MODKIT_DMR's own stub never reads its inputs. + files.collectMany { f -> + ['hp1', 'hp2'].collect { haplotype -> [ meta + [ id: "\${meta.id}_\${haplotype}" ], f, [], [] ] } + } + } + } + input[1] = "compress" + input[2] = true + input[3] = [] + """ + } + } + } + + test("[bedmethyl_a, tbi], [bedmethyl_b, tbi], regions_bed, fasta") { + + when { + params { + module_args = '--base C' + } + process { + """ + input[0] = HTSLIB_BGZIPTABIX.out.output + .join(HTSLIB_BGZIPTABIX.out.index) + .filter { meta, bedgz, tbi -> meta.id.endsWith('_hp1') } + .map { meta, bedgz, tbi -> [ [ id: 'test' ], bedgz, tbi ] } + input[1] = HTSLIB_BGZIPTABIX.out.output + .join(HTSLIB_BGZIPTABIX.out.index) + .filter { meta, bedgz, tbi -> meta.id.endsWith('_hp2') } + .map { meta, bedgz, tbi -> [ [ id: 'test' ], bedgz, tbi ] } + input[2] = Channel.of('chr22\t0\t40001') + .collectFile(name: 'chr22.bed', newLine: true) + .map { file -> [ [ id:'chr22' ], file ] } + input[3] = [ + [ id: 'test_ref' ], + 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) + ] + """ + } + } + + then { + assertAll ( + { assert process.success }, + // log-file content includes a timestamp/runtime that varies between runs -- + // assert existence only, and snapshot everything else. + { assert process.out.log.size() == 1 }, + { assert snapshot(process.out.bed, process.out.findAll { key, val -> key.startsWith('versions') } + ).match() } + ) + } + + } + + test("[bedmethyl_a, tbi], [bedmethyl_b, tbi], [], fasta") { + + when { + params { + module_args = '--base C' + } + process { + """ + input[0] = HTSLIB_BGZIPTABIX.out.output + .join(HTSLIB_BGZIPTABIX.out.index) + .filter { meta, bedgz, tbi -> meta.id.endsWith('_hp1') } + .map { meta, bedgz, tbi -> [ [ id: 'test' ], bedgz, tbi ] } + input[1] = HTSLIB_BGZIPTABIX.out.output + .join(HTSLIB_BGZIPTABIX.out.index) + .filter { meta, bedgz, tbi -> meta.id.endsWith('_hp2') } + .map { meta, bedgz, tbi -> [ [ id: 'test' ], bedgz, tbi ] } + input[2] = [ [ id:'regions' ], [] ] + input[3] = [ + [ id: 'test_ref' ], + 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) + ] + """ + } + } + + then { + assertAll ( + { assert process.success }, + // log-file content includes a timestamp/runtime that varies between runs -- + // assert existence only, and snapshot everything else. + { assert process.out.log.size() == 1 }, + { assert snapshot(process.out.bed, process.out.findAll { key, val -> key.startsWith('versions') } + ).match() } + ) + } + + } + + test("[bedmethyl_a, tbi], [bedmethyl_b, tbi], regions_bed, fasta - stub") { + + options "-stub" + + when { + params { + module_args = '--base C' + } + process { + """ + input[0] = HTSLIB_BGZIPTABIX.out.output + .join(HTSLIB_BGZIPTABIX.out.index) + .filter { meta, bedgz, tbi -> meta.id.endsWith('_hp1') } + .map { meta, bedgz, tbi -> [ [ id: 'test' ], bedgz, tbi ] } + input[1] = HTSLIB_BGZIPTABIX.out.output + .join(HTSLIB_BGZIPTABIX.out.index) + .filter { meta, bedgz, tbi -> meta.id.endsWith('_hp2') } + .map { meta, bedgz, tbi -> [ [ id: 'test' ], bedgz, tbi ] } + input[2] = Channel.of('chr22\t0\t40001') + .collectFile(name: 'chr22.bed', newLine: true) + .map { file -> [ [ id:'chr22' ], file ] } + input[3] = [ + [ id: 'test_ref' ], + 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) + ] + """ + } + } + + then { + assertAll ( + { assert process.success }, + // log-file content includes a timestamp/runtime that varies between runs -- + // assert existence only, and snapshot everything else. + { assert process.out.log.size() == 1 }, + { assert snapshot(process.out.bed, process.out.findAll { key, val -> key.startsWith('versions') } +).match() } + ) + } + + } + +} diff --git a/modules/nf-core/modkit/dmr/tests/main.nf.test.snap b/modules/nf-core/modkit/dmr/tests/main.nf.test.snap new file mode 100644 index 00000000..a7039d7d --- /dev/null +++ b/modules/nf-core/modkit/dmr/tests/main.nf.test.snap @@ -0,0 +1,80 @@ +{ + "[bedmethyl_a, tbi], [bedmethyl_b, tbi], regions_bed, fasta": { + "content": [ + [ + [ + { + "id": "test" + }, + "test.bed:md5,c57c0b35ce975dbafa19f3b4b5604dec" + ] + ], + { + "versions_modkit": [ + [ + "MODKIT_DMR", + "modkit", + "0.6.1" + ] + ] + } + ], + "meta": { + "nf-test": "0.9.0", + "nextflow": "25.10.2" + }, + "timestamp": "2026-09-29T11:00:29.834449326" + }, + "[bedmethyl_a, tbi], [bedmethyl_b, tbi], [], fasta": { + "content": [ + [ + [ + { + "id": "test" + }, + "test.bed:md5,e1ed04436fd4aee1c05fa4a612846417" + ] + ], + { + "versions_modkit": [ + [ + "MODKIT_DMR", + "modkit", + "0.6.1" + ] + ] + } + ], + "meta": { + "nf-test": "0.9.0", + "nextflow": "25.10.2" + }, + "timestamp": "2026-09-29T11:00:38.70209544" + }, + "[bedmethyl_a, tbi], [bedmethyl_b, tbi], regions_bed, fasta - stub": { + "content": [ + [ + [ + { + "id": "test" + }, + "test.bed:md5,d41d8cd98f00b204e9800998ecf8427e" + ] + ], + { + "versions_modkit": [ + [ + "MODKIT_DMR", + "modkit", + "0.6.1" + ] + ] + } + ], + "meta": { + "nf-test": "0.9.0", + "nextflow": "25.10.2" + }, + "timestamp": "2026-09-29T11:00:46.52021922" + } +} \ No newline at end of file diff --git a/modules/nf-core/modkit/dmr/tests/nextflow.config b/modules/nf-core/modkit/dmr/tests/nextflow.config new file mode 100644 index 00000000..a6b69568 --- /dev/null +++ b/modules/nf-core/modkit/dmr/tests/nextflow.config @@ -0,0 +1,14 @@ +process { + withName: 'MODKIT_PILEUP' { + // --modified-bases 5mC 5hmC (matching modkit/pileup's own phased test): if this BAM's + // basecalling calls both, requesting only one leaves the other counted as "other" in + // modkit's bedMethyl output, which MODKIT_DMR's internal consistency check rejects. + ext.args = '--phased --modified-bases 5mC 5hmC' + } + withName: 'HTSLIB_BGZIPTABIX' { + ext.args2 = '-p bed' + } + withName: 'MODKIT_DMR' { + ext.args = params.module_args + } +} diff --git a/nextflow.config b/nextflow.config index 696d5161..66dc85aa 100644 --- a/nextflow.config +++ b/nextflow.config @@ -29,6 +29,14 @@ params { modkit_args = '--cpg --modified-bases 5mC' modkit_phased = false + // DMR options -- only take effect when modkit_phased is also true. dmr_cpg_islands_bed/ + // dmr_gencode_gene_bed default to null here; --genome GRCh38/CHM13 fill them in from + // igenomes.config unless explicitly overridden, and any other genome must supply both + // explicitly or DMR is skipped with a warning (see workflows/lrsomatic.nf). + skip_dmr = false + dmr_cpg_islands_bed = null + dmr_gencode_gene_bed = null + // PON Options clairsto_pon_vcfs = null clairsto_pon_flags = null diff --git a/nextflow_schema.json b/nextflow_schema.json index 713e2baa..44e04a3e 100644 --- a/nextflow_schema.json +++ b/nextflow_schema.json @@ -170,6 +170,20 @@ "description": "Verdict CNA resource directory for ClairS-TO, holding loci_files/, allele_files/ and one GC_*.txt. Read only with --skip_ascat.", "help_text": "Read only on `--skip_ascat` runs; otherwise Verdict's germline tagging comes from the pipeline's ASCAT run and this is ignored. The image ships the GRCh38 set and `--genome CHM13 --skip_ascat` builds a CHM13 one, so this is only needed for another assembly or as an override. The directory must be self-contained, with no links reaching outside it. Add an `RT_.txt` to correct LogR for replication timing as well; without it, correction is GC-only.", "fa_icon": "fas fa-folder-open" + }, + "dmr_cpg_islands_bed": { + "type": "string", + "format": "file-path", + "description": "CpG-islands BED to restrict haplotype DMR calling to. Has a built-in default for --genome GRCh38/CHM13; required as an explicit override for any other genome when --modkit_phased is set and --skip_dmr is not.", + "help_text": "Only used by the DMR subworkflow (see --skip_dmr). Without a value -- built-in or supplied -- DMR is skipped with a warning rather than the run failing.", + "fa_icon": "fas fa-map-marker-alt" + }, + "dmr_gencode_gene_bed": { + "type": "string", + "format": "file-path", + "description": "Gene-model BED (used for nearest-gene annotation of DMR calls) to use for haplotype DMR calling. Has a built-in default for --genome GRCh38/CHM13; required as an explicit override for any other genome when --modkit_phased is set and --skip_dmr is not.", + "help_text": "Only used by the DMR subworkflow (see --skip_dmr). Without a value -- built-in or supplied -- DMR is skipped with a warning rather than the run failing.", + "fa_icon": "fas fa-dna" } } }, @@ -547,6 +561,10 @@ "type": "boolean", "description": "Skip SAVANA (SV + copy-number calling)" }, + "skip_dmr": { + "type": "boolean", + "description": "Skip differential methylation region (DMR) calling between haplotypes. Only takes effect when --modkit_phased is also set; has no effect otherwise." + }, "skip_m6a": { "type": "boolean", "description": "Skip m6a calling by Fibertools" diff --git a/subworkflows/local/dmr.nf b/subworkflows/local/dmr.nf new file mode 100644 index 00000000..1d7d4fa5 --- /dev/null +++ b/subworkflows/local/dmr.nf @@ -0,0 +1,135 @@ +// IMPORT MODULES +include { HTSLIB_BGZIPTABIX } from '../../modules/nf-core/htslib/bgziptabix/main' +include { MODKIT_DMR } from '../../modules/nf-core/modkit/dmr/main' +include { DMR_HAPLOTYPE_REGIONS } from '../../modules/local/dmr/haplotype_regions/main' +include { DMR_NEAREST_GENE } from '../../modules/local/dmr/nearest_gene/main' + +workflow DMR { + + take: + modkit_bedgz // [meta, bed.gz(es)] -- MODKIT_PILEUP.out.bedgz; only meaningful when modkit_phased is true + fasta // [[:], fasta] + fai // [[:], fai] + cpg_islands_bed // path + gencode_gene_bed // path + + main: + ch_versions = channel.empty() + + // Pick the hp1/hp2 pair out of MODKIT_PILEUP's glob-collected output list (which also + // includes a _combined file when --modkit_phased is set). An unphased run's single + // unsuffixed file matches neither pattern, so it's dropped by the filter below -- + // defence in depth alongside the modkit_phased gate on the caller's side. + modkit_bedgz + .map { meta, files -> + def flist = files instanceof List ? files : [files] + def hp1 = flist.find { it.name.endsWith('_hp1.bed.gz') } + def hp2 = flist.find { it.name.endsWith('_hp2.bed.gz') } + return [meta, hp1, hp2] + } + .filter { _meta, hp1, hp2 -> hp1 && hp2 } + .set { haplotype_bedmethyl } + // haplotype_bedmethyl: [meta, hp1_bedgz, hp2_bedgz] + + // + // MODULE: HTSLIB_BGZIPTABIX (label: process_low) + // Tabix-index each haplotype's bedMethyl -- modkit dmr pair requires a .tbi alongside each + // bgzip input. The bedMethyl is already bgzip-compressed (modkit pileup --bgzf), so this + // only adds the index. hp1/hp2 are tagged onto meta and mixed into one call, then split + // back apart below. + // + haplotype_bedmethyl + .flatMap { meta, hp1, hp2 -> + return [ + [meta + [haplotype: 'hp1'], hp1, [], []], + [meta + [haplotype: 'hp2'], hp2, [], []] + ] + } + .set { bgziptabix_input } + // bgziptabix_input: [meta+haplotype, bedmethyl, [], []] + + HTSLIB_BGZIPTABIX ( + bgziptabix_input, + 'compress', + true, + 'bed' + ) + + HTSLIB_BGZIPTABIX.out.output + .join(HTSLIB_BGZIPTABIX.out.index) + .map { meta, bedgz, tbi -> + def haplotype = meta.haplotype + def sample_meta = meta.findAll { it.key != 'haplotype' } + return [sample_meta, haplotype, bedgz, tbi] + } + .branch { _meta, haplotype, _bedgz, _tbi -> + hp1: haplotype == 'hp1' + hp2: haplotype == 'hp2' + } + .set { indexed_branched } + + indexed_branched.hp1 + .map { meta, _haplotype, bedgz, tbi -> [meta, bedgz, tbi] } + .set { hp1_indexed } + indexed_branched.hp2 + .map { meta, _haplotype, bedgz, tbi -> [meta, bedgz, tbi] } + .set { hp2_indexed } + + hp1_indexed + .join(hp2_indexed) + .set { dmr_haplotype_input } + // dmr_haplotype_input: [meta, hp1_bedgz, hp1_tbi, hp2_bedgz, hp2_tbi] + + // + // MODULE: DMR_HAPLOTYPE_REGIONS (label: process_single) + // Restrict cpg_islands_bed to islands covered in both haplotypes. + // + DMR_HAPLOTYPE_REGIONS ( + dmr_haplotype_input, + cpg_islands_bed, + fai.first() + ) + ch_versions = ch_versions.mix(DMR_HAPLOTYPE_REGIONS.out.versions) + + // + // MODULE: MODKIT_DMR (label: process_medium) + // Compare methylation between the two haplotypes over the restricted regions. + // + dmr_haplotype_input + .join(DMR_HAPLOTYPE_REGIONS.out.regions_bed) + .multiMap { meta, hp1_bedgz, hp1_tbi, hp2_bedgz, hp2_tbi, regions_bed -> + hp1: [meta, hp1_bedgz, hp1_tbi] + hp2: [meta, hp2_bedgz, hp2_tbi] + regions: [meta, regions_bed] + } + .set { modkit_dmr_input } + + // MODKIT_DMR's fai input only avoids rebuilding the index on every task -- see the + // module's own script comment. fasta/fai are each a single [[:], path] value channel + // (take: comments above), so .first() on each stays correct after the join. + fasta.first() + .combine(fai.first()) + .map { meta, fasta_file, _meta2, fai_file -> [meta, fasta_file, fai_file] } + .set { fasta_with_fai } + + MODKIT_DMR ( + modkit_dmr_input.hp1, + modkit_dmr_input.hp2, + modkit_dmr_input.regions, + fasta_with_fai + ) + + // + // MODULE: DMR_NEAREST_GENE (label: process_single) + // Annotate each DMR with its nearest gene. + // + DMR_NEAREST_GENE ( + MODKIT_DMR.out.bed, + gencode_gene_bed + ) + ch_versions = ch_versions.mix(DMR_NEAREST_GENE.out.versions) + + emit: + nearest_gene_bed = DMR_NEAREST_GENE.out.nearest_gene_bed // [meta, bed] -- DMRs annotated with nearest gene + versions = ch_versions +} diff --git a/subworkflows/local/tests/dmr.nf.test b/subworkflows/local/tests/dmr.nf.test new file mode 100644 index 00000000..c940e906 --- /dev/null +++ b/subworkflows/local/tests/dmr.nf.test @@ -0,0 +1,98 @@ +nextflow_workflow { + + name "Test Workflow DMR" + script "../dmr.nf" + workflow "DMR" + config "./nextflow.config" + + tag "subworkflows" + tag "subworkflows_local" + tag "dmr" + // "small" is what .github/workflows/nf-test.yml selects on for pull_request + tag "small" + + // Real (not stub), public, non-confidential coverage of the full haplotype-DMR chain: + // MODKIT_PILEUP (--phased) -> HTSLIB_BGZIPTABIX -> DMR_HAPLOTYPE_REGIONS -> MODKIT_DMR -> + // DMR_NEAREST_GENE. Uses the same public, phased nf-core/test-datasets BAM already used by + // modkit/pileup's and modkit/dmr's own module tests -- no new fixture needed. cpg_islands_bed + // reuses the exact chr22:0-40001 region modkit/dmr's own module test already proves has real + // dual-haplotype coverage; gencode_gene_bed is a small synthetic gene entry over the same span + // for DMR_NEAREST_GENE to annotate against. + setup { + run("MODKIT_PILEUP") { + script "../../../modules/nf-core/modkit/pileup/main.nf" + process { + """ + // lrsomatic's vendored modkit/pileup takes fasta/fai as separate 2-tuples plus + // a bed-region tuple -- 4 inputs, not the 3-tuple/3-input shape upstream + // nf-core/modules now uses. Matches its real call in workflows/lrsomatic.nf. + input[0] = [ + [ id: 'test', type: 'tumor' ], + 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: 'test_ref' ], + file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.fasta', checkIfExists: true) + ] + input[2] = [ + [ id: 'test_ref' ], + file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.fasta.fai', checkIfExists: true) + ] + input[3] = [[:],[]] + """ + } + } + } + + test("real phased public data -> haplotype DMR, annotated with nearest gene") { + + when { + workflow { + """ + def cpg_islands_bed = Channel.of('chr22\\t0\\t40001') + .collectFile(name: 'cpg_islands.bed', newLine: true) + + def gencode_gene_bed = Channel.of('chr22\\t0\\t40001\\tTEST_GENE\\t.\\t+') + .collectFile(name: 'gencode_gene.bed', newLine: true) + + // DMR's internal fasta.first()/fai.first() need real single-emission value + // channels -- a bare list literal here is a 2-item queue channel (meta, then + // path, as separate emissions), and .first() would silently take just the meta. + input[0] = MODKIT_PILEUP.out.bedgz + input[1] = channel.value([ + [ id: 'test_ref' ], + file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.fasta', checkIfExists: true) + ]) + input[2] = channel.value([ + [ id: 'test_ref' ], + file(params.modules_testdata_base_path + 'genomics/homo_sapiens/genome/genome.fasta.fai', checkIfExists: true) + ]) + input[3] = cpg_islands_bed + input[4] = gencode_gene_bed + """ + } + } + + then { + def processes = workflow.trace.tasks().collect { task -> task.name } + + assertAll( + { assert workflow.success }, + { assert processes.any { name -> name.contains('HTSLIB_BGZIPTABIX') } }, + { assert processes.any { name -> name.contains('DMR_HAPLOTYPE_REGIONS') } }, + { assert processes.any { name -> name.contains('MODKIT_DMR') } }, + { assert processes.any { name -> name.contains('DMR_NEAREST_GENE') } }, + { assert workflow.out.nearest_gene_bed.size() == 1 }, + // Real content, not just a file existing: at least one DMR call, annotated + // with the synthetic gene above. + { + def bed = workflow.out.nearest_gene_bed[0][1] + def lines = file(bed).readLines() + assert lines.size() > 0 + assert lines.every { line -> line.contains('TEST_GENE') } + } + ) + } + } +} diff --git a/subworkflows/local/tests/nextflow.config b/subworkflows/local/tests/nextflow.config new file mode 100644 index 00000000..a1704d19 --- /dev/null +++ b/subworkflows/local/tests/nextflow.config @@ -0,0 +1,19 @@ +// A subworkflow-level nf-test session does not reliably pick up conf/modules.config's +// DMR-scoped withName blocks (observed directly: HTSLIB_BGZIPTABIX's ext.prefix/ext.args2 and +// MODKIT_DMR's ext.args did not apply, even though the DMR: prefix was present in the process +// name). Set everything this real (non-stub) chain needs explicitly here instead of relying on +// the pipeline's shared config, so the test is self-contained and its behaviour is guaranteed. +process { + withName: 'MODKIT_PILEUP' { + ext.args = '--phased --modified-bases 5mC 5hmC' + } + withName: '.*HTSLIB_BGZIPTABIX' { + // Without meta.haplotype in the prefix, hp1/hp2 both output ".bed.gz" and collide + // once staged together into DMR_HAPLOTYPE_REGIONS. + ext.prefix = { "${meta.id}_${meta.haplotype}" } + ext.args2 = '-p bed' + } + withName: '.*MODKIT_DMR' { + ext.args = '--base C' + } +} diff --git a/workflows/lrsomatic.nf b/workflows/lrsomatic.nf index 37bd1c4c..b668e4cf 100644 --- a/workflows/lrsomatic.nf +++ b/workflows/lrsomatic.nf @@ -60,6 +60,7 @@ include { PAIRED_SMALLVAR_GERMLINE } from '../subworkflows/local/paired/p include { PHASING_HAPLOTYPING } from '../subworkflows/local/phasing_haplotyping' include { TUMORONLY_SAVANA } from '../subworkflows/local/tumor_only/tumoronly_savana' include { PAIRED_SAVANA } from '../subworkflows/local/paired/paired_savana' +include { DMR } from '../subworkflows/local/dmr' @@ -111,6 +112,13 @@ workflow LRSOMATIC { params.bed_file = getGenomeAttribute('bed_file') params.savana_contigs = getGenomeAttribute('savana_contigs') params.savana_g1000_vcf = getGenomeAttribute('savana_g1000_vcf') + // An explicit --dmr_cpg_islands_bed/--dmr_gencode_gene_bed wins over the per-genome + // default (getGenomeAttribute alone would silently clobber it, same as every other + // getGenomeAttribute-assigned param above -- these two are the only ones a user is + // expected to override directly, since GRCh38/CHM13 are the only genomes with a built-in + // default and any other genome needs its own values). + params.dmr_cpg_islands_bed = params.dmr_cpg_islands_bed ?: getGenomeAttribute('dmr_cpg_islands_bed') + params.dmr_gencode_gene_bed = params.dmr_gencode_gene_bed ?: getGenomeAttribute('dmr_gencode_gene_bed') params.vep_genome = getGenomeAttribute('vep_genome') params.vep_species = getGenomeAttribute('vep_species') params.sigprofiler_genome = getGenomeAttribute('sigprofiler_genome') @@ -776,6 +784,36 @@ workflow LRSOMATIC { : ch_index_minimap // ch_modkit_input: [meta, bam, bai] -- BAM to pile up; meta.type selects the publish directory MODKIT_PILEUP(ch_modkit_input, ch_fasta, ch_fai, [[:],[]]) + + // + // SUBWORKFLOW: DMR (label: process_medium) + // Differential methylation between a sample's two haplotypes. Only possible when + // --modkit_phased produced hp1/hp2 bedMethyl in the first place; --skip_dmr additionally + // turns it off on top of that. dmr_cpg_islands_bed/dmr_gencode_gene_bed have a built-in + // default only for --genome GRCh38/CHM13 (see getGenomeAttribute() above) -- any other + // genome needs both supplied explicitly, and DMR is skipped with a warning rather than + // the run failing when they are not (e.g. a custom --fasta with no --genome, or a + // --genome this pipeline doesn't ship DMR references for). + if (params.modkit_phased && !params.skip_dmr) { + if (params.dmr_cpg_islands_bed && params.dmr_gencode_gene_bed) { + DMR ( + MODKIT_PILEUP.out.bedgz, + ch_fasta, + ch_fai, + file(params.dmr_cpg_islands_bed, checkIfExists: true), + file(params.dmr_gencode_gene_bed, checkIfExists: true) + ) + ch_versions = ch_versions.mix(DMR.out.versions) + } + else { + log.warn( + "--modkit_phased is set but DMR is being skipped: --dmr_cpg_islands_bed/" + + "--dmr_gencode_gene_bed have no default for --genome '${params.genome}' " + + "(only GRCh38/CHM13 ship one). Pass both explicitly to enable DMR calling, " + + "or set --skip_dmr to silence this warning." + ) + } + } } // Prepare phased VCFs for VEP: add empty 'extra' list required by ENSEMBLVEP_VEP