Skip to content

Commit

Permalink
Vignette updated.
Browse files Browse the repository at this point in the history
  • Loading branch information
@sebastianbossert
sebastianbossert committed Dec 18, 2023
1 parent 01a1d24 commit f852df6
Show file tree
Hide file tree
Showing 2 changed files with 90 additions and 49 deletions.
4 changes: 3 additions & 1 deletion R/modelling.R
View file
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
#' @title getModelFits
#'
#' @description Fits dose-response curves for the specified dose-repsonse models, based on the posterior distributions.
#' @description Fits dose-response curves for the specified dose-response models, based on the posterior distributions.
#' For the simplified fit, multivariate normal distributions will be approximated and reduced by one-dimensional normal distributions.
#' For the default case, the Nelder-Mead algorithm is used. Will be further updated and links to publication as well as references will be added.
#'
Expand All @@ -17,6 +17,7 @@ getModelFits <- function (
models,
dose_levels,
posterior,
#avg_fit=FALSE, if possible we should add the average fit directly here
simple = FALSE

) {
Expand All @@ -32,6 +33,7 @@ getModelFits <- function (
attr(model_fits, "posterior") <- posterior
class(model_fits) <- "modelFits"


return (model_fits)

}
Expand Down
135 changes: 87 additions & 48 deletions vignettes/analysis_normal.Rmd
View file
Original file line number Diff line number Diff line change
Expand Up @@ -28,12 +28,12 @@ In this vignette we will show the use of the Bayesian MCPMod package for continu
Hereby the focus is on the utilization of an informative prior and the BayesianMCPMod evaluation of a single trial.
We will use data that is included in the clinDR package.
More specifically trial results for BRINTELLIX will be used to illustrate the specification of an informative prior and the usage of such a prior for the bayesian evaluation of a (hypothetical) new trial.
More information around BRINTELLIX will be added...
BRINTELLIX is a medication used to treat major depressive disorder. Various clinical trials with different dose groups (including a control group) were conducted.


# Calculation of a MAP prior
In a first step a meta analytic prior will be calculated using historical data from 4 trials (with main endpoint CfB in MADRS score after 8 weeks).
Please note that only information from the control group will be integrated (leading to an informative multi-component prior for the control group), while for the active groups a non-informative prior will be specified.
In a first step a meta analytic prior will be calculated using historical data from 4 trials (with main endpoint Change from baseline in MADRS score after 8 weeks).
Please note that only information from the control group will be integrated (leading to an informative mixture prior for the control group), while for the active groups a non-informative prior will be specified.
```{r Historical Data for Control Arm}
data("metaData")
dataset <- filter(as.data.frame(metaData), bname == "BRINTELLIX")
Expand Down Expand Up @@ -113,10 +113,13 @@ dose_levels <- c(0, 2.5, 5, 10)
prior_list <- getPriorList(
hist_data = hist_data,
dose_levels = dose_levels)
dose_levels = dose_levels,
robust_weight = 0.3)
getESS(prior_list)
```
# Specifications new trial
We will use the trial with ct.gov number NCT00635219 as exemplary new trial.
We will use the trial with ct.gov number NCT00735709 as exemplary new trial.
To be able to apply the bayesian MCPMod approach, candidate models need to be specified.
Since there are only 3 active dose levels we will limit the set of candidate models to a linear, an exponential, and an emax model.
Note that the linear candidate model does not require a guesstimate.
Expand Down Expand Up @@ -144,66 +147,76 @@ new_trial <- filter(dataset,
indication == "MAJOR DEPRESSIVE DISORDER",
protid == 5)
n_patients <- c(150, 142, 147, 149)
n_patients <- c(128, 124, 129, 122)
```

# Combination of prior information and trial results

As outlined in citePaper, in a first step the posterior is calculated combining the prior information with the estimated results of the new trial.
As outlined in Fleischer et al. [Bayesian MCPMod. Pharmaceutical Statistics. 2022; 21(3): 654-670.], in a first step the posterior is calculated combining the prior information with the estimated results of the new trial. Via the summary function it is possible to print out the summary information of the posterior distributions.

The simulation takes the variability of the historical data into account.

```{r Trial results}
#Simulation of new trial
##Note: This step should not be required, as we provide summary measures directly from the new trial
data_emax <- simulateData(
n_patients = n_patients,
dose_levels = dose_levels,
sd = with(hist_data, sum(sd * n) / sum(n)),
mods = mods,
true_model = "emax")
posterior <- getPosterior(prior = prior_list,
data = data_emax)
posterior <- getPosterior(prior = prior_list,
mu_hat = new_trial$rslt,
se_hat = new_trial$se)
se_hat = new_trial$se,
calc_ess = TRUE)
summary(posterior)
```

# Execution of Bayesian MCPMod Test step
For the execution of the testing step of Bayesian MCPMod a critical value (on the probability scale) will be determined for a given alpha level. In addition a contrast matrix is created. Please note that there are different possibilities how to generate contrasts.
This information is then used as input for the Bayesian MCP testing function.
For the execution of the testing step of Bayesian MCPMod a critical value (on the probability scale) will be determined for a given alpha level. This critical value is calculated using (the re-estimated) contrasts for the frequentist MCPMod to ensure that when using non-informative priors the actual error level for falsely declaring a significant trial in the Bayesian MCPMod is controlled (by the specified alpha level).

A pseudo-optimal contrast matrix is generated based on the variability of the posterior distribution (see Fleischer et al. 2022 for more details). Please note that there are different possibilities how to establish the contrasts.



```{r Execution of Bayesian MCPMod Test step}
```{r Preparation of input for Bayesian MCPMod Test step}
crit_pval <- getCritProb(
mods = mods,
dose_levels = dose_levels,
dose_weights = n_patients,
alpha_crit_val = 0.1)
mods = mods,
dose_levels = dose_levels,
se_new_trial = new_trial$se,
alpha_crit_val = 0.05)
contr_mat<- getContr(
mods = mods,
dose_levels = dose_levels,
sd_posterior = summary(posterior)[,2])
#Alternative ways to come up with a contrast would be to use
#i) the frequentist contrast
#contr_mat_prior <- getContr(
# mods = mods,
# dose_levels = dose_levels,
# dose_weights = n_patients,
# prior_list = prior_list)
# ii) re-estimated frequentist contrasts
# contr_mat_prior <- getContr(
# mods = mods,
# dose_levels = dose_levels,
# se_new_trial = new_trial$se)
# iii) Bayesian approach using number of patients for new trial and prior distribution
#contr_mat_prior <- getContr(
# mods = mods,
# dose_levels = dose_levels,
# dose_weights = n_patients,
# prior_list = prior_list)
#This would be the most reasonable output, but since it is not exported it is currently not usable.
#BMCP_result<-BayesMCPi(posterior_i = posterior,
# contr_mat = contr_mat_prior,
# crit_prob = crit_pval)
contr_mat_prior <- getContr(
mods = mods,
dose_levels = dose_levels,
dose_weights = n_patients,
prior_list = prior_list)
#This would be the most reasonable output, but since it is not exported it is currently not usable.
#BMCP_result<-BayesMCPi(posterior_i = posterior,
# contr_mat = contr_mat_prior,
# crit_prob = crit_pval)
BMCP_result <- performBayesianMCP(
posterior_list = posterior,
contr = contr_mat_prior,
crit_prob_adj = crit_pval)
BMCP_result
#BMCP_result2 <- performBayesianMCPMod(
# posteriors_list = posterior_emax,
Expand All @@ -212,13 +225,25 @@ BMCP_result
#BMCP_result2
```
The testing step is significant indicating a non-flat dose-response shape. In detail the p-values for the emax model is the most significant one.
The bayesian MCP testing step is then executed based on the posterior information, the provided contrasts and the (multiplicity adjusted) critical value.
```{r Executionr of Bayesian MCPMod Test step}
BMCP_result <- performBayesianMCP(
posterior_list = posterior,
contr = contr_mat,
crit_prob_adj = crit_pval)
BMCP_result
```

The testing step is significant indicating a non-flat dose-response shape. In detail, all 3 models are significant and the p-value for the emax model is the most significant one.

# Model fitting and visualization of results
In the model fitting step the posterior distribution is used as basis. Both simplified and full fitting are performed.
Please note that all models are fitted to the data even though they were not significant in the testing step.
For the plotting bootstrap based credible intervals (for 50% and 95%) are shown as well.

In the model fitting step the posterior distribution is used as basis. Both simplified and full fitting are performed.
For the simplified fit the multivariate normal distribution of the control group is approximated and reduced by a one-dimensional normal distribution.The actual fit (on this approximated posterior distribution) is then performed using generalized least squares criterion.
For the full fitting the (multi-dimensional) non-linear optimization problem is adressed via the Nelder Mead algorithm.

The output of the fit includes information about the predicted effects for the included dose levels the generalized AIC (and the corresponding weights). For the considered case the simplified and the full fit are very similar.
```{r Model fitting}
#Model fit
post_observed <- posterior
Expand All @@ -237,21 +262,35 @@ fit <- getModelFits(
dose_levels = dose_levels,
posterior = post_observed,
simple = FALSE)
```

Via the predict function one can also receive estimates for dose levels that were not included in the trial.

It is possible to plot the fitted dose response models and an AIC based average model (black lines). To assess the uncertainty one can in addition visualize credible intervals (yellow shaded areas, the default is set to 50% and 95%). These credible intervals are calculated as follows.
Samples from the posterior distribution are drawn and for every sample the simplified fitting step and a prediction is performed. These fits are then used to identify and visualize the specified quantiles.
The bootstrap based quantiles can also directly be calculated and displayed via the gotbootstrapQuantiles function.


```{r Visualization and interpretation of the results}
predict(fit, doses = c(0, 2.5,4, 5,7, 10))
plot(fit, cr_bands = TRUE)
plot(fit_simple, cr_bands = TRUE)
predict(fit, doses = c(4, 7))
getBootstrapQuantiles(fit, quantiles = c(0.025, 0.975), doses = c(4, 7))
getBootstrapQuantiles(fit, quantiles = c(0.025,0.5, 0.975), doses = c(0, 2.5,4, 5,7, 10))
```

Technical note: The median quantile of the bootstrap based procedure is not necessary similar to the main model fit, as they are derived via different procedures (main fit, i.e. black lines in the plot, show the best fit of a certain model based on minimizing the residuals for the posterior distribution while bootstrap based 50% quantile indicates the median fit of the random sampling and fitting procedure).

# Additional notes
It is also possible to perform the testing and modelling step in combined fashion via the performBayesianMCPMod function
It is also possible to perform the testing and modelling step in combined fashion via the performBayesianMCPMod function.
```{r}
performBayesianMCPMod(
posterior_list = posterior,
contr = contr_mat_prior,
crit_prob_adj = crit_pval)
#performBayesianMCPMod(
# posterior_list = posterior,
# contr = contr_mat,
# crit_prob_adj = crit_pval)
```


Expand Down

0 comments on commit f852df6

Please sign in to comment.