feat: enable gzip output for bedtools wrappers

bio/bedtools/bamtobed/environment.yaml

@@ -4,3 +4,4 @@ channels:
   - nodefaults
 dependencies:
   - bedtools =2.31.1
+  - htslib =1.21

bio/bedtools/bamtobed/test/Snakefile

@@ -1,3 +1,4 @@
+
 rule bamtobed:
     input:
         "{sample}.bam",
@@ -9,3 +10,15 @@ rule bamtobed:
         extra="-bedpe",  # optional parameters
     wrapper:
         "master/bio/bedtools/bamtobed"
+
+rule bamtobed_gz:
+    input:
+        "{sample}.bam",
+    output:
+        "{sample}.bed.gz",
+    log:
+        "logs/bamtobed/{sample}.gz.log",
+    params:
+        extra="-bedpe",  # optional parameters
+    wrapper:
+        "master/bio/bedtools/bamtobed"

bio/bedtools/bamtobed/wrapper.py

@@ -2,18 +2,19 @@
 __copyright__ = "Copyright 2022, Filipe G. Vieira"
 __license__ = "MIT"
 
-
-from snakemake.shell import shell
-
+import snakemake  # type: ignore
+from snakemake.shell import shell  # type: ignore
 
 log = snakemake.log_fmt_shell(stdout=False, stderr=True)
 extra = snakemake.params.get("extra", "")
 
+compress = "| bgzip" if snakemake.output[0].endswith(".gz") else ""
 
 shell(
     "(bamToBed"
     " {extra}"
     " -i {snakemake.input[0]}"
+    " {compress}"
     " > {snakemake.output[0]}"
     ") {log}"
 )

bio/bedtools/complement/environment.yaml

@@ -4,3 +4,4 @@ channels:
   - nodefaults
 dependencies:
   - bedtools =2.31.1
+  - htslib =1.21

bio/bedtools/complement/test/Snakefile

@@ -12,6 +12,20 @@ rule bedtools_complement_bed:
     wrapper:
         "master/bio/bedtools/complement"
 
+rule bedtools_complement_bed_gz:
+    input:
+        in_file="a.bed",
+        genome="dummy.genome"
+    output:
+        "results/bed-complement/a.complement.bed.gz"
+    params:
+        ## Add optional parameters
+        extra="-L"
+    log:
+        "logs/a.complement.bed.gz.log"
+    wrapper:
+        "master/bio/bedtools/complement"
+
 rule bedtools_complement_vcf:
     input:
         in_file="a.vcf",

bio/bedtools/complement/wrapper.py

@@ -7,13 +7,16 @@
 
 extra = snakemake.params.get("extra", "")
 
+compress = "| bgzip" if snakemake.output[0].endswith(".gz") else ""
+
 log = snakemake.log_fmt_shell(stdout=True, stderr=True)
 
 shell(
     "(bedtools complement"
     " {extra}"
     " -i {snakemake.input.in_file}"
     " -g {snakemake.input.genome}"
+    " {compress}"
     " > {snakemake.output[0]})"
     " {log}"
 )

bio/bedtools/coveragebed/environment.yaml

@@ -4,3 +4,4 @@ channels:
   - nodefaults
 dependencies:
   - bedtools =2.31.1
+  - htslib =1.21

bio/bedtools/coveragebed/wrapper.py

@@ -3,7 +3,6 @@
 __email__ = "[email protected]"
 __license__ = "MIT"
 
-
 from snakemake.shell import shell
 
 shell.executable("bash")
@@ -17,6 +16,8 @@
 
 output_file = snakemake.output[0]
 
+compress = "| bgzip" if snakemake.output[0].endswith(".gz") else ""
+
 if not isinstance(output_file, str) and len(snakemake.output) != 1:
     raise ValueError("Output should be one file: " + str(output_file) + "!")
 
@@ -25,6 +26,7 @@
     " -a {input_a}"
     " -b {input_b}"
     " {extra_params}"
+    " {compress}"
     " > {output_file}"
     " {log}"
 )

bio/bedtools/genomecov/environment.yaml

@@ -4,3 +4,4 @@ channels:
   - nodefaults
 dependencies:
   - bedtools =2.31.1
+  - htslib =1.21

bio/bedtools/genomecov/wrapper.py

@@ -19,10 +19,13 @@
         input_file = "-i " + snakemake.input.get("bed")
         genome = "-g " + snakemake.input.get("ref")
 
+compress = "| bgzip" if snakemake.output[0].endswith(".gz") else ""
+
 shell(
     "(genomeCoverageBed"
     " {snakemake.params}"
     " {input_file}"
     " {genome}"
+    " {compress}"
     " > {snakemake.output[0]}) {log}"
 )

bio/bedtools/intersect/environment.yaml

@@ -4,3 +4,4 @@ channels:
   - nodefaults
 dependencies:
   - bedtools =2.31.1
+  - htslib =1.21

bio/bedtools/intersect/test/Snakefile

@@ -1,4 +1,4 @@
-rule bedtools_merge:
+rule bedtools_intersect:
     input:
         left="A.bed",
         right="B.bed"
@@ -11,3 +11,18 @@ rule bedtools_merge:
         "logs/intersect/A_B.log"
     wrapper:
         "master/bio/bedtools/intersect"
+
+
+rule bedtools_intersect_gz:
+    input:
+        left="A.bed",
+        right="B.bed"
+    output:
+        "A_B.intersected.bed.gz"
+    params:
+        ## Add optional parameters 
+        extra="-wa -wb" ## In this example, we want to write original entries in A and B for each overlap.
+    log:
+        "logs/intersect/A_B_gz.log"
+    wrapper:
+        "master/bio/bedtools/intersect"

bio/bedtools/intersect/wrapper.py

@@ -3,17 +3,21 @@
 __email__ = "[email protected]"
 __license__ = "MIT"
 
-from snakemake.shell import shell
+import snakemake  # type: ignore
+from snakemake.shell import shell  # type: ignore
 
 ## Extract arguments
 extra = snakemake.params.get("extra", "")
 log = snakemake.log_fmt_shell(stdout=True, stderr=True)
 
+compress = "| bgzip" if snakemake.output[0].endswith(".gz") else ""
+
 shell(
     "(bedtools intersect"
     " {extra}"
     " -a {snakemake.input.left}"
     " -b {snakemake.input.right}"
+    " {compress}"
     " > {snakemake.output})"
     " {log}"
 )

bio/bedtools/merge/environment.yaml

@@ -4,3 +4,4 @@ channels:
   - nodefaults
 dependencies:
   - bedtools =2.31.1
+  - htslib =1.21

bio/bedtools/merge/wrapper.py

@@ -9,6 +9,9 @@
 extra = snakemake.params.get("extra", "")
 
 log = snakemake.log_fmt_shell(stdout=True, stderr=True)
+
+compress = "| bgzip" if snakemake.output[0].endswith(".gz") else ""
+
 if len(snakemake.input) > 1:
     if all(f.endswith(".gz") for f in snakemake.input):
         cat = "zcat"
@@ -20,14 +23,17 @@
         "({cat} {snakemake.input} | "
         "sort -k1,1 -k2,2n | "
         "bedtools merge {extra} "
-        "-i stdin > {snakemake.output}) "
+        " -i stdin "
+        " {compress}"
+        "> {snakemake.output[0]}) "
         " {log}"
     )
 else:
     shell(
         "( bedtools merge"
         " {extra}"
         " -i {snakemake.input}"
-        " > {snakemake.output})"
+        " {compress}"
+        " > {snakemake.output[0]})"
         " {log}"
     )

bio/bedtools/slop/environment.yaml

@@ -4,3 +4,4 @@ channels:
   - nodefaults
 dependencies:
   - bedtools =2.31.1
+  - htslib =1.21

bio/bedtools/slop/wrapper.py

@@ -10,11 +10,14 @@
 
 log = snakemake.log_fmt_shell(stdout=True, stderr=True)
 
+compress = "| bgzip" if snakemake.output[0].endswith(".gz") else ""
+
 shell(
     "(bedtools slop"
     " {extra}"
     " -i {snakemake.input[0]}"
     " -g {snakemake.params.genome}"
+    " {compress}"
     " > {snakemake.output})"
     " {log}"
 )

bio/bedtools/sort/environment.yaml

@@ -4,3 +4,4 @@ channels:
   - nodefaults
 dependencies:
   - bedtools =2.31.1
+  - htslib =1.21

bio/bedtools/sort/wrapper.py

@@ -11,6 +11,8 @@
 
 log = snakemake.log_fmt_shell(stdout=True, stderr=True)
 
+compress = "| bgzip" if snakemake.output[0].endswith(".gz") else ""
+
 if genome:
     extra += " -g {}".format(genome)
 elif faidx:
@@ -20,6 +22,7 @@
     "(bedtools sort"
     " {extra}"
     " -i {snakemake.input.in_file}"
+    " {compress}"
     " > {snakemake.output[0]})"
     " {log}"
 )

bio/bedtools/split/environment.yaml

@@ -4,3 +4,4 @@ channels:
   - nodefaults
 dependencies:
   - bedtools =2.31.1
+  - htslib =1.21

bio/bedtools/split/wrapper.py

@@ -10,6 +10,8 @@
 
 n_subfiles = len(snakemake.output)
 
+compress = "| bgzip" if snakemake.output[0].endswith(".gz") else ""
+
 with tempfile.TemporaryDirectory() as tmpdir:
     shell(
         "bedtools split"
@@ -23,4 +25,4 @@
     for i in range(n_subfiles):
         out_tmp = f"{tmpdir}/out.{i+1:05d}.bed"
         out = snakemake.output[i]
-        shell("cat {out_tmp} > {out}")
+        shell("cat {out_tmp} {compress} > {out}")

test_wrappers.py

@@ -1909,7 +1909,7 @@ def test_bcftools_view_uncompressed_bcf(run):
 def test_bedtools_bamtobed(run):
     run(
         "bio/bedtools/bamtobed",
-        ["snakemake", "--cores", "1", "a.bed", "--use-conda", "-F"],
+        ["snakemake", "--cores", "1", "a.bed", "a.bed.gz", "--use-conda", "-F"],
     )
 
 
@@ -1936,10 +1936,16 @@ def test_bedtools_complement(run):
             "--cores",
             "1",
             "results/bed-complement/a.complement.bed",
+            "results/bed-complement/a.complement.bed.gz",
             "results/vcf-complement/a.complement.vcf",
             "--use-conda",
             "-F",
         ],
+        compare_results_with_expected={
+            "results/bed-complement/a.complement.bed": "expected/bed-complement/a.complement.bed",
+            "results/bed-complement/a.complement.bed.gz": "expected/bed-complement/a.complement.bed.gz",
+            "results/vcf-complement/a.complement.vcf": "expected/vcf-complement/a.complement.vcf"
+        }
     )
 
 
@@ -1954,6 +1960,7 @@ def test_bedtools_sort(run):
             "results/bed-sorted/a.sorted_by_file.bed",
             "results/vcf-sorted/a.sorted_by_file.vcf",
             "--use-conda",
+
             "-F",
         ],
     )
@@ -1977,7 +1984,19 @@ def test_bedtools_split(run):
 def test_bedtools_intersect(run):
     run(
         "bio/bedtools/intersect",
-        ["snakemake", "--cores", "1", "A_B.intersected.bed", "--use-conda", "-F"],
+        [ 
+            "snakemake", 
+            "--cores", 
+            "1", 
+            "A_B.intersected.bed",
+            "A_B.intersected.bed.gz",
+            "--use-conda",
+             "-F"
+        ],
+         compare_results_with_expected={
+            "A_B.intersected.bed": "expected/A_B.intersected.bed",
+            "A_B.intersected.bed.gz": "expected/A_B.intersected.bed.gz",
+        }
     )