I want to compare the result between limma + voom with DESeq2. But I don't know how to make the design matrix for my date when I am using Limma.
I have two groups: WT, ABX. Each group has 6 samples: X1~X6(WT), X7~X12(ABX). It looks like:
head(exprSet, n=1) X.1 X.2 X.3 X.4 X.5 X.6 X.7 X.8 ENSMUSG00000051951 88 113 83 101 57 45 45 50 X.9 X.10 X.11 X.12 ENSMUSG00000051951 31 51 59 43
The row name is the geneID, the column is the count for each sample.
My code is below:
## read count matrix from file exprSet=read.table("ex_matrix_g1g2_h.txt", sep="\t", header=T, row.names = 1, stringsAsFactors = F) group_list <- factor(c(rep("WT",6), rep("ABX",6))) design <- model.matrix(~0+group_list) ## Is this correct? colnames(design) <- levels(group_list) rownames(design) <- colnames(exprSet) v <- voom(exprSet,design,normalize="quantile", plot = T) fit <- lmFit(v,design) fit2 <- eBayes(fit) tempOutput = topTable(fit2, n=Inf, adjust.method = 'BH', coef=2) ## if I want to compare ABX - WT, how to set coef?
Unfortunately, the above code generated very different results compared with DESeq2. In limma's result, almost all genes are significantly different.
I think that maybe I used the wrong design matrix. But I am not familiar with limma. Could some help me to figure out what the matter of my code?
Thank you so much!