Please help me build a design that accounts for both the matched case-control design (there were 1 case [disease=1] and 2 controls [disease=0] in each pair [having the same matchid], matched for age and gender) and the correlation between multiple measurements at different timepoints (time 1, 2, 3) in a given subject. My objective is to compare the microbial differential abundance between cases and controls over all timepoints (not treating each timepoint separately). Thanks a lot!

My dataset is as following:

If I use edgeR, is the design correct?

design <- model.matrix(~matchid+personid+disease+time+disease:time)

(There are paired examples and repeated measurement examples in the edgeRUsersGuide. However, in repeated measurement examples, repeated measurements at different timepoints were obtained from independent subjects. By contrast, in my study, repeated measurements at different timepoints were obtained from the same subjects.)   

Unfortunately, setting up the design to answer your question is much more complicated than you might think, because personid is nested within disease. This means that you cannot estimate the log-fold change between disease and healthy individuals and block on the person effect at the same time (which would be bad if there are strong patient-to-patient effects, which are probably present in any kind of human patient data).

There are two solutions to this, but one of them involves subsetting the dataset by timepoint and is not appropriate here given that you want to use information from all timepoints at once. So this leaves us with the only other solution, which is to switch to limma and voom. Basically, set up your design matrix like so:

personid <- factor(rep(1:9, each=3))
matchid <- factor(rep(1:3, each=9))
time <- rep(1:3, 3)
disease <- factor(rep(rep(c(1,0,0), each=3), 3))
design <- model.matrix(~matchid + time + disease)

... and run an analysis using duplicateCorrelation with the blocking factor set to personid. This will use information from all timepoints, account for any person-to-person correlations and allow you to compare directly between disease states.

• Note the factor, which is important to ensure that the numbers are not treated as real-valued covariates.
• You're assuming that the effect of time on microbial abundance is linear.
• personid is nested within matchid, and conventionally (i.e., when blocking on all factors in the design matrix) you wouldn't need the latter. In this case, though, I've used both of them, to avoid violating homogeneity assumptions in duplicateCorrelation when there are systematic differences between matchid levels.