diff --git a/.github/workflows/benchmark.yml b/.github/workflows/benchmark.yml new file mode 100644 index 0000000..cf0d5eb --- /dev/null +++ b/.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/* + diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..53c185d --- /dev/null +++ b/.gitignore @@ -0,0 +1,2 @@ +.snakemake +output \ No newline at end of file diff --git a/README.md b/README.md index 999afec..c588353 100644 --- a/README.md +++ b/README.md @@ -9,6 +9,7 @@ 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, @@ -16,8 +17,10 @@ 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 ``` diff --git a/config/config.yaml b/config/config.yaml new file mode 100644 index 0000000..a96b26a --- /dev/null +++ b/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 diff --git a/data/SRR9218144.888888.fastq.gz b/data/SRR9218144.888888.fastq.gz new file mode 100644 index 0000000..b545ac2 Binary files /dev/null and b/data/SRR9218144.888888.fastq.gz differ diff --git a/plots/benchmarks.pdf b/plots/benchmarks.pdf index 891d36b..03eec90 100644 Binary files a/plots/benchmarks.pdf and b/plots/benchmarks.pdf differ diff --git a/plots/benchmarks.png b/plots/benchmarks.png index 701a793..1f25d46 100644 Binary files a/plots/benchmarks.png and b/plots/benchmarks.png differ diff --git a/plots/grades.pdf b/plots/grades.pdf index 4dcb769..e57d8f5 100644 Binary files a/plots/grades.pdf and b/plots/grades.pdf differ diff --git a/plots/grades.png b/plots/grades.png index 0e0b9af..9d381d7 100644 Binary files a/plots/grades.png and b/plots/grades.png differ diff --git a/workflow/Snakefile b/workflow/Snakefile index abf7da7..61f1af7 100644 --- a/workflow/Snakefile +++ b/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} """ + diff --git a/workflow/envs/adapterremoval.yaml b/workflow/envs/adapterremoval.yaml index a463b41..4a92cdd 100644 --- a/workflow/envs/adapterremoval.yaml +++ b/workflow/envs/adapterremoval.yaml @@ -2,6 +2,5 @@ name: adapterremoval channels: - conda-forge - bioconda - - default dependencies: - - adapterremoval=2.3.3 + - adapterremoval=2.3.4 diff --git a/workflow/envs/cutadapt.yaml b/workflow/envs/cutadapt.yaml index f5d02ec..513ed3c 100644 --- a/workflow/envs/cutadapt.yaml +++ b/workflow/envs/cutadapt.yaml @@ -2,6 +2,5 @@ name: cutadapt channels: - conda-forge - bioconda - - default dependencies: - - cutadapt=4.7 + - cutadapt=5.1 diff --git a/workflow/envs/fasten.yaml b/workflow/envs/fasten.yaml new file mode 100644 index 0000000..55797d3 --- /dev/null +++ b/workflow/envs/fasten.yaml @@ -0,0 +1,6 @@ +name: fasten +channels: + - conda-forge + - bioconda +dependencies: + - fasten=0.8.7 diff --git a/workflow/envs/fastp.yaml b/workflow/envs/fastp.yaml index f805274..e1863b6 100644 --- a/workflow/envs/fastp.yaml +++ b/workflow/envs/fastp.yaml @@ -2,6 +2,5 @@ name: fastp channels: - conda-forge - bioconda - - default dependencies: - - fastp=0.23.4 + - fastp=0.25.0 diff --git a/workflow/envs/github.yaml b/workflow/envs/github.yaml new file mode 100644 index 0000000..59e89b7 --- /dev/null +++ b/workflow/envs/github.yaml @@ -0,0 +1,6 @@ +name: smk +channels: + - conda-forge + - bioconda +dependencies: + - snakemake-minimal diff --git a/workflow/envs/main.yaml b/workflow/envs/main.yaml index 6724a6c..12fd386 100644 --- a/workflow/envs/main.yaml +++ b/workflow/envs/main.yaml @@ -2,7 +2,6 @@ name: adaptertest channels: - conda-forge - bioconda - - default dependencies: - wget=1.20.1 - bbmap=38.62 diff --git a/workflow/envs/trimmomatic.yaml b/workflow/envs/trimmomatic.yaml index 1f3e90d..178e10c 100644 --- a/workflow/envs/trimmomatic.yaml +++ b/workflow/envs/trimmomatic.yaml @@ -2,6 +2,5 @@ name: trimmomatic channels: - conda-forge - bioconda - - default dependencies: - trimmomatic=0.39