Resfinder improvements, added support for Shigella spp., added XDR Shigella prediction

tasks/gene_typing/task_resfinder.wdl

@@ -6,11 +6,13 @@ task resfinder {
     String samplename
     String? organism # Species in the sample, species should be entered with their full scientific names (e.g. "escherichia coli"), using quotation marks
     Boolean acquired = true # Run resfinder for acquired resistance genes
-    Float min_cov = 0.6 # Minimum (breadth-of) coverage of ResFinder
-    Float threshold = 0.9 # Threshold for identity of ResFinder
-    Boolean point = false # Run pointfinder for chromosomal mutations
+    Float min_cov = 0.5 # Minimum (breadth-of) coverage of ResFinder
+    Float min_id = 0.9 # Threshold for identity of ResFinder
+    Boolean call_pointfinder = false # Run pointfinder for chromosomal mutations
     String docker = "us-docker.pkg.dev/general-theiagen/staphb/resfinder:4.1.11"
     Int disk_size = 100
+    Int cpu = 2
+    Int memory = 8
   }
   command <<<
     date | tee DATE
@@ -28,7 +30,8 @@ task resfinder {
       resfinder_organism="enterococcus faecalis"
     elif [[ "~{organism}" == *"Enterococcus"*"faecium"* ]]; then 
       resfinder_organism="enterococcus faecium"
-    elif [[ "~{organism}" == *"Escherichia"*"coli"* ]]; then 
+    # to allow for both E. coli and any Shigella species to be processed via PointFinder as E. coli
+    elif [[ "~{organism}" == *"Escherichia"*"coli"* ]] || [[ "~{organism}" == *"Shigella"* ]]; then 
       resfinder_organism="escherichia coli"
     elif [[ "~{organism}" == *"Klebsiella"* ]]; then 
       resfinder_organism="klebsiella"
@@ -45,6 +48,7 @@ task resfinder {
     else 
       echo "Either Gambit predicted taxon is not supported by resfinder or the user did not supply an organism as input."
       echo "Skipping the use of resfinder --species optional parameter."
+      echo "WARNING: This will disable PointFinder due to the requirement of --species flag."
     fi
 
     # if resfinder_organism variable is set, use --species flag, otherwise do not use --species flag
@@ -55,52 +59,141 @@ task resfinder {
         --species "${resfinder_organism}" \
         ~{true="--acquired" false="" acquired} \
         ~{'--min_cov ' + min_cov} \
-        ~{'--threshold ' + threshold} \
-        ~{true="--point" false="" point} 
+        ~{'--threshold ' + min_id} \
+        ~{true="--point" false="" call_pointfinder} 
     else 
+      # pointfinder requires the use of the --species flag, so if resfinder_organism is not set, do not run pointfinder
       run_resfinder.py \
         --inputfasta ~{assembly} \
         --outputPath . \
         --species "other" \
         ~{true="--acquired" false="" acquired} \
         ~{'--min_cov ' + min_cov} \
-        ~{'--threshold ' + threshold} 
+        ~{'--threshold ' + min_id}
     fi
 
     # replace space in resfinder_organism with underscore
     resfinder_organism="${resfinder_organism// /_}"
 
     # rename files
-    mv pheno_table.txt ~{samplename}_pheno_table.txt
-    if [ -f pheno_table_${resfinder_organism}.txt ]; then
-      mv pheno_table_${resfinder_organism}.txt ~{samplename}_pheno_table_species.txt
+    mv -v pheno_table.txt ~{samplename}_pheno_table.tsv
+    if [ -f "pheno_table_${resfinder_organism}.txt" ]; then
+      # rename file to have proper extension & samplename included
+      mv -v "pheno_table_${resfinder_organism}.txt" ~{samplename}_pheno_table_species.tsv
     fi
-    mv ResFinder_Hit_in_genome_seq.fsa ~{samplename}_ResFinder_Hit_in_genome_seq.fsa
-    mv ResFinder_Resistance_gene_seq.fsa ~{samplename}_ResFinder_Resistance_gene_seq.fsa
-    mv ResFinder_results_tab.txt ~{samplename}_ResFinder_results_tab.txt
+    mv -v ResFinder_Hit_in_genome_seq.fsa ~{samplename}_ResFinder_Hit_in_genome_seq.fsa
+    mv -v ResFinder_Resistance_gene_seq.fsa ~{samplename}_ResFinder_Resistance_gene_seq.fsa
+    mv -v ResFinder_results_tab.txt ~{samplename}_ResFinder_results_tab.tsv
+
+    # if pointfinder was run, rename files
     if [ -f PointFinder_prediction.txt ]; then
-      mv PointFinder_prediction.txt ~{samplename}_PointFinder_prediction.txt
-      mv PointFinder_results.txt ~{samplename}_PointFinder_results.txt
+      mv -v PointFinder_prediction.txt ~{samplename}_PointFinder_prediction.tsv
+      mv -v PointFinder_results.txt ~{samplename}_PointFinder_results.tsv
+    fi
+
+    # parse ~{samplename}_pheno_table.tsv for predicted phenotypes and genes associated
+    # strip off 18 lines from top of file (18th line is the header with the columns: antimicrobial, class, WGS-predicted phenotype, Match, Genetic Background)
+    tail +18 ~{samplename}_pheno_table.tsv > ~{samplename}_pheno_table.headerless.tsv
+
+    # convert all letters in first column (antibiotic) to uppercase. For readability of output string
+    awk -F '\t' 'BEGIN{OFS="\t"} { $1=toupper($1) } 1' ~{samplename}_pheno_table.headerless.tsv > ~{samplename}_pheno_table.headerless.uppercase.tsv
+
+    # if column 3 shows 'Resistant', then print list of drugs followed by the genes/point mutations responsible
+    awk -F '\t' 'BEGIN{OFS=":"; ORS="; "} { if($3 == "Resistant") {print $1,$5}}' ~{samplename}_pheno_table.headerless.uppercase.tsv \
+    | sed 's/..$//' > RESFINDER_PREDICTED_PHENO_RESISTANCE.txt
+
+    # check for XDR Shigella status, based on CDC definition here: https://emergency.cdc.gov/han/2023/han00486.asp
+    # requirements:
+    # organism input (i.e. gambit_predicted_taxon) must contain "Shigella"
+    # predicted resistance to antimicrobials must include ALL: "ceftriaxone", "azithromycin", "ciprofloxacin", "trimethoprim", "sulfamethoxazole", and "ampicillin"
+    if [[ "~{organism}" != *"Shigella"* ]]; then
+      echo 'Either the input gambit_predicted_taxon does not contain "Shigella" or the user did not supply the organism as an input string to the workflow.'
+      echo "Skipping XDR Shigella check."
+      echo "Not Shigella based on gambit_predicted_taxon or user input" | tee RESFINDER_PREDICTED_XDR_SHIGELLA.txt
+    # if organism input string DOES contain the word "Shigella", check for resistance predictions to 6 drugs in XDR definition
+    elif [[ "~{organism}" == *"Shigella"* ]]; then
+      # nested if: if grep finds the all drugs, set output to XDR shigella, but if not, set it to "Not XDR Shigella"
+      if grep -qi "ceftriaxone" RESFINDER_PREDICTED_PHENO_RESISTANCE.txt && \
+      grep -qi "azithromycin" RESFINDER_PREDICTED_PHENO_RESISTANCE.txt && \
+      grep -qi "ciprofloxacin" RESFINDER_PREDICTED_PHENO_RESISTANCE.txt && \
+      grep -qi "trimethoprim" RESFINDER_PREDICTED_PHENO_RESISTANCE.txt && \
+      grep -qi "sulfamethoxazole" RESFINDER_PREDICTED_PHENO_RESISTANCE.txt && \
+      grep -qi "ampicillin" RESFINDER_PREDICTED_PHENO_RESISTANCE.txt; then
+        echo "XDR Shigella based on predicted resistance to ceftriaxone, azithromycin, ciprofloxacin, trimethoprim, sulfamethoxazole, and ampicillin. Please verify by reviewing ~{samplename}_pheno_table.tsv and ~{samplename}_ResFinder_results_tab.tsv"
+        echo "XDR Shigella" > RESFINDER_PREDICTED_XDR_SHIGELLA.txt
+      # ~{organism} does contain the word "Shigella", but one of the greps failed, meaning not all drug resistances' were predicted
+      else
+        echo "Not XDR Shigella" | tee RESFINDER_PREDICTED_XDR_SHIGELLA.txt
+      fi
+    fi
+
+    # set output strings for Resistance or no Resistance predicted for each of 6 drugs
+    # if grep finds the drug in the RESFINDER_PREDICTED_PHENO_RESISTANCE.txt file, then set the output string to "Resistance"
+    # if grep does not find the drug in the RESFINDER_PREDICTED_PHENO_RESISTANCE.txt file, then set the output string to "No resistance predicted"
+    # ampicillin
+    if grep -qi "ampicillin" RESFINDER_PREDICTED_PHENO_RESISTANCE.txt; then
+      echo "Resistance" > RESFINDER_PREDICTED_RESISTANCE_AMP.txt
+    else
+      echo "No resistance predicted" > RESFINDER_PREDICTED_RESISTANCE_AMP.txt
+    fi
+    # azithromycin
+    if grep -qi "azithromycin" RESFINDER_PREDICTED_PHENO_RESISTANCE.txt; then
+      echo "Resistance" > RESFINDER_PREDICTED_RESISTANCE_AZM.txt
+    else
+      echo "No resistance predicted" > RESFINDER_PREDICTED_RESISTANCE_AZM.txt
+    fi
+    # ceftriaxone
+    if grep -qi "ceftriaxone" RESFINDER_PREDICTED_PHENO_RESISTANCE.txt; then
+      echo "Resistance" > RESFINDER_PREDICTED_RESISTANCE_AXO.txt
+    else
+      echo "No resistance predicted" > RESFINDER_PREDICTED_RESISTANCE_AXO.txt
+    fi
+    # ciprofloxacin
+    if grep -qi "ciprofloxacin" RESFINDER_PREDICTED_PHENO_RESISTANCE.txt; then
+      echo "Resistance" > RESFINDER_PREDICTED_RESISTANCE_CIP.txt
+    else
+      echo "No resistance predicted" > RESFINDER_PREDICTED_RESISTANCE_CIP.txt
+    fi
+    # sulfamethoxazole
+    if grep -qi "sulfamethoxazole" RESFINDER_PREDICTED_PHENO_RESISTANCE.txt; then
+      echo "Resistance" > RESFINDER_PREDICTED_RESISTANCE_SMX.txt
+    else
+      echo "No resistance predicted" > RESFINDER_PREDICTED_RESISTANCE_SMX.txt
+    fi
+    # trimethoprim
+    if grep -qi "trimethoprim" RESFINDER_PREDICTED_PHENO_RESISTANCE.txt; then
+      echo "Resistance" > RESFINDER_PREDICTED_RESISTANCE_TMP.txt
+    else
+      echo "No resistance predicted" > RESFINDER_PREDICTED_RESISTANCE_TMP.txt
     fi
 
   >>>
   output {
-    File resfinder_pheno_table = "~{samplename}_pheno_table.txt"
-    File? resfinder_pheno_table_species = "~{samplename}_pheno_table_species.txt"
+    File resfinder_pheno_table = "~{samplename}_pheno_table.tsv"
+    File? resfinder_pheno_table_species = "~{samplename}_pheno_table_species.tsv"
     File resfinder_hit_in_genome_seq = "~{samplename}_ResFinder_Hit_in_genome_seq.fsa"
     File resfinder_resistance_gene_seq = "~{samplename}_ResFinder_Resistance_gene_seq.fsa"
-    File resfinder_results_tab = "~{samplename}_ResFinder_results_tab.txt"
-    File? pointfinder_pheno_table = "~{samplename}_PointFinder_prediction.txt"
-    File? pointfinder_results = "~{samplename}_PointFinder_results.txt"
+    File resfinder_results_tab = "~{samplename}_ResFinder_results_tab.tsv"
+    File? pointfinder_pheno_table = "~{samplename}_PointFinder_prediction.tsv"
+    File? pointfinder_results = "~{samplename}_PointFinder_results.tsv"
+    String resfinder_predicted_pheno_resistance = read_string("RESFINDER_PREDICTED_PHENO_RESISTANCE.txt")
+    String resfinder_predicted_xdr_shigella = read_string("RESFINDER_PREDICTED_XDR_SHIGELLA.txt")
+    String resfinder_predicted_resistance_Amp = read_string("RESFINDER_PREDICTED_RESISTANCE_AMP.txt")
+    String resfinder_predicted_resistance_Azm = read_string("RESFINDER_PREDICTED_RESISTANCE_AZM.txt")
+    String resfinder_predicted_resistance_Axo = read_string("RESFINDER_PREDICTED_RESISTANCE_AXO.txt")
+    String resfinder_predicted_resistance_Cip = read_string("RESFINDER_PREDICTED_RESISTANCE_CIP.txt")
+    String resfinder_predicted_resistance_Smx = read_string("RESFINDER_PREDICTED_RESISTANCE_SMX.txt")
+    String resfinder_predicted_resistance_Tmp = read_string("RESFINDER_PREDICTED_RESISTANCE_TMP.txt")
     String resfinder_docker = "~{docker}"
     String resfinder_version = read_string("RESFINDER_VERSION")
     String resfinder_db_version = read_string("RESFINDER_DB_VERSION")
   }
   runtime {
-    memory: "8 GB"
-    cpu: 4
+    memory: "~{memory} GB"
+    cpu: cpu
     docker: docker
     disks: "local-disk " + disk_size + " SSD"
     disk: disk_size + " GB"
-    maxRetries: 3  }
+    maxRetries: 3  
+  }
 }

tasks/utilities/task_broad_terra_tools.wdl

@@ -117,6 +117,14 @@ task export_taxon_tables {
     File? resfinder_results 
     File? resfinder_pointfinder_pheno_table 
     File? resfinder_pointfinder_results 
+    String? resfinder_predicted_pheno_resistance
+    String? resfinder_predicted_xdr_shigella
+    String? resfinder_predicted_resistance_Amp
+    String? resfinder_predicted_resistance_Azm
+    String? resfinder_predicted_resistance_Axo
+    String? resfinder_predicted_resistance_Cip
+    String? resfinder_predicted_resistance_Smx
+    String? resfinder_predicted_resistance_Tmp
     String? resfinder_db_version 
     String? resfinder_docker 
     File? ts_mlst_results
@@ -609,6 +617,14 @@ task export_taxon_tables {
       "resfinder_results": "~{resfinder_results}",
       "resfinder_pointfinder_pheno_table": "~{resfinder_pointfinder_pheno_table}",
       "resfinder_pointfinder_results": "~{resfinder_pointfinder_results}",
+      "resfinder_predicted_pheno_resistance": "~{resfinder_predicted_pheno_resistance}",
+      "resfinder_predicted_xdr_shigella": "~{resfinder_predicted_xdr_shigella}",
+      "resfinder_predicted_resistance_Amp": "~{resfinder_predicted_resistance_Amp}",
+      "resfinder_predicted_resistance_Azm": "~{resfinder_predicted_resistance_Azm}",
+      "resfinder_predicted_resistance_Axo": "~{resfinder_predicted_resistance_Axo}",
+      "resfinder_predicted_resistance_Cip": "~{resfinder_predicted_resistance_Cip}",
+      "resfinder_predicted_resistance_Smx": "~{resfinder_predicted_resistance_Smx}",
+      "resfinder_predicted_resistance_Tmp": "~{resfinder_predicted_resistance_Tmp}",
       "resfinder_db_version": "~{resfinder_db_version}",
       "resfinder_docker": "~{resfinder_docker}",
       "prokka_gff": "~{prokka_gff}",

tests/workflows/theiaprok/test_wf_theiaprok_illumina_pe.yml

@@ -560,7 +560,7 @@
     - path: miniwdl_run/wdl/tasks/gene_typing/task_prokka.wdl
       md5sum: 8e8423b99d5cd3e92bb89b3c3211d432
     - path: miniwdl_run/wdl/tasks/gene_typing/task_resfinder.wdl
-      md5sum: 18d0dd6560f88ed2541ad342e12412e8
+      md5sum: 27528633723303b462d095b642649453
     - path: miniwdl_run/wdl/tasks/gene_typing/task_snippy_variants.wdl
       md5sum: bc048390626650b6a939fd740c1aaf13
     - path: miniwdl_run/wdl/tasks/quality_control/task_bbduk.wdl
@@ -632,9 +632,9 @@
     - path: miniwdl_run/wdl/tasks/taxon_id/task_midas.wdl
       md5sum: faacd87946ee3fbdf70f3a15b79ce547
     - path: miniwdl_run/wdl/tasks/utilities/task_broad_terra_tools.wdl
-      md5sum: 26ac9a20c8043a28d373bfe0ca361cc6
+      md5sum: 5dc54b8446b6a430fc7375ae364908f0
     - path: miniwdl_run/wdl/workflows/theiaprok/wf_theiaprok_illumina_pe.wdl
-      md5sum: ac4971ad992c3e8abee5d3817928a8f2
+      md5sum: 3acf4dcddbb44d547b69f597761cc048
     - path: miniwdl_run/wdl/workflows/utilities/wf_merlin_magic.wdl
       md5sum: 00bd2489b2a7aa5b88340a940961a857
     - path: miniwdl_run/wdl/workflows/utilities/wf_read_QC_trim_pe.wdl

tests/workflows/theiaprok/test_wf_theiaprok_illumina_se.yml

@@ -528,7 +528,7 @@
     - path: miniwdl_run/wdl/tasks/gene_typing/task_prokka.wdl
       md5sum: 8e8423b99d5cd3e92bb89b3c3211d432
     - path: miniwdl_run/wdl/tasks/gene_typing/task_resfinder.wdl
-      md5sum: 18d0dd6560f88ed2541ad342e12412e8
+      md5sum: 27528633723303b462d095b642649453
     - path: miniwdl_run/wdl/tasks/gene_typing/task_snippy_variants.wdl
       md5sum: bc048390626650b6a939fd740c1aaf13
     - path: miniwdl_run/wdl/tasks/quality_control/task_bbduk.wdl
@@ -598,9 +598,9 @@
     - path: miniwdl_run/wdl/tasks/taxon_id/task_midas.wdl
       md5sum: faacd87946ee3fbdf70f3a15b79ce547
     - path: miniwdl_run/wdl/tasks/utilities/task_broad_terra_tools.wdl
-      md5sum: 26ac9a20c8043a28d373bfe0ca361cc6
+      md5sum: 5dc54b8446b6a430fc7375ae364908f0
     - path: miniwdl_run/wdl/workflows/theiaprok/wf_theiaprok_illumina_se.wdl
-      md5sum: 7303247badeb82119ccef528f7367f89
+      md5sum: 3e19938fc8a624c7948b57867865561a
     - path: miniwdl_run/wdl/workflows/utilities/wf_merlin_magic.wdl
       md5sum: 00bd2489b2a7aa5b88340a940961a857
     - path: miniwdl_run/wdl/workflows/utilities/wf_read_QC_trim_se.wdl

workflows/theiaprok/wf_theiaprok_fasta.wdl

@@ -200,6 +200,14 @@ workflow theiaprok_fasta {
         resfinder_results = resfinder_task.resfinder_results_tab,
         resfinder_pointfinder_pheno_table = resfinder_task.pointfinder_pheno_table,
         resfinder_pointfinder_results = resfinder_task.pointfinder_results,
+        resfinder_predicted_pheno_resistance = resfinder_task.resfinder_predicted_pheno_resistance,
+        resfinder_predicted_xdr_shigella = resfinder_task.resfinder_predicted_xdr_shigella,
+        resfinder_predicted_resistance_Amp = resfinder_task.resfinder_predicted_resistance_Amp,
+        resfinder_predicted_resistance_Azm = resfinder_task.resfinder_predicted_resistance_Azm,
+        resfinder_predicted_resistance_Axo = resfinder_task.resfinder_predicted_resistance_Axo,
+        resfinder_predicted_resistance_Cip = resfinder_task.resfinder_predicted_resistance_Cip,
+        resfinder_predicted_resistance_Smx = resfinder_task.resfinder_predicted_resistance_Smx,
+        resfinder_predicted_resistance_Tmp = resfinder_task.resfinder_predicted_resistance_Tmp,
         resfinder_db_version = resfinder_task.resfinder_db_version,
         resfinder_docker = resfinder_task.resfinder_docker,
         ts_mlst_results = ts_mlst.ts_mlst_results,
@@ -456,6 +464,14 @@ workflow theiaprok_fasta {
     File? resfinder_results = resfinder_task.resfinder_results_tab
     File? resfinder_pointfinder_pheno_table = resfinder_task.pointfinder_pheno_table
     File? resfinder_pointfinder_results = resfinder_task.pointfinder_results
+    String? resfinder_predicted_pheno_resistance = resfinder_task.resfinder_predicted_pheno_resistance
+    String? resfinder_predicted_xdr_shigella = resfinder_task.resfinder_predicted_xdr_shigella
+    String? resfinder_predicted_resistance_Amp = resfinder_task.resfinder_predicted_resistance_Amp
+    String? resfinder_predicted_resistance_Azm = resfinder_task.resfinder_predicted_resistance_Azm
+    String? resfinder_predicted_resistance_Axo = resfinder_task.resfinder_predicted_resistance_Axo
+    String? resfinder_predicted_resistance_Cip = resfinder_task.resfinder_predicted_resistance_Cip
+    String? resfinder_predicted_resistance_Smx = resfinder_task.resfinder_predicted_resistance_Smx
+    String? resfinder_predicted_resistance_Tmp = resfinder_task.resfinder_predicted_resistance_Tmp
     String? resfinder_db_version = resfinder_task.resfinder_db_version
     String? resfinder_docker = resfinder_task.resfinder_docker
     # MLST Typing

workflows/theiaprok/wf_theiaprok_illumina_pe.wdl

@@ -325,6 +325,14 @@ workflow theiaprok_illumina_pe {
             resfinder_results = resfinder_task.resfinder_results_tab,
             resfinder_pointfinder_pheno_table = resfinder_task.pointfinder_pheno_table,
             resfinder_pointfinder_results = resfinder_task.pointfinder_results,
+            resfinder_predicted_pheno_resistance = resfinder_task.resfinder_predicted_pheno_resistance,
+            resfinder_predicted_xdr_shigella = resfinder_task.resfinder_predicted_xdr_shigella,
+            resfinder_predicted_resistance_Amp = resfinder_task.resfinder_predicted_resistance_Amp,
+            resfinder_predicted_resistance_Azm = resfinder_task.resfinder_predicted_resistance_Azm,
+            resfinder_predicted_resistance_Axo = resfinder_task.resfinder_predicted_resistance_Axo,
+            resfinder_predicted_resistance_Cip = resfinder_task.resfinder_predicted_resistance_Cip,
+            resfinder_predicted_resistance_Smx = resfinder_task.resfinder_predicted_resistance_Smx,
+            resfinder_predicted_resistance_Tmp = resfinder_task.resfinder_predicted_resistance_Tmp,
             resfinder_db_version = resfinder_task.resfinder_db_version,
             resfinder_docker = resfinder_task.resfinder_docker,
             ts_mlst_results = ts_mlst.ts_mlst_results,
@@ -650,6 +658,14 @@ workflow theiaprok_illumina_pe {
     File? resfinder_results = resfinder_task.resfinder_results_tab
     File? resfinder_pointfinder_pheno_table = resfinder_task.pointfinder_pheno_table
     File? resfinder_pointfinder_results = resfinder_task.pointfinder_results
+    String? resfinder_predicted_pheno_resistance = resfinder_task.resfinder_predicted_pheno_resistance
+    String? resfinder_predicted_xdr_shigella = resfinder_task.resfinder_predicted_xdr_shigella
+    String? resfinder_predicted_resistance_Amp = resfinder_task.resfinder_predicted_resistance_Amp
+    String? resfinder_predicted_resistance_Azm = resfinder_task.resfinder_predicted_resistance_Azm
+    String? resfinder_predicted_resistance_Axo = resfinder_task.resfinder_predicted_resistance_Axo
+    String? resfinder_predicted_resistance_Cip = resfinder_task.resfinder_predicted_resistance_Cip
+    String? resfinder_predicted_resistance_Smx = resfinder_task.resfinder_predicted_resistance_Smx
+    String? resfinder_predicted_resistance_Tmp = resfinder_task.resfinder_predicted_resistance_Tmp
     String? resfinder_db_version = resfinder_task.resfinder_db_version
     String? resfinder_docker = resfinder_task.resfinder_docker
     # MLST Typing