Hello, 

I have a fairly complicated design of experiment and I would like some help/feedback on designing the model.matrix. The data is coming from an experiment for which there are two groups of mice (young/old), the cells of which have been used for sorting populations with two markers (sort1/sort2) and each sort has positive cells and negative cells. The problem is that the mice from which the cells are coming are nested within both the sorts and the age groups and that some groups have 3 some 4 and some 5 mice. To explain a bit better, my samples matrix looks like this:

Initially I thought of splitting the data in two (sort1 and sort2 groups) and then using a nested design within that:

  design <- model.matrix( ~ cell + age +  age:cell + age:mouse_nest)

which works for sort2 but not for sort1 because the groups of mice are different for the two groups of young and old (5 vs 4 samples per group). As far as I understand the way to resolve this is either to remove a pair of samples so that I have 4 samples in each group or to remove the age:mouse_nest term. However, neither of the two solutions sounds great to me because a) don't like removing samples and b) there seem to be differences according to the mice. How do people go about choosing which is best, looking at the dispersion estimates? Any other ways to resolve this?

Also I would like to be able to compare between the positive cells of one marker (sort1) with the positive cells of the other marker (sort2) so I would like to put all the samples together but then the problem with the nesting becomes even greater because of the differences in group sizes. Is the best way to just put these samples together (sort1+ve vs sort2+ve) for young and old and forget about the nesting all together, using a design matrix like the below:

  design <- model.matrix( ~ age + sort + age:sort)

(or the equivalent form of combining them into one factor and using contrasts)

Thank you in advance for all your help!

Best wishes,

Emma

The only factor that is nested inside another factor is mouse. So I think your best bet is to use the limma-voom with duplicateCorrelation model. Assuming you want the complete 3-way interaction between age, sort, and cell, I would create a group variable and use that in the design:

library(limma)
group <- interaction(sort, cell, age, sep=".", drop=TRUE)
design <- model.matrix(~0 + group)

and then proceed with the analysis as described here[1], using mouse as the block argument to duplicateCorrelation. If you don't want the 3-way interaction, use whatever design you like involving those three variables, but leave out mouse, since that is handled as a random effect by duplicateCorrelation.

I don't think the unequal group sizes are an issue in this design, as long as you are properly modelling all the variables involved.

[1]: A: using duplicateCorrelation with limma+voom for RNA-seq data