feat: add wrapper for metaDMG

bio/metadmg/compressbam/environment.yaml

@@ -0,0 +1,6 @@
+channels:
+  - conda-forge
+  - bioconda
+  - nodefaults
+dependencies:
+  - metadmg =0.4

bio/metadmg/compressbam/meta.yaml

@@ -0,0 +1,12 @@
+name: CompressBam
+url: https://github.com/metaDMG-dev/metaDMG-cpp
+description: metaDMG-cpp is a fast and efficient method for estimating mutation and damage rates in ancient DNA data
+authors:
+  - Filipe G. Vieira
+input:
+  - aln: SAM/BAM/CRAM file
+  - ref: reference file in FASTA format (mandatory if CRAM input)
+output:
+  - output BAM file
+params:
+  - extra: additional program arguments.

bio/metadmg/compressbam/test/Snakefile

@@ -0,0 +1,15 @@
+
+rule compressbam:
+    input:
+        aln="{sample}.bam",
+    output:
+        "results/compressbam/{sample}.bam",
+    log:
+        "logs/compressbam/{sample}.log",
+    params:
+        extra="",  # optional
+    threads: 4
+    resources:
+        mem_mb=1024,
+    wrapper:
+        "master/bio/metadmg/compressbam"

bio/metadmg/compressbam/wrapper.py

@@ -0,0 +1,22 @@
+__author__ = "Filipe G. Vieira"
+__copyright__ = "Copyright 2023, Filipe G. Vieira"
+__license__ = "MIT"
+
+
+import tempfile
+from pathlib import Path
+from snakemake.shell import shell
+
+
+extra = snakemake.params.get("extra", "")
+log = snakemake.log_fmt_shell(stdout=True, stderr=True)
+
+
+ref = snakemake.input.get("ref", "")
+if ref:
+    ref = f"--ref {ref}"
+
+
+shell(
+    "compressbam --threads {snakemake.threads} --input {snakemake.input.aln} {ref} {extra} --output {snakemake.output[0]} {log}"
+)

bio/metadmg/dfit/environment.yaml

@@ -0,0 +1,6 @@
+channels:
+  - conda-forge
+  - bioconda
+  - nodefaults
+dependencies:
+  - metadmg =0.4

bio/metadmg/dfit/meta.yaml

@@ -0,0 +1,19 @@
+name: metaDMG dfit
+url: https://github.com/metaDMG-dev/metaDMG-cpp
+description: metaDMG-cpp is a fast and efficient method for estimating mutation and damage rates in ancient DNA data
+authors:
+  - Filipe G. Vieira
+input:
+  - aln: SAM/BAM/CRAM file
+  - names: taxonomy file "names.dmp"
+  - nodes: taxonomy file "nodes.dmp"
+  - acc2taxid: TSV with correspondence between accesions and taxa IDs
+output:
+  - dmg: path to TSV file containing counts of mismatchs conditional on strand and cycle.
+  - lca: path to TSV file with LCA results.
+  - stat: path to TSV file with general stats.
+params:
+  - extra: additional program arguments.
+notes: |
+  * Input BAM file has to be sorted by read name.
+  * More information about output formats in https://github.com/metaDMG-dev/metaDMG-cpp/blob/master/doc/formats.pdf

bio/metadmg/dfit/test/Snakefile

@@ -0,0 +1,20 @@
+
+rule metadmg_dfit:
+    input:
+        dmg="{sample}.bam",
+        names="names.dmp.gz",
+        nodes="nodes.dmp.gz",
+        lca_stats="{sample}.stat",
+    output:
+        dfit="results/dfit/{sample}.out.gz",
+        stats_dfit="stats/dfit/{sample}.dfit.tsv.gz",
+        stats_boot="stats/dfit/{sample}.boot.tsv.gz",
+    log:
+        "logs/dfit/{sample}.log",
+    params:
+        extra="--nopt 10 --doboot 1 --nbootstrap 20 --showfits 2 --lib ds",
+    threads: 1
+    resources:
+        mem_mb=1024,
+    wrapper:
+        "master/bio/metadmg/dfit"

bio/metadmg/dfit/wrapper.py

@@ -0,0 +1,36 @@
+__author__ = "Filipe G. Vieira"
+__copyright__ = "Copyright 2023, Filipe G. Vieira"
+__license__ = "MIT"
+
+
+import tempfile
+from snakemake.shell import shell
+
+
+extra = snakemake.params.get("extra", "")
+log = snakemake.log_fmt_shell(stdout=True, stderr=True)
+
+
+names = snakemake.input.get("names", "")
+if names:
+    names = f"--names {names}"
+
+nodes = snakemake.input.get("nodes", "")
+if nodes:
+    nodes = f"--nodes {nodes}"
+
+lca_stats = snakemake.input.get("lca_stats", "")
+if lca_stats:
+    lca_stats = f"--lcastat {lca_stats}"
+
+
+with tempfile.TemporaryDirectory() as tmpdir:
+    shell(
+        "metaDMG-cpp dfit {snakemake.input.dmg} {names} {nodes} {lca_stats} {extra} --out_prefix {tmpdir}/out {log}"
+    )
+
+    for output in snakemake.output:
+        for ext in [".dfit.txt.gz", ".dfit.stat.txt.gz", ".boot.stat.txt.gz"]:
+            if output.endswith(ext):
+                shell("cat {tmpdir}/out{ext} > {output}")
+                continue

bio/metadmg/getdamage/environment.yaml

@@ -0,0 +1,6 @@
+channels:
+  - conda-forge
+  - bioconda
+  - nodefaults
+dependencies:
+  - metadmg =0.4

bio/metadmg/getdamage/meta.yaml

@@ -0,0 +1,17 @@
+name: metaDMG getdamage
+url: https://github.com/metaDMG-dev/metaDMG-cpp
+description: metaDMG-cpp is a fast and efficient method for estimating mutation and damage rates in ancient DNA data
+authors:
+  - Filipe G. Vieira
+input:
+  - aln: SAM/BAM/CRAM file
+  - ref: reference file in FASTA format (mandatory if CRAM input)
+output:
+  - dmg: path to TSV file containing counts of mismatchs conditional on strand and cycle.
+  - res: path to TSV file with estimates of damage.
+  - stat: path to TSV file with general stats.
+params:
+  - extra: additional program arguments.
+notes: |
+  * Input BAM file has to be sorted by read name.
+  * More information about output formats in https://github.com/metaDMG-dev/metaDMG-cpp/blob/master/doc/formats.pdf

bio/metadmg/getdamage/test/Snakefile

@@ -0,0 +1,17 @@
+
+rule metadmg_getdamage:
+    input:
+        aln="{sample}.bam",
+    output:
+        res="results/getdamage/{sample}.out.gz",
+        dmg="results/getdamage/{sample}.dmg.gz",
+        stats="stats/getdamage/{sample}.tsv",
+    log:
+        "logs/getdamage/{sample}.log",
+    params:
+        extra="--min_length 30 --print_length 30 --run_mode 1",
+    threads: 4
+    resources:
+        mem_mb=1024,
+    wrapper:
+        "master/bio/metadmg/getdamage"

bio/metadmg/getdamage/wrapper.py

@@ -0,0 +1,28 @@
+__author__ = "Filipe G. Vieira"
+__copyright__ = "Copyright 2023, Filipe G. Vieira"
+__license__ = "MIT"
+
+
+import tempfile
+from snakemake.shell import shell
+
+
+extra = snakemake.params.get("extra", "")
+log = snakemake.log_fmt_shell(stdout=True, stderr=True)
+
+
+ref = snakemake.input.get("ref", "")
+if ref:
+    ref = f"--fasta {ref}"
+
+
+with tempfile.TemporaryDirectory() as tmpdir:
+    shell(
+        "metaDMG-cpp getdamage --threads {snakemake.threads} {ref} {extra} --out_prefix {tmpdir}/out {snakemake.input.aln} {log}"
+    )
+
+    for output in snakemake.output:
+        for ext in [".bdamage.gz", ".res.gz", ".stat"]:
+            if output.endswith(ext):
+                shell("cat {tmpdir}/out{ext} > {output}")
+                continue

bio/metadmg/lca/environment.yaml

@@ -0,0 +1,6 @@
+channels:
+  - conda-forge
+  - bioconda
+  - nodefaults
+dependencies:
+  - metadmg =0.4

bio/metadmg/lca/meta.yaml

@@ -0,0 +1,19 @@
+name: metaDMG lca
+url: https://github.com/metaDMG-dev/metaDMG-cpp
+description: metaDMG-cpp is a fast and efficient method for estimating mutation and damage rates in ancient DNA data
+authors:
+  - Filipe G. Vieira
+input:
+  - aln: SAM/BAM/CRAM file
+  - names: taxonomy file "names.dmp"
+  - nodes: taxonomy file "nodes.dmp"
+  - acc2taxid: TSV with correspondence between accesions and taxa IDs
+output:
+  - dmg: path to TSV file containing counts of mismatchs conditional on strand and cycle.
+  - lca: path to TSV file with LCA results.
+  - stat: path to TSV file with general stats.
+params:
+  - extra: additional program arguments.
+notes: |
+  * Input BAM file has to be sorted by read name.
+  * More information about output formats in https://github.com/metaDMG-dev/metaDMG-cpp/blob/master/doc/formats.pdf

bio/metadmg/lca/test/Snakefile

@@ -0,0 +1,20 @@
+
+rule metadmg_lca:
+    input:
+        aln="{sample}.bam",
+        names="names.dmp.gz",
+        nodes="nodes.dmp.gz",
+        acc2taxid="acc2taxid.tsv",
+    output:
+        res="results/lca/{sample}.out.gz",
+        lca="results/lca/{sample}.lca.gz",
+        stats="stats/lca/{sample}.tsv",
+    log:
+        "logs/lca/{sample}.log",
+    params:
+        extra="--sim_score_low 0.95 --sim_score_high 1.0 --min_mapq 30 --how_many 30 --lca_rank genus --fix_ncbi 0",
+    threads: 4
+    resources:
+        mem_mb=1024,
+    wrapper:
+        "master/bio/metadmg/lca"

bio/metadmg/lca/test/acc2taxid.tsv

@@ -0,0 +1,2 @@
+accession	accession.version	taxid	gi
+NC_023100	NC_023100.1	1425170	1

bio/metadmg/lca/wrapper.py

@@ -0,0 +1,23 @@
+__author__ = "Filipe G. Vieira"
+__copyright__ = "Copyright 2023, Filipe G. Vieira"
+__license__ = "MIT"
+
+
+import tempfile
+from snakemake.shell import shell
+
+
+extra = snakemake.params.get("extra", "")
+log = snakemake.log_fmt_shell(stdout=True, stderr=True)
+
+
+with tempfile.TemporaryDirectory() as tmpdir:
+    shell(
+        "metaDMG-cpp lca --threads {snakemake.threads} --bam {snakemake.input.aln} --names {snakemake.input.names} --nodes {snakemake.input.nodes} --acc2tax {snakemake.input.acc2taxid} --temp {tmpdir} {extra} --out_prefix {tmpdir}/out {log}"
+    )
+
+    for output in snakemake.output:
+        for ext in [".bdamage.gz", ".lca.gz", ".stat"]:
+            if output.endswith(ext):
+                shell("cat {tmpdir}/out{ext} > {output}")
+                continue

test.py

@@ -139,6 +139,72 @@ def run(wrapper, cmd, check_log=None):
             os.chdir(origdir)
 
 
+@skip_if_not_modified
+def test_metadmg_getdamage():
+    run(
+        "bio/metadmg/getdamage",
+        [
+            "snakemake",
+            "--cores",
+            "1",
+            "--use-conda",
+            "-F",
+            "results/getdamage/a.out.gz",
+            "results/getdamage/a.dmg.gz",
+            "stats/getdamage/a.tsv",
+        ],
+    )
+
+
+@skip_if_not_modified
+def test_metadmg_lca():
+    run(
+        "bio/metadmg/lca",
+        [
+            "snakemake",
+            "--cores",
+            "1",
+            "--use-conda",
+            "-F",
+            "results/lca/a.out.gz",
+            "results/lca/a.lca.gz",
+            "stats/lca/a.tsv",
+        ],
+    )
+
+
+@skip_if_not_modified
+def test_metadmg_dfit():
+    run(
+        "bio/metadmg/dfit",
+        [
+            "snakemake",
+            "--cores",
+            "1",
+            "--use-conda",
+            "-F",
+            "results/dfit/a.out.gz",
+            "stats/dfit/a.dfit.tsv.gz",
+            "stats/dfit/a.boot.tsv.gz",
+        ],
+    )
+
+
+@skip_if_not_modified
+def test_metadmg_compressbam():
+    run(
+        "bio/metadmg/compressbam",
+        [
+            "snakemake",
+            "--cores",
+            "1",
+            "--use-conda",
+            "-F",
+            "results/compressbam/a.bam",
+        ],
+    )
+
+
 @skip_if_not_modified
 def test_galah():
     run(