4
$\begingroup$

Imagine that I have a treatment that reduces the likelihood of response to a stimulus. This could be anything you like, but the simplest example is of a treatment (e.g., hand washing, mask wearing, etc.) that prevents disease when exposed.

For simplicity, I built an (overly complicated) model of the risk of responding to the stimulus over time under the two treatments (all of this is done in R using the tidyverse):

day_prob <-
  list(
    A = c(0.06, 0.06, pmax(dchisq(1:13, 3)*1, 0.06))
    , B = c(0.06, 0.06, pmax(dchisq(seq(1, 25, 2), 5)*6, 0.06))
  )

This can be plotted to show the risk:

day_prob %>%
  lapply(function(x){
    tibble(
      Day = 1:length(x)
      , Prob = x
    )
  }) %>%
  bind_rows(.id = "Group") %>%
  ggplot(aes(x = Day
             , y = Prob
             , col = Group)) +
  geom_line() +
  geom_point() +
  scale_color_brewer(palette = "Dark2")

Plot:

enter image description here

Note that there is a similar baseline risk for both group (0.06), some latency to develop after exposure (risk rises on day 3), and that risk falls over time back to baseline.

Now, assuming that I randomize individuals into the two treatments, what is my best approach to identify this effect? I can just analyze each day separately, though that raises some repeated-measures questions since (unlike the sample data below) it is likely that a positive individual is more (or even less) likely to test positive on subsequent days.

I've searched through a number of alternative questions, but nothing seems to quite capture this issue. I could try repeated-measures analyses, but the latency and return-to-baseline should both swamp my effect if they are sufficiently long. Further, it would be ideal to actually know when the effect is observed.

A sample data set:

make_obs <- function(group, day){
  sapply(1:length(group), function(idx){
    rbinom(1, 1, day_prob[[group[idx]]][day[idx]] )
  })
}

set.seed(12345)
example_data <-
  tibble(
    Group = rep(c("A", "B"), each = 50)
    , Ind = 1:100
    , Day = 1
  ) %>%
  complete(nesting(Group, Ind), Day = 1:15) %>%
  mutate(
    Obs = make_obs(Group, Day)
  )

Plotting to show the observed data:

example_data %>%
  group_by(Group, Day) %>%
  summarise(Prop_pos = mean(Obs)) %>%
  ungroup() %>%
  ggplot(aes(x = Day
             , y = Prop_pos
             , col = Group)) +
  geom_line() +
  scale_color_brewer(palette = "Dark2")

enter image description here

$\endgroup$
  • $\begingroup$ The chi-squared call you use in the simulated data - does that represent an actual model of the outcome, or is it just a function that produces good-looking results? If the latter is the case, would you be able to describe the response to treatment as a parametric curve and/or differential equation? Vaguely reminds me of PK/PD (pharmacokinetic/pharmacodynamic) models - is there a connection? $\endgroup$ – Martin Modrák Jun 3 at 19:58
  • $\begingroup$ The chi-square model was simply handy to produce roughly the pattern that I was looking for (after lots of tinkering with df etc.) -- it does not describe any particular underlying result. I am fine using any underlying model; the output data are probably the more relevant, but I thought that it would be instructive to include the generative model. There is no direct connection of PK/PD to the model I have in mind, though it could explain the pattern in similar cases. $\endgroup$ – Mark Peterson Jun 3 at 20:52

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.