In maomlab/BayesPharma: Tools for Bayesian Analysis of Non-Linear Pharmacology Models

#| echo=FALSE

options(rmarkdown.html_vignette.check_title = FALSE)

# inspired by https://www.jumpingrivers.com/blog/knitr-default-options-settings-hooks/
knitr::opts_chunk$set(
  cache = TRUE,
  echo = FALSE,
  fig.path = "implement_tidymodels_files/",
  fig.retina = 2, # Control using dpi
  fig.width = 6,  # generated images
  fig.height = 5, # generated images
  fig.pos = "t",  # pdf mode
  fig.align = "center",
  dpi = if (knitr::is_latex_output()) 72 else 300, 
  out.width = "100%")

library(here)
suppressWarnings(suppressMessages(library(tidymodels)))
library(BayesPharma)
suppressWarnings(suppressMessages(library(bayesplot)))
library(bayesian)

bayesplot::bayesplot_theme_set(
  new = ggplot2::theme_bw())

set.seed(0)

Using BayesPharma models with Tidymodels

Tidymodels is a platform for running and analyzing models in R. This vignette demonstrates how to use the BayesPharma models with Tidy Models.

#| echo=TRUE
load(here::here("data", "kor_antag.rda"))
kor_antag <- kor_antag |> dplyr::transmute(
  substance_id = Drug,
  cell_id = cell_number,
  log_dose = log_dose,
  response = normalized_measurement)

Define a recipe and assign roles to different columns

#| echo=TRUE

kor_recipe <- kor_antag |>
  recipes::recipe() |>
  recipes::update_role(response, new_role = "outcome") |>
  recipes::update_role(log_dose, new_role = "predictor") |>
  recipes::update_role(substance_id, new_role = "predictor")

Define the model using the \code{bayesian} wrapper package for Tidymodels.

#| echo=TRUE
kor_model <- bayesian::bayesian(
  cores = 4,
  chains = 4,
  iter = 2000,
  control = list(adapt_delta = 0.99),
  stanvars = BayesPharma::sigmoid_stanvar(),
  init = BayesPharma::sigmoid_antagonist_init(top = 100),
  prior = BayesPharma::sigmoid_antagonist_prior(top = 100)) |>
  parsnip::set_engine("brms") |>
  parsnip::set_mode("regression") |>
  recipes::update(
    formula.override = bayesian::bayesian_formula(
      ..y ~ sigmoid(ic50, hill, top, bottom, log_dose),
      ic50 + hill + top + bottom ~ 0 + substance_id,
      nl = TRUE))

Build and run the workflow by adding the recipe and model and then running fit.

#| echo=TRUE
kor_workflow <- workflows::workflow() |>
  workflows::add_recipe(kor_recipe) |>
  workflows::add_model(spec = kor_model) |>
  parsnip::fit(data = kor_antag)

kor_workflow

maomlab/BayesPharma documentation built on Aug. 24, 2024, 8:45 a.m.