RNA-seq DE analysis: repeated measures design (options available, mixed models?)

Question

[bio.infanatics]

I am performing DE analysis on RNA-seq data. The experimental design is as follows:

Animals were ranked/assigned based on the observed phenotype of interest. Factor variable with 2 levels: T or N (binary). T group has n=20 animals and the N group has n=22.

The same animals were sampled on 2 separate days, (hot day=day1, cool day=day2). Day collected was recorded for each sample, factor variable with 2 levels: day1 or day2 (binary).

I have other pieces of data collected, but I created another variable rankday, which was coded as: T on day 1=0, N on day 1=1, T on day 2=2, N on day 2=3.

So, the design is repeated measures with a total of 42 animals and 84 samples. I know I need to do a linear mixed model because of the repeated measures, and treat rank as a random effect.

I coded all of this and ran the DE using the dream function in the variancePartition package: (https://bioconductor.org/packages/release/bioc/vignettes/variancePartition/inst/doc/dream.html) similar to this but fitting my data.

I did DE contrast comparisons on T1 vs T2 , T1 vs N1, N1 vs N2, and T2 vs N2. We are interested in the differences in the T and N group and how they changed between day 1 and day 2.

What I'm wondering is, are there other ways/software/packages that can fit a linear mixed model for repeated measures when running a DE analysis? Or is variancePartition:dream my best option for this?

Edit: also is it ok to use the dream function for all of these contrasts, even if T1 vs N1 and T2 vs N2 do not technically have repeated measures when comparing these? Thanks in advance

Answer 1

[ATpoint]

You can use limma-voom and take the repeat measures into account by combining with its duplicateCorrelation function. This has been recommended before many times at the Bioconductor support site to deal with repeated measures.

You could use the voomLmFit function from the edgeR package and provide the repeated measure information to the block argument which then runs duplicateCorrelation as well as voom and lmFit interally. It lives in the edgeR rather than limma package because it uses edgeR functionality to do some other things under the hood which are not relevant here though.