-
Notifications
You must be signed in to change notification settings - Fork 25
/
sensitivity_eval.R
114 lines (96 loc) · 2.67 KB
/
sensitivity_eval.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
#' Run Sensitivity Analysis
#'
#' @param model An evaluated Markov model.
#' @param dsa An object returned by
#' [define_dsa()].
#' @return A `data.frame` with one row per model and
#' parameter value.
#' @export
#'
#' @example inst/examples/example_run_dsa.R
run_dsa <- function(model, dsa) {
if (! all(c(".cost", ".effect") %in% names(get_model_results(model)))) {
stop("No cost and/or effect defined, sensitivity analysis unavailable.")
}
init <- get_uneval_init(model)
cycles <- get_cycles(model)
method <- get_method(model)
strategy_names <- get_strategy_names(model)
n_par <- length(dsa$variables)
pos_par <- cumsum(c(1, rep(c(n_par, n_par+1), n_par)))
pos_par <- pos_par[-length(pos_par)]
list_res <- list()
e_newdata <- list()
for (n in strategy_names) {
message(sprintf(
"Running DSA on strategy '%s'...", n
))
tab <- eval_strategy_newdata(
model,
strategy = n,
newdata = dsa$dsa
)
res <- tab %>%
dplyr::mutate_if(
names(tab) %in% dsa$variables,
to_text_dots,
name = FALSE
)
list_res <- c(
list_res,
list(res)
)
e_newdata <- c(
e_newdata,
list(unlist(lapply(
tab$.mod,
function(x) x$complete_parameters[1, dsa$variables]))[pos_par]))
names(e_newdata)[length(e_newdata)] <- n
}
for (i in seq_along(strategy_names)) {
list_res[[i]]$.strategy_names <- strategy_names[i]
}
res <-
dplyr::bind_rows(list_res) %>%
reshape_long(
key_col = ".par_names", value_col = ".par_value",
gather_cols = dsa$variables, na.rm = TRUE) %>%
dplyr::rowwise()
e_newdata <- lapply(e_newdata, function(x){
split(x, ceiling(seq_along(x)/2))
}) %>%
tibble::as_tibble() %>%
apply(1, FUN = c) %>%
unlist()
res <- res %>%
dplyr::do(get_total_state_values(.data$.mod)) %>%
dplyr::bind_cols(res %>% dplyr::select(-.data$.mod)) %>%
dplyr::ungroup() %>%
dplyr::mutate(
.par_value_eval = unlist(e_newdata)) %>%
dplyr::mutate(
!!! compat_lazy_dots(get_ce(model)))
structure(
list(
dsa = res,
variables = dsa$variables,
model = model
),
class = c("dsa", "list")
)
}
get_model.dsa <- function(x) {
x$model
}
digits_at_diff <- function(x, y, addl_digits = 1){
stopifnot(length(x) == length(y))
diff <- abs(x - y)
num_digits <- -floor(log(diff, 10)) + addl_digits
round_x <-
sapply(seq(along = x),
function(i){round(x[i], num_digits[i])})
round_y <-
sapply(seq(along = y),
function(i){round(y[i], num_digits[i])})
list(x = round_x, y = round_y, nd = num_digits)
}