Confidence intervals and hypothesis tests for binary choice data with repeated measures

I have data with two between-subjects factors (each with two levels; bs1 has levels ‘hi’ and ‘lo’, and bs2 has levels ‘happy’ and ‘sad’) and one within-subjects factor (ws1, with levels ‘good’, ‘neutral’, and ‘bad’). Participants complete two trials in each within-subjects condition. That is, participants are split into four groups and then complete six trials across the three within-subjects conditions. I measure a binary outcome variable that can be interpreted as success or failure. I thus have data that looks like this:

id <- rep(1:4, each = 6)
trial <- rep(1:6, 4)
bs1 <- rep(c('hi', 'lo'), each = 6, times = 2)
bs2 <- rep(c('happy', 'sad'), each = 6, times = 2)
ws1 <- rep(c('good', 'neutral', 'bad'), each = 2, times = 4)
choice <- rbinom(24, 1, 0.5)
data <- as.data.frame(cbind(id, trial, bs1, bs2, ws1, choice))


Question 1

I want to plot the proportion of successes by condition with 95% confidence intervals, but I am not sure of the appropriate way to generate these confidence intervals. If I were only observing participants once in each condition, I would use binom.test() in R; however (at least as far as I can tell), this package does not allow me to specify that measurements are repeated. Can anyone suggest a method for calculating these intervals?

Question 2

I also want to test a series of hypotheses on this dataset. To do this I specify a series of generalized estimating equations models using the R package ‘geepack’. This approach takes the nested structure of my data into account. I am using GEE rather than a mixed effects model because my primary interest is in population average effects rather than individual-specific effects (as per When to use generalized estimating equations vs. mixed effects models?).

Hypothesis 1 tests whether there is an effect of the within-subjects manipulation, collapsing across the between-subjects conditions. I hypothesise that the proportion of successes will increase as we move from the ‘bad’ to the ‘neutral’ to the ‘good’ condition (i.e., the probability of success by condition is ‘bad’ < ‘neutral’ < ‘good’).

I thus specify the model

geeM1 <- geeglm(choice ~ ws1, data = data, id = id, family = binomial,
corstr = "exchangeable")


with the condition ‘neutral’ as the reference and examine whether the coefficients are significantly positive for the ‘good’ term and significantly negative fro the ‘bad’ term.

Hypothesis 2 tests for an interaction between bs1 and ws1 levels ‘neutral’ and ‘good’. I hypothesise that the difference in the proportion of successes between the ‘neutral’ and ‘good’ conditions will be greater in the bs1 ‘hi’ condition than in the bs1 ‘lo’ condition. That is,

(p(success) | bs1 ‘hi’ * ws1 ‘good’) - (p(success) | bs1 ‘hi’ * ws1 ‘neutral’) > (p(success) | bs1 ‘lo’ * ws1 ‘good’) - (p(success) | bs1 ‘lo’ * ws1 ‘neutral’).

I thus specify the model

geeM2 <- geeglm(choice ~ ws1 * bs1, data = data, id = id, family =
binomial, corstr = "exchangeable")


(with ws1 condition ‘neutral’ and bs1 condition ‘lo’ as the references) and examine whether the coefficient for the interaction term ws1 ‘good’ * bs1 ‘hi’ is positive and significantly different from zero.

Hypothesis 3 tests for an interaction between bs1 and bs2 within ws1 level ‘good’. I hypothesise that, for those in the bs1 ‘lo’ condition, being in the bs2 ‘happy’ group will increase the proportion of successes relative to those in the bs2 ‘sad’ group. For those in the bs1 ‘hi’ condition, the proportion of success will be similar across the two bs2 groups. That is,

(p(success) | bs1 ‘lo’ * bs2 ‘happy’ * ws1 ‘good’) - (p(success) | bs1 ‘lo’ * bs2 ‘happy’ * ws1 ‘good’) > (p(success) | bs1 ‘hi’ * bs2 ‘happy’ * ws1 ‘good’) - (p(success) | bs1 ‘hi’ * bs2 ‘happy’ * ws1 ‘good’).

I thus take a subset of my data

dataSub <- subset(data, ws1 == ‘good’)


and specify the model

geeM3 <- geeglm(choice ~ bs1 * bs2, data = dataSub, id = id, family =
binomial, corstr = "exchangeable")


(with bs1 condition ‘lo’ and bs2 condition ‘sad’ as the references) and examine whether the coefficient for the interaction term bs1 ‘hi’ * bs2 ‘happy’ is positive and significantly different from zero.

It strikes me that there are other ways to test these hypotheses. As such, I’m wondering whether the GEE approach that I have taken is appropriate.