because I don't have a statistician who could help me verify my design matrix. So, I am looking for someone who has time to help me on this.
Hi, I am pretty new to the RNA seq data analysis and I need some advice on the design matrix formulation:
This is my sample file with a column group id which I made to group all the variables in one to distinguish one sample from another and then use it in the design matrix:
Sample Cell Type Concentration time GroupID
1 AB1 Control 0 5 AB1_Control_0_5
2 AB1 Control 0 5 AB1_Control_0_5
3 AB1 Control 0 5 AB1_Control_0_5
4 AB1 Treatment 5 5 AB1_Treatment_5_5
5 AB1 Treatment 5 5 AB1_Treatment_5_5
6 AB1 Treatment 5 5 AB1_Treatment_5_5
7 ST1 Control 0 5 ST1_Control_0_5
8 ST1 Control 0 5 ST1_Control_0_5
9 ST1 Control 0 5 ST1_Control_0_5
10 ST1 Treatment 5 5 ST1_Treatment_5_5
11 ST1 Treatment 5 5 ST1_Treatment_5_5
12 ST1 Treatment 5 5 ST1_Treatment_5_5
13 AB1 Control 0 8 AB1_Control_0_8
14 AB1 Control 0 8 AB1_Control_0_8
15 AB1 Control 0 8 AB1_Control_0_8
16 AB1 Treatment 8 8 AB1_Treatment_8_8
17 AB1 Treatment 8 8 AB1_Treatment_8_8
18 AB1 Treatment 8 8 AB1_Treatment_8_8
19 ST1 Control 0 8 ST1_Control_0_8
20 ST1 Control 0 8 ST1_Control_0_8
21 ST1 Control 0 8 ST1_Control_0_8
22 ST1 Treatment 8 8 ST1_Treatment_8_8
23 ST1 Treatment 8 8 ST1_Treatment_8_8
24 ST1 Treatment 8 8 ST1_Treatment_8_8
Contrasts I need
for my contrast matrix I want to compare:
AB1.Treatmentat concentration 5 and time point 5 vsAB1.controlat concentration 0 and time point 5AB1.Treatmentat concentration 8 and time point 8 vsAB1.controlat concentration 0 and time point 8ST1.Treatmentat concentration 5 and time point 5 vsST1.controlat concentration 0 and time point 5ST1.Treatmentat concentration 8 and time point 8 vsST1.controlat concentration 0 and time point 8
design matrix I used:
dds = DESeqDataSetFromMatrix(countData = Count_data,
colData = Meta_data,
design = ~ GroupID)
results(dds,contrast=c("GroupID","AB1_Treatment_5_5","AB1_Control_0_5"))
results(dds,contrast=c("GroupID","ST1_Treatment_5_5","ST1_Control_0_5"))
results(dds,contrast=c("GroupID","AB1_Treatment_8_8","AB1_Control_0_8"))
results(dds,contrast=c("GroupID","ST1_Treatment_8_8","ST1_Control_0_8"))
and when I run resultnames(dds) I see some contrast I don't need and not the ones I need. For example:
GroupIDAB1_Treatment_5_5vsAB1_Control_0_5
GroupIDST1_Treatment_5_5vsAB1_Control_0_5 ( this is not what I want)
but I want GroupIDST1_Treatment_5_5vsST1_Control_0_5.
Also I get the results from this but I see some only 1 as adj values for all the genes in some contrasts. So, is my design matrix right? Could somebody help me on this?
Please tell me if I am doing wrong somewhere
Thanks in advance!
1 answer
You already got answers on the bioconductor support site, why are you asking the same question here?
I'd recommend creating 4 design matrices or using only a subset of the samples for each DE analysis instead of using them all at once.
For example, for your first DE query, use just the first 6 samples.
EDIT: See swbarnes2's comment below.
I think in generally, it's preferred not to do this; it's preferred to keep all the samples in the dds object to get better estimate of the gene dispersion. In this particular case, you might consider separating the two types of cells out, if they are very different; keeping them together might make the library normalization wonky, or if one cell type has a much higher dispersion than the other, using the higher estimates on the other one will hurt the p-value unnecessarily.
thanks for the tip @swbarnes2.
Log in to answer this question.
Hello riya!
It appears that your post has been cross-posted to another site: https://support.bioconductor.org/p/128057/
This is typically not recommended as it runs the risk of annoying people in both communities.
Thanks for the pointer, swbarnes2