Hello, If someone can help I will be very grateful. My question is about which design should have my data to run deseq2 and see some comparisons.
First, I have the following data: 2genotypes, 2 tissue, and a treatment evaluated in 4 different times (control, t1, t2 and t3). In a matrix it is like this:
- library genotype tissue treatment replicate
- 1 A tissue1 Control 1
- 2 A tissue1 Time1 1
- 3 A tissue1 Time2 1
- 4 A tissue1 Time3 1
- 5 A tissue2 Control 1
- 6 A tissue2 Time1 1
- 7 A tissue2 Time2 1
- 8 A tissue2 Time3 1
- 9 A tissue1 Control 2
- 10 A tissue1 Time1 2
- 11 A tissue1 Time2 2
- 12 A tissue1 Time3 2
- 13 A tissue2 Control 2
- 14 A tissue2 Time1 2
- 15 A tissue2 Time2 2
- 16 A tissue2 Time3 2
- 17 A tissue1 Control 3
- 18 A tissue1 Time1 3
- 19 A tissue1 Time2 3
- 20 A tissue1 Time3 3
- 21 A tissue2 Control 3
- 22 A tissue2 Time1 3
- 23 A tissue2 Time2 3
- 24 A tissue2 Time3 3
I have the same for genotype B. So in total I have 48 samples.
I would like to know how gene expression changes over control, time1, time2, and time3 in one genotype and specific tissue, for example in the genotype B, in the tissue1.
So, I was applying this design:
dds <- DESeqDataSetFromHTSeqCount(sampleTable = sampleTable, directory = directory,
design = ~genotype + tissue + treatment)
dds$group <- factor(paste0(dds$genotype,dds$tissue,dds$treatment))
design(dds) <- ~ group
dds <- DESeq(dds)
then: resultsNames(dds):
[1] "Intercept" "group_ATissue1T1_vs_ATissue1Control"
[3] "group_ATissue1T2_vs_ATissue1Control" "group_ATissue1T3_vs_ATissue1Control"
[5] "group_ATissue2Control_vs_ATissue1Control" "group_ATissue2T1_vs_ATissue1Control"
[7] "group_ATissue2T2_vs_ATissue1Control" "group_ATissue2T3_vs_ATissue1Control"
[9] "group_BTissue1Control_vs_ATissue1Control" "group_BTissue1T1_vs_ATissue1Control"
[11] "group_BTissue1T2_vs_ATissue1Control" "group_BTissue1T3_vs_ATissue1Control"
[13] "group_BTissue2Control_vs_ATissue1Control" "group_BTissue2T1_vs_ATissue1Control"
[15] "group_BTissue2T2_vs_ATissue1Control" "group_BTissue2T3_vs_ATissue1Control"
There I don´t have "group_BTissue1T1_vs_BTissue1Control", to compare the fold change between the T1 and the control in the tissue1 in the genotype B.
However, I do not know why, but I can run this:
res <- results(dds, alpha = 0.05, contrast = c("group", "BTissue1T1","BTissue1Control"))
But to generate MA plot, I need to do the shrinkage (lfcShrink), and for that, I need the coef :
coef = c("group_BTissue1T1_vs_BTissue1Control") which I don't have.
Please how should be my design to see the fold change across the different treatment, in a specific genotype and tissue?. I was employing the default wald test, but because treatment is evaluated in 4 different times (control, t1, t2 and t3), maybe I should have to employ the LRT test.
Thanks in advance for your time.
thanks!, but after doing the relevel then I can´t see other combinations
You need a design without an intercept to make all possible contrasts. I am not familiar enough with DESeq2 for hands-on help, but in general you need
design = ~ 0 + (...). I am sure there is online help for this. In edgeR it is easy with themakeContrastsfunction, pretty sure there is something similar in DESeq2....and please comment to answers with
ADD COMMENT, not the answer field, to keep things logically organized.