| Title: | Decision Analysis Modelling Package with Parameters Estimation Ability from Individual Patient Level Data |
|---|---|
| Description: | A collection of functions to construct Markov model for model-based cost-effectiveness analysis. This includes creating Markov model (both time homogenous and time dependent models), decision analysis, sensitivity analysis (deterministic and probabilistic). The package allows estimation of parameters for the Markov model from a given individual patient level data, provided the data file follows some standard data entry rules. |
| Authors: | Sheeja Manchira Krishnan [aut, cre] |
| Maintainer: | Sheeja Manchira Krishnan <[email protected]> |
| License: | GPL-3 |
| Version: | 1.1.0 |
| Built: | 2024-09-18 05:59:30 UTC |
| Source: | https://github.com/cranhaven/cranhaven.r-universe.dev |
Function to get sum of entries of resource per individual at diff timepoints if same cateogry has listed multiple time for same id of participant , this method comes in handy to get the sum
add_entries_sameuse_timepoint( use_data, timepointcol, timepointval, idcolumn, result_col )add_entries_sameuse_timepoint( use_data, timepointcol, timepointval, idcolumn, result_col )
use_data |
the data where the observations are held |
timepointcol |
columnname in the data where the timepoints are noted |
timepointval |
which time point is considered now at which the descriptive analysis is done |
idcolumn |
id for each participant |
result_col |
name of the column where the sum of entries to be saved |
the data with added sum of resource use
eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 7, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"), time = c(1,1,2,2), id = c(1, 1, 1, 2))) add_entries_sameuse_timepoint(eg_data, "time", 1, "id", ("mark_at_2"))eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 7, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"), time = c(1,1,2,2), id = c(1, 1, 1, 2))) add_entries_sameuse_timepoint(eg_data, "time", 1, "id", ("mark_at_2"))
adl_scoring table
adl_scoring.dfadl_scoring.df
A 41 by 3 dataframe
created on Jan 15, 2020
Function to assign the values of nested parameters from the parameter list
assign_parameters(param_list)assign_parameters(param_list)
param_list |
list of parameters, can be nested, or can be used the list returned from define_parameters() |
The parameter list should be a list of parameters in the form name value pairs. If the name value pairs is given as a string it throws error as in assign_parameters(c("cost_A = 100", "a = 10")) even if you use assign_parameters(define_parameters(c("cost_A = 100","a = 10"))) but this will be ok if you use the below forms assign_list2 <- c(a = 10, cost_A = "a + 100", cost_B = 10) assign_parameters(assign_list2) OR param_list <- define_parameters(a = 10, cost_A = "a + 100", cost_B = 10) assign_list <- assign_parameters(param_list) Also for nested parameters, remember to give the parameters in order so that at run time, the parameters can be evaluated for example, assign_list = define_parameters(cost_A="a+100", a=10) assign_parameters(assign_list) will throw an error, while assign_list = define_parameters( a = 10, cost_A = "a + 100") assign_parameters(assign_list) will successfully assign parameters as the parameters 'a' is visible before the calculation of 'cost_A' Another thing to note is that while using define_parameters, just enumerate them, no need to create as a list by using c() or list function
list of assigned parameters
param_list <- define_parameters( cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, cost_zido = 2456, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido" ) assign_parameters(param_list)param_list <- define_parameters( cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, cost_zido = 2456, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido" ) assign_parameters(param_list)
Parameter table created
blank.dfblank.df
A 2 column 1 observation
created on September 5, 2020
Estimation of ICER and NMB
calculate_icer_nmb(list_markov, threshold, comparator = NULL)calculate_icer_nmb(list_markov, threshold, comparator = NULL)
list_markov |
list of Markov model objects with their Markov trace, cost matrix and utility matrix |
threshold |
threshold value of WTP |
comparator |
the strategy to be compared with |
ICER and NMB for all the strategies compared to comparator
well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "control") this_markov <- markov_model(this.strategy, 24, c(1000, 0, 0), c(0,0)) well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 10, utility = 0.5) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.4, 0.4, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "intervention") sec_markov <- markov_model(this.strategy, 24, c(1000, 0, 0), c(0,0)) list_markov <- combine_markov(this_markov, sec_markov) calculate_icer_nmb(list_markov, 20000, comparator = "control")well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "control") this_markov <- markov_model(this.strategy, 24, c(1000, 0, 0), c(0,0)) well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 10, utility = 0.5) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.4, 0.4, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "intervention") sec_markov <- markov_model(this.strategy, 24, c(1000, 0, 0), c(0,0)) list_markov <- combine_markov(this_markov, sec_markov) calculate_icer_nmb(list_markov, 20000, comparator = "control")
Function to check the equality of column contents between two data sets with omitting NA
check_equal_columncontents_NAomitted( data1, data2, samecol, column_data1, column_data2 )check_equal_columncontents_NAomitted( data1, data2, samecol, column_data1, column_data2 )
data1 |
first data set |
data2 |
second data set |
samecol |
a unique col in both datasets, like tno or id |
column_data1 |
column name in data 1 to be compared |
column_data2 |
column name in data 2 to be compared |
0 if they are equal else error message
eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 7, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"), time = c(1,1,2,2), id = c(1, 1, 1, 2))) eg_data2 <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 27, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"), time = c(1,1,2,2), id = c(1, 1, 1, 2))) check_equal_columncontents_NAomitted(eg_data, eg_data2, "no", "mark_at_1","mark_at_2")eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 7, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"), time = c(1,1,2,2), id = c(1, 1, 1, 2))) eg_data2 <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 27, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"), time = c(1,1,2,2), id = c(1, 1, 1, 2))) check_equal_columncontents_NAomitted(eg_data, eg_data2, "no", "mark_at_1","mark_at_2")
Function to check the sum of column contents between two data sets with omitting NA
check_equal_sumcolumncontents_NAomitted( data1, data2, samecol, column_data1, column_data2 )check_equal_sumcolumncontents_NAomitted( data1, data2, samecol, column_data1, column_data2 )
data1 |
first data set |
data2 |
second data set |
samecol |
a unique col in both datasets, like tno or id |
column_data1 |
column name in data 1 to be compared |
column_data2 |
column name in data 2 to be compared |
0 if they are equal else error message
eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 7, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"), time = c(1,1,2,2), id = c(1, 1, 1, 2))) eg_data2 <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 27, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"), time = c(1,1,2,2), id = c(1, 1, 1, 2))) check_equal_sumcolumncontents_NAomitted(eg_data, eg_data2, "no", "mark_at_1","mark_at_2")eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 7, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"), time = c(1,1,2,2), id = c(1, 1, 1, 2))) eg_data2 <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 27, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"), time = c(1,1,2,2), id = c(1, 1, 1, 2))) check_equal_sumcolumncontents_NAomitted(eg_data, eg_data2, "no", "mark_at_1","mark_at_2")
Function to find the keyword for family of distribution in glm
check_link_glm(family, link)check_link_glm(family, link)
family |
family of distribution |
link |
function to be used |
Check and get the link function for the method glm
the link if they can be accepted else error
check_link_glm("gaussian", "identity")check_link_glm("gaussian", "identity")
check the list of Markov models
check_list_markov_models(list_markov)check_list_markov_models(list_markov)
list_markov |
list of Markov model objects with their Markov trace, cost matrix and utility matrix |
0 if success else error
well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "example") this_markov <- markov_model(this.strategy, 24, c(1000, 0, 0), c(0,0)) well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 10, utility = 0.5) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.4, 0.4, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "example_two") sec_markov <- markov_model(this.strategy, 24, c(1000, 0, 0), c(0, 0)) list_markov <- combine_markov(this_markov, sec_markov) check_list_markov_models(list_markov)well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "example") this_markov <- markov_model(this.strategy, 24, c(1000, 0, 0), c(0,0)) well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 10, utility = 0.5) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.4, 0.4, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "example_two") sec_markov <- markov_model(this.strategy, 24, c(1000, 0, 0), c(0, 0)) list_markov <- combine_markov(this_markov, sec_markov) check_list_markov_models(list_markov)
Function to check the variable null or NA
check_null_na(variable)check_null_na(variable)
variable |
name of variable or list of variable to check |
-1 or -2 as error, else return 0 as success
var = c("a") check_null_na(var)var = c("a") check_null_na(var)
Check the transition probabilities for numeric values and unity row sum
check_trans_prob(trans_mat)check_trans_prob(trans_mat)
trans_mat |
transition matrix |
checking for rowsum - checks for the class of transition matrix, value of rowsum (to be 1) and numeric values
0 if they add to 1 else error
tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- populate_transition_matrix(2, tmat, list_prob = c(0.5, 0.5, 0, 1)) check_trans_prob(tm)tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- populate_transition_matrix(2, tmat, list_prob = c(0.5, 0.5, 0, 1)) check_trans_prob(tm)
Function to return treatment arm
check_treatment_arm(arm)check_treatment_arm(arm)
arm |
the arm of the trial |
0, if success -1, if failure
check_treatment_arm("control")check_treatment_arm("control")
Check if the values of health states are provided
check_values_states(health_states)check_values_states(health_states)
health_states |
list of health_state objects |
This is to check if the values are numeric during the run time, else to throw an error
true or false
well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1)) health_states <- combine_state(well, disabled, dead) check_values_states(health_states)well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1)) health_states <- combine_state(well, disabled, dead) check_values_states(health_states)
Checks the input to run the Markov cycles and picks correct method
checks_markov_pick_method( current_strategy, initial_state, discount, method, half_cycle_correction, startup_cost, startup_util, state_cost_only_prevalent, state_util_only_prevalent )checks_markov_pick_method( current_strategy, initial_state, discount, method, half_cycle_correction, startup_cost, startup_util, state_cost_only_prevalent, state_util_only_prevalent )
current_strategy |
strategy object |
initial_state |
value of states initially |
discount |
rate of discount for costs and qalys |
method |
what type of half cycle correction needed |
half_cycle_correction |
boolean to indicate half cycle correction |
startup_cost |
cost of states initially |
startup_util |
utility of states initially if any |
state_cost_only_prevalent |
boolean parameter to indicate if the costs for state occupancy is only for those in the state excluding those that transitioned new. This is relevant when the transition cost is provided for eg. in a state with dialysis the cost of previous dialysis is different from the newly dialysis cases.Then the state_cost_only_prevalent should be TRUE |
state_util_only_prevalent |
boolean parameter to indicate if the utilities for state occupancy is only for those in the state excluding those that transitioned new. |
changed method name
tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- populate_transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0, 0, TRUE) health_states <- combine_state(a, b) this.strategy <- strategy(tm, health_states, "intervention") checks_markov_pick_method(this.strategy, c(1, 0), c(0, 0), "half cycle correction", TRUE,NULL,NULL)tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- populate_transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0, 0, TRUE) health_states <- combine_state(a, b) this.strategy <- strategy(tm, health_states, "intervention") checks_markov_pick_method(this.strategy, c(1, 0), c(0, 0), "half cycle correction", TRUE,NULL,NULL)
Function to do some checks before plotting sensitivity analysis results
checks_plot_dsa( result_dsa_control, plotfor, type, result_dsa_treat, threshold, comparator )checks_plot_dsa( result_dsa_control, plotfor, type, result_dsa_treat, threshold, comparator )
result_dsa_control |
result from deterministic sensitivity analysis for first or control model |
plotfor |
the variable to plotfor e.g. cost, utility NMB etc |
type |
type of analysis, range or difference |
result_dsa_treat |
result from deterministic sensitivity analysis for the comparative Markov model |
threshold |
threshold value of WTP |
comparator |
the strategy to be compared with |
the plot variable
Join Markov model objects
combine_markov(markov1, ...)combine_markov(markov1, ...)
markov1 |
object 1 of class markov_model |
... |
any additional objects |
Combining Markov models for easiness of comparison
joined objects of type markov_model
well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "example") this_markov <- markov_model(this.strategy, 24, c(1000, 0, 0)) well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 10, utility = 0.5) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.4, 0.4, 0.2, 0.6, 0.4, 1)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "example") sec_markov <- markov_model(this.strategy, 24, c(1000, 0, 0)) list_markov <- combine_markov(this_markov, sec_markov)well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "example") this_markov <- markov_model(this.strategy, 24, c(1000, 0, 0)) well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 10, utility = 0.5) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.4, 0.4, 0.2, 0.6, 0.4, 1)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "example") sec_markov <- markov_model(this.strategy, 24, c(1000, 0, 0)) list_markov <- combine_markov(this_markov, sec_markov)
Join health states
combine_state(...)combine_state(...)
... |
any additional objects |
checking each state is a health state and join them
joined health states
a <- health_state("IT", 100, 0.4, 0, FALSE) b <- health_state("PT", 100, 0.4, 0, FALSE) combine_state(a, b)a <- health_state("IT", 100, 0.4, 0, FALSE) b <- health_state("PT", 100, 0.4, 0, FALSE) combine_state(a, b)
Convert frequency medication to given basis
convert_freq_diff_basis(freq_given, basis = "day")convert_freq_diff_basis(freq_given, basis = "day")
freq_given |
given frequency |
basis |
given basis, default is daily |
converted frequency
convert_freq_diff_basis("once daily") convert_freq_diff_basis("bd", "week") convert_freq_diff_basis("Every 4 days", "day")convert_freq_diff_basis("once daily") convert_freq_diff_basis("bd", "week") convert_freq_diff_basis("Every 4 days", "day")
Convert period to given basis
convert_to_given_timeperiod(given_time, basis_time = "day")convert_to_given_timeperiod(given_time, basis_time = "day")
given_time |
given time |
basis_time |
given basis, default is "day" |
converted unit
convert_to_given_timeperiod("4 weeks") convert_to_given_timeperiod("a month") convert_to_given_timeperiod("1 week")convert_to_given_timeperiod("4 weeks") convert_to_given_timeperiod("a month") convert_to_given_timeperiod("1 week")
Convert volume to given basis
convert_volume_basis(given_unit, basis = "ml")convert_volume_basis(given_unit, basis = "ml")
given_unit |
given unit |
basis |
given basis, default is "ml" |
converted unit
convert_volume_basis("ml", "liter")convert_volume_basis("ml", "liter")
Convert unit strength to given basis
convert_weight_diff_basis(given_unit, basis = "mg")convert_weight_diff_basis(given_unit, basis = "mg")
given_unit |
given unit |
basis |
given basis, default is "mg" |
converted unit
convert_weight_diff_basis("mg") convert_weight_diff_basis("kilogram", "micro gram")convert_weight_diff_basis("mg") convert_weight_diff_basis("kilogram", "micro gram")
Convert weight per time to given basis
convert_wtpertimediff_basis(given_unit, basis = "mcg/hour")convert_wtpertimediff_basis(given_unit, basis = "mcg/hour")
given_unit |
given unit |
basis |
given basis, default is "mg" |
converted unit
convert_wtpertimediff_basis("mg/day") convert_wtpertimediff_basis("mcg/day") convert_wtpertimediff_basis("mg/hour")convert_wtpertimediff_basis("mg/day") convert_wtpertimediff_basis("mcg/day") convert_wtpertimediff_basis("mg/hour")
Convert wt per unit volume to given basis
convert_wtpervoldiff_basis(given_unit, basis = "mg/ml")convert_wtpervoldiff_basis(given_unit, basis = "mg/ml")
given_unit |
given unit |
basis |
given basis, default is "mg/ml" |
converted unit
convert_wtpervoldiff_basis("g/ml")convert_wtpervoldiff_basis("g/ml")
cost matrix
cost_data.dfcost_data.df
A 11 by 2 dataframe
created on Nov 26, 2019 from tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, FALSE) b <- health_state("Dead", 1, 0, TRUE) health_states <- combine_state(a, b) this.strategy <- strategy(tm, health_states, "intervention")
Function to estimate the cost of inpatient admission but taken from GP records where code or description known
costing_AandE_admission( ind_part_data, code_ae, descrip_ae, number_use_ae, type_admit_ae, unit_cost_data, code_col, type_admit_col, description_col, unit_cost_col, cost_calculated_in = "attendance", sheet = NULL )costing_AandE_admission( ind_part_data, code_ae, descrip_ae, number_use_ae, type_admit_ae, unit_cost_data, code_col, type_admit_col, description_col, unit_cost_col, cost_calculated_in = "attendance", sheet = NULL )
ind_part_data |
IPD |
code_ae |
column name of code (for inpatient admission) |
descrip_ae |
column name of description for inpatient admission |
number_use_ae |
the number of days spent in each admission if that is a criteria to be included. Otherwise each admission will be costed |
type_admit_ae |
term indicating admission and type of attendance |
unit_cost_data |
unit cost data file with code/descriptions and unit costs are listed for inpatient admission |
code_col |
code column name in unit cost data |
type_admit_col |
colname that describes type of the attendance and |
description_col |
column name of description of inpatient admission in the unit cost data |
unit_cost_col |
column name of unit cost in unit_cost_data |
cost_calculated_in |
name of unit where the cost is calculated assumed to be per admission |
sheet |
sheet where the unit costs are listed in the unit costs data file |
the calculated cost of inpatient admission long with original data
costs_file <- system.file("extdata", "National_schedule_of_NHS_costs_2019_AandE.csv", package = "packDAMipd") datafile <- system.file("extdata", "resource_use_ae_ip.csv", package = "packDAMipd") ind_part_data <- packDAMipd::load_trial_data(datafile) unit_cost_data <- packDAMipd::load_trial_data(costs_file) result <- costing_AandE_admission(ind_part_data = ind_part_data, code_ae = "code", descrip_ae = NULL, number_use_ae = "number_use", type_admit_ae = "type_admit", unit_cost_data = unit_cost_data, code_col = "Currency_Code", type_admit_col = "Service_Code", description_col = NULL, unit_cost_col ="National_Average_Unit_Cost", cost_calculated_in = "attendance")costs_file <- system.file("extdata", "National_schedule_of_NHS_costs_2019_AandE.csv", package = "packDAMipd") datafile <- system.file("extdata", "resource_use_ae_ip.csv", package = "packDAMipd") ind_part_data <- packDAMipd::load_trial_data(datafile) unit_cost_data <- packDAMipd::load_trial_data(costs_file) result <- costing_AandE_admission(ind_part_data = ind_part_data, code_ae = "code", descrip_ae = NULL, number_use_ae = "number_use", type_admit_ae = "type_admit", unit_cost_data = unit_cost_data, code_col = "Currency_Code", type_admit_col = "Service_Code", description_col = NULL, unit_cost_col ="National_Average_Unit_Cost", cost_calculated_in = "attendance")
Function to estimate the cost of inpatient admission but taken from GP records where HRG code or description known
costing_inpatient_daycase_admission( ind_part_data, hrg_code_ip_admi, descrip_ip_admi, number_use_ip_admi, elective_col, unit_cost_data, hrg_code_col, description_col, unit_cost_col, cost_calculated_in = "admission" )costing_inpatient_daycase_admission( ind_part_data, hrg_code_ip_admi, descrip_ip_admi, number_use_ip_admi, elective_col, unit_cost_data, hrg_code_col, description_col, unit_cost_col, cost_calculated_in = "admission" )
ind_part_data |
IPD |
hrg_code_ip_admi |
column name of hrg code (for inpatient admission) |
descrip_ip_admi |
column name of description for inpatient admission |
number_use_ip_admi |
the number of days spent in each admission if that is a criteria to be included. Otherwise each admission will be costed |
elective_col |
colname to say whether it is an elective admission or non elective admission |
unit_cost_data |
unit cost data file with hrg code/descriptions and unit costs are listed for inpatient admission |
hrg_code_col |
hrg code column name in unit cost data |
description_col |
column name of description of inpatient admission in the unit cost data |
unit_cost_col |
column name of unit cost in unit_cost_data |
cost_calculated_in |
name of unit where the cost is calculated assumed to be per admission |
the calculated cost of inpatient admission long with original data
costs_file <- system.file("extdata", "National_schedule_of_NHS_costs_2019.csv", package = "packDAMipd") datafile <- system.file("extdata", "resource_use_hc_ip.csv", package = "packDAMipd") ind_part_data <- packDAMipd::load_trial_data(datafile) unit_cost_data <- packDAMipd::load_trial_data(costs_file) result <- costing_inpatient_daycase_admission(ind_part_data, hrg_code_ip_admi = "HRGcode", descrip_ip_admi = NULL, number_use_ip_admi = "number_use", elective_col = "EL", unit_cost_data, hrg_code_col = "Currency_Code", description_col = NULL, unit_cost_col ="National_Average_Unit_Cost", cost_calculated_in = "admission")costs_file <- system.file("extdata", "National_schedule_of_NHS_costs_2019.csv", package = "packDAMipd") datafile <- system.file("extdata", "resource_use_hc_ip.csv", package = "packDAMipd") ind_part_data <- packDAMipd::load_trial_data(datafile) unit_cost_data <- packDAMipd::load_trial_data(costs_file) result <- costing_inpatient_daycase_admission(ind_part_data, hrg_code_ip_admi = "HRGcode", descrip_ip_admi = NULL, number_use_ip_admi = "number_use", elective_col = "EL", unit_cost_data, hrg_code_col = "Currency_Code", description_col = NULL, unit_cost_col ="National_Average_Unit_Cost", cost_calculated_in = "admission")
Function to estimate the cost of liquids when IPD is in long format
costing_opioid_liquids_averageMED_long( the_columns, ind_part_data_long, name_med, brand_med = NULL, bottle_size, bottle_size_unit = NULL, bottle_lasts, bottle_lasts_unit = NULL, preparation_dose, preparation_unit = NULL, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = NULL, basis_strength_unit = NULL )costing_opioid_liquids_averageMED_long( the_columns, ind_part_data_long, name_med, brand_med = NULL, bottle_size, bottle_size_unit = NULL, bottle_lasts, bottle_lasts_unit = NULL, preparation_dose, preparation_unit = NULL, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = NULL, basis_strength_unit = NULL )
the_columns |
columns that are to be used to convert the data from long to wide |
ind_part_data_long |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
bottle_size |
size of the bottle used |
bottle_size_unit |
unit of bottle volume |
bottle_lasts |
how long the bottle lasted |
bottle_lasts_unit |
time unit of how long the bottle lasted |
preparation_dose |
dose if preparation is given |
preparation_unit |
unit of preparatio dose |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit where the cost is calculated |
strength_column |
column column name that has strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
list_of_code_bottle_size_unit |
list of bottle size units and codes |
list_of_code_bottle_lasts_unit |
list of time of bottle lasts and codes |
list_preparation_dose_unit |
list of preparation dose units and codes |
eqdose_covtab |
table to get the conversion factor for equivalent doses, optional |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_liq_brand_empty.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) res <- costing_opioid_liquids_averageMED_wide( ind_part_data = ind_part_data, name_med = "liq_name", brand_med = "liq_brand", bottle_size = "liq_bottle_size", bottle_size_unit = NULL, bottle_lasts = "liq_lasts", bottle_lasts_unit = NULL, preparation_dose = "liq_strength", preparation_unit = NULL, timeperiod = "1 day", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = table, basis_strength_unit = NULL)med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_liq_brand_empty.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) res <- costing_opioid_liquids_averageMED_wide( ind_part_data = ind_part_data, name_med = "liq_name", brand_med = "liq_brand", bottle_size = "liq_bottle_size", bottle_size_unit = NULL, bottle_lasts = "liq_lasts", bottle_lasts_unit = NULL, preparation_dose = "liq_strength", preparation_unit = NULL, timeperiod = "1 day", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = table, basis_strength_unit = NULL)
Function to estimate the cost of liquids taken (from IPD)
costing_opioid_liquids_averageMED_wide( ind_part_data, name_med, brand_med = NULL, bottle_size, bottle_size_unit = NULL, bottle_lasts, bottle_lasts_unit = NULL, preparation_dose, preparation_unit = NULL, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = NULL, basis_strength_unit = NULL )costing_opioid_liquids_averageMED_wide( ind_part_data, name_med, brand_med = NULL, bottle_size, bottle_size_unit = NULL, bottle_lasts, bottle_lasts_unit = NULL, preparation_dose, preparation_unit = NULL, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = NULL, basis_strength_unit = NULL )
ind_part_data |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
bottle_size |
size of the bottle used |
bottle_size_unit |
unit of bottle volume |
bottle_lasts |
how long the bottle lasted |
bottle_lasts_unit |
time unit of how long the bottle lasted |
preparation_dose |
dose if preparation is given |
preparation_unit |
unit of preparation dose |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit where the cost is calculated |
strength_column |
column column name that has strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
list_of_code_bottle_size_unit |
list of bottle size units and codes |
list_of_code_bottle_lasts_unit |
list of time of bottle lasts and codes |
list_preparation_dose_unit |
list of preparation dose units and codes |
eqdose_covtab |
table to get the conversion factor for equivalent doses, optional, but the column names have to unique Similar to c("Drug", "form", "unit", "factor") or c("Drug", "form", "unit", "conversion") |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_liq.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx",package = "packDAMipd") table <- load_trial_data(conv_file) res <- microcosting_liquids_wide( ind_part_data = ind_part_data, name_med = "liq_name", brand_med = NULL, dose_med = "liq_strength", unit_med = NULL, bottle_size = "liq_bottle_size", bottle_size_unit = NULL, bottle_lasts = "liq_lasts", bottle_lasts_unit = NULL, preparation_dose = NULL, preparation_unit = NULL, timeperiod = "4 months", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_dose_unit = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = table, basis_strength_unit = NULL)med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_liq.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx",package = "packDAMipd") table <- load_trial_data(conv_file) res <- microcosting_liquids_wide( ind_part_data = ind_part_data, name_med = "liq_name", brand_med = NULL, dose_med = "liq_strength", unit_med = NULL, bottle_size = "liq_bottle_size", bottle_size_unit = NULL, bottle_lasts = "liq_lasts", bottle_lasts_unit = NULL, preparation_dose = NULL, preparation_unit = NULL, timeperiod = "4 months", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_dose_unit = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = table, basis_strength_unit = NULL)
#'########################################################################### Function to estimate the cost of patches when IPD is in long format using a IPD data of long format
costing_opioid_patches_averageMED_long( the_columns, ind_part_data_long, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )costing_opioid_patches_averageMED_long( the_columns, ind_part_data_long, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )
the_columns |
columns that are to be used to convert the data from long to wide |
ind_part_data_long |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
dose_med |
dose of medication used |
unit_med |
unit of medication ; use null if its along with the dose |
no_taken |
how many taken |
freq_taken |
frequency of medication |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit in the cost is calculated |
strength_column |
column column name that contain strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_freq |
if frequency is coded, give the code:frequency pairs, optional |
list_of_code_dose_unit |
if unit is coded, give the code:unit pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
eqdose_cov_tab |
table to get the conversion factor for equivalent doses, optional |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx",package = "packDAMipd") table <- load_trial_data(conv_file) names <- colnames(ind_part_data) ending <- length(names) ind_part_data_long <- tidyr::gather(ind_part_data, measurement, value, names[2]:names[ending], factor_key = TRUE) the_columns <- c("measurement", "value") res <- costing_opioid_patches_averageMED_long(the_columns, ind_part_data_long = ind_part_data_long, name_med = "patch_name", brand_med = "patch_brand", dose_med = "patch_strength",unit_med = NULL, no_taken = "patch_no_taken", freq_taken = "patch_frequency", timeperiod = "4 months", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = table, basis_strength_unit = "mcg/hr")med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx",package = "packDAMipd") table <- load_trial_data(conv_file) names <- colnames(ind_part_data) ending <- length(names) ind_part_data_long <- tidyr::gather(ind_part_data, measurement, value, names[2]:names[ending], factor_key = TRUE) the_columns <- c("measurement", "value") res <- costing_opioid_patches_averageMED_long(the_columns, ind_part_data_long = ind_part_data_long, name_med = "patch_name", brand_med = "patch_brand", dose_med = "patch_strength",unit_med = NULL, no_taken = "patch_no_taken", freq_taken = "patch_frequency", timeperiod = "4 months", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = table, basis_strength_unit = "mcg/hr")
Function to estimate the cost of patches taken (from IPD)
costing_opioid_patches_averageMED_wide( ind_part_data, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )costing_opioid_patches_averageMED_wide( ind_part_data, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )
ind_part_data |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
dose_med |
dose of medication used |
unit_med |
unit of medication ; use null if its along with the dose |
no_taken |
how many taken |
freq_taken |
frequency of medication |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit where the cost is calculated |
strength_column |
column column name that contain strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_freq |
if frequency is coded, give the code:frequency pairs, optional |
list_of_code_dose_unit |
if unit is coded, give the code:unit pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
eqdose_cov_tab |
table to get the conversion factor for equivalent doses, optional, but the column names have to unique Similar to c("Drug", "form", "unit", "factor") or c("Drug", "form", "unit", "conversion") |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
Assumes individual level data has name of medication, dose, dose unit, number taken, frequency taken, and basis time Assumes unit cost data contains the name of medication, form/type, strength, unit of strength (or the unit in which the cost calculated), preparation, unit cost, size and size unit (in which name, forms, size, size unit, and preparation are not passed on) @importFrom dplyr %>%
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) res <- costing_opioid_patches_averageMED_wide( ind_part_data = ind_part_data, name_med = "patch_name", brand_med = "patch_brand", dose_med = "patch_strength", unit_med = NULL, no_taken = "patch_no_taken", freq_taken = "patch_frequency", timeperiod = "4 months", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = table, basis_strength_unit = "mcg/hr")med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) res <- costing_opioid_patches_averageMED_wide( ind_part_data = ind_part_data, name_med = "patch_name", brand_med = "patch_brand", dose_med = "patch_strength", unit_med = NULL, no_taken = "patch_no_taken", freq_taken = "patch_frequency", timeperiod = "4 months", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = table, basis_strength_unit = "mcg/hr")
Function to estimate the cost of tablets when IPD is in long format
costing_opioid_tablets_averageMED_long( the_columns, ind_part_data_long, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )costing_opioid_tablets_averageMED_long( the_columns, ind_part_data_long, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )
the_columns |
columns that are to be used to convert the data from long to wide |
ind_part_data_long |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
dose_med |
dose of medication used |
unit_med |
unit of medication ; use null if its along with the dose |
no_taken |
how many taken |
freq_taken |
frequency of medication |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit where the cost is calculated |
strength_column |
column column name that contain strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_freq |
if frequency is coded, give the code:frequency pairs, optional |
list_of_code_dose_unit |
if unit is coded, give the code:unit pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
eqdose_cov_tab |
table to get the conversion factor for equivalent doses, optional |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_all.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) names <- colnames(ind_part_data) ending <- length(names) ind_part_data_long <- tidyr::gather(ind_part_data, measurement, value, names[2]:names[ending], factor_key = TRUE) the_columns <- c("measurement", "value") res <- microcosting_tablets_long(the_columns, ind_part_data_long = ind_part_data_long, name_med = "tab_name", brand_med = "tab_brand", dose_med = "tab_strength", unit_med = "tab_str_unit", no_taken = "tab_no_taken", freq_taken = "tab_frequency", timeperiod = "2 months",unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL,list_of_code_dose_unit = NULL, eqdose_cov_tab = table, basis_strength_unit = "mg")med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_all.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) names <- colnames(ind_part_data) ending <- length(names) ind_part_data_long <- tidyr::gather(ind_part_data, measurement, value, names[2]:names[ending], factor_key = TRUE) the_columns <- c("measurement", "value") res <- microcosting_tablets_long(the_columns, ind_part_data_long = ind_part_data_long, name_med = "tab_name", brand_med = "tab_brand", dose_med = "tab_strength", unit_med = "tab_str_unit", no_taken = "tab_no_taken", freq_taken = "tab_frequency", timeperiod = "2 months",unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL,list_of_code_dose_unit = NULL, eqdose_cov_tab = table, basis_strength_unit = "mg")
Function to estimate the cost of tablets taken as an average cost per equivalent dose in the opioid scenario is the morphine equivalent dose (from IPD)
costing_opioid_tablets_averageMED_wide( ind_part_data, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )costing_opioid_tablets_averageMED_wide( ind_part_data, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )
ind_part_data |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
dose_med |
dose of medication used |
unit_med |
unit of medication ; use null if its along with the dose |
no_taken |
how many taken |
freq_taken |
frequency of medication |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit where the cost is calculated |
strength_column |
column column name that contain strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_freq |
if frequency is coded, give the code:frequency pairs, optional |
list_of_code_dose_unit |
if unit is coded, give the code:unit pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
eqdose_cov_tab |
table to get the conversion factor for equivalent doses, optional, but the column names have to unique Similar to c("Drug", "form", "unit", "factor") or c("Drug", "form", "unit", "conversion") |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
Assumes individual level data has name of medication, dose, dose unit, number taken, frequency taken, and basis time Assumes unit cost data contains the name of medication, form/type, strength, unit of strength (or the unit in which the cost calculated), preparation, unit cost, size and size unit (in which name, forms, size, size unit, and preparation are not passed on) @importFrom dplyr %>%
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_all_brandNull.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) res <- costing_opioid_tablets_MED_wide(ind_part_data = ind_part_data, name_med = "tab_name", brand_med = "tab_brand", dose_med = "tab_str", unit_med = "tab_unit", no_taken = "tab_no_taken", freq_taken = "tab_frequency", timeperiod = "one day", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, eqdose_cov_tab = table, basis_strength_unit = "mg")med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_all_brandNull.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) res <- costing_opioid_tablets_MED_wide(ind_part_data = ind_part_data, name_med = "tab_name", brand_med = "tab_brand", dose_med = "tab_str", unit_med = "tab_unit", no_taken = "tab_no_taken", freq_taken = "tab_frequency", timeperiod = "one day", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, eqdose_cov_tab = table, basis_strength_unit = "mg")
Function to estimate the cost of tablets taken (from IPD)
costing_opioid_tablets_MED_wide( ind_part_data, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )costing_opioid_tablets_MED_wide( ind_part_data, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )
ind_part_data |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
dose_med |
dose of medication used |
unit_med |
unit of medication ; use null if its along with the dose |
no_taken |
how many taken |
freq_taken |
frequency of medication |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit where the cost is calculated |
strength_column |
column column name that contain strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_freq |
if frequency is coded, give the code:frequency pairs, optional |
list_of_code_dose_unit |
if unit is coded, give the code:unit pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
eqdose_cov_tab |
table to get the conversion factor for equivalent doses, optional, but the column names have to unique Similar to c("Drug", "form", "unit", "factor") or c("Drug", "form", "unit", "conversion") |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
Assumes individual level data has name of medication, dose, dose unit, number taken, frequency taken, and basis time Assumes unit cost data contains the name of medication, form/type, strength, unit of strength (or the unit in which the cost calculated), preparation, unit cost, size and size unit (in which name, forms, size, size unit, and preparation are not passed on) @importFrom dplyr %>%
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_all_brandNull.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) res <- costing_opioid_tablets_MED_wide(ind_part_data = ind_part_data, name_med = "tab_name", brand_med = "tab_brand", dose_med = "tab_str", unit_med = "tab_unit", no_taken = "tab_no_taken", freq_taken = "tab_frequency", timeperiod = "one day", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, eqdose_cov_tab = table, basis_strength_unit = "mg")med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_all_brandNull.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) res <- costing_opioid_tablets_MED_wide(ind_part_data = ind_part_data, name_med = "tab_name", brand_med = "tab_brand", dose_med = "tab_str", unit_med = "tab_unit", no_taken = "tab_no_taken", freq_taken = "tab_frequency", timeperiod = "one day", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, eqdose_cov_tab = table, basis_strength_unit = "mg")
Function to estimate the cost of resource use taken (from IPD)
costing_resource_use( ind_part_data, name_use_col, each_length_num_use = NULL, each_use_provider_indicator = NULL, unit_length_use = "day", unit_cost_data, name_use_unit_cost, unit_cost_column, cost_calculated_in, list_code_use_indicator = NULL, list_code_provider_indicator = NULL )costing_resource_use( ind_part_data, name_use_col, each_length_num_use = NULL, each_use_provider_indicator = NULL, unit_length_use = "day", unit_cost_data, name_use_unit_cost, unit_cost_column, cost_calculated_in, list_code_use_indicator = NULL, list_code_provider_indicator = NULL )
ind_part_data |
IPD |
name_use_col |
name of the column containing resource use |
each_length_num_use |
list of column names that shows length/number of repeated use eg. hospital admission |
each_use_provider_indicator |
list of column names that shows the bool indicators for the use of resource if this is to be included for the particular provider, say an nhs hospital use |
unit_length_use |
the column name that contains how many or how long used |
unit_cost_data |
unit costs data where the assumption is that the unit cost for resources such as hospital use, gp visit are listed in column resource/resource use with unit costs in another column and the units calculated as in another column |
name_use_unit_cost |
name of resource use (the column name in the unit cost data is assumed to be name/resource/type etc) in unit cost data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_in |
column name of unit where the cost is calculated |
list_code_use_indicator |
if the column name_use_col shows codes to indicate the resource use provide the list of codes and resource use for eg., list(c("yes", "no", c(1,2))) |
list_code_provider_indicator |
column each_use_provider_indicator shows codes to indicate the resource use provide the list of codes and resource use for eg., list(c("yes", "no", c(1,2))) |
the calculated cost of resource uses along with original data
costs_file <- system.file("extdata", "costs_resource_use.csv", package = "packDAMipd") datafile <- system.file("extdata", "resource_use_hc_2.csv", package = "packDAMipd") ind_part_data <- load_trial_data(datafile) unit_cost_data <- load_trial_data(costs_file) res <- costing_resource_use( ind_part_data[1, ], "hospital_admission_1", list("length_1", "length_2"), list("nhs_1", "nhs_2"), "day", unit_cost_data, "Inpatient hospital admissions", "UnitCost", "UnitUsed", NULL, NULL )costs_file <- system.file("extdata", "costs_resource_use.csv", package = "packDAMipd") datafile <- system.file("extdata", "resource_use_hc_2.csv", package = "packDAMipd") ind_part_data <- load_trial_data(datafile) unit_cost_data <- load_trial_data(costs_file) res <- costing_resource_use( ind_part_data[1, ], "hospital_admission_1", list("length_1", "length_2"), list("nhs_1", "nhs_2"), "day", unit_cost_data, "Inpatient hospital admissions", "UnitCost", "UnitUsed", NULL, NULL )
create new dataset while keeping cox regression results and returned coefficients
create_new_dataset(var, covar, dataset, categorical)create_new_dataset(var, covar, dataset, categorical)
var |
variable for which the levels have to be identified usually indep variable |
covar |
the other covariates |
dataset |
the dataset where these variables contain |
categorical |
are these variables categorical? True of false |
new data frame
dataset <- survival::lung new = create_new_dataset("status", c("age"), dataset, c(FALSE))dataset <- survival::lung new = create_new_dataset("status", c("age"), dataset, c(FALSE))
Create a table to compare the descriptive analysis (short) from gtsummary of two groups, but at different timepoints
create_shorttable_from_gtsummary_compare_twogroups_timpoints( variables, gtsummary, name_use, timepoints )create_shorttable_from_gtsummary_compare_twogroups_timpoints( variables, gtsummary, name_use, timepoints )
variables |
variables that interested |
gtsummary |
a gtsummary object that contains summary parameters |
name_use |
name of the variable or category |
timepoints |
the timepoints at which the descriptive analysis is done |
the table
eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 34, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 23, 45))) outcome_summary <- IPDFileCheck::get_summary_gtsummary(eg_data, c("gender", "mark_at_1", "mark_at_2"), byvar = "gender") variables <- "Mark" k <- create_shorttable_from_gtsummary_compare_twogroups_timpoints(variables, outcome_summary, "Category", c(1, 2))eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 34, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 23, 45))) outcome_summary <- IPDFileCheck::get_summary_gtsummary(eg_data, c("gender", "mark_at_1", "mark_at_2"), byvar = "gender") variables <- "Mark" k <- create_shorttable_from_gtsummary_compare_twogroups_timpoints(variables, outcome_summary, "Category", c(1, 2))
Create a table to compare the descriptive analysis from gtsummary of two groups
create_table_from_gtsummary_compare_twogroups(variables, gtsummary, name_use)create_table_from_gtsummary_compare_twogroups(variables, gtsummary, name_use)
variables |
variables that interested |
gtsummary |
a gtsummary object that contains summary parameters |
name_use |
name of the variable or category |
the table
eg_data <- as.data.frame(list(no = c(1,2,3,4), mark = c(12,34,23,45), gender = c("M", "F", "M", "F"))) outcome_summary <- IPDFileCheck::get_summary_gtsummary(eg_data, c("gender", "mark"), byvar = "gender") variables <- "Mark" create_table_from_gtsummary_compare_twogroups(variables, outcome_summary, "Category")eg_data <- as.data.frame(list(no = c(1,2,3,4), mark = c(12,34,23,45), gender = c("M", "F", "M", "F"))) outcome_summary <- IPDFileCheck::get_summary_gtsummary(eg_data, c("gender", "mark"), byvar = "gender") variables <- "Mark" create_table_from_gtsummary_compare_twogroups(variables, outcome_summary, "Category")
Function to return a list of parameters given
define_parameters(...)define_parameters(...)
... |
any parameters set of name value pairs expected |
To return a list of parameters For using with assign_parameters() just list or enumerate the parameters, do not use c() or list() to create a data type list
a list of parameters
define_parameters(rr = 1)define_parameters(rr = 1)
Define parameter lists for deterministic sensitivity analysis
define_parameters_psa(base_param_list, sample_list)define_parameters_psa(base_param_list, sample_list)
base_param_list |
list of parameters that used to define Markov model |
sample_list |
list of parameter values with their sampling distributions |
table for probability sensitivity analysis
param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A+ cost_comm_care_A", cost_health_B = "cost_direct_med_B+ cost_comm_care_B", cost_health_C = "cost_direct_med_C+ cost_comm_care_C", cost_drug = "cost_zido" ) sample_list <- define_parameters(cost_zido = "gamma(mean = 2756, sd = sqrt(2756))") param_table <- define_parameters_psa(param_list, sample_list)param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A+ cost_comm_care_A", cost_health_B = "cost_direct_med_B+ cost_comm_care_B", cost_health_C = "cost_direct_med_C+ cost_comm_care_C", cost_drug = "cost_zido" ) sample_list <- define_parameters(cost_zido = "gamma(mean = 2756, sd = sqrt(2756))") param_table <- define_parameters_psa(param_list, sample_list)
ISLR flextable huxtable nlme tm Define parameter lists for deterministic sensitivity analysis
define_parameters_sens_anal(param_list, low_values, upp_values)define_parameters_sens_anal(param_list, low_values, upp_values)
param_list |
list of parameters that used to define Markov model |
low_values |
list of lower values of those parameters for whom the sensitivity is to be estimated |
upp_values |
list of upper values of those parameters for whom the sensitivity is to be estimated |
Get the parameter list, min and maximum values of the parameters. The min and max values should have same entries, but they should be contained in param_list too. Copy the exact values of parameters that are in param list but not in min and max values
table for sensitivity analysis
param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido" ) low_values <- define_parameters(cost_direct_med_B = 177.4, cost_comm_care_C = 205.9) upp_values <- define_parameters(cost_direct_med_B = 17740, cost_comm_care_C = 20590) param_table <- define_parameters_sens_anal(param_list, low_values, upp_values)param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido" ) low_values <- define_parameters(cost_direct_med_B = 177.4, cost_comm_care_C = 205.9) upp_values <- define_parameters(cost_direct_med_B = 17740, cost_comm_care_C = 20590) param_table <- define_parameters_sens_anal(param_list, low_values, upp_values)
Define the table for transition
define_transition_table(tmat)define_transition_table(tmat)
tmat |
transition matrix in the format as in package 'mstate' |
Generating a table for transition matrix for efficient understanding and checking The transition matrix in the format as per 'mstate' package is transformed to a table. if tmat is not a square matrix, it gives error else it spells out the transition number, probability name and from state to state
the transition table with the probabilities
tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") define_transition_table(tmat)tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") define_transition_table(tmat)
Function to do probabilistic sensitivity analysis
do_psa(this_markov, psa_table, num_rep)do_psa(this_markov, psa_table, num_rep)
this_markov |
Markov model object |
psa_table |
table object from define_parameters_psa |
num_rep |
number of repetitions |
result after sensitivity analysis
param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A+ cost_comm_care_A", cost_health_B = "cost_direct_med_B+ cost_comm_care_B", cost_health_C = "cost_direct_med_C+ cost_comm_care_C", cost_drug = "cost_zido" ) A <- health_state("A", cost = "cost_health_A+ cost_drug ", utility = 1) B <- health_state("B", cost = "cost_health_B + cost_drug", utility = 1) C <- health_state("C", cost = "cost_health_C + cost_drug", utility = 1) D <- health_state("D", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3, 4), c(NA, 5, 6, 7), c(NA, NA, 8, 9), c(NA, NA, NA, 10)) colnames(tmat) <- rownames(tmat) <- c("A", "B", "C", "D") tm <- populate_transition_matrix(4, tmat, c( "tpAtoA", "tpAtoB", "tpAtoC", "tpAtoD", "tpBtoB", "tpBtoC", "tpBtoD", "tpCtoC", "tpCtoD", "tpDtoD" ), colnames(tmat)) health_states <- combine_state(A, B, C, D) mono_strategy <- strategy(tm, health_states, "mono") mono_markov <- markov_model(mono_strategy, 20, discount = c(0.06, 0), initial_state =c(1,0,0,0),param_list) sample_list <- define_parameters(cost_zido = "gamma(mean = 2756, sd = sqrt(2756))") param_table <- define_parameters_psa(param_list, sample_list) result <- do_psa(mono_markov, param_table, 10)param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A+ cost_comm_care_A", cost_health_B = "cost_direct_med_B+ cost_comm_care_B", cost_health_C = "cost_direct_med_C+ cost_comm_care_C", cost_drug = "cost_zido" ) A <- health_state("A", cost = "cost_health_A+ cost_drug ", utility = 1) B <- health_state("B", cost = "cost_health_B + cost_drug", utility = 1) C <- health_state("C", cost = "cost_health_C + cost_drug", utility = 1) D <- health_state("D", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3, 4), c(NA, 5, 6, 7), c(NA, NA, 8, 9), c(NA, NA, NA, 10)) colnames(tmat) <- rownames(tmat) <- c("A", "B", "C", "D") tm <- populate_transition_matrix(4, tmat, c( "tpAtoA", "tpAtoB", "tpAtoC", "tpAtoD", "tpBtoB", "tpBtoC", "tpBtoD", "tpCtoC", "tpCtoD", "tpDtoD" ), colnames(tmat)) health_states <- combine_state(A, B, C, D) mono_strategy <- strategy(tm, health_states, "mono") mono_markov <- markov_model(mono_strategy, 20, discount = c(0.06, 0), initial_state =c(1,0,0,0),param_list) sample_list <- define_parameters(cost_zido = "gamma(mean = 2756, sd = sqrt(2756))") param_table <- define_parameters_psa(param_list, sample_list) result <- do_psa(mono_markov, param_table, 10)
Function to do deterministic sensitivity analysis
do_sensitivity_analysis(this_markov, param_table)do_sensitivity_analysis(this_markov, param_table)
this_markov |
Markov model object |
param_table |
table object from define_parameters_sens_anal() with parameters (base case value, lower and upper) |
result after sensitivity analysis
param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido") low_values <- define_parameters(cost_direct_med_B = 177.4, cost_comm_care_C = 205.9) upp_values <- define_parameters(cost_direct_med_B = 17740, cost_comm_care_C = 20590) A <- health_state("A", cost = "cost_health_A + cost_drug ", utility = 1) B <- health_state("B", cost = "cost_health_B + cost_drug", utility = 1) C <- health_state("C", cost = "cost_health_C + cost_drug", utility = 1) D <- health_state("D", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3, 4), c(NA, 5, 6, 7), c(NA, NA, 8, 9), c(NA, NA, NA, 10)) colnames(tmat) <- rownames(tmat) <- c("A", "B", "C", "D") tm <- populate_transition_matrix(4, tmat, c("tpAtoA", "tpAtoB", "tpAtoC", "tpAtoD","tpBtoB", "tpBtoC", "tpBtoD", "tpCtoC", "tpCtoD", "tpDtoD"), colnames(tmat)) health_states <- combine_state(A, B, C, D) mono_strategy <- strategy(tm, health_states, "mono") mono_markov <- markov_model(mono_strategy, 20, c(1, 0, 0, 0), discount = c(0.06, 0), param_list) param_table <- define_parameters_sens_anal(param_list, low_values, upp_values) result <- do_sensitivity_analysis(mono_markov, param_table)param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido") low_values <- define_parameters(cost_direct_med_B = 177.4, cost_comm_care_C = 205.9) upp_values <- define_parameters(cost_direct_med_B = 17740, cost_comm_care_C = 20590) A <- health_state("A", cost = "cost_health_A + cost_drug ", utility = 1) B <- health_state("B", cost = "cost_health_B + cost_drug", utility = 1) C <- health_state("C", cost = "cost_health_C + cost_drug", utility = 1) D <- health_state("D", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3, 4), c(NA, 5, 6, 7), c(NA, NA, 8, 9), c(NA, NA, NA, 10)) colnames(tmat) <- rownames(tmat) <- c("A", "B", "C", "D") tm <- populate_transition_matrix(4, tmat, c("tpAtoA", "tpAtoB", "tpAtoC", "tpAtoD","tpBtoB", "tpBtoC", "tpBtoD", "tpCtoC", "tpCtoD", "tpDtoD"), colnames(tmat)) health_states <- combine_state(A, B, C, D) mono_strategy <- strategy(tm, health_states, "mono") mono_markov <- markov_model(mono_strategy, 20, c(1, 0, 0, 0), discount = c(0.06, 0), param_list) param_table <- define_parameters_sens_anal(param_list, low_values, upp_values) result <- do_sensitivity_analysis(mono_markov, param_table)
Function to get the codes and the corresponding entries
encode_codes_data(list_code_values, data_column_nos, the_data)encode_codes_data(list_code_values, data_column_nos, the_data)
list_code_values |
list of codes and values, given as list of lists |
data_column_nos |
the column numbers of data to look for the entries |
the_data |
the data where to look for |
weight and vol units
data_file <- system.file("extdata", "medication_liq_codes.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) data_column_nos = c(2,12) list_of_code_names = list(c(1, 2),c("Morphine", "Oxycodone")) encode_codes_data(list_of_code_names, data_column_nos, ind_part_data)data_file <- system.file("extdata", "medication_liq_codes.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) data_column_nos = c(2,12) list_of_code_names = list(c(1, 2),c("Morphine", "Oxycodone")) encode_codes_data(list_of_code_names, data_column_nos, ind_part_data)
Attribute parameters to probabilities of transition matrix
eval_assign_trans_prob(tm, parameter_values)eval_assign_trans_prob(tm, parameter_values)
tm |
A transition matrix in the format from the package 'mstate' |
parameter_values |
name value pairs of parameter values in the probability matrix |
Once the transition matrix is populated, the probabilities in transition matrix gets evaluated and assigned in this function call If the entry in transition matrix is NA, replaces it with zero similarly to evaluate and assign health states, the parameter values is excepted to be a list from assign_parameter() and define_parameter(). The exception is that if the parameters are defined directly and no nested calculation is required. For eg. assign_list = c(p1 = 0.2, p2 = 0.3, p3 = 0.4, p4 = 0.5) prob <- eval_assign_trans_prob(tmat, assign_list) will work For those with nested calculations, this has to be defined as below assign_list<-assign_parameters(define_parameters(p1 = 0.2, p2 = 0.3, p3 = 0.4, p4 = 0.5)) prob <- eval_assign_trans_prob(tmat, assign_list) The below will give error assign_list <- c(p1=0.1, p2 = "p1 + 0.2", p3=0, p4=0.3) prob <- eval_assign_trans_prob(tmat, assign_list)
the transition table with the probabilities
tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tmat <- populate_transition_matrix(2, tmat, list_prob = c("p1", "p2", "p3", "p4")) tmat_assigned <- eval_assign_trans_prob(tmat, c(p1 = 0.2, p2 = 0.3, p3 = 0.4, p4 = 0.5))tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tmat <- populate_transition_matrix(2, tmat, list_prob = c("p1", "p2", "p3", "p4")) tmat_assigned <- eval_assign_trans_prob(tmat, c(p1 = 0.2, p2 = 0.3, p3 = 0.4, p4 = 0.5))
Attribute values in health states
eval_assign_values_states(health_states, assigned_param)eval_assign_values_states(health_states, assigned_param)
health_states |
list of health_state objects |
assigned_param |
name value pairs of parameter values in the probability matrix expected created using function assign_parameters() |
Assigning the param is done for the cost and utility if the param is not numeric, check if it can be evaluated at the run time if yes, assign the evaluated numeric value if not get the parameters between operators, and assign the values to each individual parameters and then evaluate. only works for two levels. For the example shown the cost is sum of cost_A and cost_B which will only get added in the call eval_assign_values_states While initialising the state "well" it will be only saved as expression(cost_A + cost_B) assigned_param (a list) can be expected to be created using assign_parameters() the exception is if parameter is directly assigned with no nested calculation and no missing parameters. For example assigned_param = c(cost_a = 10, cost_b=10) will be ok but not assigned_param = c(a=10, cost_A = "a+100", cost_B =10) as it requires a nested calculation then use define_parameters() with assign_parameters() as in param_list <- define_parameters(a = 10, cost_A = "a + 100", cost_B = 10) assign_list <- assign_parameters(param_list)
health states with assigned values
well <- health_state("well", cost = "cost_A + cost_B", utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1)) health_states <- combine_state(well, disabled, dead) eval_assign_values_states(health_states, c(cost_A = 100, cost_B = 11))well <- health_state("well", cost = "cost_A + cost_B", utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1)) health_states <- combine_state(well, disabled, dead) eval_assign_values_states(health_states, c(cost_A = 100, cost_B = 11))
Function to find the keyword for family of distribution in glm
find_glm_distribution(text)find_glm_distribution(text)
text |
distribution |
Find the family for glm method
the keyword - the name of distribution
find_glm_distribution("gamma")find_glm_distribution("gamma")
Function to find the keyword for generating random numbers the distribution
find_keyword_rand_generation(text)find_keyword_rand_generation(text)
text |
name of the probability distribution |
This function returns the keyword for generating random number using a keyword provided that is generally used for prob distribution (but R might require a different keyword)
the keyword that should be used in R for generating random numbers
find_keyword_rand_generation("gamma")find_keyword_rand_generation("gamma")
Function to find the keyword for regression methods
find_keyword_regression_method(text, additional_info = NA)find_keyword_regression_method(text, additional_info = NA)
text |
regression method |
additional_info |
additional information required |
This function returns the keyword to use in regression methods. For example linear regression requires lm in R some regression methods require additional info and it has to be provided
the keyword that should be used in R for regression analysis
find_keyword_regression_method("linear")find_keyword_regression_method("linear")
Function to return parameters with in a expression containing operators
find_parameters_btn_operators(expr)find_parameters_btn_operators(expr)
expr |
an expression |
This function returns the parameters between the operators if the state value or probabilities are defined as expressions, we need to extract the parameters and then assign First the position of all operators are found and then return the parameters separated by those operators This happens only for one level find_parameters_btn_operators("a+b") provides a and b but for find_parameters_btn_operators("mean(a,b)+b") provides mean(a,b) and b
parameters in the expression expr
find_parameters_btn_operators("a+b")find_parameters_btn_operators("a+b")
Function to find the parameters that determine the probability distribution
find_required_parameter_combs(name_distri)find_required_parameter_combs(name_distri)
name_distri |
name of the probability distribution |
For each of the probability distribution we require certain parameters and this function provides that required list of parameters.
the parameters that determine the distribution
find_required_parameter_combs("gamma")find_required_parameter_combs("gamma")
Function to get sum of multiple columns of observations
find_rowwise_sum_multiplecol(the_data, colnames, sumcolname)find_rowwise_sum_multiplecol(the_data, colnames, sumcolname)
the_data |
the data where the observations are held |
colnames |
columnname in the data whose sum to be obtianed |
sumcolname |
name of the new column where sum to be saved |
the data with added sum
eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 7, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"), time = c(1,1,2,2), id = c(1, 1, 1, 2))) find_rowwise_sum_multiplecol(eg_data, c("mark_at_1","mark_at_2"), "totalmarks")eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 7, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"), time = c(1,1,2,2), id = c(1, 1, 1, 2))) find_rowwise_sum_multiplecol(eg_data, c("mark_at_1","mark_at_2"), "totalmarks")
Function to find the keyword for survreg distribution
find_survreg_distribution(text)find_survreg_distribution(text)
text |
distribution |
For surveg method, find the distribution
the keyword - the name of distribution
find_survreg_distribution("weibull")find_survreg_distribution("weibull")
Form expression to use with glm()
form_expression_glm( param_to_be_estimated, indep_var, family, covariates, interaction, naaction, link )form_expression_glm( param_to_be_estimated, indep_var, family, covariates, interaction, naaction, link )
param_to_be_estimated |
parameter of interest |
indep_var |
the independent variable (column name in data file) |
family |
distribution name eg. for logistic regression -binomial |
covariates |
list of covariates |
interaction |
boolean value to indicate interaction in the case of generalised linear models, false by default |
naaction |
action to be taken with the missing values |
link |
link function if not the default for each family |
Form expression for the method glm
the formula for glm
formula <- form_expression_glm("admit", indep_var = "gre", family = "binomial", covariates = c("gpa", "rank"), interaction = FALSE, naaction = "na.omit", link = NA)formula <- form_expression_glm("admit", indep_var = "gre", family = "binomial", covariates = c("gpa", "rank"), interaction = FALSE, naaction = "na.omit", link = NA)
Form expression to use with lm()
form_expression_lm(param_to_be_estimated, indep_var, covariates, interaction)form_expression_lm(param_to_be_estimated, indep_var, covariates, interaction)
param_to_be_estimated |
parameter of interest |
indep_var |
the independent variable (column name in data file) |
covariates |
list of covariates |
interaction |
boolean value to indicate interaction in the case of linear regression, false by default |
This function helps to create the expression for liner regression model it takes care of covariates and interaction
the formula for lm
formula <- form_expression_lm("gre", indep_var = "gpa", covariates = NA, interaction = FALSE)formula <- form_expression_lm("gre", indep_var = "gpa", covariates = NA, interaction = FALSE)
Form expression to use with mixed models
form_expression_mixed_model_lme4( param_to_be_estimated, dataset, fix_eff, fix_eff_interact_vars, random_intercept_vars, nested_intercept_vars_pairs, cross_intercept_vars_pairs, uncorrel_slope_intercept_pairs, random_slope_intercept_pairs, family, link )form_expression_mixed_model_lme4( param_to_be_estimated, dataset, fix_eff, fix_eff_interact_vars, random_intercept_vars, nested_intercept_vars_pairs, cross_intercept_vars_pairs, uncorrel_slope_intercept_pairs, random_slope_intercept_pairs, family, link )
param_to_be_estimated |
column name of dependent variable |
dataset |
a dataframe |
fix_eff |
names of variables as fixed effect predictors |
fix_eff_interact_vars |
if interaction -true |
random_intercept_vars |
names of variables for random intercept |
nested_intercept_vars_pairs |
those of the random intercept variables with nested effect |
cross_intercept_vars_pairs |
those of the random intercept variables with crossed effect |
uncorrel_slope_intercept_pairs |
variables with correlated intercepts |
random_slope_intercept_pairs |
random slopes intercept pairs - this is a list of paired variables |
family |
family of distribution for non gaussian distribution of predicted variable |
link |
link function for the variance |
Form the expression for mixed model
result regression result with plot if success and -1, if failure
datafile <- system.file("extdata", "data_linear_mixed_model.csv", package = "packDAMipd") dt = utils::read.csv(datafile, header = TRUE) formula <- form_expression_mixed_model_lme4("extro", dataset = dt, fix_eff = c("open", "agree", "social"), fix_eff_interact_vars = NULL, random_intercept_vars = c("school", "class"), nested_intercept_vars_pairs = list(c("school", "class")), cross_intercept_vars_pairs = NULL, uncorrel_slope_intercept_pairs = NULL, random_slope_intercept_pairs = NULL, family = "binomial", link = NA )datafile <- system.file("extdata", "data_linear_mixed_model.csv", package = "packDAMipd") dt = utils::read.csv(datafile, header = TRUE) formula <- form_expression_mixed_model_lme4("extro", dataset = dt, fix_eff = c("open", "agree", "social"), fix_eff_interact_vars = NULL, random_intercept_vars = c("school", "class"), nested_intercept_vars_pairs = list(c("school", "class")), cross_intercept_vars_pairs = NULL, uncorrel_slope_intercept_pairs = NULL, random_slope_intercept_pairs = NULL, family = "binomial", link = NA )
Function to get the weight and time units
generate_wt_time_units()generate_wt_time_units()
weight and time units
ans <- generate_wt_time_units()ans <- generate_wt_time_units()
Function to get the weight and volume units
generate_wt_vol_units()generate_wt_vol_units()
weight and vol units
ans <- generate_wt_vol_units()ans <- generate_wt_vol_units()
Function to get the details of the age column
get_age_details(trialdata)get_age_details(trialdata)
trialdata |
data containing individual level trial data |
expecting the data contains the information on age preferably column names "age", "dob" or "yob" or "date of birth". "year of birth", "birth year" If multiple column names match these, then first match will be chosen.
the name of the variable related to age and the unique contents if success, else error
get_age_details(data.frame("Age" = c(21, 15), "arm" = c("control", "intervention")))get_age_details(data.frame("Age" = c(21, 15), "arm" = c("control", "intervention")))
Function to keep the column name, coded values and non response code into a dataframe
get_colnames_codedvalues(variable, name, code, nrcode = NA)get_colnames_codedvalues(variable, name, code, nrcode = NA)
variable |
name of the variable in the column |
name |
column name |
code |
coded values |
nrcode |
code for non response |
data frame with all the above information
get_colnames_codedvalues("arm", "pat_trial_arm", c("Y", "N"))get_colnames_codedvalues("arm", "pat_trial_arm", c("Y", "N"))
Function to extract the unit hospital inpatient admission by matching code
get_cost_AandE_code( code, type_admit, ref_cost_data_file, col_name_code, unit_cost_col, type_admit_col, sheet = NULL )get_cost_AandE_code( code, type_admit, ref_cost_data_file, col_name_code, unit_cost_col, type_admit_col, sheet = NULL )
code |
code for AE attendance |
type_admit |
term indicating admission and type of attendance |
ref_cost_data_file |
file that has unit cost |
col_name_code |
name of the column that has the code |
unit_cost_col |
name of the column with the unit cost |
type_admit_col |
colname that describes type of the attendance and that indicates admitted or not |
sheet |
sheet if excel file is given |
unit cost the unit cost matching the code
ref_cost_data_file <- system.file("extdata", "National_schedule_of_NHS_costs_2019_AandE.csv", package = "packDAMipd") re = get_cost_AandE_code("VB02Z", "T01A", ref_cost_data_file, "Currency_Code","National_Average_Unit_Cost", "Service_Code")ref_cost_data_file <- system.file("extdata", "National_schedule_of_NHS_costs_2019_AandE.csv", package = "packDAMipd") re = get_cost_AandE_code("VB02Z", "T01A", ref_cost_data_file, "Currency_Code","National_Average_Unit_Cost", "Service_Code")
Function to extract the unit cost by descirption of AandE att matching description
get_cost_AandE_description( description, type_admit, ref_cost_data_file, col_name_description, unit_cost_col, type_admit_col, sheet = NULL )get_cost_AandE_description( description, type_admit, ref_cost_data_file, col_name_description, unit_cost_col, type_admit_col, sheet = NULL )
description |
description of the AE attendance |
type_admit |
term indicating admission and type of attendance |
ref_cost_data_file |
file that has unit cost |
col_name_description |
name of the column that has the description |
unit_cost_col |
name of the column with the unit cost |
type_admit_col |
colname that descirbes type of the attendance and |
sheet |
sheet if excel file is given |
unit cost the unit cost matching the hrg code
ref_cost_data_file <- system.file("extdata", "National_schedule_of_NHS_costs_2019_AandE.csv", package = "packDAMipd") re = get_cost_AandE_description("Emergency Medicine", "T01A", ref_cost_data_file, "Currency_Description", "National_Average_Unit_Cost", "Service_Code")ref_cost_data_file <- system.file("extdata", "National_schedule_of_NHS_costs_2019_AandE.csv", package = "packDAMipd") re = get_cost_AandE_description("Emergency Medicine", "T01A", ref_cost_data_file, "Currency_Description", "National_Average_Unit_Cost", "Service_Code")
Function to extract the unit hospital inpatient admission by matching description
get_cost_ip_dc_description( description, ref_cost_data_file, col_name_description, unit_cost_col, sheet = NULL )get_cost_ip_dc_description( description, ref_cost_data_file, col_name_description, unit_cost_col, sheet = NULL )
description |
description corresponding to the inpatient admission |
ref_cost_data_file |
file that has unit cost |
col_name_description |
name of the column that has the description |
unit_cost_col |
name of the column with the unit cost |
sheet |
sheet if excel file is given |
unit cost the unit cost matching the hrg code
ref_cost_data_file <- system.file("extdata", "National_schedule_of_NHS_costs_2019.csv", package = "packDAMipd") result <- get_cost_ip_dc_description("Cerebrovascular Accident", ref_cost_data_file, "Currency_Description", "National_Average_Unit_Cost")ref_cost_data_file <- system.file("extdata", "National_schedule_of_NHS_costs_2019.csv", package = "packDAMipd") result <- get_cost_ip_dc_description("Cerebrovascular Accident", ref_cost_data_file, "Currency_Description", "National_Average_Unit_Cost")
Function to extract the unit hospital inpatient admission by matching HRG code
get_cost_ip_dc_hrg( hrg, ref_cost_data_file, col_name_hrg_code, unit_cost_col, sheet = NULL )get_cost_ip_dc_hrg( hrg, ref_cost_data_file, col_name_hrg_code, unit_cost_col, sheet = NULL )
hrg |
hrg code corresponding to the inpatient admission |
ref_cost_data_file |
file that has unit cost |
col_name_hrg_code |
name of the column that has the hrg code |
unit_cost_col |
name of the colum with the unit cost |
sheet |
sheet if excel file is given |
unit cost the unit cost matching the hrg code
ref_cost_data_file <- system.file("extdata", "National_schedule_of_NHS_costs_2019.csv", package = "packDAMipd") get_cost_ip_dc_hrg("AA22C", ref_cost_data_file, "Currency_Code", "National_Average_Unit_Cost")ref_cost_data_file <- system.file("extdata", "National_schedule_of_NHS_costs_2019.csv", package = "packDAMipd") get_cost_ip_dc_hrg("AA22C", ref_cost_data_file, "Currency_Code", "National_Average_Unit_Cost")
Convert the combined dose to its individual component numerical value or can be unit/unit
get_doses_combination(the_string, separator = "/")get_doses_combination(the_string, separator = "/")
the_string |
given combined unit |
separator |
given character for separation , default is "/" |
separated texts
get_doses_combination("g/ml")get_doses_combination("g/ml")
Convert the combined dose to its individual component numerical value and units
get_doses_combination_units(the_string, separator = "/")get_doses_combination_units(the_string, separator = "/")
the_string |
given combined unit |
separator |
given character for separation , default is "/" |
separated numerical value and its units
get_doses_combination_units("10g/2ml")get_doses_combination_units("10g/2ml")
Function to get the details of the EQ5D column
get_eq5d_details(trialdata)get_eq5d_details(trialdata)
trialdata |
data containing individual level trial data |
Specific to the EQ5D data - the column names are given as certain sets, Tried to give 15 sets as the column names
the name of the variable related to EQ5D and the unique contents if success, else error
get_eq5d_details(data.frame( "MO" = c(1, 2), "SC" = c(1, 2), "UA" = c(1, 2), "PD" = c(1, 2), "AD" = c(1, 2) ))get_eq5d_details(data.frame( "MO" = c(1, 2), "SC" = c(1, 2), "UA" = c(1, 2), "PD" = c(1, 2), "AD" = c(1, 2) ))
Function to get extension of a file name
get_extension_file(filename)get_extension_file(filename)
filename |
name of a file |
if there is no "." character returns error else returns last characters those after string split using "."
the extension
get_extension_file("data.txt")get_extension_file("data.txt")
Function to get the details of the gender column
get_gender_details(trialdata)get_gender_details(trialdata)
trialdata |
data containing individual level trial data |
expecting the data contains the information on gender preferably column names "gender", "sex" or "male" or "female". If multiple column names match these, then first match will be chosen.
the name of the variable related to gender and the unique contents if success, else error
get_gender_details(data.frame("Age" = c(21, 15), "sex" = c("m", "f")))get_gender_details(data.frame("Age" = c(21, 15), "sex" = c("m", "f")))
Function to return mean age from a data frame
get_mean_sd_age(this_data, age_nrcode)get_mean_sd_age(this_data, age_nrcode)
this_data |
the data containing column with age |
age_nrcode |
non response code |
Age data is complete with the nr code given and get the mean and sd
mean and sd, if success -1, if failure
this_data <- as.data.frame(cbind(num = c(1, 2, 3, 4), age = c(14, 25, 26, 30))) get_mean_sd_age(this_data, NA)this_data <- as.data.frame(cbind(num = c(1, 2, 3, 4), age = c(14, 25, 26, 30))) get_mean_sd_age(this_data, NA)
Get the mortality rate values from reading a file
get_mortality_from_file(paramfile, age, mortality_colname, gender = NULL)get_mortality_from_file(paramfile, age, mortality_colname, gender = NULL)
paramfile |
parameter file to get the mortality eg.national life table data |
age |
age to get the age specific data |
mortality_colname |
column name with the mortality rates if it is not gender specific |
gender |
gender details to get the gender specific mortality data |
Provides the mortality rates as age and gender dependent Assumes the data contains mortality rate for single year and once it extracted per gender will retrieve single value Age column can consists of range of values, or a particular value also assumes that the mortality rate for each gender is listed under the gender column for gender specific values. if the mortality is not gender specific, the column name should be passed on to the function if gender is not null, mortality_name will be ignored
the paramvalue
paramfile <- system.file("extdata", "LifeTable_USA_Mx_2015.csv", package = "packDAMipd" ) a <- get_mortality_from_file(paramfile, age = 10, mortality_colname = "total", gender = NULL)paramfile <- system.file("extdata", "LifeTable_USA_Mx_2015.csv", package = "packDAMipd" ) a <- get_mortality_from_file(paramfile, age = 10, mortality_colname = "total", gender = NULL)
Function to return the two parameters from a given expression separated by comma,
get_name_value_probdistrb_def(expr)get_name_value_probdistrb_def(expr)
expr |
an expression |
It will return the parameters of the distribution separated by commas and given in usual notation as brackets. It will identify those in between first occurrence of "( "and last occurrence of ")" and from the characters in between search for comma to indicate different parameters then it will extract (from those extracted parameters separated by commas) that on the left side of "equal" sign get_name_value_probdistrb_def("gamma(mean = sqrt(2), b =17)") will be ok but get_name_value_probdistrb_def("gamma(shape, scale")) and get_name_value_probdistrb_def("gamma(shape =1 & scale =1")) will show error
parameters in the expression expr
get_name_value_probdistrb_def("gamma(mean = 10, sd =1)")get_name_value_probdistrb_def("gamma(mean = 10, sd =1)")
Function to get the details of the outcome column
get_outcome_details(trialdata, name, related_words, multiple = FALSE)get_outcome_details(trialdata, name, related_words, multiple = FALSE)
trialdata |
data containing individual level trial data |
name |
name of the variable |
related_words |
probable column names |
multiple |
indicates true if there are multiple columns |
if the words related to outcome is given, the function will get the columns and the codes used for the outcome, the difference here is that certain outcomes can be distributed in multiple columns
the name of the variable related to health outcome (any) and the unique contents if success, else error
get_outcome_details( data.frame("qol.MO" = c(1, 2), "qol.PD" = c(1, 2), "qol.AD" = c(1, 2)), "eq5d", "qol", TRUE )get_outcome_details( data.frame("qol.MO" = c(1, 2), "qol.PD" = c(1, 2), "qol.AD" = c(1, 2)), "eq5d", "qol", TRUE )
Get the definition of given parameter distribution defined in a file
get_parameter_def_distribution( parameter, paramfile, colnames_paramdistr, strategycol = NA, strategyname = NA )get_parameter_def_distribution( parameter, paramfile, colnames_paramdistr, strategycol = NA, strategyname = NA )
parameter |
parameter of interest |
paramfile |
data file to be provided |
colnames_paramdistr |
list of column names for the parameters that define the distribution |
strategycol |
treatment strategy column name |
strategyname |
treatment strategy name in the column strategycol |
This function reads the parameter distribution from a file and return the parameter obtained This assumes that the file contains parameter, distribution colnames for parameter values for the distribution are passed on to the function assumes the name of each parameter and value are given in the consecutive columns. Once the expression is created using the parameters given in the file, it gets checked for correctness of specifying the distribution in R context using the function check_estimate_substitute_proper_params and then evaluated.
the definition of parameter from the given distribution
paramfile <- system.file("extdata", "table_param.csv", package = "packDAMipd") a <- get_parameter_def_distribution("rr", paramfile, c("Param1_name", "Param1_value"))paramfile <- system.file("extdata", "table_param.csv", package = "packDAMipd") a <- get_parameter_def_distribution("rr", paramfile, c("Param1_name", "Param1_value"))
Get the parameter values from reading a file
get_parameter_direct(parameter, paramvalue)get_parameter_direct(parameter, paramvalue)
parameter |
parameter of interest |
paramvalue |
parameter value to be assigned |
Basic function to assign a parameter directly
the paramvalue
a <- get_parameter_direct("cost_IT", paramvalue = 100)a <- get_parameter_direct("cost_IT", paramvalue = 100)
Get the parameter values using the provided statistical regression methods
get_parameter_estimated_regression( param_to_be_estimated, data, method, indep_var, info_get_method = NA, info_distribution = NA, covariates = NA, timevar_survival = NA, interaction = FALSE, fix_eff = NA, fix_eff_interact_vars = NA, random_intercept_vars = NA, nested_intercept_vars_pairs = NA, cross_intercept_vars_pairs = NA, uncorrel_slope_intercept_pairs = NA, random_slope_intercept_pairs = NA, naaction = "stats::na.omit", param2_to_be_estimated = NA, covariates2 = NA, interaction2 = FALSE, link = NA, cluster_var = NA, package_mixed_model = NA )get_parameter_estimated_regression( param_to_be_estimated, data, method, indep_var, info_get_method = NA, info_distribution = NA, covariates = NA, timevar_survival = NA, interaction = FALSE, fix_eff = NA, fix_eff_interact_vars = NA, random_intercept_vars = NA, nested_intercept_vars_pairs = NA, cross_intercept_vars_pairs = NA, uncorrel_slope_intercept_pairs = NA, random_slope_intercept_pairs = NA, naaction = "stats::na.omit", param2_to_be_estimated = NA, covariates2 = NA, interaction2 = FALSE, link = NA, cluster_var = NA, package_mixed_model = NA )
param_to_be_estimated |
parameter of interest |
data |
data to be provided or the data file containing dataset |
method |
method of estimation (for example, linear, logistic regression etc) |
indep_var |
the independent variable (column name in data file) |
info_get_method |
additional information on methods e.g Kaplan-Meier ot hazard |
info_distribution |
distribution name eg. for logistic regression -binomial |
covariates |
list of covariates-calculations to be done before passing |
timevar_survival |
time variable for survival analysis |
interaction |
boolean value to indicate interaction in the case of linear regression |
fix_eff |
boolean value to indicate interaction in the case of linear regression |
fix_eff_interact_vars |
boolean value to indicate interaction in the case of linear regression |
random_intercept_vars |
boolean value to indicate interaction in the case of linear regression |
nested_intercept_vars_pairs |
boolean value to indicate interaction in the case of linear regression |
cross_intercept_vars_pairs |
boolean value to indicate interaction in the case of linear regression |
uncorrel_slope_intercept_pairs |
boolean value to indicate interaction in the case of linear regression |
random_slope_intercept_pairs |
boolean value to indicate interaction in the case of linear regression |
naaction |
what action to be taken for the missing values, default is a missing value. |
param2_to_be_estimated |
parameter of interest for equation 2 in bivariate regression |
covariates2 |
list of covariates - for equation 2 in bivariate regression |
interaction2 |
boolean value to indicate interaction for equation 2 in bivariate regression |
link |
link function to be provided if not using the default link for each of the info_distribution |
cluster_var |
cluster variable if any |
package_mixed_model |
package to be used for mixed model ie nlme or lme4 |
This function is the top in the layer of functions used for regression analysis Thus it contains many parameters to be passed on The required ones are parameter to be estimated, data that contains the observation, the method of regression to be used, the independent variable and the information for the distribution and method. if the data is given as a file name. it will load the data in that file Then it calls the appropriate functions depending on the regression method that specified. The methods that are considered : Survival analysis, linear regression, logistic regression,generalised linear model, linear multilevel or mixed model, and seemingly unrelated regression
results the results of the regression analysis
result <- get_parameter_estimated_regression( param_to_be_estimated = "Direction", data = ISLR::Smarket, method = "logistic", indep_var = "Lag1", info_get_method = NA, info_distribution = "binomial", covariates = c("Lag2", "Lag3"), interaction = FALSE, naaction = "na.omit", link = NA)result <- get_parameter_estimated_regression( param_to_be_estimated = "Direction", data = ISLR::Smarket, method = "logistic", indep_var = "Lag1", info_get_method = NA, info_distribution = "binomial", covariates = c("Lag2", "Lag3"), interaction = FALSE, naaction = "na.omit", link = NA)
Get the parameter values from reading a file
get_parameter_read(parameter, paramfile, strategycol = NA, strategyname = NA)get_parameter_read(parameter, paramfile, strategycol = NA, strategyname = NA)
parameter |
parameter of interest |
paramfile |
parameter file to be provided |
strategycol |
treatment strategy |
strategyname |
treatment strategy name in the column strategycol |
This function read the parameter from a file given that the file has these column names (at least) Parameter and Value Strategy col and name are optional. Check if the data file contains column names parameter and value and then get the results.
the paramvalue
a <- get_parameter_read("cost_IT", paramfile = system.file("extdata", "table_param.csv", package = "packDAMipd"))a <- get_parameter_read("cost_IT", paramfile = system.file("extdata", "table_param.csv", package = "packDAMipd"))
Function to get cols for the pattern given
get_single_col_multiple_pattern(pattern, the_data)get_single_col_multiple_pattern(pattern, the_data)
pattern |
the pattern to look for |
the_data |
data where to look at |
zero or -1
the_data <- as.data.frame(cbind(c("one", "two"), c("a", "b"), c("aa", "bb"))) colnames(the_data) <- c("name", "brand_one", "two") get_single_col_multiple_pattern(c("brand", "trade"), the_data)the_data <- as.data.frame(cbind(c("one", "two"), c("a", "b"), c("aa", "bb"))) colnames(the_data) <- c("name", "brand_one", "two") get_single_col_multiple_pattern(c("brand", "trade"), the_data)
help function to keep slope and intercept portion ready in mixed model expression
get_slope_intercept( expression, random_intercept_vars, random_slope_intercept_pairs, uncorrel_slope_intercept_pairs )get_slope_intercept( expression, random_intercept_vars, random_slope_intercept_pairs, uncorrel_slope_intercept_pairs )
expression |
expression created so far |
random_intercept_vars |
names of variables for random intercept |
random_slope_intercept_pairs |
random slopes intercept pairs this is a list of paired variables |
uncorrel_slope_intercept_pairs |
variables with correlated intercepts |
expression expression created
help function to keep slope and intercept portion ready in mixed model expression
get_slope_intercept_cross( expression, random_intercept_vars, intercept_vars_pairs, random_slope_intercept_pairs, uncorrel_slope_intercept_pairs )get_slope_intercept_cross( expression, random_intercept_vars, intercept_vars_pairs, random_slope_intercept_pairs, uncorrel_slope_intercept_pairs )
expression |
expression created so far |
random_intercept_vars |
names of variables for random intercept |
intercept_vars_pairs |
those of the random intercept variables with nested effect |
random_slope_intercept_pairs |
random slopes intercept pairs this is a list of paired variables |
uncorrel_slope_intercept_pairs |
variables with correlated intercepts |
expression expression created
help function to keep slope and intercept portion ready in mixed model expression
get_slope_intercept_nested( expression, random_intercept_vars, intercept_vars_pairs, random_slope_intercept_pairs, uncorrel_slope_intercept_pairs )get_slope_intercept_nested( expression, random_intercept_vars, intercept_vars_pairs, random_slope_intercept_pairs, uncorrel_slope_intercept_pairs )
expression |
expression created so far |
random_intercept_vars |
names of variables for random intercept |
intercept_vars_pairs |
those of the random intercept variables with nested effect |
random_slope_intercept_pairs |
random slopes intercept pairs this is a list of paired variables |
uncorrel_slope_intercept_pairs |
variables with correlated intercepts |
expression expression created
Function to get the details of the time point column
get_timepoint_details(trialdata)get_timepoint_details(trialdata)
trialdata |
data containing individual level trial data |
expecting the data contains the information on timepoints preferably column names "time point", "times" or "time" or "timepoint". If multiple column names match these, then first match will be chosen.
the name of the variable related to time point and the unique contents if success, else error
get_timepoint_details(data.frame("time" = c(21, 15), "arm" = c("control", "intervention")))get_timepoint_details(data.frame("time" = c(21, 15), "arm" = c("control", "intervention")))
Function to get the details of the trial arm
get_trial_arm_details(trialdata)get_trial_arm_details(trialdata)
trialdata |
data containing individual level trial data |
expecting the data contains the information on trial arm preferably column names "arm", "trial" or "trial arm". If multiple column names match these, then first match will be chosen.
the name of the variable related to trial arm and the unique contents if success, else error
get_trial_arm_details(data.frame( "Age" = c(21, 15), "arm" = c("control", "intervention") ))get_trial_arm_details(data.frame( "Age" = c(21, 15), "arm" = c("control", "intervention") ))
Get the attribute for the health state
get_var_state(state, var)get_var_state(state, var)
state |
object of class health state |
var |
attribute of the health state |
After checking the given state is a health state and given variable is defined in the health state, the value of the variable is returned
modified health state
get_var_state(health_state("IT", 100, 0.4, 0, FALSE), "cost")get_var_state(health_state("IT", 100, 0.4, 0, FALSE), "cost")
Definition of health state class or health state constructor
health_state(name, cost, utility, state_time = 0, absorb = FALSE)health_state(name, cost, utility, state_time = 0, absorb = FALSE)
name |
name of the health state |
cost |
value or expression that represents cost of the health state |
utility |
value or expression that represents utility of the health state |
state_time |
time denoting how long in the state |
absorb |
boolean indicating health state absorbing or not |
Initialising the name, cost, utility and time spent for the health state name is the name of the health state cost/utility can be defined as characters e.g. "cost_A" if they are characters, the value is assigned after parsing the text. state_time is integer and absorb is boolean
value of the state
st <- health_state("IT", 100, 0.4, 0, FALSE) st <- health_state("IT", "cost_A", 0.4, 0, FALSE)st <- health_state("IT", 100, 0.4, 0, FALSE) st <- health_state("IT", "cost_A", 0.4, 0, FALSE)
Define an all zero trace matrix
init_trace(health_states, cycles)init_trace(health_states, cycles)
health_states |
health states |
cycles |
no of cycles |
Initialise the trace matrix with all zeros trace matrix will be with no_cycles+1 by no_states matrix
trace matrix -all zero
a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0.5, 0, FALSE) health_states <- combine_state(a, b) init_trace(health_states, 10)a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0.5, 0, FALSE) health_states <- combine_state(a, b) init_trace(health_states, 10)
Define an all zero trace matrix
init_trace_sjtime(health_states, cycles)init_trace_sjtime(health_states, cycles)
health_states |
health states |
cycles |
no of cycles |
Initialise the trace matrix with all zeros trace matrix will be with no_cycles+1 by no_states matrix
trace matrix -all zero
a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0.5, 0, FALSE) health_states <- combine_state(a, b) init_trace(health_states, 10)a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0.5, 0, FALSE) health_states <- combine_state(a, b) init_trace(health_states, 10)
Function to do some checks before plotting sensitivity analysis results
keep_results_plot_dsa( result_dsa_control, plotfor, result_dsa_treat, plot_variable, threshold, comparator )keep_results_plot_dsa( result_dsa_control, plotfor, result_dsa_treat, plot_variable, threshold, comparator )
result_dsa_control |
result from deterministic sensitivity analysis for first or control model |
plotfor |
the variable to plotfor e.g. cost, utility NMB etc |
result_dsa_treat |
result from deterministic sensitivity analysis for the comparative Markov model |
plot_variable |
variable for plotting |
threshold |
threshold value of WTP |
comparator |
the strategy to be compared with |
results to plot dsa
Function to list probabilistic sensitivity analysis results parameterwise
list_paramwise_psa_result( result_psa_params_control, result_psa_params_treat, threshold, comparator )list_paramwise_psa_result( result_psa_params_control, result_psa_params_treat, threshold, comparator )
result_psa_params_control |
result from probabilistic sensitivity analysis for first or control model |
result_psa_params_treat |
result from probabilistic sensitivity analysis for the comparative Markov model |
threshold |
threshold value of WTP |
comparator |
the strategy to be compared with |
plot of sensitivity analysis
param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido" ) A <- health_state("A", cost = "cost_health_A + cost_drug ", utility = 1) B <- health_state("B", cost = "cost_health_B + cost_drug", utility = 1) C <- health_state("C", cost = "cost_health_C + cost_drug", utility = 1) D <- health_state("D", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3, 4), c(NA, 5, 6, 7), c(NA, NA, 8, 9), c(NA, NA, NA, 10)) colnames(tmat) <- rownames(tmat) <- c("A", "B", "C", "D") tm <- populate_transition_matrix(4, tmat, c( "tpAtoA", "tpAtoB", "tpAtoC", "tpAtoD", "tpBtoB", "tpBtoC", "tpBtoD", "tpCtoC", "tpCtoD", "tpDtoD" ), colnames(tmat)) health_states <- combine_state(A, B, C, D) mono_strategy <- strategy(tm, health_states, "mono") mono_markov <- markov_model(mono_strategy, 20, initial_state = c(1,0,0,0), discount = c(0.06, 0),param_list) sample_list <- define_parameters(cost_zido = "gamma(mean = 2756, sd = sqrt(2756))") param_table <- define_parameters_psa(param_list, sample_list) result <- do_psa(mono_markov, param_table, 10) list_paramwise_psa_result(result, NULL, NULL, NULL)param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido" ) A <- health_state("A", cost = "cost_health_A + cost_drug ", utility = 1) B <- health_state("B", cost = "cost_health_B + cost_drug", utility = 1) C <- health_state("C", cost = "cost_health_C + cost_drug", utility = 1) D <- health_state("D", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3, 4), c(NA, 5, 6, 7), c(NA, NA, 8, 9), c(NA, NA, NA, 10)) colnames(tmat) <- rownames(tmat) <- c("A", "B", "C", "D") tm <- populate_transition_matrix(4, tmat, c( "tpAtoA", "tpAtoB", "tpAtoC", "tpAtoD", "tpBtoB", "tpBtoC", "tpBtoD", "tpCtoC", "tpCtoD", "tpDtoD" ), colnames(tmat)) health_states <- combine_state(A, B, C, D) mono_strategy <- strategy(tm, health_states, "mono") mono_markov <- markov_model(mono_strategy, 20, initial_state = c(1,0,0,0), discount = c(0.06, 0),param_list) sample_list <- define_parameters(cost_zido = "gamma(mean = 2756, sd = sqrt(2756))") param_table <- define_parameters_psa(param_list, sample_list) result <- do_psa(mono_markov, param_table, 10) list_paramwise_psa_result(result, NULL, NULL, NULL)
Function to load the file containing trial data and return it
load_trial_data(file = NULL, sheet = NULL)load_trial_data(file = NULL, sheet = NULL)
file |
name of the file in full |
sheet |
name of the sheet if excel work book is given |
trial data if success, else -1
load_trial_data(system.file("extdata", "trial_data.csv", package = "packDAMipd" ))load_trial_data(system.file("extdata", "trial_data.csv", package = "packDAMipd" ))
Function to map EQ5D5L scores to EQ5D3L scores and then add to IPD data
map_eq5d5Lto3L_VanHout(ind_part_data, eq5d_nrcode)map_eq5d5Lto3L_VanHout(ind_part_data, eq5d_nrcode)
ind_part_data |
a data frame |
eq5d_nrcode |
non response code for EQ5D5L, default is NA |
qaly included modified data, if success -1, if failure
http://eprints.whiterose.ac.uk/121473/1/Devlin_et_al-2017-Health_Economics.pdf
library(valueEQ5D) datafile <- system.file("extdata", "trial_data.csv", package = "packDAMipd") trial_data <- load_trial_data(datafile) map_eq5d5Lto3L_VanHout(trial_data, NA)library(valueEQ5D) datafile <- system.file("extdata", "trial_data.csv", package = "packDAMipd") trial_data <- load_trial_data(datafile) map_eq5d5Lto3L_VanHout(trial_data, NA)
Function to to read the text form of resource use and replace it with standard texts of resoure use ie. some one can describe GP visit as GP surgery visit, surgery visit or general practioners visit etc. Here all these texts should be given in a excel or csv file and then corresponidng standard form will be read from the file and will be replaced.
map_resource_use_categories( the_data, service_actual, new_column, mapped_data, mapped_use, analysis, replace_only, relevant_column = NULL, check_value_relevant = NULL, nhs_use_column = NULL, check_value_nhs_use = NULL )map_resource_use_categories( the_data, service_actual, new_column, mapped_data, mapped_use, analysis, replace_only, relevant_column = NULL, check_value_relevant = NULL, nhs_use_column = NULL, check_value_nhs_use = NULL )
the_data |
the data where the observations are held |
service_actual |
columna name of the actual service use |
new_column |
the name of the column where the mapped resource use to be |
mapped_data |
data where the service name and mapped service name has been stored |
mapped_use |
columan name of mapped resource use in mapped_data |
analysis |
base case or secondary |
replace_only |
if we want to replace only certain resource use |
relevant_column |
the name of the column where the mapped resource use is indicated as relevant or not |
check_value_relevant |
how is the mapped resource is indicated as relevant by a value |
nhs_use_column |
the name of the column where the mapped resource use comes under NHS or not |
check_value_nhs_use |
value that is used to indicated the nhs use |
the data with added sum
Definition of Markov model and trace
markov_model( current_strategy, cycles, initial_state, discount = c(0, 0), parameter_values = NULL, half_cycle_correction = TRUE, state_cost_only_prevalent = FALSE, state_util_only_prevalent = FALSE, method = "half cycle correction", startup_cost = NULL, startup_util = NULL )markov_model( current_strategy, cycles, initial_state, discount = c(0, 0), parameter_values = NULL, half_cycle_correction = TRUE, state_cost_only_prevalent = FALSE, state_util_only_prevalent = FALSE, method = "half cycle correction", startup_cost = NULL, startup_util = NULL )
current_strategy |
strategy object |
cycles |
no of cycles |
initial_state |
value of states initially |
discount |
rate of discount for costs and qalys |
parameter_values |
parameters for assigning health states and probabilities |
half_cycle_correction |
boolean to indicate half cycle correction |
state_cost_only_prevalent |
boolean parameter to indicate if the costs for state occupancy is only for those in the state excluding those that transitioned new. This is relevant when the transition cost is provided for eg. in a state with dialysis the cost of previous dialysis is different from the newly dialysis cases. Then the state_cost_only_prevalent should be TRUE |
state_util_only_prevalent |
boolean parameter to indicate if the utilities for state occupancy is only for those in the state excluding those that transitioned new. |
method |
what type of half cycle correction needed |
startup_cost |
cost of states initially |
startup_util |
utility of states initially if any |
Use the strategy, cycles, initial state values creating the markov model and trace. As many probabilities /cost/utility value depend on age/time the evaluation and assignment happens during each cycle. At the heart it does a matrix multiplication using the previous row of the trace matrix and the columns of the transition matrix. Also checks for population loss, calculates cumulative costs and qalys (accounts for discounting and half cycle correction)
Markov trace
tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- populate_transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0, 0, TRUE) health_states <- combine_state(a, b) this.strategy <- strategy(tm, health_states, "intervention") markov_model(this.strategy, 10, c(1, 0))tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- populate_transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0, 0, TRUE) health_states <- combine_state(a, b) this.strategy <- strategy(tm, health_states, "intervention") markov_model(this.strategy, 10, c(1, 0))
Definition of Markov model and trace
markov_model_sojourntime( current_strategy, nCycles, initial_state, discount = c(0, 0), parameter_values = NULL, half_cycle_correction = TRUE, state_cost_only_prevalent = FALSE, state_util_only_prevalent = FALSE, method = "half cycle correction", startup_cost = NULL, startup_util = NULL )markov_model_sojourntime( current_strategy, nCycles, initial_state, discount = c(0, 0), parameter_values = NULL, half_cycle_correction = TRUE, state_cost_only_prevalent = FALSE, state_util_only_prevalent = FALSE, method = "half cycle correction", startup_cost = NULL, startup_util = NULL )
current_strategy |
strategy object |
nCycles |
no of cycles |
initial_state |
value of states initially |
discount |
rate of discount for costs and qalys |
parameter_values |
parameters for assigning health states and probabilities |
half_cycle_correction |
boolean to indicate half cycle correction |
state_cost_only_prevalent |
boolean parameter to indicate if the costs for state occupancy is only for those in the state excluding those that transitioned new. This is relevant when the transition cost is provided for eg. in a state with dialysis the cost of previous dialysis is different from the newly dialysis cases. Then the state_cost_only_prevalent should be TRUE |
state_util_only_prevalent |
boolean parameter to indicate if the utilities for state occupancy is only for those in the state excluding those that transitioned new. |
method |
what type of half cycle correction needed |
startup_cost |
cost of states initially |
startup_util |
utility of states initially if any |
Use the strategy, cycles, initial state values creating the markov model and trace. As many probabilities /cost/utility value depend on age/time the evaluation and assignment happens during each cycle. At the heart it does a matrix multiplication using the previous row of the trace matrix and the columns of the transition matrix. Also checks for population loss, calculates cumulative costs and qalys (accounts for discounting and half cycle correction)
Markov trace
tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- populate_transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0, 0, TRUE) health_states <- combine_state(a, b) this.strategy <- strategy(tm, health_states, "intervention") markov_model(this.strategy, 10, c(1, 0))tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- populate_transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0, 0, TRUE) health_states <- combine_state(a, b) this.strategy <- strategy(tm, health_states, "intervention") markov_model(this.strategy, 10, c(1, 0))
Function to estimate the cost of liquids when IPD is in long format
microcosting_liquids_long( the_columns, ind_part_data_long, name_med, brand_med = NULL, dose_med, unit_med = NULL, bottle_size, bottle_size_unit = NULL, bottle_lasts, bottle_lasts_unit = NULL, preparation_dose = NULL, preparation_unit = NULL, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_dose_unit = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = NULL, basis_strength_unit = NULL )microcosting_liquids_long( the_columns, ind_part_data_long, name_med, brand_med = NULL, dose_med, unit_med = NULL, bottle_size, bottle_size_unit = NULL, bottle_lasts, bottle_lasts_unit = NULL, preparation_dose = NULL, preparation_unit = NULL, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_dose_unit = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = NULL, basis_strength_unit = NULL )
the_columns |
columns that are to be used to convert the data from long to wide |
ind_part_data_long |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
dose_med |
dose of medication used |
unit_med |
unit of medication ; use null if its along with the dose |
bottle_size |
size of the bottle used |
bottle_size_unit |
unit of bottle volume |
bottle_lasts |
how long the bottle lasted |
bottle_lasts_unit |
time unit of how long the bottle lasted |
preparation_dose |
dose if preparation is given |
preparation_unit |
unit of preparatio dose |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit where the cost is calculated |
strength_column |
column column name that has strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
list_of_code_dose_unit |
if unit is coded, give the code:unit pairs, optional |
list_of_code_bottle_size_unit |
list of bottle size units and codes |
list_of_code_bottle_lasts_unit |
list of time of bottle lasts and codes |
list_preparation_dose_unit |
list of preparation dose units and codes |
eqdose_covtab |
table to get the conversion factor for equivalent doses, optional |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_liq.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) names <- colnames(ind_part_data) ending <- length(names) ind_part_data_long <- tidyr::gather(ind_part_data, measurement, value, names[2]:names[ending], factor_key = TRUE) the_columns <- c("measurement", "value") res <- microcosting_liquids_long(the_columns, ind_part_data_long = ind_part_data_long, name_med = "liq_name", brand_med = NULL, dose_med = "liq_strength", unit_med = NULL, bottle_size = "liq_bottle_size",bottle_size_unit = NULL, bottle_lasts = "liq_lasts",bottle_lasts_unit = NULL,preparation_dose = NULL, preparation_unit = NULL,timeperiod = "4 months",unit_cost_data = med_costs, unit_cost_column = "UnitCost",cost_calculated_per = "Basis", strength_column = "Strength",list_of_code_names = NULL, list_of_code_brand = NULL,list_of_code_dose_unit = NULL, list_of_code_bottle_size_unit = NULL,list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL,eqdose_covtab = table, basis_strength_unit = NULL)med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_liq.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) names <- colnames(ind_part_data) ending <- length(names) ind_part_data_long <- tidyr::gather(ind_part_data, measurement, value, names[2]:names[ending], factor_key = TRUE) the_columns <- c("measurement", "value") res <- microcosting_liquids_long(the_columns, ind_part_data_long = ind_part_data_long, name_med = "liq_name", brand_med = NULL, dose_med = "liq_strength", unit_med = NULL, bottle_size = "liq_bottle_size",bottle_size_unit = NULL, bottle_lasts = "liq_lasts",bottle_lasts_unit = NULL,preparation_dose = NULL, preparation_unit = NULL,timeperiod = "4 months",unit_cost_data = med_costs, unit_cost_column = "UnitCost",cost_calculated_per = "Basis", strength_column = "Strength",list_of_code_names = NULL, list_of_code_brand = NULL,list_of_code_dose_unit = NULL, list_of_code_bottle_size_unit = NULL,list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL,eqdose_covtab = table, basis_strength_unit = NULL)
Function to estimate the cost of liquids taken (from IPD)
microcosting_liquids_wide( ind_part_data, name_med, brand_med = NULL, dose_med, unit_med = NULL, bottle_size, bottle_size_unit = NULL, bottle_lasts, bottle_lasts_unit = NULL, preparation_dose = NULL, preparation_unit = NULL, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_dose_unit = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = NULL, basis_strength_unit = NULL )microcosting_liquids_wide( ind_part_data, name_med, brand_med = NULL, dose_med, unit_med = NULL, bottle_size, bottle_size_unit = NULL, bottle_lasts, bottle_lasts_unit = NULL, preparation_dose = NULL, preparation_unit = NULL, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_dose_unit = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = NULL, basis_strength_unit = NULL )
ind_part_data |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
dose_med |
dose of medication used |
unit_med |
unit of medication ; use null if its along with the dose |
bottle_size |
size of the bottle used |
bottle_size_unit |
unit of bottle volume |
bottle_lasts |
how long the bottle lasted |
bottle_lasts_unit |
time unit of how long the bottle lasted |
preparation_dose |
dose if preparation is given |
preparation_unit |
unit of preparatio dose |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit where the cost is calculated |
strength_column |
column column name that has strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
list_of_code_dose_unit |
if unit is coded, give the code:unit pairs, optional |
list_of_code_bottle_size_unit |
list of bottle size units and codes |
list_of_code_bottle_lasts_unit |
list of time of bottle lasts and codes |
list_preparation_dose_unit |
list of preparation dose units and codes |
eqdose_covtab |
table to get the conversion factor for equivalent doses, optional, but the column names have to unique Similar to c("Drug", "form", "unit", "factor") or c("Drug", "form", "unit", "conversion") |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_liq.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx",package = "packDAMipd") table <- load_trial_data(conv_file) res <- microcosting_liquids_wide( ind_part_data = ind_part_data, name_med = "liq_name", brand_med = NULL, dose_med = "liq_strength", unit_med = NULL, bottle_size = "liq_bottle_size", bottle_size_unit = NULL, bottle_lasts = "liq_lasts", bottle_lasts_unit = NULL, preparation_dose = NULL, preparation_unit = NULL, timeperiod = "4 months", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_dose_unit = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = table, basis_strength_unit = NULL)med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_liq.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx",package = "packDAMipd") table <- load_trial_data(conv_file) res <- microcosting_liquids_wide( ind_part_data = ind_part_data, name_med = "liq_name", brand_med = NULL, dose_med = "liq_strength", unit_med = NULL, bottle_size = "liq_bottle_size", bottle_size_unit = NULL, bottle_lasts = "liq_lasts", bottle_lasts_unit = NULL, preparation_dose = NULL, preparation_unit = NULL, timeperiod = "4 months", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_brand = NULL, list_of_code_dose_unit = NULL, list_of_code_bottle_size_unit = NULL, list_of_code_bottle_lasts_unit = NULL, list_preparation_dose_unit = NULL, eqdose_covtab = table, basis_strength_unit = NULL)
#'########################################################################### Function to estimate the cost of patches when IPD is in long format using a IPD data of long format
microcosting_patches_long( the_columns, ind_part_data_long, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )microcosting_patches_long( the_columns, ind_part_data_long, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )
the_columns |
columns that are to be used to convert the data from long to wide |
ind_part_data_long |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
dose_med |
dose of medication used |
unit_med |
unit of medication ; use null if its along with the dose |
no_taken |
how many taken |
freq_taken |
frequency of medication |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit in the cost is calculated |
strength_column |
column column name that contain strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_freq |
if frequency is coded, give the code:frequency pairs, optional |
list_of_code_dose_unit |
if unit is coded, give the code:unit pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
eqdose_cov_tab |
table to get the conversion factor for equivalent doses, optional |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx",package = "packDAMipd") table <- load_trial_data(conv_file) names <- colnames(ind_part_data) ending <- length(names) ind_part_data_long <- tidyr::gather(ind_part_data, measurement, value, names[2]:names[ending], factor_key = TRUE) the_columns <- c("measurement", "value") res <- microcosting_patches_long(the_columns, ind_part_data_long = ind_part_data_long, name_med = "patch_name", brand_med = "patch_brand", dose_med = "patch_strength",unit_med = NULL, no_taken = "patch_no_taken", freq_taken = "patch_frequency", timeperiod = "4 months", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = table, basis_strength_unit = "mcg/hr")med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx",package = "packDAMipd") table <- load_trial_data(conv_file) names <- colnames(ind_part_data) ending <- length(names) ind_part_data_long <- tidyr::gather(ind_part_data, measurement, value, names[2]:names[ending], factor_key = TRUE) the_columns <- c("measurement", "value") res <- microcosting_patches_long(the_columns, ind_part_data_long = ind_part_data_long, name_med = "patch_name", brand_med = "patch_brand", dose_med = "patch_strength",unit_med = NULL, no_taken = "patch_no_taken", freq_taken = "patch_frequency", timeperiod = "4 months", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = table, basis_strength_unit = "mcg/hr")
Function to estimate the cost of patches taken (from IPD)
microcosting_patches_wide( ind_part_data, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )microcosting_patches_wide( ind_part_data, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )
ind_part_data |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
dose_med |
dose of medication used |
unit_med |
unit of medication ; use null if its along with the dose |
no_taken |
how many taken |
freq_taken |
frequency of medication |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit where the cost is calculated |
strength_column |
column column name that contain strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_freq |
if frequency is coded, give the code:frequency pairs, optional |
list_of_code_dose_unit |
if unit is coded, give the code:unit pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
eqdose_cov_tab |
table to get the conversion factor for equivalent doses, optional, but the column names have to unique Similar to c("Drug", "form", "unit", "factor") or c("Drug", "form", "unit", "conversion") |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
Assumes individual level data has name of medication, dose, dose unit, number taken, frequency taken, and basis time Assumes unit cost data contains the name of medication, form/type, strength, unit of strength (or the unit in which the cost calculated), preparation, unit cost, size and size unit (in which name, forms, size, size unit, and preparation are not passed on) @importFrom dplyr %>% a patient use 1 mg/hr patches 5 patches once a week that patch comes in a pack of 4 with cost £2.50 we want to estimate the cost for 3 months that means amount of medication 3 months = 21 weeks number of patches taken = 215 = 105 patches packs = (105/4) almost 27 packs cost = 272.50
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) res <- microcosting_patches_wide( ind_part_data = ind_part_data, name_med = "patch_name", brand_med = "patch_brand", dose_med = "patch_strength", unit_med = NULL, no_taken = "patch_no_taken", freq_taken = "patch_frequency", timeperiod = "4 months", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = table, basis_strength_unit = "mcg/hr")med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) res <- microcosting_patches_wide( ind_part_data = ind_part_data, name_med = "patch_name", brand_med = "patch_brand", dose_med = "patch_strength", unit_med = NULL, no_taken = "patch_no_taken", freq_taken = "patch_frequency", timeperiod = "4 months", unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = table, basis_strength_unit = "mcg/hr")
Function to estimate the cost of tablets when IPD is in long format
microcosting_tablets_long( the_columns, ind_part_data_long, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )microcosting_tablets_long( the_columns, ind_part_data_long, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )
the_columns |
columns that are to be used to convert the data from long to wide |
ind_part_data_long |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
dose_med |
dose of medication used |
unit_med |
unit of medication ; use null if its along with the dose |
no_taken |
how many taken |
freq_taken |
frequency of medication |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit where the cost is calculated |
strength_column |
column column name that contain strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_freq |
if frequency is coded, give the code:frequency pairs, optional |
list_of_code_dose_unit |
if unit is coded, give the code:unit pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
eqdose_cov_tab |
table to get the conversion factor for equivalent doses, optional |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_all.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) names <- colnames(ind_part_data) ending <- length(names) ind_part_data_long <- tidyr::gather(ind_part_data, measurement, value, names[2]:names[ending], factor_key = TRUE) the_columns <- c("measurement", "value") res <- microcosting_tablets_long(the_columns, ind_part_data_long = ind_part_data_long, name_med = "tab_name", brand_med = "tab_brand", dose_med = "tab_strength", unit_med = "tab_str_unit", no_taken = "tab_no_taken", freq_taken = "tab_frequency", timeperiod = "2 months",unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL,list_of_code_dose_unit = NULL, eqdose_cov_tab = table, basis_strength_unit = "mg")med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_all.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx", package = "packDAMipd") table <- load_trial_data(conv_file) names <- colnames(ind_part_data) ending <- length(names) ind_part_data_long <- tidyr::gather(ind_part_data, measurement, value, names[2]:names[ending], factor_key = TRUE) the_columns <- c("measurement", "value") res <- microcosting_tablets_long(the_columns, ind_part_data_long = ind_part_data_long, name_med = "tab_name", brand_med = "tab_brand", dose_med = "tab_strength", unit_med = "tab_str_unit", no_taken = "tab_no_taken", freq_taken = "tab_frequency", timeperiod = "2 months",unit_cost_data = med_costs, unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL,list_of_code_dose_unit = NULL, eqdose_cov_tab = table, basis_strength_unit = "mg")
Function to estimate the cost of tablets taken (from IPD)
microcosting_tablets_wide( ind_part_data, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )microcosting_tablets_wide( ind_part_data, name_med, brand_med = NULL, dose_med, unit_med = NULL, no_taken, freq_taken, timeperiod, unit_cost_data, unit_cost_column, cost_calculated_per, strength_column, list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, list_of_code_brand = NULL, eqdose_cov_tab = NULL, basis_strength_unit = NULL )
ind_part_data |
IPD |
name_med |
name of medication |
brand_med |
brand name of medication if revealed |
dose_med |
dose of medication used |
unit_med |
unit of medication ; use null if its along with the dose |
no_taken |
how many taken |
freq_taken |
frequency of medication |
timeperiod |
time period for cost calculation |
unit_cost_data |
unit costs data |
unit_cost_column |
column name of unit cost in unit_cost_data |
cost_calculated_per |
column name of unit where the cost is calculated |
strength_column |
column column name that contain strength of medication |
list_of_code_names |
if names is coded, give the code:name pairs, optional |
list_of_code_freq |
if frequency is coded, give the code:frequency pairs, optional |
list_of_code_dose_unit |
if unit is coded, give the code:unit pairs, optional |
list_of_code_brand |
if brand names are coded, give the code:brand pairs, optional |
eqdose_cov_tab |
table to get the conversion factor for equivalent doses, optional, but the column names have to unique Similar to c("Drug", "form", "unit", "factor") or c("Drug", "form", "unit", "conversion") |
basis_strength_unit |
strength unit to be taken as basis required for total medication calculations |
Assumes individual level data has name of medication, dose, dose unit, number taken, frequency taken, and basis time Assumes unit cost data contains the name of medication, form/type, strength, unit of strength (or the unit in which the cost calculated), preparation, unit cost, size and size unit (in which name, forms, size, size unit, and preparation are not passed on) @importFrom dplyr %>%
the calculated cost of tablets along with original data
med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_all.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx",package = "packDAMipd") table <- load_trial_data(conv_file) res <- microcosting_tablets_wide(ind_part_data = ind_part_data, name_med = "tab_name", brand_med = "tab_brand", dose_med = "tab_strength", unit_med = "tab_str_unit", no_taken = "tab_no_taken", freq_taken = "tab_frequency",timeperiod = "2 months", unit_cost_data = med_costs,unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, eqdose_cov_tab = table, basis_strength_unit = "mg")med_costs_file <- system.file("extdata", "medicaton_costs_all.xlsx", package = "packDAMipd") data_file <- system.file("extdata", "medication_all.xlsx", package = "packDAMipd") ind_part_data <- load_trial_data(data_file) med_costs <- load_trial_data(med_costs_file) conv_file <- system.file("extdata", "Med_calc.xlsx",package = "packDAMipd") table <- load_trial_data(conv_file) res <- microcosting_tablets_wide(ind_part_data = ind_part_data, name_med = "tab_name", brand_med = "tab_brand", dose_med = "tab_strength", unit_med = "tab_str_unit", no_taken = "tab_no_taken", freq_taken = "tab_frequency",timeperiod = "2 months", unit_cost_data = med_costs,unit_cost_column = "UnitCost", cost_calculated_per = "Basis", strength_column = "Strength", list_of_code_names = NULL, list_of_code_freq = NULL, list_of_code_dose_unit = NULL, eqdose_cov_tab = table, basis_strength_unit = "mg")
Function to plot CEAC
plot_ceac_psa( control_markov, trt_markov, psa_table_ctrl, psa_table_trt, thresholds, num_rep, comparator )plot_ceac_psa( control_markov, trt_markov, psa_table_ctrl, psa_table_trt, thresholds, num_rep, comparator )
control_markov |
markov model for control |
trt_markov |
markov model for treatment |
psa_table_ctrl |
param table for psa for control |
psa_table_trt |
param table for psa for treatment |
thresholds |
threshold values of WTP |
num_rep |
number of repetitions |
comparator |
the strategy to be compared with |
plot of ceac
Function to plot results of sensitivity analysis do_sensitivity_analysis()
plot_dsa( result_dsa_control, plotfor, type = "range", result_dsa_treat = NULL, threshold = NULL, comparator = NULL )plot_dsa( result_dsa_control, plotfor, type = "range", result_dsa_treat = NULL, threshold = NULL, comparator = NULL )
result_dsa_control |
result from deterministic sensitivity analysis for first or control model |
plotfor |
the variable to plotfor e.g. cost, utility NMB etc |
type |
type of analysis, range or difference |
result_dsa_treat |
result from deterministic sensitivity analysis for the comparative Markov model |
threshold |
threshold value of WTP |
comparator |
the strategy to be compared with |
plot of sensitivity analysis
param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido") low_values <- define_parameters(cost_direct_med_B = 177.4, cost_comm_care_C = 205.9) upp_values <- define_parameters(cost_direct_med_B = 17740, cost_comm_care_C = 20590) A <- health_state("A", cost = "cost_health_A + cost_drug ", utility = 1) B <- health_state("B", cost = "cost_health_B + cost_drug", utility = 1) C <- health_state("C", cost = "cost_health_C + cost_drug", utility = 1) D <- health_state("D", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3, 4), c(NA, 5, 6, 7), c(NA, NA, 8, 9), c(NA, NA, NA, 10)) colnames(tmat) <- rownames(tmat) <- c("A", "B", "C", "D") tm <- populate_transition_matrix(4, tmat, c("tpAtoA", "tpAtoB", "tpAtoC", "tpAtoD","tpBtoB", "tpBtoC", "tpBtoD", "tpCtoC", "tpCtoD", "tpDtoD"), colnames(tmat)) health_states <- combine_state(A, B, C, D) mono_strategy <- strategy(tm, health_states, "mono") mono_markov <- markov_model(mono_strategy, 20, c(1, 0, 0, 0), discount = c(0.06, 0), param_list) param_table <- define_parameters_sens_anal(param_list, low_values, upp_values) result <- do_sensitivity_analysis(mono_markov, param_table) param_list_treat <- define_parameters( cost_zido = 3000, cost_direct_med_A = 890, cost_comm_care_A = 8976, cost_direct_med_B = 2345, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido") treat_strategy <- strategy(tm, health_states, "treat") treat_markov <- markov_model(treat_strategy, 20, c(1, 0, 0, 0), discount = c(0.06, 0), param_list_treat) treat_low_values <- define_parameters(cost_direct_med_B = 234.5, cost_comm_care_C = 694.8) treat_upp_values <- define_parameters(cost_direct_med_B = 23450, cost_comm_care_C = 69480) param_table_treat <- define_parameters_sens_anal(param_list_treat, treat_low_values,treat_upp_values) result_treat <- do_sensitivity_analysis(treat_markov, param_table) plot_dsa(result,"NMB","range",result_treat, 20000, "treat")param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido") low_values <- define_parameters(cost_direct_med_B = 177.4, cost_comm_care_C = 205.9) upp_values <- define_parameters(cost_direct_med_B = 17740, cost_comm_care_C = 20590) A <- health_state("A", cost = "cost_health_A + cost_drug ", utility = 1) B <- health_state("B", cost = "cost_health_B + cost_drug", utility = 1) C <- health_state("C", cost = "cost_health_C + cost_drug", utility = 1) D <- health_state("D", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3, 4), c(NA, 5, 6, 7), c(NA, NA, 8, 9), c(NA, NA, NA, 10)) colnames(tmat) <- rownames(tmat) <- c("A", "B", "C", "D") tm <- populate_transition_matrix(4, tmat, c("tpAtoA", "tpAtoB", "tpAtoC", "tpAtoD","tpBtoB", "tpBtoC", "tpBtoD", "tpCtoC", "tpCtoD", "tpDtoD"), colnames(tmat)) health_states <- combine_state(A, B, C, D) mono_strategy <- strategy(tm, health_states, "mono") mono_markov <- markov_model(mono_strategy, 20, c(1, 0, 0, 0), discount = c(0.06, 0), param_list) param_table <- define_parameters_sens_anal(param_list, low_values, upp_values) result <- do_sensitivity_analysis(mono_markov, param_table) param_list_treat <- define_parameters( cost_zido = 3000, cost_direct_med_A = 890, cost_comm_care_A = 8976, cost_direct_med_B = 2345, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido") treat_strategy <- strategy(tm, health_states, "treat") treat_markov <- markov_model(treat_strategy, 20, c(1, 0, 0, 0), discount = c(0.06, 0), param_list_treat) treat_low_values <- define_parameters(cost_direct_med_B = 234.5, cost_comm_care_C = 694.8) treat_upp_values <- define_parameters(cost_direct_med_B = 23450, cost_comm_care_C = 69480) param_table_treat <- define_parameters_sens_anal(param_list_treat, treat_low_values,treat_upp_values) result_treat <- do_sensitivity_analysis(treat_markov, param_table) plot_dsa(result,"NMB","range",result_treat, 20000, "treat")
Function to do some checks before plotting sensitivity analysis results
plot_dsa_difference(ob_results, plotfor, plot_var)plot_dsa_difference(ob_results, plotfor, plot_var)
ob_results |
results from deterministic sensitivity analysis |
plotfor |
the quantity plotting |
plot_var |
the variable |
plot
Function to do some checks before plotting sensitivity analysis results
plot_dsa_icer_range(ob_results, plot_var)plot_dsa_icer_range(ob_results, plot_var)
ob_results |
results from deterministic sensitivity analysis |
plot_var |
the variable |
plot
Function to do some checks before plotting sensitivity analysis results
plot_dsa_nmb_range(ob_results, plot_var)plot_dsa_nmb_range(ob_results, plot_var)
ob_results |
results from deterministic sensitivity analysis |
plot_var |
the variable |
plot
Function to do some checks before plotting sensitivity analysis results
plot_dsa_others_range(ob_results, plot_var)plot_dsa_others_range(ob_results, plot_var)
ob_results |
results from deterministic sensitivity analysis |
plot_var |
the variable |
plot
Plot efficiency frontier
plot_efficiency_frontier(results_calculate_icer_nmb, threshold)plot_efficiency_frontier(results_calculate_icer_nmb, threshold)
results_calculate_icer_nmb |
results from cea (from calculate_icer_nmb method) |
threshold |
threshold value |
plot well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "control") this_markov <- markov_model(this.strategy, 24, c(1000, 0, 0), c(0,0)) well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 10, utility = 0.5) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.4, 0.4, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "intervention") sec_markov <- markov_model(this.strategy, 24, c(1000, 0, 0), c(0,0)) list_markov <- combine_markov(this_markov, sec_markov) results_cea <- calculate_icer_nmb(list_markov, 20000, comparator = "control") plot_efficiency_frontier(results_cea, c(1000, 2000))
Function to plot mean and SE for longitudinal observations for twogroups compared
plot_histogram_onetimepoint_twogroups( thedata, colname, timepointstring, xstring, ylimits = NULL, nbins = NULL )plot_histogram_onetimepoint_twogroups( thedata, colname, timepointstring, xstring, ylimits = NULL, nbins = NULL )
thedata |
the data where the observations are held |
colname |
columnname in the data where the intersted observations |
timepointstring |
the text that correspond to timepoints at which the descriptive analysis is done |
xstring |
xlable text |
ylimits |
on hsitogram plot |
nbins |
nbins the number of bins to plot histogram |
the historgram plot at a timepoint for two groups
eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 7, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"))) plot_histogram_onetimepoint_twogroups(eg_data, c("mark_at_2"), c("1","2"), "mark")eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 7, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"))) plot_histogram_onetimepoint_twogroups(eg_data, c("mark_at_2"), c("1","2"), "mark")
Function to plot mean and SE for longitudinal observations for twogroups compared
plot_meanSE_longitudinal_twogroups( thedata, columnnames, timepoints, observation )plot_meanSE_longitudinal_twogroups( thedata, columnnames, timepoints, observation )
thedata |
the data where the observations are held |
columnnames |
columnnames in the data where the intersted observations |
timepoints |
the timepoints at which the descriptive analysis is done |
observation |
name of the observations |
the plot that shows mean and SE
eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 7, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"))) plot_meanSE_longitudinal_twogroups(eg_data, c("mark_at_1", "mark_at_2"), c("1","2"), "mark")eg_data <- as.data.frame(list(no = c(1, 2, 3, 4), mark_at_1 = c(12, 7, 23, 45), gender = c("M", "F", "M", "F"), mark_at_2 = c(12, 34, 89, 45), trialarm = c("1","1","2","2"))) plot_meanSE_longitudinal_twogroups(eg_data, c("mark_at_1", "mark_at_2"), c("1","2"), "mark")
E1. Plot a Markov model
plot_model(markov)plot_model(markov)
markov |
markov_model object |
plots
tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- populate_transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0, 0, TRUE) health_states <- combine_state(a, b) this.strategy <- strategy(tm, health_states, "intervention") this_markov <- markov_model(this.strategy, 10, c(1, 0), c(0, 0)) p <- plot_model(this_markov)tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- populate_transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0, 0, TRUE) health_states <- combine_state(a, b) this.strategy <- strategy(tm, health_states, "intervention") this_markov <- markov_model(this.strategy, 10, c(1, 0), c(0, 0)) p <- plot_model(this_markov)
Plot cost effectiveness acceptability curve
plot_nmb_lambda(list_markov, threshold_values, comparator, currency = "GBP")plot_nmb_lambda(list_markov, threshold_values, comparator, currency = "GBP")
list_markov |
markov_model objects |
threshold_values |
list of threshold values |
comparator |
the comparator |
currency |
currency |
plots
well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "control") this_markov <- markov_model(this.strategy, 24, c(1000, 0, 0),c(0, 0)) well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 10, utility = 0.5) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.4, 0.4, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "intervention") sec_markov <- markov_model(this.strategy, 24, c(1000, 0, 0), c(0, 0)) list_markov <- combine_markov(this_markov, sec_markov) plot_nmb_lambda(list_markov, c(1000, 2000, 3000), comparator = "control")well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 100, utility = 1) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.6, 0.2, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "control") this_markov <- markov_model(this.strategy, 24, c(1000, 0, 0),c(0, 0)) well <- health_state("well", cost = 0, utility = 1) disabled <- health_state("disabled", cost = 10, utility = 0.5) dead <- health_state("dead", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("well", "disabled", "dead") tm <- populate_transition_matrix(3, tmat, c(0.4, 0.4, 0.2, 0.6, 0.4, 1), colnames(tmat)) health_states <- combine_state(well, disabled, dead) this.strategy <- strategy(tm, health_states, "intervention") sec_markov <- markov_model(this.strategy, 24, c(1000, 0, 0), c(0, 0)) list_markov <- combine_markov(this_markov, sec_markov) plot_nmb_lambda(list_markov, c(1000, 2000, 3000), comparator = "control")
Plot the predicted survival curves for covariates keeping the others fixed
plot_prediction_parametric_survival( param_to_be_estimated, indep_var, covariates, dataset, fit, timevar_survival )plot_prediction_parametric_survival( param_to_be_estimated, indep_var, covariates, dataset, fit, timevar_survival )
param_to_be_estimated |
parameter to be estimated |
indep_var |
variable for which the levels have to be identified |
covariates |
the covariates |
dataset |
the dataset where these variables contain |
fit |
the fit result survreg |
timevar_survival |
time variable from the dataset |
plot
data_for_survival <- survival::lung surv_estimated <- use_parametric_survival("status", data_for_survival, "sex", info_distribution = "weibull",covariates = c("ph.ecog"),"time") plot_prediction_parametric_survival("status", "sex", covariates = c("ph.ecog"),data_for_survival, surv_estimated$fit, "time")data_for_survival <- survival::lung surv_estimated <- use_parametric_survival("status", data_for_survival, "sex", info_distribution = "weibull",covariates = c("ph.ecog"),"time") plot_prediction_parametric_survival("status", "sex", covariates = c("ph.ecog"),data_for_survival, surv_estimated$fit, "time")
Plotting and return the residuals after cox proportional hazard model
plot_return_residual_cox( param_to_be_estimated, indep_var, covariates, fit, dataset )plot_return_residual_cox( param_to_be_estimated, indep_var, covariates, fit, dataset )
param_to_be_estimated |
parameter to be estimated |
indep_var |
independent variable |
covariates |
covariates |
fit |
fit object from coxph method |
dataset |
data used for cox ph model |
plot and the residuals
data_for_survival <- survival::lung surv_estimated <- use_coxph_survival("status", data_for_survival, "sex", covariates = c("ph.ecog"), "time") plot_return_residual_cox("status", "sex", covariates = c("ph.ecog"), surv_estimated$fit,data_for_survival )data_for_survival <- survival::lung surv_estimated <- use_coxph_survival("status", data_for_survival, "sex", covariates = c("ph.ecog"), "time") plot_return_residual_cox("status", "sex", covariates = c("ph.ecog"), surv_estimated$fit,data_for_survival )
Plotting and return the residuals after survival model
plot_return_residual_survival( param_to_be_estimated, indep_var, covariates, fit )plot_return_residual_survival( param_to_be_estimated, indep_var, covariates, fit )
param_to_be_estimated |
parameter to be estimated |
indep_var |
independent variable |
covariates |
covariates |
fit |
fit object from survreg method |
plot and the residuals
data_for_survival <- survival::lung surv_estimated <- use_parametric_survival("status", data_for_survival, "sex", info_distribution = "weibull",covariates = c("ph.ecog"), "time") plot_return_residual_survival("status", "sex", covariates = c("ph.ecog"),surv_estimated$fit)data_for_survival <- survival::lung surv_estimated <- use_parametric_survival("status", data_for_survival, "sex", info_distribution = "weibull",covariates = c("ph.ecog"), "time") plot_return_residual_survival("status", "sex", covariates = c("ph.ecog"),surv_estimated$fit)
Plotting survival function for all covariates using survfit
plot_return_survival_curve( param_to_be_estimated, dataset, indep_var, covariates, timevar_survival )plot_return_survival_curve( param_to_be_estimated, dataset, indep_var, covariates, timevar_survival )
param_to_be_estimated |
parameter to be estimated |
dataset |
param describing the methods |
indep_var |
independent variable |
covariates |
covariates |
timevar_survival |
time variable for survival analysis |
plot and the survival function values
data_for_survival <- survival::lung plot_return_survival_curve(param_to_be_estimated = "status", dataset = data_for_survival,indep_var = "sex",covariates = c("ph.ecog"), timevar_survival = "time")data_for_survival <- survival::lung plot_return_survival_curve(param_to_be_estimated = "status", dataset = data_for_survival,indep_var = "sex",covariates = c("ph.ecog"), timevar_survival = "time")
Plotting survival function for all covariates calculated from cox regression results and returned coefficients
plot_survival_cox_covariates( coxfit, dataset, param_to_be_estimated, covariates, indep_var )plot_survival_cox_covariates( coxfit, dataset, param_to_be_estimated, covariates, indep_var )
coxfit |
cox regression fit result |
dataset |
param describing the methods |
param_to_be_estimated |
parameter to be estimated |
covariates |
covariates |
indep_var |
independent variable |
plot and the survival function values
data_for_survival <- survival::lung surv_estimated <- use_coxph_survival("status", data_for_survival, "sex", covariates = c("ph.ecog"), "time") plot_survival_cox_covariates(surv_estimated$fit,data_for_survival, "status", covariates = c("ph.ecog"), "sex")data_for_survival <- survival::lung surv_estimated <- use_coxph_survival("status", data_for_survival, "sex", covariates = c("ph.ecog"), "time") plot_survival_cox_covariates(surv_estimated$fit,data_for_survival, "status", covariates = c("ph.ecog"), "sex")
Populate transition matrix
populate_transition_matrix(no_states, tmat, list_prob, name_states = NULL)populate_transition_matrix(no_states, tmat, list_prob, name_states = NULL)
no_states |
number of the health states |
tmat |
A transition matrix in the format from the package 'mstate' |
list_prob |
list of probabilities as in the order of transitions (row wise) |
name_states |
names of the health states |
If the state names are null, they are replaced with numbers starting from 1 First find those missing probabilities, and fill a list from the given list of probabilities and fill those are not NA in the matrix Note that the probabilities need not be numeric here and no checks are needed for sum
value of the transition matrix
tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") populate_transition_matrix(2, tmat, list_prob = c(0.2, 0.5, 0, 0.3))tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") populate_transition_matrix(2, tmat, list_prob = c(0.2, 0.5, 0, 0.3))
Predict risk/hazard function for cox ph regression
predict_coxph( coxfit, dataset, param_to_be_estimated, covariates, indep_var, timevar_survival )predict_coxph( coxfit, dataset, param_to_be_estimated, covariates, indep_var, timevar_survival )
coxfit |
cox regression fit result |
dataset |
param describing the methods |
param_to_be_estimated |
parameter to be estimated |
covariates |
covariates |
indep_var |
independent variable |
timevar_survival |
time variable |
"risk" option for "type" returns the hazard ratio relative to mean e.g given below For lung data with data_for_survival <- survival::lung fit <- use_coxph_survival("status", data_for_survival, "sex", covariates = c("ph.ecog"), "time") coeffit = fit$coefficients r1234 <- exp(coeffit("sex")lung$sex+ coeffit("ph.ecog")lung$ph.ecog) rMean <- exp(sum(coef(fit) * fit$means, na.rm=TRUE)) rr <- r1234/rMean
plot and the survival function values
data_for_survival <- survival::lung surv_estimated <- use_coxph_survival("status", data_for_survival, "sex",covariates = c("ph.ecog"), "time") predict_coxph(surv_estimated$fit,data_for_survival, "status","sex", covariates = c("ph.ecog"), "time")data_for_survival <- survival::lung surv_estimated <- use_coxph_survival("status", data_for_survival, "sex",covariates = c("ph.ecog"), "time") predict_coxph(surv_estimated$fit,data_for_survival, "status","sex", covariates = c("ph.ecog"), "time")
promis 3a scoring table
promis3a_scoring.dfpromis3a_scoring.df
A 14 by 3 dataframe
created on April 08, 2021
Function to report deterministic sensitivity analysis
report_sensitivity_analysis( result_dsa_control, result_dsa_treat = NULL, threshold = NULL, comparator = NULL )report_sensitivity_analysis( result_dsa_control, result_dsa_treat = NULL, threshold = NULL, comparator = NULL )
result_dsa_control |
result from deterministic sensitivity analysis for first or control model |
result_dsa_treat |
result from deterministic sensitivity analysis for the comparative Markov model |
threshold |
threshold value of WTP |
comparator |
the strategy to be compared with |
report in the form of a table
param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido") low_values <- define_parameters(cost_direct_med_B = 177.4, cost_comm_care_C = 205.9) upp_values <- define_parameters(cost_direct_med_B = 17740, cost_comm_care_C = 20590) A <- health_state("A", cost = "cost_health_A + cost_drug ", utility = 1) B <- health_state("B", cost = "cost_health_B + cost_drug", utility = 1) C <- health_state("C", cost = "cost_health_C + cost_drug", utility = 1) D <- health_state("D", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3, 4), c(NA, 5, 6, 7), c(NA, NA, 8, 9), c(NA, NA, NA, 10)) colnames(tmat) <- rownames(tmat) <- c("A", "B", "C", "D") tm <- populate_transition_matrix(4, tmat, c("tpAtoA", "tpAtoB", "tpAtoC", "tpAtoD","tpBtoB", "tpBtoC", "tpBtoD", "tpCtoC", "tpCtoD", "tpDtoD"), colnames(tmat)) health_states <- combine_state(A, B, C, D) mono_strategy <- strategy(tm, health_states, "mono") mono_markov <- markov_model(mono_strategy, 20, c(1, 0, 0, 0), discount = c(0.06, 0), param_list) param_table <- define_parameters_sens_anal(param_list, low_values, upp_values) result <- do_sensitivity_analysis(mono_markov, param_table) reporting <- report_sensitivity_analysis(result)param_list <- define_parameters( cost_zido = 2278, cost_direct_med_A = 1701, cost_comm_care_A = 1055, cost_direct_med_B = 1774, cost_comm_care_B = 1278, cost_direct_med_C = 6948, cost_comm_care_C = 2059, tpAtoA = 1251 / (1251 + 483), tpAtoB = 350 / (350 + 1384), tpAtoC = 116 / (116 + 1618), tpAtoD = 17 / (17 + 1717), tpBtoB = 731 / (731 + 527), tpBtoC = 512 / (512 + 746), tpBtoD = 15 / (15 + 1243), tpCtoC = 1312 / (1312 + 437), tpCtoD = 437 / (437 + 1312), tpDtoD = 1, cost_health_A = "cost_direct_med_A + cost_comm_care_A", cost_health_B = "cost_direct_med_B + cost_comm_care_B", cost_health_C = "cost_direct_med_C + cost_comm_care_C", cost_drug = "cost_zido") low_values <- define_parameters(cost_direct_med_B = 177.4, cost_comm_care_C = 205.9) upp_values <- define_parameters(cost_direct_med_B = 17740, cost_comm_care_C = 20590) A <- health_state("A", cost = "cost_health_A + cost_drug ", utility = 1) B <- health_state("B", cost = "cost_health_B + cost_drug", utility = 1) C <- health_state("C", cost = "cost_health_C + cost_drug", utility = 1) D <- health_state("D", cost = 0, utility = 0) tmat <- rbind(c(1, 2, 3, 4), c(NA, 5, 6, 7), c(NA, NA, 8, 9), c(NA, NA, NA, 10)) colnames(tmat) <- rownames(tmat) <- c("A", "B", "C", "D") tm <- populate_transition_matrix(4, tmat, c("tpAtoA", "tpAtoB", "tpAtoC", "tpAtoD","tpBtoB", "tpBtoC", "tpBtoD", "tpCtoC", "tpCtoD", "tpDtoD"), colnames(tmat)) health_states <- combine_state(A, B, C, D) mono_strategy <- strategy(tm, health_states, "mono") mono_markov <- markov_model(mono_strategy, 20, c(1, 0, 0, 0), discount = c(0.06, 0), param_list) param_table <- define_parameters_sens_anal(param_list, low_values, upp_values) result <- do_sensitivity_analysis(mono_markov, param_table) reporting <- report_sensitivity_analysis(result)
Function to get the subset of data compared to a string after trimming the white spaces
return_equal_liststring_col(col, the_data, list_str)return_equal_liststring_col(col, the_data, list_str)
col |
the form of medication either tablet or patch |
the_data |
the data to be get the subset from |
list_str |
list of strings to be compared |
the subset data
the_data <- as.data.frame(cbind(c("one", "two"), c("a", "b"))) colnames(the_data) <- c("name", "brand") ans <- return_equal_liststring_col(2, the_data, c("a", "cc"))the_data <- as.data.frame(cbind(c("one", "two"), c("a", "b"))) colnames(the_data) <- c("name", "brand") ans <- return_equal_liststring_col(2, the_data, c("a", "cc"))
Function to get the subset of data compared to a string after trimming the white spaces
return_equal_liststring_listcol(col, the_data, list_str)return_equal_liststring_listcol(col, the_data, list_str)
col |
the form of medication either tablet or patch |
the_data |
the data to be get the subset from |
list_str |
list of strings to be compared |
the subset data
the_data <- as.data.frame(cbind(c("one", "two"), c("tablet", "tablets"), c("aa", "bb"))) colnames(the_data) <- c("name", "brand_a", "xx") ans <- return_equal_liststring_listcol(2, the_data, c("tablet", "tablets"))the_data <- as.data.frame(cbind(c("one", "two"), c("tablet", "tablets"), c("aa", "bb"))) colnames(the_data) <- c("name", "brand_a", "xx") ans <- return_equal_liststring_listcol(2, the_data, c("tablet", "tablets"))
Function to get the subset of data compared to a string after trimming the white spaces
return_equal_str_col(col, the_data, the_str)return_equal_str_col(col, the_data, the_str)
col |
the form of medication either tablet or patch |
the_data |
the data to be get the subset from |
the_str |
the string to be compared |
the subset data
the_data <- as.data.frame(cbind(c("one", "two"), c("a", "b"))) colnames(the_data) <- c("name", "brand") ans <- return_equal_str_col(2, the_data, "a")the_data <- as.data.frame(cbind(c("one", "two"), c("a", "b"))) colnames(the_data) <- c("name", "brand") ans <- return_equal_str_col(2, the_data, "a")
Function to return 0 if the param is not null or NA trimming the white spaces
return0_if_not_null_na(param)return0_if_not_null_na(param)
param |
the form of medication either tablet or patch |
zero or -1
parame = NULL ans <- return0_if_not_null_na(parame) parame = 1 ans <- return0_if_not_null_na(parame)parame = NULL ans <- return0_if_not_null_na(parame) parame = 1 ans <- return0_if_not_null_na(parame)
Set the attribute for the health state
set_var_state(state, var, new_value)set_var_state(state, var, new_value)
state |
object of class health state |
var |
attribute of the health state |
new_value |
new value to be assigned |
After checking the given state is a health state the value of the variable is set if the value is not numeric, it is being parsed to form an expression
modified health state
set_var_state(health_state("IT", 100, 0.4, 0, FALSE), "cost", 1)set_var_state(health_state("IT", 100, 0.4, 0, FALSE), "cost", 1)
Definition of strategy - or arm
strategy(trans_mat, states, name, trans_cost = NULL, trans_util = NULL)strategy(trans_mat, states, name, trans_cost = NULL, trans_util = NULL)
trans_mat |
transition matrix |
states |
health states |
name |
name of the strategy |
trans_cost |
values of costs if these are attached to transitions |
trans_util |
values of utility if these are attached to transitions |
Defining strategy keeping all transition matrix, states and names together to use in defining Markov model
object strategy
tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- populate_transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0.5, 0, FALSE) states <- combine_state(a, b) strategy(tm, states, "intervention")tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- populate_transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, 0, FALSE) b <- health_state("Dead", 1, 0.5, 0, FALSE) states <- combine_state(a, b) strategy(tm, states, "intervention")
Function to summarise and plot probabilistic sensitivity analysis
summary_plot_psa( result_psa_params_control, result_psa_params_treat = NULL, threshold = NULL, comparator = NULL )summary_plot_psa( result_psa_params_control, result_psa_params_treat = NULL, threshold = NULL, comparator = NULL )
result_psa_params_control |
result from probabilistic sensitivity analysis for first or control model |
result_psa_params_treat |
result from probabilistic sensitivity analysis for the comparative Markov model |
threshold |
threshold value of WTP |
comparator |
the strategy to be compared with |
plot of sensitivity analysis
param_list <- define_parameters( cost_direct_med_A = 1701, cost_direct_med_B = 1774, tpAtoA = 0.2, tpAtoB = 0.5, tpAtoC = 0.3, tpBtoB = 0.3, tpBtoC = 0.7, tpCtoC = 1,cost_health_A = "cost_direct_med_A", cost_health_B = "cost_direct_med_B") sample_list <- define_parameters(cost_direct_med_A = "gamma(mean = 1701, sd = sqrt(1701))") A <- health_state("A", cost = "cost_health_A ", utility = 1) B <- health_state("B", cost = "cost_health_B", utility = 1) C <- health_state("C", cost = 0, utility = 0, absorb = "TRUE") tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("A", "B", "C") tm <- populate_transition_matrix(3, tmat, c( "tpAtoA", "tpAtoB", "tpAtoC", "tpBtoB", "tpBtoC", "tpCtoC"), colnames(tmat)) health_states <- combine_state(A, B, C) mono_strategy <- strategy(tm, health_states, "mono") mono_markov <- markov_model(mono_strategy, 20, initial_state =c(1,0,0), discount = c(0.06, 0),param_list) param_table <- define_parameters_psa(param_list, sample_list) result <- do_psa(mono_markov, param_table, 3) result_plot <- summary_plot_psa(result, NULL, NULL, NULL) param_list_comb <- define_parameters( cost_direct_med_A = 1800, cost_direct_med_B = 1774, tpAtoA = 0.6, tpAtoB = 0.1, tpAtoC = 0.3,tpBtoB = 0.3, tpBtoC = 0.7,tpCtoC = 1, cost_health_A = "cost_direct_med_A",cost_health_B = "cost_direct_med_B") comb_strategy <- strategy(tm, health_states, "comb") comb_markov <- markov_model(comb_strategy, 20, c(1, 0, 0), discount = c(0.06, 0), param_list) param_table_comb <- define_parameters_psa(param_list_comb, sample_list) result_comb <- do_psa(comb_markov, param_table_comb, 3) summary_plot_psa(result, result_comb, 2000, "mono")param_list <- define_parameters( cost_direct_med_A = 1701, cost_direct_med_B = 1774, tpAtoA = 0.2, tpAtoB = 0.5, tpAtoC = 0.3, tpBtoB = 0.3, tpBtoC = 0.7, tpCtoC = 1,cost_health_A = "cost_direct_med_A", cost_health_B = "cost_direct_med_B") sample_list <- define_parameters(cost_direct_med_A = "gamma(mean = 1701, sd = sqrt(1701))") A <- health_state("A", cost = "cost_health_A ", utility = 1) B <- health_state("B", cost = "cost_health_B", utility = 1) C <- health_state("C", cost = 0, utility = 0, absorb = "TRUE") tmat <- rbind(c(1, 2, 3), c(NA, 4, 5), c(NA, NA, 6)) colnames(tmat) <- rownames(tmat) <- c("A", "B", "C") tm <- populate_transition_matrix(3, tmat, c( "tpAtoA", "tpAtoB", "tpAtoC", "tpBtoB", "tpBtoC", "tpCtoC"), colnames(tmat)) health_states <- combine_state(A, B, C) mono_strategy <- strategy(tm, health_states, "mono") mono_markov <- markov_model(mono_strategy, 20, initial_state =c(1,0,0), discount = c(0.06, 0),param_list) param_table <- define_parameters_psa(param_list, sample_list) result <- do_psa(mono_markov, param_table, 3) result_plot <- summary_plot_psa(result, NULL, NULL, NULL) param_list_comb <- define_parameters( cost_direct_med_A = 1800, cost_direct_med_B = 1774, tpAtoA = 0.6, tpAtoB = 0.1, tpAtoC = 0.3,tpBtoB = 0.3, tpBtoC = 0.7,tpCtoC = 1, cost_health_A = "cost_direct_med_A",cost_health_B = "cost_direct_med_B") comb_strategy <- strategy(tm, health_states, "comb") comb_markov <- markov_model(comb_strategy, 20, c(1, 0, 0), discount = c(0.06, 0), param_list) param_table_comb <- define_parameters_psa(param_list_comb, sample_list) result_comb <- do_psa(comb_markov, param_table_comb, 3) summary_plot_psa(result, result_comb, 2000, "mono")
Parameter table created
table_param.dftable_param.df
A 11 by 2 dataframe
created on Jan 15, 2020
Trace matrix
trace_data.dftrace_data.df
A 11 by 2 dataframe
created on Nov 26, 2019 from tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, FALSE) b <- health_state("Dead", 1, 0, TRUE) health_states <- combine_state(a, b) this.strategy <- strategy(tm, health_states, "intervention")
Create the the values of cost and utility while transition
transition_cost_util( no_states, tmat_cost_util, list_values, name_states = NULL )transition_cost_util( no_states, tmat_cost_util, list_values, name_states = NULL )
no_states |
number of the health states |
tmat_cost_util |
A transition matrix for the cost/utility values in the format from the package 'mstate' use NA to indicate if the value is zero |
list_values |
list of probabilities as in the order of transitions (row wise) |
name_states |
names of the health states |
Similar to transition matrix but for denoting one time change during transitions
value of the transition matrix
tmat_cost <- rbind(c(NA, 1), c(NA, NA)) colnames(tmat_cost) <- rownames(tmat_cost) <- c("Healthy", "Dead") transition_cost_util(2, tmat_cost, list_values = c(500))tmat_cost <- rbind(c(NA, 1), c(NA, NA)) colnames(tmat_cost) <- rownames(tmat_cost) <- c("Healthy", "Dead") transition_cost_util(2, tmat_cost, list_values = c(500))
Example trial data
trial_data.dftrial_data.df
A 31 by 33 dataframe
created on Jan 15, 2020
########################################################################### Get the parameter values using the survival analysis using cox proportional hazard
use_coxph_survival( param_to_be_estimated, dataset, indep_var, covariates, timevar_survival )use_coxph_survival( param_to_be_estimated, dataset, indep_var, covariates, timevar_survival )
param_to_be_estimated |
parameter of interest |
dataset |
data set to be provided |
indep_var |
the independent variable (column name in data file) |
covariates |
list of covariates - calculations to be done before passing |
timevar_survival |
time variable for survival analysis, default is NA false by default |
plots baseline cumulative hazard function, survival function for each covariate while keeping the other fixed at the mean value (using plot_survival_cox_covariates), survival function for each combination of covariate using survfit (using plot_return_survival_curve) and test for cox regression results It also returns risk relative to mean (predicted at mean value of each covariate) along with the fit results coefficients, SE of coefficients, summary, and analysis of deviance
the results of the regression analysis
data_for_survival <- survival::aml surv_estimated <- use_coxph_survival("status", data_for_survival, "x", covariates = NA, "time") data_for_survival <- survival::lung surv_estimated <- use_coxph_survival("status", data_for_survival, "sex", covariates = c("ph.ecog"), "time")data_for_survival <- survival::aml surv_estimated <- use_coxph_survival("status", data_for_survival, "x", covariates = NA, "time") data_for_survival <- survival::lung surv_estimated <- use_coxph_survival("status", data_for_survival, "sex", covariates = c("ph.ecog"), "time")
########################################################################### Get the parameter values using the survival analysis method FH
use_fh_survival( param_to_be_estimated, dataset, indep_var, covariates, timevar_survival )use_fh_survival( param_to_be_estimated, dataset, indep_var, covariates, timevar_survival )
param_to_be_estimated |
parameter of interest |
dataset |
data set to be provided |
indep_var |
the independent variable (column name in data file) |
covariates |
list of covariates |
timevar_survival |
time variable for survival analysis, default is NA |
This function is for survival analysis using FH. This plots the cumulative survival function for each combination of covariate If the covariate is numeric, R takes it as different levels. The plot uses the returned list of survfit and extracts the time and the strata from summary of the fit (implemented in plot_return_survival_curve function)
the results of the regression analysis
data_for_survival <- survival::aml surv_estimated <- use_fh_survival("status", data_for_survival, "x", covariates = NA, "time" )data_for_survival <- survival::aml surv_estimated <- use_fh_survival("status", data_for_survival, "x", covariates = NA, "time" )
############################################################################ Get the parameter values using the survival analysis using FH2 method
use_fh2_survival( param_to_be_estimated, dataset, indep_var, covariates, timevar_survival )use_fh2_survival( param_to_be_estimated, dataset, indep_var, covariates, timevar_survival )
param_to_be_estimated |
parameter of interest |
dataset |
data set to be provided |
indep_var |
the independent variable (column name in data file) |
covariates |
list of covariates |
timevar_survival |
time variable for survival analysis, default is NA false by default |
This function is for survival analysis using FH2. This plots the cumulative survival function for each combination of covariate If the covariate is numeric, R takes it as different levels. The plot uses the returned list of survfit and extracts the time and the strata from summary of the fit (implemented in plot_return_survival_curve function)
the results of the regression analysis
data_for_survival <- survival::aml surv_estimated <- use_fh2_survival("status", data_for_survival, "x", covariates = NA, "time")data_for_survival <- survival::aml surv_estimated <- use_fh2_survival("status", data_for_survival, "x", covariates = NA, "time")
Function for generalised linear mixed model
use_generalised_linear_mixed_model( param_to_be_estimated, dataset, fix_eff, fix_eff_interact_vars, random_intercept_vars, nested_intercept_vars_pairs, cross_intercept_vars_pairs, uncorrel_slope_intercept_pairs, random_slope_intercept_pairs, family, link, package_mixed_model )use_generalised_linear_mixed_model( param_to_be_estimated, dataset, fix_eff, fix_eff_interact_vars, random_intercept_vars, nested_intercept_vars_pairs, cross_intercept_vars_pairs, uncorrel_slope_intercept_pairs, random_slope_intercept_pairs, family, link, package_mixed_model )
param_to_be_estimated |
column name of dependent variable |
dataset |
a dataframe |
fix_eff |
names of variables as fixed effect predictors |
fix_eff_interact_vars |
those of the fixed effect predictors that show interaction |
random_intercept_vars |
names of variables for random intercept |
nested_intercept_vars_pairs |
those of the random intercept variables with nested effect |
cross_intercept_vars_pairs |
those of the random intercept variables with crossed effect |
uncorrel_slope_intercept_pairs |
variables with no correlated intercepts |
random_slope_intercept_pairs |
random slopes intercept pairs - this is a list of paired variables |
family |
family of distributions for the response variable |
link |
link function for the variances |
package_mixed_model |
package to be used for mixed model |
result regression result with plot if success and -1, if failure
datafile <- system.file("extdata", "culcita_data.csv", package = "packDAMipd") dataset <- read.csv(datafile) results1 = use_generalised_linear_mixed_model("predation", dataset = datafile,fix_eff = c("ttt"), family = "binomial", fix_eff_interact_vars = NULL, random_intercept_vars = c("block"), nested_intercept_vars_pairs = NULL, cross_intercept_vars_pairs = NULL, uncorrel_slope_intercept_pairs = NULL, random_slope_intercept_pairs = NULL, link = NA, package_mixed_model = NA)datafile <- system.file("extdata", "culcita_data.csv", package = "packDAMipd") dataset <- read.csv(datafile) results1 = use_generalised_linear_mixed_model("predation", dataset = datafile,fix_eff = c("ttt"), family = "binomial", fix_eff_interact_vars = NULL, random_intercept_vars = c("block"), nested_intercept_vars_pairs = NULL, cross_intercept_vars_pairs = NULL, uncorrel_slope_intercept_pairs = NULL, random_slope_intercept_pairs = NULL, link = NA, package_mixed_model = NA)
############################################################################ Get the parameter values using logistic regression
use_generalised_linear_model( param_to_be_estimated, dataset, indep_var, family, covariates, interaction, naaction, link = NA )use_generalised_linear_model( param_to_be_estimated, dataset, indep_var, family, covariates, interaction, naaction, link = NA )
param_to_be_estimated |
parameter of interest |
dataset |
data set to be provided |
indep_var |
the independent variable (column name in data file) |
family |
distribution name eg. for logistic regression -binomial |
covariates |
list of covariates-calculations to be done before passing |
interaction |
boolean value to indicate interaction in the case of linear regression |
naaction |
action to be taken with the missing values |
link |
link function if not the default for each family |
This function returns the results and plots after doing linear regression Requires param to be estimated, dataset, independent variables and information on covariates, and interaction variables if there are Uses form_expression_glm to create the expression as per R standard for e.g glm(y ~ x ). Returns the fit result,s summary results as returned by summary(), confidence interval for fit coefficients (ci_coeff), variance covariance matrix, cholesky decomposition matrix, results from correlation test, plot of diagnostic tests and model fit assumptions, plot of model prediction diagnostic include AIC, R2, and BIC. The results of the prediction ie predicted values for fixed other variables will be returned in prediction matrix
the results of the regression analysis
gm_result <- use_generalised_linear_model( param_to_be_estimated = "Direction", dataset = ISLR::Smarket, indep_var = "Lag1", family = "binomial", covariates = c("Lag2", "Lag3"), interaction = FALSE, naaction = "na.omit", link = NA)gm_result <- use_generalised_linear_model( param_to_be_estimated = "Direction", dataset = ISLR::Smarket, indep_var = "Lag1", family = "binomial", covariates = c("Lag2", "Lag3"), interaction = FALSE, naaction = "na.omit", link = NA)
########################################################################### Get the parameter values using the Kaplan-Meier survival analysis
use_km_survival( param_to_be_estimated, dataset, indep_var, covariates, timevar_survival )use_km_survival( param_to_be_estimated, dataset, indep_var, covariates, timevar_survival )
param_to_be_estimated |
parameter of interest |
dataset |
data set to be provided |
indep_var |
the independent variable (column name in data file) |
covariates |
list of covariates |
timevar_survival |
time variable for survival analysis, default is NA |
This function is for survival analysis using Kaplan Meier. This plots the cumulative survival function for each combination of covariate If the covariate is numeric, R takes it as different levels. The plot uses the returned list of survfit and extracts the time and the strata from summary of the fit (implemented in plot_return_survival_curve function)
the results of the regression analysis, fit results, summary and plot
data_for_survival <- survival::aml surv_estimated <- use_km_survival("status", data_for_survival, "x", covariates = NA, "time") data_for_survival <- survival::lung surv_estimated <- use_km_survival("status", data_for_survival, "sex", covariates = c("ph.ecog"), "time")data_for_survival <- survival::aml surv_estimated <- use_km_survival("status", data_for_survival, "x", covariates = NA, "time") data_for_survival <- survival::lung surv_estimated <- use_km_survival("status", data_for_survival, "sex", covariates = c("ph.ecog"), "time")
Function for mixed effect regression
use_linear_mixed_model( param_to_be_estimated, dataset, fix_eff, fix_eff_interact_vars, random_intercept_vars, nested_intercept_vars_pairs, cross_intercept_vars_pairs, uncorrel_slope_intercept_pairs, random_slope_intercept_pairs, package_mixed_model )use_linear_mixed_model( param_to_be_estimated, dataset, fix_eff, fix_eff_interact_vars, random_intercept_vars, nested_intercept_vars_pairs, cross_intercept_vars_pairs, uncorrel_slope_intercept_pairs, random_slope_intercept_pairs, package_mixed_model )
param_to_be_estimated |
column name of dependent variable |
dataset |
a dataframe |
fix_eff |
names of variables as fixed effect predictors |
fix_eff_interact_vars |
those of the fixed effect predictors that show interaction |
random_intercept_vars |
names of variables for random intercept |
nested_intercept_vars_pairs |
those of the random intercept variables with nested effect |
cross_intercept_vars_pairs |
those of the random intercept variables with crossed effect |
uncorrel_slope_intercept_pairs |
variables with no correlated intercepts |
random_slope_intercept_pairs |
random slopes intercept pairs - this is a list of paired variables |
package_mixed_model |
package to be used for mixed model |
result regression result with plot if success and -1, if failure
datafile <- system.file("extdata", "data_linear_mixed_model.csv", package = "packDAMipd") dataset = utils::read.table(datafile, header = TRUE, sep = ",", na.strings = "NA", dec = ".", strip.white = TRUE) result <- use_linear_mixed_model("extro", dataset = dataset, fix_eff = c("open", "agree", "social"), fix_eff_interact_vars = NULL, random_intercept_vars = c("school", "class"), nested_intercept_vars_pairs = list(c("school", "class")), cross_intercept_vars_pairs = NULL, uncorrel_slope_intercept_pairs = NULL, random_slope_intercept_pairs = NULL, package_mixed_model = NA)datafile <- system.file("extdata", "data_linear_mixed_model.csv", package = "packDAMipd") dataset = utils::read.table(datafile, header = TRUE, sep = ",", na.strings = "NA", dec = ".", strip.white = TRUE) result <- use_linear_mixed_model("extro", dataset = dataset, fix_eff = c("open", "agree", "social"), fix_eff_interact_vars = NULL, random_intercept_vars = c("school", "class"), nested_intercept_vars_pairs = list(c("school", "class")), cross_intercept_vars_pairs = NULL, uncorrel_slope_intercept_pairs = NULL, random_slope_intercept_pairs = NULL, package_mixed_model = NA)
########################################################################### Get the parameter values using the linear regression
use_linear_regression( param_to_be_estimated, dataset, indep_var, covariates, interaction )use_linear_regression( param_to_be_estimated, dataset, indep_var, covariates, interaction )
param_to_be_estimated |
parameter of interest |
dataset |
data set to be provided |
indep_var |
the independent variable (column name in data file) |
covariates |
list of covariates-calculations to be done before passing |
interaction |
boolean value to indicate interaction in the case of linear regression, false by default |
This function returns the results and plots after doing linear regression Requires param to be estimated, dataset, independent variables and information on covariates, and interaction variables if there are Uses form_expression_lm to create the expression as per R standard for e.g lm(y ~ x ). Returns the fit result,s summary results as returned by summary(), confidence interval for fit coefficients (ci_coeff), variance covariance matrix, cholesky decomposition matrix, Results from correlation test, plot of diagnostic tests and model fit assumptions, plot of model prediction diagnostic include AIC, R2, and BIC. The results of the prediction ie predicted values when each of covariate is fixed will be returned in prediction matrix predicted values will provide the mean value of param_to_to_estimated as calculated by the linear regression formula. ref:https://www.statmethods.net/stats/regression.html
the results of the regression analysis
results_lm <- use_linear_regression("dist", dataset = cars, indep_var = "speed", covariates = NA, interaction = FALSE) library(car) results_lm <- use_linear_regression("mpg", dataset = mtcars, indep_var = "disp", covariates = c("hp", "wt", "drat"), interaction = FALSE)results_lm <- use_linear_regression("dist", dataset = cars, indep_var = "speed", covariates = NA, interaction = FALSE) library(car) results_lm <- use_linear_regression("mpg", dataset = mtcars, indep_var = "disp", covariates = c("hp", "wt", "drat"), interaction = FALSE)
########################################################################### Get the parameter values using the survival analysis parametric survival
use_parametric_survival( param_to_be_estimated, dataset, indep_var, info_distribution, covariates, timevar_survival, cluster_var = NA )use_parametric_survival( param_to_be_estimated, dataset, indep_var, info_distribution, covariates, timevar_survival, cluster_var = NA )
param_to_be_estimated |
parameter of interest |
dataset |
data set to be provided |
indep_var |
the independent variable (column name in data file) |
info_distribution |
distribution name eg. for logistic regression -binomial |
covariates |
list of covariates |
timevar_survival |
time variable for survival analysis, default is NA |
cluster_var |
cluster variable for survival analysis |
This function is the last in the layer of function for parametric survival analysis. This then returns the parameters of interest, plots the results etc if the distribution is weibull it uses the package SurvRegCensCov for easy interpretation of results Returns the fit result, summary of regression, variance-covariance matrix of coeff, cholesky decomposition, the parameters that define the assumed distribution and the plot of model prediction Using survfit from survival package to plot the survival curve R's weibull distribution is defined as std weibull in terms of a and b as (a/b) (x/b)^ (a-1) exp((-x/b)^a) where a is the shape and b is the scale In HE the weibull distribution is parameterised as bit different it is like gamma.lambda. t^(gamma-1) .exp(-lambda*t^gamma) where gamma is the shape and lambda is the scale. The relationship is as below. HE_shape = rweibull_shape HE_scale = rweibull_scale ^(-rweibull_shape) The survreg shape and scale are again bit different and they are rweibull's shape and scale as below. rweibull_shape = 1/fit$scale rweibull_scale = exp(fit intercept)= exp(fit$coefficients) remember to use 1st of coefficients This has been utilised in SurvRegCensCov::ConvertWeibull predict() for survreg object with type =quantile will provide the failure times as survival function is 1-CDF of failure time.
the results of the regression analysis
data_for_survival <- survival::lung surv_estimated <- use_parametric_survival("status", data_for_survival, "sex", info_distribution = "weibull", covariates = c("ph.ecog"), "time")data_for_survival <- survival::lung surv_estimated <- use_parametric_survival("status", data_for_survival, "sex", info_distribution = "weibull", covariates = c("ph.ecog"), "time")
############################################################################ Get the parameter values using the survival analysis
use_survival_analysis( param_to_be_estimated, dataset, indep_var, info_get_method, info_distribution, covariates, timevar_survival, cluster_var = NA )use_survival_analysis( param_to_be_estimated, dataset, indep_var, info_get_method, info_distribution, covariates, timevar_survival, cluster_var = NA )
param_to_be_estimated |
parameter of interest |
dataset |
data set to be provided |
indep_var |
the independent variable (column name in data file) |
info_get_method |
additional information on methods e.g Kaplan-Meier ot hazard |
info_distribution |
distribution name eg. for logistic regression -binomial |
covariates |
list of covariates - calculations to be done before passing |
timevar_survival |
time variable for survival analysis, default is NA |
cluster_var |
cluster variable for survival analysis |
This function helps to get the parameter values after the survival analysis Takes into account many different methods like KM.FH, Cox proportional etc. and then calls appropriate functions to do the survival analysis
the results of the regression analysis
data_for_survival <- survival::aml surv_estimated_aml <- use_survival_analysis("status", data_for_survival, "x", info_get_method = "parametric", info_distribution = "weibull", covariates = NA, "time")data_for_survival <- survival::aml surv_estimated_aml <- use_survival_analysis("status", data_for_survival, "x", info_get_method = "parametric", info_distribution = "weibull", covariates = NA, "time")
utility matrix
utility_data.dfutility_data.df
A 11 by 2 dataframe
created on Nov 26, 2019 from tmat <- rbind(c(1, 2), c(3, 4)) colnames(tmat) <- rownames(tmat) <- c("Healthy", "Dead") tm <- transition_matrix(2, tmat, c(0.5, 0.5, 0, 1)) a <- health_state("Healthy", 1, 1, FALSE) b <- health_state("Dead", 1, 0, TRUE) health_states <- combine_state(a, b) this.strategy <- strategy(tm, health_states, "intervention")
Function to convert ADL scores to a T score
value_ADL_scores_IPD( ind_part_data, adl_related_words, adl_nrcode, adl_scoring_table = NULL )value_ADL_scores_IPD( ind_part_data, adl_related_words, adl_nrcode, adl_scoring_table = NULL )
ind_part_data |
a data frame containing IPD data |
adl_related_words |
related words to find out which columns contain adl data |
adl_nrcode |
non response code for ADL |
adl_scoring_table |
ADL scoring table, if given as NULL use the default one |
ADL scores converted to T score included modified data, if success -1, if failure
datafile <- system.file("extdata", "trial_data.csv", package = "packDAMipd") trial_data <- load_trial_data(datafile) value_ADL_scores_IPD(trial_data,c("tpi"), NA, adl_scoring_table = NULL)datafile <- system.file("extdata", "trial_data.csv", package = "packDAMipd") trial_data <- load_trial_data(datafile) value_ADL_scores_IPD(trial_data,c("tpi"), NA, adl_scoring_table = NULL)
Function to add EQ5D3L scores to IPD data
value_eq5d3L_IPD(ind_part_data, eq5d_nrcode)value_eq5d3L_IPD(ind_part_data, eq5d_nrcode)
ind_part_data |
a dataframe |
eq5d_nrcode |
non response code for EQ5D3L, default is NA |
qaly included modified data, if success -1, if failure
http://eprints.whiterose.ac.uk/121473/1/Devlin_et_al-2017-Health_Economics.pdf
datafile <- system.file("extdata", "trial_data.csv", package = "packDAMipd") trial_data <- load_trial_data(datafile) value_eq5d5L_IPD(trial_data, NA)datafile <- system.file("extdata", "trial_data.csv", package = "packDAMipd") trial_data <- load_trial_data(datafile) value_eq5d5L_IPD(trial_data, NA)
Function to add EQ5D5L scores to IPD data
value_eq5d5L_IPD(ind_part_data, eq5d_nrcode)value_eq5d5L_IPD(ind_part_data, eq5d_nrcode)
ind_part_data |
a dataframe |
eq5d_nrcode |
non response code for EQ5D5L, default is NA |
qaly included modified data, if success -1, if failure
http://eprints.whiterose.ac.uk/121473/1/Devlin_et_al-2017-Health_Economics.pdf
datafile <- system.file("extdata", "trial_data.csv", package = "packDAMipd") trial_data <- load_trial_data(datafile) value_eq5d5L_IPD(trial_data, NA)datafile <- system.file("extdata", "trial_data.csv", package = "packDAMipd") trial_data <- load_trial_data(datafile) value_eq5d5L_IPD(trial_data, NA)
Function to convert promis3a scores to a T score
value_promis3a_scores_IPD( ind_part_data, promis3a_related_words, promis3a_nrcode, promis3a_scoring_table = NULL )value_promis3a_scores_IPD( ind_part_data, promis3a_related_words, promis3a_nrcode, promis3a_scoring_table = NULL )
ind_part_data |
a data frame containing IPD data |
promis3a_related_words |
related words to find out which columns contain promis3a data |
promis3a_nrcode |
non response code for promis3a |
promis3a_scoring_table |
promis3a scoring table, if given as NULL use the default one |
promis3a scores converted to T score included modified data, if success -1, if failure
trial_data <- data.frame("tpi.q1" = c(1, 2), "tpi.q2" = c(1, 2), "tpi.q3" = c(1, 2)) results <- value_promis3a_scores_IPD(trial_data, c("tpi"), NA, NULL)trial_data <- data.frame("tpi.q1" = c(1, 2), "tpi.q2" = c(1, 2), "tpi.q3" = c(1, 2)) results <- value_promis3a_scores_IPD(trial_data, c("tpi"), NA, NULL)
Function to estimate the cost of tablets taken (from IPD)
value_Shows_IPD(ind_part_data, shows_related_words, shows_nrcode)value_Shows_IPD(ind_part_data, shows_related_words, shows_nrcode)
ind_part_data |
a data frame containing IPD |
shows_related_words |
a data frame containing IPD |
shows_nrcode |
non response code for ADL, default is NA |
sum of scores, if success -1, if failure
datafile <- system.file("extdata", "trial_data.csv", package = "packDAMipd") trial_data <- load_trial_data(datafile) value_Shows_IPD(trial_data, "qsy", NA)datafile <- system.file("extdata", "trial_data.csv", package = "packDAMipd") trial_data <- load_trial_data(datafile) value_Shows_IPD(trial_data, "qsy", NA)
Function to check the variable null or NA
word2num(word)word2num(word)
word |
word for the number |
https://stackoverflow.com/questions/18332463/convert-written-number-to-number-in-r
return the number
answer <- word2num("one forty one") answer <- word2num("forty one and five hundred") answer <- word2num("five thousand two hundred and eight")answer <- word2num("one forty one") answer <- word2num("forty one and five hundred") answer <- word2num("five thousand two hundred and eight")