In skgallagher/loewesadditivity: Loewe's Additivity

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.align = "center"
)
library(loewesadditivity)

Installation

To install loewesadditivity, run the following code.

if(!require(loewesadditivity)){
  library(devtools)
  devtools::install_github("shannong19/loewesadditivity")
}
library(loewesadditivity)

# hidden load currently
#library(devtools)
devtools::load_all()

Quickstart

library(ggplot2)
library(dplyr)
library(knitr)
library(kableExtra)

results <- rh5_ama1ron2 %>% 
  rename(dose_A = "RH5", dose_B = "AMA1RON2")  %>% 
  fortify_gia_data() %>%
  estimate_params() 

results$params_est

Look at the raw data

This is a raw set of GIA data that has been collected and uploaded to this package. We measure the GIA% for different experiments, plates (a or b) and repetitions. We are using the compounds AMA1RON2 and RH5.

The variables iRBC and uRBC are optional rows which measure the minimal and maximum GIA respectively.

data(rh5_ama1ron2)

rh5_ama1ron2 %>% head() %>% knitr::kable() %>%
  kableExtra::kable_styling(full_width = T, font_size = 7)

Fortify data

The standard output we use is that of the long, data.frame with the columns

1. dose_A
2. dose_B
3. GIA

where dose_A and dose_B are the unscaled doses (e.g. mg/mL) and GIA is the percent GIA.

We provide the function fortify_gia_data() which works when the data is in the above format with column headings of exp{x}{y}rep{z} along with the two doses (must have names dose_A and dose_B.

rh5_ama1ron2_f <-
  rh5_ama1ron2 %>% rename(dose_A = 'RH5', dose_B = 'AMA1RON2') %>%
  fortify_gia_data()

  rh5_ama1ron2_f %>% head() %>%  kable(align = "c") %>%
  kableExtra::kable_styling(full_width = F)

The fortification has averaged over the repetitions but has kept the plates and the experiments.

Estimate parameters

You can use the function estimate_params() to estimate the parameters. The parameters and their meaning can be found in the documentation of ?base_GIA.

## Set up initial guesses and parameters
model_params <- c(
  "beta_A" = .5,
  "beta_B" = .5,
  "gamma_A" = .5,
  "gamma_B" = .5,
  "tau_1" = 0,
  "tau_2" = 0
  )
  n_boot <- 10
  fn_list <- NULL
  alpha <- .05
  verbose <- TRUE


out <- estimate_params(
    data = rh5_ama1ron2_f,
    init_params = model_params,
    n_boot = n_boot,
    alpha = alpha,
    verbose = verbose
  )

Hewlett's S can be found below.

out$S_est  %>%  kable(align = "c", digits = 3) %>%
  kableExtra::kable_styling(full_width = F)

The parameter estimates from the model are shown here:

out$params_est %>%  kable(align = "c", digits = 3) %>%
  kableExtra::kable_styling(full_width = F)

Finally, the (first six) individual observation estimates and 95% CIs can be found here:

out$GIA_est %>% head()  %>%  kable(align = "c", digits = 3) %>%
  kableExtra::kable_styling(full_width = F)

Plot results

We can plot the estimated surface, individual curves, and the isobologram.

# Note xlab is dose_A and ylab is dose_B
g1 <- plot_surface(out, xlab = "RH5", ylab = "AMA1RON2")

g2 <- plot_curves(out, dose_A = "RH5", dose_B = "AMA1RON2" )

g3 <- plot_isobologram(out, dose_A = "RH5", dose_B = "AMA1RON2" )

Simulation of coverage

We provide the function simulate_coverage() as a way to get an estimate of the power of the experiment to detect synergy or antagonism. For given model parameters, an experimental grid of dose combinations, and an assumed noise structure for GIA, we can determine 1) percent of times we expect 0 to in the 95% CI of the interaction parameter $\tau_1$ and 2) percent of times we expect the true given value of the model parameters to be in the 95% CI.

WARNING It is advised to use at least 100 bootstraps for each of 10 simulations but is not shown here.

  model_params <- c("beta_A" = .250, "beta_B" = .146,
                    "gamma_A" = .527, "gamma_B" = .921,
                    "tau_1" = -.058, "tau_2" = -.289)
  experimental_grid <- make_grid(par = model_params, 
                                 n = 5)
  n_boot <- 3
  n_sims <- 5
  GIA_fn <- base_GIA
  S_fn <- calc_S_base
  fn_list <- NULL
  alpha <- .05
  verbose <- TRUE
  out <- simulate_coverage(n_sims = n_sims,
                         n_boot = n_boot,
                         verbose = FALSE,
                         experimental_grid = experimental_grid,
                         model_par = model_params,
                         alpha = .05,
                         noise_par = c("a0" = 3, "a1" = .01),
                         GIA_fn = base_GIA,
                         fn_list = NULL)

  out

Design an experiment

A unique feature of loewesadditivity is the ability to generate code to use to create your own coverage simulation in conjunction with a (hopefully) intuitive shiny App. Run the below code, copy and paste the results into R and see what happens!

Ex.

 out <- design_experiment(n_rep = 2)

skgallagher/loewesadditivity documentation built on April 6, 2021, 1:49 a.m.