From a5cd79382e4737dfb9bc2a71f77bb47c2c85126a Mon Sep 17 00:00:00 2001 From: Andersen Date: Mon, 9 Oct 2023 16:50:06 +0200 Subject: [PATCH] styled vignette! --- vignettes/analysis_normal.Rmd | 108 +++++++++++++++++++--------------- 1 file changed, 62 insertions(+), 46 deletions(-) diff --git a/vignettes/analysis_normal.Rmd b/vignettes/analysis_normal.Rmd index 60069f9..ccf3f91 100644 --- a/vignettes/analysis_normal.Rmd +++ b/vignettes/analysis_normal.Rmd @@ -26,16 +26,17 @@ In a first step a meta analytic prior will be calculated using historical data. ```{r} data("metaData") testdata <- as.data.frame(metaData) -dataset<-filter(testdata, bname == "VIAGRA") -histcontrol<-filter(dataset, dose==0, primtime==12,indication=="ERECTILE DYSFUNCTION") +dataset <- filter(testdata, bname == "VIAGRA") +histcontrol <- filter(dataset, dose == 0, primtime == 12, indication == "ERECTILE DYSFUNCTION") -##Create MAP Prior +## Create MAP Prior hist_data <- data.frame( - trial=c(1,2,3,4), - est = histcontrol$rslt, - se = histcontrol$se, - sd = histcontrol$sd, - n = histcontrol$sampsize) + trial = c(1, 2, 3, 4), + est = histcontrol$rslt, + se = histcontrol$se, + sd = histcontrol$sd, + n = histcontrol$sampsize +) sd_tot <- with(hist_data, sum(sd * n) / sum(n)) @@ -47,39 +48,46 @@ gmap <- RBesT::gMAP( family = gaussian, beta.prior = cbind(0, 100 * sd_tot), tau.dist = "HalfNormal", - tau.prior = cbind(0, sd_tot / 4)) + tau.prior = cbind(0, sd_tot / 4) +) -prior_ctr <- RBesT::robustify( +prior_ctr <- RBesT::robustify( priormix = RBesT::automixfit(gmap), weight = 0.5, mean = 1.4, - sigma = sd_tot) + sigma = sd_tot +) -#RBesT::ess(prior_ctr) +# RBesT::ess(prior_ctr) ## derive prior for treatment ## weak informative using same parameters as for robustify component prior_trt <- RBesT::mixnorm( comp1 = c(w = 1, m = 1.4, n = 1), sigma = sd_tot, - param = "mn") + param = "mn" +) dose_levels <- c(0, 50, 100, 200) ## combine priors in list prior_list <- c(list(prior_ctr), rep(list(prior_trt), times = length(dose_levels[-1]))) ``` # ii) Pre-Specification of candidate models ```{r} -#Pre-Specification (B)MCPMod +# Pre-Specification (B)MCPMod ## candidate models for MCPMod # linear function - no guestimates needed -exp <- DoseFinding::guesst(d = 50, - p = c(0.2), - model = "exponential", - Maxd = max(dose_levels)) -emax <- DoseFinding::guesst(d = 100, - p = c(0.9), - model = "emax") +exp <- DoseFinding::guesst( + d = 50, + p = c(0.2), + model = "exponential", + Maxd = max(dose_levels) +) +emax <- DoseFinding::guesst( + d = 100, + p = c(0.9), + model = "emax" +) mods <- DoseFinding::Mods( @@ -88,46 +96,50 @@ mods <- DoseFinding::Mods( exponential = exp, doses = dose_levels, maxEff = 10, - placEff = 1.4) + placEff = 1.4 +) ``` # iii) Trial results ```{r} -#Simulation of new trial -##Note: This part will be simplified and direct results from one trial will be used +# Simulation of new trial +## Note: This part will be simplified and direct results from one trial will be used mods_sim <- DoseFinding::Mods( emax = emax, doses = dose_levels, maxEff = 12, - placEff = 1.4) + placEff = 1.4 +) -n_patients<-c(50,50,50,50) +n_patients <- c(50, 50, 50, 50) data <- simulateData( n_patients = n_patients, dose_levels = dose_levels, sd = sd_tot, mods = mods_sim, - n_sim = 1) + n_sim = 1 +) -data_emax <- data[, c("simulation", "dose", "emax")] +data_emax <- data[, c("simulation", "dose", "emax")] names(data_emax)[3] <- "response" posterior_emax <- getPosterior( data = data_emax, - prior_list = prior_list) - + prior_list = prior_list +) ``` # iv) Execution of Bayesian MCPMod Test step ```{r} -#Evaluation of Bayesian MCPMod +# Evaluation of Bayesian MCPMod contr_mat <- DoseFinding::optContr( models = mods, doses = dose_levels, - w = n_patients) -##Calculation of critical value can be done with critVal function -crit_val_equal<-DoseFinding:::critVal(contr_mat$corMat, alpha=0.05, df=0,alternative = "one.sided") -crit_pval <-pnorm(crit_val_equal) + w = n_patients +) +## Calculation of critical value can be done with critVal function +crit_val_equal <- DoseFinding:::critVal(contr_mat$corMat, alpha = 0.05, df = 0, alternative = "one.sided") +crit_pval <- pnorm(crit_val_equal) ess_prior <- round(unlist(lapply(prior_list, RBesT::ess))) @@ -135,12 +147,14 @@ ess_prior <- round(unlist(lapply(prior_list, RBesT::ess))) contr_mat_prior <- DoseFinding::optContr( models = mods, doses = dose_levels, - w = n_patients + ess_prior) + w = n_patients + ess_prior +) BMCP_result <- BMCPMod( posteriors_list = list(posterior_emax), - contr_mat = contr_mat_prior, - crit_prob = crit_pval) + contr_mat = contr_mat_prior, + crit_prob = crit_pval +) BMCP_result ``` @@ -148,10 +162,10 @@ BMCP_result # v) Model fitting and Plotting In the model fitting step the posterior distribution is used as basis. ```{r} -#Model fit -#This part is currently not working +# Model fit +# This part is currently not working post_observed <- posterior_emax -model_shapes <- c("linear", "emax", "exponential") +model_shapes <- c("linear", "emax", "exponential") # Option a) Simplified approach by using frequentist idea @@ -159,18 +173,20 @@ fit_simple <- getModelFits( models = model_shapes, dose_levels = dose_levels, posterior = post_observed, - simple = TRUE) + simple = TRUE +) # Option b) Making use of the complete posterior distribution -fit<- getModelFits( +fit <- getModelFits( models = model_shapes, dose_levels = dose_levels, posterior = post_observed, - simple = FALSE) + simple = FALSE +) ``` ```{r} -plot_modelFits(fit,CrI= TRUE) -plot_modelFits(fit_simple, CrI= TRUE) +plot_modelFits(fit, CrI = TRUE) +plot_modelFits(fit_simple, CrI = TRUE) ``` For the plotting the credible intervals are shown as well.