Analyzing replicates of contingency tables

Question

Problem statement

We would like to measure the effect of treatment on cell differentiation, in other words whether or not specific treatments cause cells to change their cell types. So we treat cells, measure number of cells that belong to different cell types and we perform multiple biological replicates of this experiment (so they are independent), 3 to be exact. The processed data looks as following:

Cell type	Control, replicate 1	Treatment, replicate 1	Control, replicate 2	Treatment, replicate 2	...
A	200	300	500	200	...
B	500	700	300	900	...

Question

How can we combine results of individual replicates in a statistically sound way? Ideally we would like to say whether or not treatment has an effect, rather than treatment has an effect in replicate 1, but not in replicate 2.

Possible approaches

• We run Fisher's exact test on the contingency table for each replicate, then correct for multiple hypothesis testing. We get multiple p-values as a result so can't measure significance of the overall effect.

• We pool the data from multiple replicates together ( by adding appropriate columns) and run Fisher test on resulting contingency table. But that defeats the purpose of running multiple replicates in the first place: two replicates giving consistent results is more significant than just one.

• We run some sort of permutation test to get a distribution of the test statistic under the null hypothesis - treatment has no effect, and compare it to our observed statistic. I'm at a loss here, though. Permuting both treatment labels and replicate labels of the original observation is equivalent to pooling the data from all replicates. Moreover, how do we calculate one observed test statistic of the unpermuted data given that there are multiple replicates?

• We go full Bayesian and say that our observations are sampled from binomial distribution, and its parameters depend on both treatment and replicate:

  $Y \sim Binom(\theta|N)$

  $log(\theta) \sim Treatment + Replicate + \epsilon$

  That way we decompose the effect into treatment-dependent and replicate-dependent and our result will be a full posterior distribution of the treatment-dependent effect.

Am I missing something? Are there other approaches to this problem? What would be the pros and cons in your opinion?

Thanks!

Answer 1

No need to go full Bayesian, a per your last point. The model you've written down is just logistic regression, so I see no problem doing

$$ \operatorname{logit}(p) = \beta_0 + \beta_1 x$$

Here, $\beta_0$ is the log odds of cell type A, $\beta_1$ is the log odds ratio for the treatment, and $x$ is a binary indicator for treatment or not. You can very easy do this in R via the following commands

repli <- c(1,1,2,2,...)
cell_type_a <- c(200, 300, 500, 200, ...)
cell_type_b <- c(500, 700, 300, 900, ...)
txt <- c(0, 1, 0, 1, ...)

model = glm(cbind(cell_type_a, cell_type_b) ~ txt, family = binomial())

The cbind(cell_type_a, cell_type_b) will model the probability of being in cell type A. If you wanted to model the probability of being in cell type B, just do cbind(cell_type_b, cell_type_a).