Check if dropout rates are independent for an interaction of two independent variables (one with a large amount of levels) I am trying to analyse dropout rates in an experiment, but there are multiple issues which collide, and I don't know how to deal with them as a whole. First, find a list of those issues. Below, see a more elaborate description and a fictional example for illustration (which involves cakes!). I am using R for my analysis, but appreciate any theoretical input.
Issues


*

*A combination of two independent variables (stimulus $\times$ condition).


*

*... where one of the variables has a large number of levels (100 stimuli).

*... where there might be (?) dependencies: The same stimuli are used in both conditions, even though there are no repeated measures.


*In some cells, the count is zero.

*I have less participants (i.e. maximum dropout rates/cell counts) than stimuli.


What I thought of so far


*

*A $\chi^2$ test, as I am dealing with frequencies - but that doesn't cover the interaction between the two idependent variables.

*I found this tutorial using the Titanic dataset in R which seems to cover my problem and uses a log-linear model to solve it.


*

*However: The large amount of levels, zero counts and possibly (?) dependencies would still be an issue?


*I also found another tutorial heading in the same direction. Here, the loglm() function (MASS package) and the glm() function are used in R.


*

*However, I am unsure if a generalised linear model (glm() using family = poisson(link = "log") is the same as loglm() and if those are fit for my case.


*The nature of the stimuli (many of them; used in both conditions) made me think of including them as a random effect, but that seems wrong because it is my main interest to find out if the dropout rate differs across stimuli - which seems more like a fixed effect to me, but I might be on the wrong track here. In later analyses (e.g. for reaction times), my goal is to include the stimuli as random effects, but the dropout question seems to have a different nature.


I thought I'd go with a fully saturated log-linear model and then remove the interactions one by one to figure out if the dropout is independent of stimulus, condition or the interaction of those. But it appears to me that some assumptions are violated in my case.
Toy example (with code)
I created a toy data set which you can find on GitHub. Let's assume a cake experiment, where 40 participants taste 20 different cakes (stimuli) in two conditions: blindfolded or not. The experiment is between subjects: Participants each get to try 10 cakes blindfolded and 10 without a blindfold. To that end, cakes are divided into two sets (A and B) with 10 cakes each. It's always the same cakes in the same set. Each participant gets to try set A blindfolded and set B without a blindfold - or vice versa (pseudo-randomised across participants).
What I want to know
Does the number of dropouts (i.e. people who get sick and break up the experiment) differ?
a) Most importantly: Across cakes? I.e. are some cakes more likely to cause participants to quit?
b) Also: Is there an interaction between cake and condition, i.e. does the red velvet cause people to get sick in the blindfolded condition, but not without a blindfold - while it is the other way around for the chocolate cake?
Quick overview with three of the 20 cakes selected:
cakes[cakes$cake %in% c("chocolate", "banoffee", "carrot"), ]

        cake set blindfold dropout nodrop
3  chocolate   A       yes       7     33
6     carrot   A       yes       4     36
17  banoffee   B       yes      23     17
23 chocolate   A        no       1     39
26    carrot   A        no       9     31
37  banoffee   B        no      23     17

set is included in case that should play a role for my analysis. dropout is the number of participants who quit, nodrop the number of participants who did not quit, i.e. the two columns always total to 40 participants. Note that in real life, I do not only have 20 cakes, but 100 stimuli (50 in each set).
Thanks and a virtual piece of cake to everyone pointing me in the right direction!
Solution
As suggested in the answers, restructuring the dataset to single trials and doing a logistic regression was the key. I used the package brms to do it the Bayesian way. With the single-trial structure, I was not only able to include the cakes as random slope and intercept, but also the participants, which I found to explain much more variance than the ckaes - which makes sense. In the end, my model looked like this:
dropout ~ blindfold + (1 + blindfold | participant) + (1 + blindfold | cake) 

Thanks for the help!
 A: My initial reaction is to use a logistic regression.  
First, you would need to restructure your data. In the cake example, you'd create a new dataset with three variables: (1) dropout indicator, (2) cake, and (3) blindfold indicator.  Each row represents a unique person x cake x blindfold combination.  
A) Do dropout rates differ by cakes?
You can just compare the residual deviance between the null model and a model with the cake variable.  If the results are significant it means that there is at least one cake with differing dropout rates.  
mod <- glm(dropout ~ 1, data=cakes)
mod1 <- glm(dropout ~ factor(cake), data=cakes)

anova(mod, mod1, type="Chisq")

B) Differences based on interaction
Here, you could do something similar by comparing a model with just the main effects and a model with the interaction terms.
mod2 <- glm(dropout ~ factor(cake) + blindfold, data=cakes)
mod3 <- glm(dropout ~ factor(cake)*blindfold, data=cakes)

anova(mod2, mod3, type="Chisq")

If you're interested in comparing specific combinations of cake and blindfold status, you could use the multcomp package in R.  UCLA has a great tutorial.
