adding in devider for whole-segment haplotyping

Author: nrminor

modules/devider.nf

@@ -0,0 +1,30 @@
+process PHASE_READS_WITH_DEVIDER {
+
+    /* */
+
+    tag "${sample_id}"
+    // publishDir
+
+    errorStrategy { task.attempt < 3 ? 'retry' : 'ignore' }
+    maxRetries 2
+
+    cpus 4
+
+    input:
+    tuple val(sample_id), path(reads), path(vcf), path(fasta_ref)
+
+    output:
+    path "${sample_id}_devider"
+
+    script:
+    """
+    devider
+    -b ${reads} \
+    -v ${vcf} \
+    -r ${fasta_ref} \
+    -o ${sample_id}_devider \
+    -t ${task.cpus} \
+    --preset ${devider_preset}
+    """
+
+}

modules/samtools.nf

@@ -180,7 +180,7 @@ process SPLIT_SEGMENTS {
     tuple val(sample_id), path(bam), path(index)
 
     output:
-    path "${sample_id}_*.bam"
+    tuple val(sample_id), path("${sample_id}_*.bam")
 
     shell:
     '''

nextflow.config

@@ -139,6 +139,9 @@ params {
     // nextclade dataset
     nextclade_dataset = null
 
+	// devider haplotype phasing preset
+	devider_preset = "nanopore-r10" // old-long-reads, nanopore-r9, nanopore-r10, hi-fi
+
     // nextclade dataset cache
     nextclade_cache = "$launchDir/work/nextclade_datasets"
 

pixi.lock

@@ -88,6 +88,7 @@ csvtk = { version = ">=0.30.0,<0.31", channel = "bioconda"}
 snpsift = ">=5.2,<6"
 bbmap = ">=39.1,<40"
 bedtools = ">=2.31.1,<3"
+devider = ">=0.0.1,<0.0.2"
 
 [tool.uv]
 dev-dependencies = [

pyproject.toml

@@ -1,17 +1,22 @@
-#!/usr/bin/env nextflow
 
 include { SPLIT_SEGMENTS ; FASTQ_CONVERSION } from "../modules/samtools"
 include { IDENTIFY_HAPLOTYPES } from "../modules/vsearch"
 
 workflow HAPLOTYPING {
     take:
     ch_bams
+    ch_vcfs
+    ch_ref
 
     main:
     SPLIT_SEGMENTS(
         ch_bams
     )
 
+    PHASE_READS_WITH_DEVIDER(
+        SPLIT_SEGMENTS.out.join(ch_vcfs).combine(ch_ref)
+    )
+
     FASTQ_CONVERSION(
         SPLIT_SEGMENTS.out.flatten().map { bam ->
             tuple(file(bam).getSimpleName().split("_")[0..-2].join('_'), bam)

subworkflows/haplotyping.nf

@@ -50,4 +50,7 @@ workflow VARIANTS {
     MERGE_VCF_FILES(
         ANNOTATE_VCF.out.map { _label, vcf -> vcf }.collect()
     )
+
+    emit:
+    ANNOTATE_VCF.out.
 }

subworkflows/variant_calling.nf

@@ -70,7 +70,9 @@ workflow NANOPORE {
         if ( params.primer_bed && Utils.countFastaHeaders(params.refseq) == Utils.countAmplicons(params.primer_bed) ) {
 
             HAPLOTYPING (
-                ALIGNMENT.out
+                ALIGNMENT.out,
+                VARIANTS.out,
+                ch_refseq
             )
 
         }