Merge pull request #4186 from vgteam/diploid-giraffe

Diploid sampling in Giraffe
vgteam · Dec 8, 2023 · 710f00c · 710f00c · adamnovak · Dec 8, 2023
2 parents 25b02cc + f066fc8
commit 710f00c
Show file tree

Hide file tree

Showing 4 changed files with 13 additions and 11 deletions.
diff --git a/src/subcommand/gbwt_main.cpp b/src/subcommand/gbwt_main.cpp
@@ -1829,4 +1829,4 @@ void print_metadata(std::ostream& out, const GBWTHandler& gbwts) {
 //----------------------------------------------------------------------------
 
 // Register subcommand
-static vg::subcommand::Subcommand vg_gbwt("gbwt", "build and manipulate GBWTs", vg::subcommand::TOOLKIT, main_gbwt);
+static vg::subcommand::Subcommand vg_gbwt("gbwt", "build and manipulate GBWT and GBZ files", vg::subcommand::TOOLKIT, main_gbwt);
diff --git a/src/subcommand/giraffe_main.cpp b/src/subcommand/giraffe_main.cpp
@@ -1671,6 +1671,8 @@ string sample_haplotypes(const vector<pair<string, string>>& indexes, string& ba
     Haplotypes::Verbosity verbosity = (progress ? Haplotypes::verbosity_basic : Haplotypes::verbosity_silent);
     Recombinator recombinator(gbz, verbosity);
     Recombinator::Parameters parameters;
+    parameters.num_haplotypes = Recombinator::NUM_CANDIDATES;
+    parameters.diploid_sampling = true;
     parameters.include_reference = true;
     gbwt::GBWT sampled_gbwt;
     try {

diff --git a/src/subcommand/haplotypes_main.cpp b/src/subcommand/haplotypes_main.cpp
@@ -179,7 +179,7 @@ int main_haplotypes(int argc, char** argv) {
     return 0;
 }
 
-static vg::subcommand::Subcommand vg_haplotypes("haplotypes", "haplotype sampling based on kmer counts", vg::subcommand::DEVELOPMENT, main_haplotypes);
+static vg::subcommand::Subcommand vg_haplotypes("haplotypes", "haplotype sampling based on kmer counts", vg::subcommand::TOOLKIT, main_haplotypes);
 
 //----------------------------------------------------------------------------
 

diff --git a/test/t/54_vg_haplotypes.t b/test/t/54_vg_haplotypes.t
@@ -38,11 +38,18 @@ is $(vg gbwt -S -Z no_ref.gbz) 1 "1 sample"
 is $(vg gbwt -C -Z no_ref.gbz) 2 "2 contigs"
 is $(vg gbwt -H -Z no_ref.gbz) 4 "4 haplotypes"
 
+# Diploid sampling
+vg haplotypes --validate -i full.hapl -k haplotype-sampling/HG003.kff --include-reference --diploid-sampling -g diploid.gbz full.gbz
+is $? 0 "diploid sampling"
+is $(vg gbwt -S -Z diploid.gbz) 3 "1 generated + 2 reference samples"
+is $(vg gbwt -C -Z diploid.gbz) 2 "2 contigs"
+is $(vg gbwt -H -Z diploid.gbz) 4 "2 generated + 2 reference haplotypes"
+
 # Giraffe integration, guessed output name
 vg giraffe -Z full.gbz --haplotype-name full.hapl --kff-name haplotype-sampling/HG003.kff \
     -f haplotype-sampling/HG003.fq.gz > default.gam 2> /dev/null
 is $? 0 "Giraffe integration with a guessed output name"
-cmp indirect.gbz full.HG003.gbz
+cmp diploid.gbz full.HG003.gbz
 is $? 0 "the sampled graph is identical to a manually sampled one"
 
 # Giraffe integration, specified output name
@@ -53,16 +60,9 @@ is $? 0 "Giraffe integration with a specified output name"
 cmp full.HG003.gbz sampled.003HG.gbz
 is $? 0 "the sampled graphs are identical"
 
-# Diploid sampling
-vg haplotypes --validate -i full.hapl -k haplotype-sampling/HG003.kff --include-reference --diploid-sampling --num-haplotypes 8 -g diploid.gbz full.gbz
-is $? 0 "diploid sampling"
-is $(vg gbwt -S -Z diploid.gbz) 3 "1 generated + 2 reference samples"
-is $(vg gbwt -C -Z diploid.gbz) 2 "2 contigs"
-is $(vg gbwt -H -Z diploid.gbz) 4 "2 generated + 2 reference haplotypes"
-
 # Cleanup
 rm -r full.gbz full.ri full.dist full.hapl
 rm -f indirect.gbz direct.gbz no_ref.gbz
+rm -f diploid.gbz
 rm -f full.HG003.gbz full.HG003.dist full.HG003.min default.gam
 rm -f sampled.003HG.gbz sampled.003HG.dist sampled.003HG.min specified.gam
-rm -f diploid.gbz