HISAT2 ALIGN
Map reads with hisat2.
URL: http://daehwankimlab.github.io/hisat2
Example
This wrapper can be used in the following way:
rule hisat2_align:
input:
reads=["reads/{sample}_R1.fastq", "reads/{sample}_R2.fastq"],
idx="index/",
output:
"mapped/{sample}.bam",
log:
"logs/hisat2_align_{sample}.log",
params:
extra="",
threads: 2
wrapper:
"v4.6.0/bio/hisat2/align"
Note that input, output and log file paths can be chosen freely.
When running with
snakemake --use-conda
the software dependencies will be automatically deployed into an isolated environment before execution.
Notes
The -S flag must not be used since output is already directly piped to samtools for compression.
The –threads/-p flag must not be used since threads is set separately via the snakemake threads directive.
The wrapper does not yet handle SRA input accessions.
No reference index files checking is done since the actual number of files may differ depending on the reference sequence size. This is also why the index is supplied in the params directive instead of the input directive.
Software dependencies
hisat2=2.2.1
samtools=1.21
Input/Output
Input:
reads
: either 1 or 2 FASTQ files with reads
Output:
bam file with mapped reads
Params
idx
: prefix of index file path (required)extra
: additional parameters
Code
__author__ = "Wibowo Arindrarto"
__copyright__ = "Copyright 2016, Wibowo Arindrarto"
__email__ = "bow@bow.web.id"
__license__ = "BSD"
import os
from pathlib import Path
from snakemake.shell import shell
# Placeholder for optional parameters
extra = snakemake.params.get("extra", "")
# Run log
log = snakemake.log_fmt_shell()
# Input file wrangling
reads = snakemake.input.get("reads")
if isinstance(reads, str):
input_flags = "-U {0}".format(reads)
elif len(reads) == 1:
input_flags = "-U {0}".format(reads[0])
elif len(reads) == 2:
input_flags = "-1 {0} -2 {1}".format(*reads)
else:
raise RuntimeError(
"Reads parameter must contain at least 1 and at most 2" " input files."
)
ht2_files = Path(snakemake.input.idx).glob("*.ht2")
idx_prefix = os.path.commonprefix(list(ht2_files)).rstrip(".")
# Executed shell command
shell(
"(hisat2 {extra} "
"--threads {snakemake.threads} "
" -x {idx_prefix} {input_flags} "
" | samtools view -Sbh -o {snakemake.output[0]} -) "
" {log}"
)