Pipeline2: added more helpers

.ocamlinit

@@ -4,7 +4,7 @@ open Pipeline2
 let h0 =
   bppseqgen_simulation
     ~hyp:Convergence_hypothesis.(H0 (Fixed 5.))
-    ~tree:"example/trees_analyses/cyp_coding.Chrysithr_root.nhx"
+    ~tree:(`NHX "example/trees_analyses/cyp_coding.Chrysithr_root.nhx")
     ~profiles:"example/aa_fitness/263SelectedProfiles.tsv"
     ~nb_sites:10
     ~seed:42
@@ -12,37 +12,57 @@ let h0 =
 let ha =
   bppseqgen_simulation
     ~hyp:Convergence_hypothesis.(HaPC (Fixed 5.))
-    ~tree:"example/trees_analyses/cyp_coding.Chrysithr_root.nhx"
+    ~tree:(`NHX "example/trees_analyses/cyp_coding.Chrysithr_root.nhx")
     ~profiles:"example/aa_fitness/263SelectedProfiles.tsv"
     ~nb_sites:10
     ~seed:42
 
 let sim =
-  Convdet_simulation {
-    tree = `NHX "example/trees_test/tree_small_bl.nhx" ;
-    profiles = "example/aa_fitness/263SelectedProfiles.tsv" ;
-    n_h0 = 10 ;
-    n_ha = 5 ;
-    ne_s = 8. ;
-    branch_factor = 1. ;
-  }
+  convdet_simulation
+    ~tree:(`NHX "example/trees_test/tree_small_bl.nhx")
+    ~profiles:"example/aa_fitness/263SelectedProfiles.tsv"
+    ~n_h0:10
+    ~n_ha:5
+    ~ne_s:8.
+    ~seed:42
+    ()
 
 let sim2 =
-  Convdet_simulation {
-    tree = `NHX "example/trees_analyses/cyp_coding.Chrysithr_root.nhx" ;
-    profiles = "example/aa_fitness/263SelectedProfiles.tsv" ;
-    n_h0 = 40 ;
-    n_ha = 10 ;
-    ne_s = 10. ;
-    branch_factor = 3. ;
-  }
+  convdet_simulation
+    ~tree:(`NHX "example/trees_analyses/cyp_coding.Chrysithr_root.nhx")
+    ~profiles:"example/aa_fitness/263SelectedProfiles.tsv"
+    ~n_h0:40
+    ~n_ha:10
+    ~ne_s:8.
+    ~seed:428
+    ()
+
+let convdet_cyp_coding =
+  convdet_simulation
+    ~tree:(`NHX "example/trees_analyses/cyp_coding.Chrysithr_root.nhx")
+    ~profiles:"example/aa_fitness/263SelectedProfiles.tsv"
+
+let convdet_pairs ~bl n =
+  convdet_simulation
+    ~tree:(`Pair_tree (bl /. 10., bl, n))
+    ~profiles:"example/aa_fitness/263SelectedProfiles.tsv"
+
+let bppseqgen_cyp_coding =
+  bppseqgen_mixed_simulation
+    ~tree:(`NHX "example/trees_analyses/cyp_coding.Chrysithr_root.nhx")
+    ~profiles:"example/aa_fitness/263SelectedProfiles.tsv"
+
+let bppseqgen_pairs ~bl n =
+  bppseqgen_mixed_simulation
+    ~tree:(`Pair_tree (bl /. 10., bl, n))
+    ~profiles:"example/aa_fitness/263SelectedProfiles.tsv"
 
 let sim3 =
   Bppseqgen_mixed {
     tree = `NHX "example/trees_analyses/cyp_coding.Chrysithr_root.nhx" ;
     profiles = "example/aa_fitness/263SelectedProfiles.tsv" ;
-    n_h0 = 40 ;
-    n_ha = 10 ;
-    ne_s = 10. ;
+    n_h0 = 400 ;
+    n_ha = 0 ;
+    ne_s = 4. ;
     seed = 42 ;
   }

lib/pipeline2.ml

@@ -8,16 +8,18 @@ let calc_fixed_seed ~(str:string) (seed:int) : int =
   let str_hash = Hashtbl.hash str in
   Hashtbl.hash (str_hash + seed)
 
+type tree = [`NHX of string | `Pair_tree of float * float * int]
+
 type dataset =
   | Bppseqgen_simulation of {
       hypothesis : Convergence_hypothesis.t ;
-      tree : [`NHX of string] ;
+      tree : tree ;
       profiles : string ;
       nb_sites : int ;
       seed : int ;
     }
   | Bppseqgen_mixed of {
-      tree : [`NHX of string] ;
+      tree : tree ;
       profiles : string ;
       seed : int ;
       n_h0 : int ;
@@ -25,42 +27,73 @@ type dataset =
       ne_s : float ;
     }
   | Convdet_simulation of {
-      tree : [`NHX of string] ;
+      tree : tree ;
       branch_factor : float ;
       profiles : string ;
       n_h0 : int ;
       n_ha : int ;
       ne_s : float ;
+      seed : int ;
     }
 
+let bppseqgen_mixed_simulation ?(ne_s = 1.) ?(seed = 0) ~tree ~profiles ~n_h0 ~n_ha () =
+  Bppseqgen_mixed {
+    tree ;
+    profiles ;
+    seed ;
+    n_ha ;
+    n_h0 ;
+    ne_s ;
+  }
+
 let bppseqgen_simulation ~hyp ~tree ~profiles ~nb_sites ~seed =
   Bppseqgen_simulation {
     hypothesis = hyp ;
-    tree = `NHX tree ;
+    tree ;
     profiles ;
     nb_sites ;
     seed ;
   }
 
-let tree_prefix = function
-  | Bppseqgen_simulation { tree = `NHX path ; _ }
-  | Bppseqgen_mixed { tree = `NHX path ; _ }
-  | Convdet_simulation { tree = `NHX path ; _ } ->
+let convdet_simulation ?(branch_factor = 1.) ?(ne_s = 1.) ?(seed = 0) ~tree ~profiles ~n_h0 ~n_ha () =
+  Convdet_simulation {
+    tree ;
+    profiles ;
+    n_h0 ;
+    n_ha ;
+    ne_s ;
+    branch_factor ;
+    seed : int ;
+  }
+
+let prefix_of_tree = function
+  | `NHX path ->
     Caml.Filename.chop_extension path
+  | `Pair_tree (bl1, bl2, np) ->
+    sprintf "pair_tree(bl1=%g,bl2=%g,np=%d)" bl1 bl2 np
 
-let tree_dataset = function
-  | Bppseqgen_simulation { tree = `NHX tree ; _ }
-  | Bppseqgen_mixed { tree = `NHX tree ; _ }
-  | Convdet_simulation { tree = `NHX tree ; _ } as sim ->
-    Tree_dataset.prepare
-      ~descr:("simulated_data." ^ (tree_prefix sim))
-      (Workflow.input tree)
+let tree_prefix = function
+  | Bppseqgen_simulation { tree ; _ }
+  | Bppseqgen_mixed { tree ; _ }
+  | Convdet_simulation { tree ; _ } ->
+    prefix_of_tree tree
 
 let tree = function
-  | Bppseqgen_simulation { tree = `NHX path ; _ }
-  | Bppseqgen_mixed { tree = `NHX path ; _ }
-  | Convdet_simulation { tree = `NHX path ; _ } ->
-    Workflow.input path
+  | Bppseqgen_simulation { tree ; _ }
+  | Bppseqgen_mixed { tree ; _ }
+  | Convdet_simulation { tree ; _ } ->
+    match tree with
+    | `NHX path ->
+      Workflow.input path
+    | `Pair_tree (branch_length1,
+                  branch_length2,
+                  npairs) ->
+    Simulator.pair_tree ~branch_length1 ~branch_length2 ~npairs
+
+let tree_dataset sim =
+  Tree_dataset.prepare
+    ~descr:("simulated_data." ^ (tree_prefix sim))
+    (tree sim)
 
 let profile ~nb_sites ~profiles =
   Profile.profile_l_of_splitted_profile
@@ -89,10 +122,10 @@ let rec nucleotide_alignment = function
     let h0 = nucleotide_alignment (Bppseqgen_simulation { hypothesis = H0 (Fixed ne_s) ; profiles ; seed ; nb_sites = n_h0 ; tree }) in
     let ha = nucleotide_alignment (Bppseqgen_simulation { hypothesis = HaPC (Fixed ne_s) ; profiles ; seed ; nb_sites = n_ha ; tree }) in
     Utils.fasta_cappend h0 ha
-  | Convdet_simulation { n_h0 ; n_ha ; profiles ; ne_s ; branch_factor ; _ } as sim ->
+  | Convdet_simulation { n_h0 ; n_ha ; profiles ; ne_s ; branch_factor ; seed ; _ } as sim ->
     let tree = tree sim in
     let fitness_profiles = Workflow.input profiles in
-    Simulator.simulator ~branch_factor ~n_ha ~n_h0 ~ne_s ~tree ~fitness_profiles ()
+    Simulator.simulator ~branch_factor ~n_ha ~n_h0 ~ne_s ~tree ~seed ~fitness_profiles ()
 
 let amino_acid_alignment d = Bppsuite.fna2faa (nucleotide_alignment d)