diff --git a/PRO-seq/PRO-seq.snakefile b/PRO-seq/PRO-seq.snakefile index 1f166157e9dce97f537f2160b896197df6da88dd..a3aef8791812295c447b7cfe57314799a1af44bc 100644 --- a/PRO-seq/PRO-seq.snakefile +++ b/PRO-seq/PRO-seq.snakefile @@ -55,6 +55,8 @@ from libworkflows import read_htseq_counts, sum_htseq_counts from libworkflows import read_feature_counts, sum_feature_counts from smincludes import rules as irules +alignment_settings = {"bowtie2": ""}, + #TRIMMERS = ["cutadapt", "fastx_clipper"] TRIMMERS = ["cutadapt"] COUNTERS = ["feature_count"] @@ -312,6 +314,7 @@ rule map_on_genome: sam = temp(OPJ(output_dir, "{trimmer}", aligner, "mapped_C_elegans", "{lib}_{rep}_{type}_on_C_elegans.sam")), nomap_fastq = OPJ(output_dir, "{trimmer}", aligner, "not_mapped_C_elegans", "{lib}_{rep}_{type}_unmapped_on_C_elegans.fastq.gz"), params: + aligner = aligner, index = index, settings = "", message: diff --git a/small_RNA-seq/small_RNA-seq.snakefile b/small_RNA-seq/small_RNA-seq.snakefile index a37b09f033deb451e960e167cadd6c021b1de8d7..0f6ff18bff87d9ecf1dbc3a6ecb235158bcdda7e 100644 --- a/small_RNA-seq/small_RNA-seq.snakefile +++ b/small_RNA-seq/small_RNA-seq.snakefile @@ -129,6 +129,7 @@ from libworkflows import filter_combinator, sum_feature_counts, sum_htseq_counts from smincludes import rules as irules strip = str.strip +alignment_settings = {"bowtie2": "-L 6 -i S,1,0.8 -N 0"}, # Positions in small RNA sequences for which to analyse nucleotide distribution #POSITIONS = ["first", "last"] @@ -733,8 +734,9 @@ rule map_on_genome: sam = temp(OPJ(output_dir, aligner, "mapped_C_elegans", "{lib}_{rep}", "%s_on_C_elegans.sam" % size_selected)), nomap_fastq = OPJ(output_dir, aligner, "not_mapped_C_elegans", "{lib}_{rep}_%s_unmapped_on_C_elegans.fastq.gz" % size_selected), params: + aligner = aligner, index = index, - settings = "-L 6 -i S,1,0.8 -N 0", + settings = alignment_settings[aligner], message: "Mapping {wildcards.lib}_{wildcards.rep}_%s on C. elegans genome." % size_selected benchmark: diff --git a/snakemake_wrappers/map_on_genome/wrapper.py b/snakemake_wrappers/map_on_genome/wrapper.py index 13ee296a8efd9791f345625ffe91c5ae28f6a157..86e8bb2497374bbc334d96af1b27cfa8cdd923ba 100644 --- a/snakemake_wrappers/map_on_genome/wrapper.py +++ b/snakemake_wrappers/map_on_genome/wrapper.py @@ -1,11 +1,29 @@ from snakemake.shell import shell -cmd = """ -genome_dir="${{HOME}}/Genomes" -genome="C_elegans" -cmd="bowtie2 --seed 123 -t {snakemake.params.settings} --mm -x {snakemake.params.index} -U {snakemake.input.fastq} --no-unal --un-gz {snakemake.output.nomap_fastq} -S {snakemake.output.sam}" - echo ${{cmd}} > {snakemake.log.log} - eval ${{cmd}} 1>> {snakemake.log.log} 2> {snakemake.log.err} +def mapping_command(aligner): + """This function returns the shell commands to run given the *aligner*.""" + if aligner == "hisat2": + shell_commands = """ +cmd="hisat2 -p {snakemake.threads} --dta --seed 123 -t {snakemake.params.settings} --mm -x {snakemake.params.index} -U {snakemake.input.fastq} --no-unal --un-gz {snakemake.output.nomap_fastq} -S {snakemake.output.sam}" +echo ${{cmd}} 1> {log.log} +eval ${{cmd}} 1>> {log.log} 2> {log.err} """ + elif aligner == "bowtie2": + shell_commands = """ +cmd="bowtie2 -p {snakemake.threads} --seed 123 -t {snakemake.params.settings} --mm -x {snakemake.params.index} -U {snakemake.input.fastq} --no-unal --un-gz {snakemake.output.nomap_fastq} -S {snakemake.output.sam}" +echo ${{cmd}} > {snakemake.log.log} +eval ${{cmd}} 1>> {snakemake.log.log} 2> {snakemake.log.err} +""" + elif aligner == "crac": + shell_commands = """ +cmd="crac --nb-threads {snakemake.threads} --summary %s --all %s {snakemake.params.settings} -i {snakemake.params.index} -r {snakemake.input.fastq} --sam {snakemake.output.sam}" +echo ${{cmd}} > {snakemake.log.log} +eval ${{cmd}} 1>> {snakemake.log.log} 2> {snakemake.log.err} +# TODO: extract non mappers from the sam output +> {snakemake.output.nomap_fastq} +""" % (snakemake.output.nomap_fastq.split("_unmapped")[0] + "_summary.txt", snakemake.output.nomap_fastq.split("_unmapped")[0] + "_output.txt") + else: + raise NotImplementedError("%s is not yet handled." % aligner) + return shell_commands -shell(cmd) +shell(mapping_command(snakemake.params.aligner))