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 themakeContrasts
function, 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.