DESeq2 complicated design - effect of replicated samples

Question

I have RNAseq data from a relatively complicated experimental design with variables = genotype, treatment, time, and batch. I have 2 biological replicates for each genotype/condition, however unfortunately in the first iteration of the experiment one sample had poor RNA quality and is unusable. To get around this lack of a replicate, we repeated the experiment for the two biological replicates in the condition of the lost sample, and did a new library prep/sequencing from this repeated experiment. The repeated experiment obviously lead to a batch effect, and although it isn't the main driver of variance, it is a significant factor. I would like to model or remove the batch effect from the experiment using DESeq2 if possible, but after reading this similar post: can I model technical replicates in DESeq2? I am concerned about inflating the power of the test for those conditions that were replicated.

Here is a diagram of the experimental design to make it more easy to understand:

My "colData" for the DESeq2 object looks like this:

My previous thought was to simply collapse the genotype+time+treatment variables into a single combined factor "condition" and then model the batch effect using a design of "~batch + condition". However this doesn't really address the concern of inflating the power of the tests due to technical replication.

Should I Just collapse the counts of these replicated conditions even though the batch affect will affect the combined count values? Or should I use a different tool than DESeq2 ? Or is there a way for me to remove the batch effect prior to collapsing the counts for replicates?

What would you guys recommend?

Answer 1

I agree with @devon-ryan that the design matrix looks reasonable, i.e. that it's okay to treat the new samples as biological replicates of the original batch, and that the batch number should be included in that matrix (assuming DESeq2 doesn't complain about it) to filter out most of the batch effects.

DESeq2 tends to be conservative in its differential expression calculations. While it's true that more replicates in one sample will increase the number of differentially-expressed genes that are statistically significant (because the biological variation can be better modelled), I wouldn't expect this to change the overall picture.

However, whole-genome differential expression testing should not be used on its own to demonstrate biological relevance. I recommend as a first-pass sanity check that people compare any differential expression results to what is seen on a genome browser to identify any quirky things going on (e.g. expression across only a part of the gene). Ideally, results should be confirmed with targeted analyses on the same (or similar) biological samples, either at the transcript level (e.g. qPCR) or at the protein level (e.g. Western Blot). When aggregating evidence from multiple analyses, the specific way in which one of those analyses was done is less important.

Aside: It'd obviously be better [in the future] to make sure that your sample replicates are both robust to drop-out, and allow for outliers to be noticed, reducing the batch-effect problem. In practise this means a minimum of four replicates per informative condition, but because funders / researchers are well known for only choosing the minimum, I'd just state that six replicates per sample as what is needed to produce good results.