run simulation for measurement 50 times

Paper/MeasurementVariance/MeasurementVariance.R

@@ -1,3 +1,68 @@
+ida <- function() {
+  lowerBounds <- c(
+    kG = 0,
+    I0 = 0,
+    IHD = 0,
+    ID = 0
+  )
+  upperBounds <- c(
+    kG = 10^10,
+    I0 = 10^2,
+    IHD = 10^10,
+    ID = 10^10
+  )
+  additionalParameters <- c(
+    host = 4.3 * 10^-6,
+    dye = 6 * 10^-6,
+    kHD = 1.7E07
+  )
+
+  tsf::batch("ida",
+    lowerBounds, upperBounds,
+    path = "../IDA.txt",
+    additionalParameters,
+    ngen = 200,
+    num_rep = 100, num_cores = 8,
+    errorThreshold = 0.6
+  )
+}
+res <- ida()
+save(res, file = "ida_100.RData")
+stop()
+
+gda <- function() {
+  lowerBounds <- c(
+    kG = 100,
+    I0 = 0,
+    IHD = 0,
+    ID = 100
+  )
+  # NOTE: upper bounds are set due to results of preliminary runs
+  upperBounds <- c(
+    kG = 10^6,
+    I0 = 10^3,
+    IHD = 10^10,
+    ID = 10^6
+  )
+  additionalParameters <- c(
+    host = 50 * 10^-6,
+    guest = 292 * 10^-6,
+    kHD = 33000
+  )
+  tsf::batch("gda",
+    lowerBounds, upperBounds,
+    path = "../NewGDA/GDA_system_3.txt",
+    additionalParameters,
+    ngen = 100,
+    num_rep = 100, num_cores = 10,
+    errorThreshold = 1
+  )
+}
+res <- gda()
+save(res, file = "gda_100.RData")
+stop()
+
+
 dba <- function() {
   lowerBounds <- c(
     kHD = 0,
@@ -19,7 +84,7 @@ dba <- function() {
     path = "../DBA.txt",
     additionalParameters,
     ngen = 150,
-    num_rep = 5, num_cores = 5,
+    num_rep = 100, num_cores = 10,
     errorThreshold = 0.6
   )
 }
@@ -32,7 +97,7 @@ p <- Reduce(rbind, p)
 p$error <- m[, 1]
 print(p)
 
-save(res, file = "dba.RData")
+save(res, file = "dba_100Runs.RData")
 stop()
 
 

Paper/MeasurementVariance/Visualisation.R

@@ -5,16 +5,71 @@ empty_plot <- ggplot() +
   theme_void()
 size_dashes <- 0.25
 
+extract_best_runs <- function(res) {
+  states <- res$states
+  params <- res$params
+  errors <- res$metrices
+  errors <- Reduce(rbind, errors)
+  states <- Reduce(rbind, states)
+  params <- Reduce(rbind, params)
+  states <- lapply(unique(errors$dataset), function(x) {
+    states_subset <- states[states$dataset == x, ]
+    errors_subset <- errors[errors$dataset == x, ]
+    errors_subset <- errors_subset[order(errors_subset$MeanSquareError), ][1:50,]
+    states_subset[states_subset$repetition %in% errors_subset$repetition, ]
+  })
+  states <- Reduce(rbind, states)
+  params <- lapply(unique(errors$dataset), function(x) {
+    params_subset <- params[params$dataset == x, ]
+    errors_subset <- errors[errors$dataset == x, ]
+    errors_subset <- errors_subset[order(errors_subset$MeanSquareError), ][1:50,]
+    params_subset[params_subset$repetition %in% errors_subset$repetition, ]
+  })
+  params <- Reduce(rbind, params)
+  return(list(states = states, params = params))
+}
+
 param_plot <- function(path) {
   load(path)
-  res[[3]] + scale_fill_brewer(palette = "Dark2")
+  params <- extract_best_runs(res[[1]])$params
+  params <- data.frame(
+    y = c(params[, 1], params[, 2], params[, 3], params[, 4]),
+    names = rep(names(params)[1:4], each = nrow(params)),
+    dataset = rep(params[, "dataset"], 4)
+  )
+  ggplot() +
+    geom_boxplot(
+      data = params,
+      aes(
+        y = y, fill = "Entire data", x = factor(0)
+      )
+    ) +
+    geom_boxplot(
+      data = params,
+      aes(
+        x = factor(dataset), y = y,
+        group = factor(dataset),
+        fill = factor(dataset)
+      )
+    ) +
+    facet_wrap(. ~ names,
+      scales = "free_y",
+      strip.position = "left"
+    ) +
+    xlab(NULL) +
+    ylab(NULL) +
+    theme(
+      panel.spacing = unit(2, "lines"),
+      strip.background = element_blank(),
+      strip.placement = "outside"
+    ) +
+    guides(fill = guide_legend(title = "Datasets")) +
+    scale_fill_brewer(palette = "Dark2")
 }
 
 plot_fct <- function(path) {
   load(path)
-  states <- res[[1]]$states
-  states <- Reduce(rbind, states)
-
+  states <- extract_best_runs(res[[1]])$states
   df <- data.frame(
     var = rep(states[, 1], 2),
     signal = c(states[, 2], states[, 3]),
@@ -48,19 +103,19 @@ plot_fct <- function(path) {
     )
 }
 
-p_dba <- plot_fct("dba.RData")
-p_ida <- plot_fct("ida.RData")
-p_gda <- plot_fct("gda.RData")
+p_dba <- plot_fct("dba_100Runs.RData")
+p_ida <- plot_fct("ida_100.RData")
+p_gda <- plot_fct("gda_100.RData")
 p_dba <- p_dba + theme(legend.position = "none")
 p_ida <- p_ida + theme(legend.position = "none")
 p1 <- plot_grid(p_dba, p_ida, p_gda,
   nrow = 1,
   labels = c("a", "b", "c")
 )
 
-p_dba <- param_plot("dba.RData") + theme(legend.position = "none")
-p_ida <- param_plot("ida.RData") + theme(legend.position = "none")
-p_gda <- param_plot("gda.RData")
+p_dba <- param_plot("dba_100Runs.RData") + theme(legend.position = "none")
+p_ida <- param_plot("ida_100.RData") + theme(legend.position = "none")
+p_gda <- param_plot("gda_100.RData")
 
 p2 <- plot_grid(p_dba, p_ida, p_gda,
   nrow = 1,
@@ -71,7 +126,7 @@ p <- plot_grid(p1, p2, nrow = 2)
 
 ggsave(p,
   bg = "white",
-  file = "variance.png",
+  file = "variance_50reps.png",
   width = 15,
   height = 8
 )