autoincremente the number of diffsel iterations

etc/docker/diffsel/Dockerfile

@@ -34,3 +34,11 @@ RUN apt-get update && \
 
 RUN pip install --upgrade pip
 RUN pip install docopt
+
+RUN apt-get update && \
+    apt-get install --no-install-recommends -qy \
+    r-base \
+    pandoc
+
+RUN echo 'install.packages(c("rmarkdown","coda"), repos="https://mirror.ibcp.fr/pub/CRAN/", dependencies=TRUE)' > /tmp/packages.R \
+&& Rscript /tmp/packages.R

lib/diffsel.ml

@@ -3,6 +3,58 @@ open Bistro.Std
 open Bistro.EDSL
 open File_formats
 
+let assign k v =
+  seq ~sep:"=" [ string k ; v ]
+
+let bash_script args code =
+  let prelude =
+    args
+    |> List.map ~f:(fun (k, v) ->
+        assign k v
+      )
+    |> seq ~sep:"\n"
+  in
+  seq ~sep:"\n" [ prelude ; string code ]
+
+
+let diffsel_add_iterations_script ~chainname =
+  let vars = [
+    "CHAIN", chainname ;
+  ]
+  in
+  bash_script vars {|
+
+#!/bin/bash
+set -e
+
+continue=true
+
+while $continue 
+do
+
+# check convergence
+Rscript -e "rmarkdown::render(\"DiffselMCMCConvergenceAnalysis.Rmd\", params=list(set_trace1=\"../dest/myrun.trace\"))"
+new_max=`tail -n 1 new_iterations.txt | cut -f 1`
+continue=`tail -n 1 new_iterations.txt | cut -f 2`
+end_it=`tail -n 1 new_iterations.txt | cut -f 3`
+
+echo new_max=$new_max
+echo continue=$continue
+echo end_it=$end_it
+
+if $continue
+then
+    echo "diffsel $CHAIN $new_max"
+    /diffsel/_build/diffsel $CHAIN $new_max
+fi
+
+done
+
+cat new_iterations.txt > $CHAIN.iterations
+
+|}
+
+
 let diffsel ~(phy_n:nucleotide_phylip workflow) ~(tree: _ workflow) ~(w_every:int) ~(n_cycles: int) ~(id: int) ~tag : [`diffsel] directory workflow =
   let env = docker_image ~account:"vlanore" ~name:"diffsel" ~tag () in
   let tmp_tree = tmp // "myrun.tree" in
@@ -10,12 +62,16 @@ let diffsel ~(phy_n:nucleotide_phylip workflow) ~(tree: _ workflow) ~(w_every:in
   let dest_tree = dest // "myrun.tree" in
   let dest_ali = dest // "myrun.ali" in
   let chainname = dest // "myrun" in
+  let n_cycles = if (n_cycles > 200) then 200 else n_cycles in
+  let script_r = tmp // "DiffselMCMCConvergenceAnalysis.Rmd" in
   (*_build/diffsel -t data/samhd1.tree -d data/samhd1.ali -ncond 3 -x 1 10000 myrun*)
   workflow ~descr:("convergence_detection.run_diffsel." ^ tag ^ "." ^(string_of_int id))  [
     docker env (
       and_list [
         mkdir_p dest;
+        cd tmp;
         cmd "echo" [string "Run chain:"; int id];
+        cmd "cp" [ file_dump (string Scripts.diffselMCMCConvergenceAnalysis) ; script_r] ;
         cmd "cp" [dep phy_n; dest_ali]; (* required dep to link the file in the env *)
         cmd "cp" [dep tree; dest_tree]; (* required dep to link the file in the env *)
         cmd "cp" [dep phy_n; tmp_ali]; (* required dep to link the file in the env *)
@@ -27,6 +83,7 @@ let diffsel ~(phy_n:nucleotide_phylip workflow) ~(tree: _ workflow) ~(w_every:in
           opt "-x" seq [ int w_every; string " "; int n_cycles];
           ident chainname ;
         ];
+       cmd "bash" [(file_dump (diffsel_add_iterations_script ~chainname ))];
       ]
     )
   ]
@@ -36,7 +93,7 @@ let check_conv run_diffsel : text_file directory workflow =
   let script = tmp // "DiffselMCMCConvergenceAnalysis.Rmd" in
   let trace = run_diffsel / selector["myrun.trace"] in
   let out = dest // "out.html" in
-  let nb_new_iterations = dest // "new_iterations.txt.txt" in
+  let nb_new_iterations = dest // "new_iterations.txt" in
   workflow ~descr:"convergence_detection.DiffselMCMCConvergenceAnalysis" [
     docker env (
       and_list [

lib/scripts/diffselMCMCConvergenceAnalysis.Rmd

@@ -66,21 +66,44 @@ We require that estimation of the 0.025 quantile should be precise at the 1% lev
 md1 <- mcmc(d1)
 print(summary(md1))
 test = raftery.diag(md1, r=0.01)
-
 print(test)
+```
 
+```{r, echo = F}
 required_it = as.numeric(test$resmatrix[2])
 run_it = dim(md1)[1]
+next_step_its = if (required_it > (200 + run_it)) {(200 + run_it)} else {required_it} # made step of 200 iterations
 
 new_iterations = required_it - run_it
-new_iterations_str = if (new_iterations > 0) {as.character(new_iterations)} else {"0"}
-writeLines(new_iterations_str, "new_iterations.txt")
-```
 
-`r run_it ` iterations were already run.
+new_iterations_bool = if (required_it > run_it) {T} else {F}
+new_iterations_bool_str = if (required_it > run_it) {"true"} else {"false"}
+new_iterations_str=paste(as.character(next_step_its), new_iterations_bool_str, required_it, sep="\t")
 
-`r required_it ` iterations are needed.
+if (new_iterations_bool_str) {
+write(new_iterations_str, "new_iterations.txt", append=T)
+} else {
+write(paste("Burnin:", str(required_it)) , "new_iterations.txt", append=T)
+}
+```
 
-`r new_iterations` additionnal iterations must be run.
+```{r , echo = F}
+if (new_iterations_bool_str) {
+cat(run_it)
+cat(" iterations were already run.")
+print("")
+cat(required_it)
+cat(" iterations are needed.")
+print("")
+cat(new_iterations)
+cat(" additionnal iterations must be run.")
+print("")
+cat(next_step_its)
+cat(" additionnal iterations will be run.")
+} else {
+cat("Burnin: ")
+cat(required_it)
+}
+```