|
| 1 | +# Load the S7 package |
| 2 | +library(S7) |
| 3 | + |
| 4 | +#================================================================ |
| 5 | +# Class Definitions |
| 6 | +#================================================================ |
| 7 | + |
| 8 | +# Define the base 'Model' class |
| 9 | +Model <- new_class("Model") |
| 10 | + |
| 11 | +# Define the 'DecisionTree' class |
| 12 | +DecisionTree <- new_class( |
| 13 | + "DecisionTree", |
| 14 | + parent = Model, |
| 15 | + properties = list( |
| 16 | + data = class_data.frame, |
| 17 | + N = class_numeric |
| 18 | + ), |
| 19 | + constructor = function(data, N = NA_real_) { |
| 20 | + new_object( |
| 21 | + .S7_object(), |
| 22 | + data = data, |
| 23 | + N = N |
| 24 | + ) |
| 25 | + } |
| 26 | +) |
| 27 | + |
| 28 | +# Define the 'MarkovModel' class |
| 29 | +MarkovModel <- new_class( |
| 30 | + "MarkovModel", |
| 31 | + parent = Model, |
| 32 | + properties = list( |
| 33 | + init_probs = class_list, |
| 34 | + trans_matrix = class_array, |
| 35 | + cost_matrix = class_array, |
| 36 | + q_matrix = class_array, |
| 37 | + mapping = class_numeric, |
| 38 | + N = class_numeric, |
| 39 | + n_cycles = class_integer |
| 40 | + ), |
| 41 | + constructor = function(init_probs = list(), |
| 42 | + trans_matrix = array(), |
| 43 | + cost_matrix = array(), |
| 44 | + q_matrix = array(), |
| 45 | + mapping = numeric(), |
| 46 | + N = NA_real_, |
| 47 | + n_cycles = 10L) { |
| 48 | + new_object( |
| 49 | + .S7_object(), |
| 50 | + init_probs = init_probs, |
| 51 | + trans_matrix = trans_matrix, |
| 52 | + cost_matrix = cost_matrix, |
| 53 | + q_matrix = q_matrix, |
| 54 | + mapping = mapping, |
| 55 | + N = N, |
| 56 | + n_cycles = n_cycles |
| 57 | + ) |
| 58 | + } |
| 59 | +) |
| 60 | + |
| 61 | +# Define the 'CombinedModel' class |
| 62 | +CombinedModel <- new_class( |
| 63 | + "CombinedModel", |
| 64 | + parent = Model, |
| 65 | + properties = list( |
| 66 | + models = class_list |
| 67 | + ), |
| 68 | + constructor = function(...) { |
| 69 | + models <- list(...) |
| 70 | + if (length(models) < 2) { |
| 71 | + stop("At least two models must be provided.") |
| 72 | + } |
| 73 | + if (!all(sapply(models, function(m) inherits(m, "S7_object")))) { |
| 74 | + stop("All arguments must be S7 objects.") |
| 75 | + } |
| 76 | + new_object(.S7_object(), models = models) |
| 77 | + } |
| 78 | +) |
| 79 | + |
| 80 | +#================================================================ |
| 81 | +# Method Definitions |
| 82 | +#================================================================ |
| 83 | + |
| 84 | +# Define generics |
| 85 | +run_model <- new_generic("run_model", "model") |
| 86 | +update_model <- new_generic("update_model", c("model", "result")) |
| 87 | +get_costs <- new_generic("get_costs", "results") |
| 88 | +get_effects <- new_generic("get_effects", "results") |
| 89 | + |
| 90 | +# Methods for DecisionTree |
| 91 | +method(run_model, DecisionTree) <- function(model) { |
| 92 | + res <- |
| 93 | + model@data |> |
| 94 | + dplyr::group_by(decision) |> |
| 95 | + dplyr::summarise( |
| 96 | + expected_cost = sum(probability * cost), |
| 97 | + expected_eff = sum(probability * effectiveness) |
| 98 | + ) |
| 99 | + |
| 100 | + terminal_prob <- split(x = model@data$probability, |
| 101 | + f = model@data$decision) |
| 102 | + |
| 103 | + if (!is.na(model@N)) { |
| 104 | + res$expected_cost <- res$expected_cost * model@N |
| 105 | + res$expected_eff <- res$expected_eff * model@N |
| 106 | + } |
| 107 | + |
| 108 | + # Return a list that can be used by other methods |
| 109 | + list(results = res, terminal_prob = terminal_prob) |
| 110 | +} |
| 111 | + |
| 112 | +method(update_model, DecisionTree) <- function(model, result) { |
| 113 | + # This is based on the map_markov_to_decision function |
| 114 | + # from your R/functions.R file |
| 115 | + model@data <- tibble::tibble( |
| 116 | + decision = c("Treatment A", "Treatment B"), |
| 117 | + outcome = c("Success", "Failure"), |
| 118 | + probability = c(0.6, 0.4), |
| 119 | + cost = c(900, 2100), |
| 120 | + effectiveness = c(0.88, 0.45) |
| 121 | + ) |
| 122 | + model |
| 123 | +} |
| 124 | + |
| 125 | + |
| 126 | +method(get_costs, DecisionTree) <- function(results) { |
| 127 | + results$expected_cost |
| 128 | +} |
| 129 | + |
| 130 | +method(get_effects, DecisionTree) <- function(results) { |
| 131 | + results$expected_eff |
| 132 | +} |
| 133 | + |
| 134 | +# Methods for MarkovModel |
| 135 | +method(run_model, MarkovModel) <- function(model) { |
| 136 | + res <- data.frame( |
| 137 | + decision = c("Treatment A", "Treatment B"), |
| 138 | + expected_cost = c(100, 100), |
| 139 | + expected_eff = c(1, 1) |
| 140 | + ) |
| 141 | + |
| 142 | + if (!is.na(model@N)) { |
| 143 | + res$expected_cost <- res$expected_cost * model@N |
| 144 | + res$expected_eff <- res$expected_eff * model@N |
| 145 | + } |
| 146 | + # Return a simple data frame for the Markov model output |
| 147 | + res |
| 148 | +} |
| 149 | + |
| 150 | +method(update_model, MarkovModel) <- function(model, result) { |
| 151 | + # This is based on the map_decision_to_markov function |
| 152 | + # from your R/functions.R file |
| 153 | + probs <- result$terminal_prob |
| 154 | + model@init_probs <- lapply(probs, \(x) tapply(x, model@mapping, sum)) |
| 155 | + model |
| 156 | +} |
| 157 | + |
| 158 | +method(get_costs, MarkovModel) <- function(results) { |
| 159 | + results$expected_cost |
| 160 | +} |
| 161 | + |
| 162 | +method(get_effects, MarkovModel) <- function(results) { |
| 163 | + results$expected_eff |
| 164 | +} |
| 165 | + |
| 166 | +# Methods for CombinedModel |
| 167 | +method(run_model, CombinedModel) <- function(model) { |
| 168 | + result <- list() |
| 169 | + models <- model@models |
| 170 | + |
| 171 | + for (i in seq_along(models)) { |
| 172 | + current_model <- models[[i]] |
| 173 | + |
| 174 | + if (i > 1) { |
| 175 | + current_model <- update_model(current_model, result[[i - 1]]) |
| 176 | + } |
| 177 | + result[[i]] <- run_model(current_model) |
| 178 | + } |
| 179 | + result |
| 180 | +} |
| 181 | + |
| 182 | +method(get_costs, CombinedModel) <- function(results) { |
| 183 | + total_cost <- 0 |
| 184 | + for (i in seq_along(results)) { |
| 185 | + # Assuming the first element of each result is the data frame |
| 186 | + total_cost <- total_cost + get_costs(results[[i]][[1]]) |
| 187 | + } |
| 188 | + total_cost |
| 189 | +} |
| 190 | + |
| 191 | +method(get_effects, CombinedModel) <- function(results) { |
| 192 | + total_eff <- 0 |
| 193 | + for (i in seq_along(results)) { |
| 194 | + # Assuming the first element of each result is the data frame |
| 195 | + total_eff <- total_eff + get_effects(results[[i]][[1]]) |
| 196 | + } |
| 197 | + total_eff |
| 198 | +} |
| 199 | + |
| 200 | +# Load libraries |
| 201 | +library(tibble) |
| 202 | +library(dplyr) |
| 203 | +library(S7) |
| 204 | + |
| 205 | +# Source the S7 functions from the file where you saved them |
| 206 | +source("s7_functions.R") # Make sure to change this to the correct file path |
| 207 | + |
| 208 | +#================================================================ |
| 209 | +# Example Data |
| 210 | +#================================================================ |
| 211 | + |
| 212 | +decision_tree_data <- tibble( |
| 213 | + decision = rep(c("Treatment A", "Treatment B"), each = 2), |
| 214 | + outcome = c("Success", "Failure", "Success", "Failure"), |
| 215 | + probability = c(0.7, 0.3, 0.6, 0.4), |
| 216 | + cost = c(1000, 2000, 800, 2500), |
| 217 | + effectiveness = c(0.9, 0.5, 0.85, 0.4) |
| 218 | +) |
| 219 | + |
| 220 | +# Define dummy Markov model parameters |
| 221 | +trans_prob_mat <- array( |
| 222 | + c(0.9, 0.1, 0.2, 0.8, 0.9, 0.1, 0.2, 0.8), |
| 223 | + dim = c(2, 2, 2), |
| 224 | + dimnames = list(NULL, NULL, c("A", "B")) |
| 225 | +) |
| 226 | + |
| 227 | +cost_mat <- array( |
| 228 | + c(1000, 2000), |
| 229 | + dim = c(1, 2, 2), |
| 230 | + dimnames = list(NULL, NULL, c("A", "B")) |
| 231 | +) |
| 232 | + |
| 233 | +q_mat <- array( |
| 234 | + c(1, 0), |
| 235 | + dim = c(1, 2, 2), |
| 236 | + dimnames = list(NULL, NULL, c("A", "B")) |
| 237 | +) |
| 238 | + |
| 239 | +# Mapping from decision tree terminal nodes to Markov states |
| 240 | +mapping <- c(1, 2) |
| 241 | + |
| 242 | +#================================================================ |
| 243 | +# Build and Run Models |
| 244 | +#================================================================ |
| 245 | + |
| 246 | +# Build individual models |
| 247 | +dt <- DecisionTree(data = decision_tree_data) |
| 248 | +mm <- MarkovModel( |
| 249 | + trans_matrix = trans_prob_mat, |
| 250 | + cost_matrix = cost_mat, |
| 251 | + q_matrix = q_mat, |
| 252 | + mapping = mapping |
| 253 | +) |
| 254 | + |
| 255 | +# Run individual models |
| 256 | +run_model(dt) |
| 257 | +run_model(mm) |
| 258 | + |
| 259 | +# Combine models and run the chain |
| 260 | +full_model <- CombinedModel(dt, mm) |
| 261 | +final_result <- run_model(full_model) |
| 262 | + |
| 263 | +# Print the final result |
| 264 | +print(final_result) |
| 265 | + |
0 commit comments