Add fasten and some of lskatz improvements

.github/workflows/benchmark.yml

@@ -0,0 +1,36 @@
+on: push
+
+jobs:
+  benchmark:
+    name: Benchmark 
+    runs-on: ubuntu-latest
+    strategy:
+      fail-fast: true
+    defaults:
+      run:
+        shell: bash -el {0}
+    steps:
+        - uses: actions/checkout@v4
+        - uses: conda-incubator/setup-miniconda@v3
+          with:
+            python-version: 3.12
+            environment-file: workflow/envs/github.yaml
+        - name: test env
+          run: |
+            conda info
+            conda list
+            which snakemake
+            snakemake --help
+            snakemake --version
+        - name: benchmarks
+          run: |
+            sed -i 's/threads: 8/threads: 2/' config/config.yaml
+            sed -i 's/replicates: 3/replicates: 1/' config/config.yaml
+            snakemake --use-conda --scheduler greedy --cores 2
+        - uses: actions/upload-artifact@v4
+          with:
+            name: benchmarks
+            path: |
+              output/plots/*
+              output/grades/*
+

.gitignore

@@ -0,0 +1,2 @@
+.snakemake
+output

README.md

@@ -9,15 +9,18 @@ The following tools are compared:
  - [Cutadapt](https://cutadapt.readthedocs.io/en/stable/)
  - [Fastp](https://github.com/OpenGene/fastp)
  - [Trimmomatic](http://www.usadellab.org/cms/?page=trimmomatic)
+ - [Fasten](https://github.com/lskatz/fasten)
 
 ## Fake adapters, "gruseq"
 The comparison adds fake truseq adapters, "gruseq", provided by Brian Bushnell,
 downloaded from:
 http://seqanswers.com/forums/attachment.php?attachmentid=2993&d=1398383571
 
 ## Test data
-A single sample is downloaded from SRA. Feel free to replace it with whatever
-you want. 
+A subset of a single sample is from SRA is provided in the repo.
+Comment/uncomment the `download_FASTQ` setting in the `config/config.yaml` to
+download the complete sample from SRA instead. Feel free to replace the URL
+with whatever FASTQ file you want.
 
 ## Results
 These plots are produced by the workflow, one for benchmarking results
@@ -32,6 +35,6 @@ and one plot with adapter removal performance.
 # Running
 Run the benchmarking workflow with something like:
 ```bash
-snakemake --use-conda --jobs 10
+snakemake --use-conda --cores 8
 ```
 

config/config.yaml

@@ -0,0 +1,6 @@
+# Config file for Snakemake workflow adapter_benchmark
+
+download_FASTQ: False
+#download_FASTQ: "ftp://ftp.sra.ebi.ac.uk/vol1/fastq/SRR921/004/SRR9218144/SRR9218144.fastq.gz"
+replicates: 3
+threads: 8

workflow/Snakefile

@@ -6,35 +6,53 @@
 #   http://seqanswers.com/forums/attachment.php?attachmentid=2993&d=1398383571
 # Fredrik Boulund 2019
 
+import random
+import gzip
+import os
+
+configfile: "config/config.yaml"
+
 TOOLS=[
     "cutadapt",
     "trimmomatic",
     "bbduk",
     "fastp",
     "adapterremoval",
+    "fasten",
 ]
-FASTQ="ftp://ftp.sra.ebi.ac.uk/vol1/fastq/SRR921/004/SRR9218144/SRR9218144.fastq.gz"
 
 rule all:
     input:
         expand("output/grades/{tool}.grade.txt", tool=TOOLS),
         expand("output/plots/{plot}.pdf", plot=["benchmarks", "grades"]),
-
 
-rule download_fastq:
-    output:
-        "output/data.fastq.gz"
-    conda: "envs/main.yaml"
-    params:
-        fastq=FASTQ
-    shell:
-        """
-        wget --output-document {output} {params.fastq}
-        """
+
+## Download data if requested by config, else use provided example testing data.
+example_data_path = "../data/SRR9218144.888888.fastq.gz"
+input_data_path = "output/data.fastq.gz"
+if config["download_FASTQ"]:
+    rule download_fastq:
+        output:
+            input_data_path
+        conda: "envs/main.yaml"
+        params:
+            fastq=config["download_FASTQ"]
+        shell:
+            """
+            wget --output-document {output} {params.fastq}
+            """
+else:
+    if os.path.isfile(input_data_path):
+        print(f"DEBUG: download_FASTQ is set to False and '{input_data_path}' already exists, doing nothing.")
+    else:
+        print(f"INFO: Symlinking example FASTQ {example_data_path} to {input_data_path}.")
+        os.makedirs("output", exist_ok=True)
+        os.symlink(example_data_path, input_data_path)
+
 
 rule add_adapters:
     input:
-        fastq=rules.download_fastq.output,
+        fastq="output/data.fastq.gz",
         adapters="gruseq.fa",
     output:
         "output/dirty.fq"
@@ -60,9 +78,9 @@ rule cutadapt:
         "output/processed/cutadapt.fq"
     log: "output/logs/cutadapt.log"
     conda: "envs/cutadapt.yaml"
-    threads: 10
+    threads: config["threads"]
     benchmark:
-        repeat("output/benchmarks/cutadapt.benchmark.txt", 3)
+        repeat("output/benchmarks/cutadapt.benchmark.txt", config["replicates"])
     shell:
         """
         cutadapt \
@@ -105,9 +123,9 @@ rule trimmomatic:
         "output/processed/trimmomatic.fq"
     conda: "envs/trimmomatic.yaml"
     log: "output/logs/trimmomatic.log"
-    threads: 10
+    threads: config["threads"]
     benchmark:
-        repeat("output/benchmarks/trimmomatic.benchmark.txt", 3)
+        repeat("output/benchmarks/trimmomatic.benchmark.txt", config["replicates"])
     shell:
         """
         trimmomatic \
@@ -130,9 +148,9 @@ rule bbduk:
         "output/processed/bbduk.fq"
     conda: "envs/main.yaml"
     log: "output/logs/bbduk.log"
-    threads: 10
+    threads: config["threads"]
     benchmark:
-        repeat("output/benchmarks/bbduk.benchmark.txt", 3)
+        repeat("output/benchmarks/bbduk.benchmark.txt", config["replicates"])
     shell:
         """
         bbduk.sh \
@@ -157,9 +175,9 @@ rule fastp:
         json="output/processed/fastp.json",
     conda: "envs/fastp.yaml"
     log: "output/logs/fastp.log"
-    threads: 10
+    threads: config["threads"]
     benchmark:
-        repeat("output/benchmarks/fastp.benchmark.txt", 3)
+        repeat("output/benchmarks/fastp.benchmark.txt", config["replicates"])
     shell:
         """
         fastp \
@@ -181,10 +199,10 @@ rule adapterremoval:
         fq="output/processed/adapterremoval.fq",
     conda: "envs/adapterremoval.yaml"
     log: "output/logs/adapterremoval.log"
-    threads: 10
+    threads: config["threads"]
     shadow: "shallow"
     benchmark:
-        repeat("output/benchmarks/adapterremoval.benchmark.txt", 3)
+        repeat("output/benchmarks/adapterremoval.benchmark.txt", config["replicates"])
     shell:
         """
         AdapterRemoval \
@@ -195,6 +213,27 @@ rule adapterremoval:
             2> {log}
         """
 
+rule fasten:
+    input:
+        fastq=rules.add_adapters.output,
+        adapters="gruseq.fa",
+    output:
+        "output/processed/fasten.fq"
+    conda: "envs/fasten.yaml"
+    log: "output/logs/fasten.log"
+    threads: config["threads"]
+    benchmark:
+        repeat("output/benchmarks/fasten.benchmark.txt", config["replicates"])
+    shell:
+        """
+        fasten_trim \
+          --adapterseqs {input.adapters} \
+          --numcpus {threads} \
+          < {input.fastq} \
+          > {output} \
+          2> {log};
+        """
+
 rule grade:
     input:
         "output/processed/{tool}.fq"
@@ -227,3 +266,4 @@ rule plot:
             {input.grades}
         """
 
+

workflow/envs/adapterremoval.yaml

@@ -2,6 +2,5 @@ name: adapterremoval
 channels:
     - conda-forge
     - bioconda
-    - default
 dependencies:
-    - adapterremoval=2.3.3
+    - adapterremoval=2.3.4

workflow/envs/cutadapt.yaml

@@ -2,6 +2,5 @@ name: cutadapt
 channels:
     - conda-forge
     - bioconda
-    - default
 dependencies:
-    - cutadapt=4.7
+    - cutadapt=5.1

workflow/envs/fasten.yaml

@@ -0,0 +1,6 @@
+name: fasten
+channels:
+    - conda-forge
+    - bioconda
+dependencies:
+    - fasten=0.8.7

workflow/envs/fastp.yaml

@@ -2,6 +2,5 @@ name: fastp
 channels:
     - conda-forge
     - bioconda
-    - default
 dependencies:
-    - fastp=0.23.4
+    - fastp=0.25.0

workflow/envs/github.yaml

@@ -0,0 +1,6 @@
+name: smk
+channels:
+    - conda-forge
+    - bioconda
+dependencies:
+    - snakemake-minimal

workflow/envs/main.yaml

@@ -2,7 +2,6 @@ name: adaptertest
 channels:
     - conda-forge
     - bioconda
-    - default
 dependencies:
     - wget=1.20.1
     - bbmap=38.62

workflow/envs/trimmomatic.yaml

@@ -2,6 +2,5 @@ name: trimmomatic
 channels:
     - conda-forge
     - bioconda
-    - default
 dependencies:
     - trimmomatic=0.39