Hello!

I am interested in doing DE analysis in 12 samples from 3 different patients ("id"). The samples have been sequenced using bulk RNA-seq in a specific group of isolated cells. I have treated the cells by doing a knock-out experiment to remove a gene (group1 as mutated "MUT" and wild-type "WT"), and I have infused the cells in mice (group2 as "pre" and "post") before and after the knock-out. In summary I have 4 measures of gene expression for each id. My main goal is to see how the gene expression from these cells change before and after the knock-out and before and after interaction with mice. This is the data:

sample group1 group2 id
1 MUT post A
2 WT post A
3 MUT pre A
4 WT pre A
5 MUT post B
6 WT post B
7 MUT pre B
8 WT pre B
9 MUT post C
10 WT post C
11 MUT pre C
12 WT pre C

I was wondering what is the correct way to do DE with limma. I don't know if I need to perform various DE analyses or include them as one in a single design matrix. I think of doing it this way, but I want to be sure is the correct way. Many thanks!!

group_id <- factor(data$id)

group1 <- factor(data$group1)

group2 <- factor(data$group2)

design <- model.matrix(~0+group_id)

design <- cbind(design, PRE=group2=="pre" & group1=="MUT")

design <- cbind(design, POST=group2=="post" & group1=="MUT")

v.DEGList.filtered.norm <- voom(myDGEList.filtered.norm, design, plot = TRUE)

fit <- lmFit(v.DEGList.filtered.norm, design)

contrast.matrix <- makeContrasts(POST - PRE, levels=colnames(design))

fits <- contrasts.fit(fit, contrast.matrix)

ebFit <- eBayes(fits)

myTopHits <- topTable(ebFit, adjust ="BH", coef=1, number=10000, sort.by="logFC")