Skip to content

Commit 2bfb9d7

Browse files
tpshea2claude
tpshea2
and
claude
committed
Add HvpViralPipeline: geNomad, VirSorter2, CheckV, ViralSummary tasks + pipeline workflow
Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
1 parent 8199316 commit 2bfb9d7

6 files changed

Lines changed: 1001 additions & 0 deletions

File tree

‎.dockstore.yml‎

Lines changed: 19 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -79,3 +79,22 @@ workflows:
7979
- /tshea\/.*/
8080
tags:
8181
- /.*/
82+
83+
- name: HvpViralPipeline
84+
subclass: WDL
85+
primaryDescriptorPath: /wdl/pipelines/PacBio/Metagenomics/HvpViralPipeline.wdl
86+
readMePath: /README.md
87+
authors:
88+
- name: Terrance Shea
89+
orcid: 0000-0002-6586-9298
90+
affiliation: Broad Institute
91+
- name: Jonn Smith
92+
orcid: 0000-0002-4279-8220
93+
affiliation: Broad Institute
94+
filters:
95+
branches:
96+
- main
97+
- /jts_.*/
98+
- /tshea\/.*/
99+
tags:
100+
- /.*/

‎wdl/pipelines/PacBio/Metagenomics/HvpViralPipeline.wdl‎

Lines changed: 138 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,138 @@
1+
version 1.0
2+
3+
import "../../../tasks/Viral/Genomad.wdl" as GnomadTasks
4+
import "../../../tasks/Viral/VirSorter2.wdl" as VS2Tasks
5+
import "../../../tasks/Viral/CheckV.wdl" as CheckVTasks
6+
import "../../../tasks/Viral/ViralSummary.wdl" as SumTasks
7+
8+
workflow HvpViralPipeline {
9+
10+
meta {
11+
description: "Single-sample HVP viral identification and quality assessment pipeline. Accepts the merged assembly+rescue FASTA from HvpReadRescue and runs geNomad and VirSorter2 in parallel to identify viral sequences, then assesses quality of each tool's output with CheckV. Produces a per-contig summary joining all four tool outputs and aggregate scalar metrics for the Terra data table."
12+
13+
allowNestedInputs: true
14+
15+
outputs: {
16+
genomad_virus_summary_tsv: "geNomad virus_summary.tsv (per-sequence scores, topology, taxonomy)",
17+
genomad_virus_fna: "Viral sequences FASTA from geNomad; input to CheckV",
18+
vs2_score_tsv: "VirSorter2 final-viral-score.tsv (per-sequence viral scores and group assignments)",
19+
vs2_viral_combined_fa: "Viral sequences FASTA from VirSorter2; input to CheckV",
20+
genomad_checkv_quality_tsv: "CheckV quality_summary.tsv for geNomad viral sequences",
21+
genomad_checkv_viruses_fna: "CheckV viruses.fna for geNomad viral sequences (complete/HQ genomes)",
22+
genomad_checkv_proviruses_fna: "CheckV proviruses.fna for geNomad viral sequences (proviral extracts)",
23+
vs2_checkv_quality_tsv: "CheckV quality_summary.tsv for VirSorter2 viral sequences",
24+
vs2_checkv_viruses_fna: "CheckV viruses.fna for VirSorter2 viral sequences (complete/HQ genomes)",
25+
vs2_checkv_proviruses_fna: "CheckV proviruses.fna for VirSorter2 viral sequences (proviral extracts)",
26+
viral_contig_summary_tsv: "Per-contig summary TSV joining geNomad and VirSorter2 scores with CheckV quality results",
27+
viral_stats_tsv: "Single-row aggregate viral metrics TSV for the sample",
28+
n_viral_contigs: "Total viral contigs called by at least one tool",
29+
n_hq_viral: "Number of high-quality viral contigs (CheckV High-quality or Complete)",
30+
n_complete_viral: "Number of complete viral genomes (CheckV Complete)",
31+
n_both_tools: "Number of contigs called viral by both geNomad and VirSorter2",
32+
n50_viral: "N50 length of viral contigs",
33+
total_viral_bases: "Total bases across all viral contigs"
34+
}
35+
}
36+
37+
parameter_meta {
38+
merged_fa_gz: "Merged assembly+rescue FASTA (gzipped) from HvpReadRescue"
39+
sample_name: "Sample identifier used as output file prefix (e.g. bc2097)"
40+
genomad_db_tgz: "geNomad database as a single compressed archive (.tar.gz or .tar.zst)"
41+
vs2_db_tgz: "VirSorter2 database as a single compressed archive (.tar.gz or .tar.zst)"
42+
checkv_db_tgz: "CheckV database as a single compressed archive (.tar.gz or .tar.zst)"
43+
genomad_extra_args: "Additional command-line args appended verbatim to the genomad invocation"
44+
vs2_min_length: "Minimum sequence length for VirSorter2 (default 1000)"
45+
vs2_extra_args: "Additional command-line args appended verbatim to the virsorter run invocation"
46+
}
47+
48+
input {
49+
File merged_fa_gz
50+
String sample_name
51+
52+
File genomad_db_tgz
53+
File vs2_db_tgz
54+
File checkv_db_tgz
55+
56+
String genomad_extra_args = ""
57+
Int vs2_min_length = 1000
58+
String vs2_extra_args = ""
59+
}
60+
61+
# geNomad and VirSorter2 run in parallel — both take the merged assembly+rescue FASTA
62+
call GnomadTasks.Genomad as t_01_Genomad {
63+
input:
64+
merged_fa_gz = merged_fa_gz,
65+
genomad_db_tgz = genomad_db_tgz,
66+
sample_name = sample_name,
67+
extra_args = genomad_extra_args
68+
}
69+
70+
call VS2Tasks.VirSorter2 as t_02_VirSorter2 {
71+
input:
72+
merged_fa_gz = merged_fa_gz,
73+
vs2_db_tgz = vs2_db_tgz,
74+
sample_name = sample_name,
75+
min_length = vs2_min_length,
76+
extra_args = vs2_extra_args
77+
}
78+
79+
# CheckV runs once per upstream tool output
80+
call CheckVTasks.CheckV as t_03_CheckVGenomd {
81+
input:
82+
virus_fna = t_01_Genomad.virus_fna,
83+
checkv_db_tgz = checkv_db_tgz,
84+
sample_name = sample_name,
85+
tool_prefix = "genomad"
86+
}
87+
88+
call CheckVTasks.CheckV as t_04_CheckVVS2 {
89+
input:
90+
virus_fna = t_02_VirSorter2.viral_combined_fa,
91+
checkv_db_tgz = checkv_db_tgz,
92+
sample_name = sample_name,
93+
tool_prefix = "vs2"
94+
}
95+
96+
# Join all four tool outputs into a per-contig summary
97+
call SumTasks.ViralContigSummary as t_05_ViralContigSummary {
98+
input:
99+
genomad_summary_tsv = t_01_Genomad.virus_summary_tsv,
100+
genomad_checkv_tsv = t_03_CheckVGenomd.quality_summary_tsv,
101+
vs2_score_tsv = t_02_VirSorter2.score_tsv,
102+
vs2_checkv_tsv = t_04_CheckVVS2.quality_summary_tsv,
103+
sample_name = sample_name
104+
}
105+
106+
# Aggregate per-contig summary to scalar metrics for the data table
107+
call SumTasks.ViralOverallSummary as t_06_ViralOverallSummary {
108+
input:
109+
viral_contig_summary_tsv = t_05_ViralContigSummary.viral_contig_summary_tsv,
110+
sample_name = sample_name
111+
}
112+
113+
output {
114+
File genomad_virus_summary_tsv = t_01_Genomad.virus_summary_tsv
115+
File genomad_virus_fna = t_01_Genomad.virus_fna
116+
117+
File vs2_score_tsv = t_02_VirSorter2.score_tsv
118+
File vs2_viral_combined_fa = t_02_VirSorter2.viral_combined_fa
119+
120+
File genomad_checkv_quality_tsv = t_03_CheckVGenomd.quality_summary_tsv
121+
File genomad_checkv_viruses_fna = t_03_CheckVGenomd.viruses_fna
122+
File genomad_checkv_proviruses_fna = t_03_CheckVGenomd.proviruses_fna
123+
124+
File vs2_checkv_quality_tsv = t_04_CheckVVS2.quality_summary_tsv
125+
File vs2_checkv_viruses_fna = t_04_CheckVVS2.viruses_fna
126+
File vs2_checkv_proviruses_fna = t_04_CheckVVS2.proviruses_fna
127+
128+
File viral_contig_summary_tsv = t_05_ViralContigSummary.viral_contig_summary_tsv
129+
130+
File viral_stats_tsv = t_06_ViralOverallSummary.viral_stats_tsv
131+
Int n_viral_contigs = t_06_ViralOverallSummary.n_viral_contigs
132+
Int n_hq_viral = t_06_ViralOverallSummary.n_hq_viral
133+
Int n_complete_viral = t_06_ViralOverallSummary.n_complete_viral
134+
Int n_both_tools = t_06_ViralOverallSummary.n_both_tools
135+
Int n50_viral = t_06_ViralOverallSummary.n50_viral
136+
Int total_viral_bases = t_06_ViralOverallSummary.total_viral_bases
137+
}
138+
}

‎wdl/tasks/Viral/CheckV.wdl‎

Lines changed: 110 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,110 @@
1+
version 1.0
2+
3+
import "../../structs/Structs.wdl"
4+
5+
task CheckV {
6+
7+
meta {
8+
description: "Run CheckV end-to-end on a viral FASTA to assess completeness and contamination. Spaces in FASTA headers are replaced with underscores before running to prevent CheckV from producing duplicate IDs (it truncates headers at whitespace, causing downstream collisions in the summary script). Called once per upstream viral detection tool (geNomad, VirSorter2)."
9+
10+
tool: "CheckV"
11+
tool_version: "1.0.3"
12+
tool_url: "https://bitbucket.org/berkeleylab/checkv"
13+
tool_citation: "Nayfach S, et al. CheckV assesses the quality and completeness of metagenome-assembled viral genomes. Nat Biotechnol. 2021;39(5):578-585."
14+
15+
outputs: {
16+
quality_summary_tsv: "CheckV quality_summary.tsv: per-sequence completeness, contamination, provirus status, and MIUVIG quality tier",
17+
viruses_fna: "CheckV viruses.fna: sequences classified as complete or high-quality viral genomes",
18+
proviruses_fna: "CheckV proviruses.fna: proviral extracts with trimmed host flanking regions"
19+
}
20+
}
21+
22+
parameter_meta {
23+
virus_fna: "Viral sequences FASTA from an upstream detection tool (geNomad or VirSorter2)"
24+
checkv_db_tgz: "CheckV database as a single compressed archive (.tar.gz or .tar.zst); extracted at runtime"
25+
sample_name: "Sample identifier used as output file prefix"
26+
tool_prefix: "Short tool identifier prepended to output filenames (e.g. 'genomad' or 'vs2')"
27+
runtime_attr_override: "Override the default runtime attributes"
28+
}
29+
30+
input {
31+
File virus_fna
32+
File checkv_db_tgz
33+
String sample_name
34+
String tool_prefix
35+
36+
RuntimeAttr? runtime_attr_override
37+
}
38+
39+
Int disk_size = 20 + ceil(3.0 * size(checkv_db_tgz, "GB")) + ceil(5.0 * size(virus_fna, "GB"))
40+
41+
command <<<
42+
set -euxo pipefail
43+
44+
# ---- Resource detection (required preamble) ----
45+
NUM_CPUS=$(grep '^processor' /proc/cpuinfo | tail -n1 | awk '{print $NF+1}')
46+
RAM_IN_GB=$(free -g | grep "^Mem" | awk '{print $2}')
47+
48+
USABLE_RAM_GB=$((RAM_IN_GB - 1))
49+
[[ "${USABLE_RAM_GB}" -lt 1 ]] && USABLE_RAM_GB=1
50+
51+
MEM_PER_THREAD_GB=$(( USABLE_RAM_GB / NUM_CPUS ))
52+
[[ "${MEM_PER_THREAD_GB}" -lt 1 ]] && MEM_PER_THREAD_GB=1
53+
54+
JAVA_MEM_GB=${USABLE_RAM_GB}
55+
56+
echo "NUM_CPUS=${NUM_CPUS} RAM_IN_GB=${RAM_IN_GB} USABLE_RAM_GB=${USABLE_RAM_GB} MEM_PER_THREAD_GB=${MEM_PER_THREAD_GB} JAVA_MEM_GB=${JAVA_MEM_GB}"
57+
# ---- end preamble ----
58+
59+
# Extract CheckV database
60+
mkdir -p checkv_db
61+
if [[ "~{checkv_db_tgz}" == *.tar.zst ]]; then
62+
zstd -d ~{checkv_db_tgz} --stdout | tar -x -C checkv_db
63+
else
64+
tar -xzf ~{checkv_db_tgz} -C checkv_db
65+
fi
66+
DB_PATH=$(find checkv_db -name 'genome_db' -type d | head -1 | xargs dirname)
67+
[[ -z "${DB_PATH}" ]] && DB_PATH=$(ls -d checkv_db/*/ | head -1)
68+
69+
# Replace spaces with underscores in headers — CheckV truncates at whitespace,
70+
# causing duplicate IDs when the same base ID appears in multiple tool outputs
71+
sed 's/ /_/g' ~{virus_fna} > cleaned_input.fna
72+
73+
checkv end_to_end \
74+
cleaned_input.fna \
75+
checkv_out \
76+
-t "${NUM_CPUS}" \
77+
-d "${DB_PATH}"
78+
79+
cp checkv_out/quality_summary.tsv ~{sample_name}.~{tool_prefix}_checkv_quality.tsv
80+
cp checkv_out/viruses.fna ~{sample_name}.~{tool_prefix}_checkv_viruses.fna
81+
cp checkv_out/proviruses.fna ~{sample_name}.~{tool_prefix}_checkv_proviruses.fna
82+
>>>
83+
84+
output {
85+
File quality_summary_tsv = "~{sample_name}.~{tool_prefix}_checkv_quality.tsv"
86+
File viruses_fna = "~{sample_name}.~{tool_prefix}_checkv_viruses.fna"
87+
File proviruses_fna = "~{sample_name}.~{tool_prefix}_checkv_proviruses.fna"
88+
}
89+
90+
#########################
91+
RuntimeAttr default_attr = object {
92+
cpu_cores: 16,
93+
mem_gb: 16,
94+
disk_gb: disk_size,
95+
boot_disk_gb: 25,
96+
preemptible_tries: 2,
97+
max_retries: 1,
98+
docker: "ghcr.io/broadinstitute/hvp-lr/hvp-monolith:0.0.3"
99+
}
100+
RuntimeAttr runtime_attr = select_first([runtime_attr_override, default_attr])
101+
runtime {
102+
cpu: select_first([runtime_attr.cpu_cores, default_attr.cpu_cores])
103+
memory: select_first([runtime_attr.mem_gb, default_attr.mem_gb]) + " GiB"
104+
disks: "local-disk " + select_first([runtime_attr.disk_gb, default_attr.disk_gb]) + " SSD"
105+
bootDiskSizeGb: select_first([runtime_attr.boot_disk_gb, default_attr.boot_disk_gb])
106+
preemptible: select_first([runtime_attr.preemptible_tries, default_attr.preemptible_tries])
107+
maxRetries: select_first([runtime_attr.max_retries, default_attr.max_retries])
108+
docker: select_first([runtime_attr.docker, default_attr.docker])
109+
}
110+
}

‎wdl/tasks/Viral/Genomad.wdl‎

Lines changed: 113 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,113 @@
1+
version 1.0
2+
3+
import "../../structs/Structs.wdl"
4+
5+
task Genomad {
6+
7+
meta {
8+
description: "Run geNomad end-to-end on a merged assembly+rescue FASTA to identify viral sequences and plasmids. Input FASTA is symlinked to a fixed name (input_contigs.fasta.gz) so output subdirectory names are predictable regardless of the upstream filename. Outputs the virus summary TSV (scores, topology, taxonomy) and the viral sequence FASTA for downstream CheckV quality assessment."
9+
10+
tool: "geNomad"
11+
tool_version: "1.9"
12+
tool_url: "https://github.com/apcamargo/genomad"
13+
tool_citation: "Camargo AP, et al. Identification of mobile genetic elements with geNomad. Nat Biotechnol. 2023."
14+
15+
outputs: {
16+
virus_summary_tsv: "geNomad virus summary TSV (one row per viral sequence: score, FDR, hallmarks, topology, taxonomy)",
17+
virus_fna: "Viral sequences FASTA produced by geNomad; input to CheckV"
18+
}
19+
}
20+
21+
parameter_meta {
22+
merged_fa_gz: "Merged assembly+rescue FASTA (gzipped) from HvpReadRescue"
23+
genomad_db_tgz: "geNomad database as a single compressed archive (.tar.gz or .tar.zst); extracted at runtime"
24+
sample_name: "Sample identifier used as output file prefix"
25+
extra_args: "Additional command-line args appended verbatim to the genomad invocation"
26+
runtime_attr_override: "Override the default runtime attributes"
27+
}
28+
29+
input {
30+
File merged_fa_gz
31+
File genomad_db_tgz
32+
String sample_name
33+
34+
String extra_args = ""
35+
36+
RuntimeAttr? runtime_attr_override
37+
}
38+
39+
Int disk_size = 20 + ceil(5.0 * size(genomad_db_tgz, "GB")) + ceil(10.0 * size(merged_fa_gz, "GB"))
40+
41+
command <<<
42+
set -euxo pipefail
43+
44+
# ---- Resource detection (required preamble) ----
45+
NUM_CPUS=$(grep '^processor' /proc/cpuinfo | tail -n1 | awk '{print $NF+1}')
46+
RAM_IN_GB=$(free -g | grep "^Mem" | awk '{print $2}')
47+
48+
USABLE_RAM_GB=$((RAM_IN_GB - 1))
49+
[[ "${USABLE_RAM_GB}" -lt 1 ]] && USABLE_RAM_GB=1
50+
51+
MEM_PER_THREAD_GB=$(( USABLE_RAM_GB / NUM_CPUS ))
52+
[[ "${MEM_PER_THREAD_GB}" -lt 1 ]] && MEM_PER_THREAD_GB=1
53+
54+
JAVA_MEM_GB=${USABLE_RAM_GB}
55+
56+
echo "NUM_CPUS=${NUM_CPUS} RAM_IN_GB=${RAM_IN_GB} USABLE_RAM_GB=${USABLE_RAM_GB} MEM_PER_THREAD_GB=${MEM_PER_THREAD_GB} JAVA_MEM_GB=${JAVA_MEM_GB}"
57+
# ---- end preamble ----
58+
59+
# Extract geNomad database
60+
mkdir -p genomad_db
61+
if [[ "~{genomad_db_tgz}" == *.tar.zst ]]; then
62+
zstd -d ~{genomad_db_tgz} --stdout | tar -x -C genomad_db
63+
else
64+
tar -xzf ~{genomad_db_tgz} -C genomad_db
65+
fi
66+
DB_ROOT=$(dirname "$(find genomad_db -name 'genomad_db' -type d | head -1)")
67+
[[ -z "${DB_ROOT}" ]] && DB_ROOT=$(find genomad_db -maxdepth 2 -name '*.mmseqs2_db' | head -1 | xargs -I{} dirname {})
68+
DB_PATH=$(find genomad_db -maxdepth 3 -name 'genomad_db' -type d | head -1)
69+
[[ -z "${DB_PATH}" ]] && DB_PATH=$(ls -d genomad_db/*/ | head -1)
70+
71+
# Symlink input to a fixed name so output subdirectory names are predictable
72+
ln -sf ~{merged_fa_gz} input_contigs.fasta.gz
73+
74+
genomad end-to-end \
75+
--threads "${NUM_CPUS}" \
76+
~{extra_args} \
77+
input_contigs.fasta.gz \
78+
genomad_out \
79+
"${DB_PATH}"
80+
81+
# Rename outputs with sample prefix
82+
cp genomad_out/input_contigs_summary/input_contigs_virus_summary.tsv \
83+
~{sample_name}.genomad_virus_summary.tsv
84+
cp genomad_out/input_contigs_virus/input_contigs_virus.fna \
85+
~{sample_name}.genomad_virus.fna
86+
>>>
87+
88+
output {
89+
File virus_summary_tsv = "~{sample_name}.genomad_virus_summary.tsv"
90+
File virus_fna = "~{sample_name}.genomad_virus.fna"
91+
}
92+
93+
#########################
94+
RuntimeAttr default_attr = object {
95+
cpu_cores: 16,
96+
mem_gb: 32,
97+
disk_gb: disk_size,
98+
boot_disk_gb: 25,
99+
preemptible_tries: 2,
100+
max_retries: 1,
101+
docker: "ghcr.io/broadinstitute/hvp-lr/hvp-monolith:0.0.3"
102+
}
103+
RuntimeAttr runtime_attr = select_first([runtime_attr_override, default_attr])
104+
runtime {
105+
cpu: select_first([runtime_attr.cpu_cores, default_attr.cpu_cores])
106+
memory: select_first([runtime_attr.mem_gb, default_attr.mem_gb]) + " GiB"
107+
disks: "local-disk " + select_first([runtime_attr.disk_gb, default_attr.disk_gb]) + " SSD"
108+
bootDiskSizeGb: select_first([runtime_attr.boot_disk_gb, default_attr.boot_disk_gb])
109+
preemptible: select_first([runtime_attr.preemptible_tries, default_attr.preemptible_tries])
110+
maxRetries: select_first([runtime_attr.max_retries, default_attr.max_retries])
111+
docker: select_first([runtime_attr.docker, default_attr.docker])
112+
}
113+
}

0 commit comments

Comments
 (0)