Question: DESeq2: making specific comparisons between samples that have two conditions and two groups
2
5.0 years ago by
erin.gill8140
erin.gill8140 wrote:

Hello,

I'm analyzing RNA-seq data from two different cell types (ESDM and BMDM) with two treatment conditions (uninfected and infected). My problem is that I'm trying to compare two specific cell types and treatments that aren't in the results table shown after running the DESeq command. This is the code I'm using:

#read .csv file containing sample info

#make sample data frame from .csv sample sheet
coldata = with(samples,
data.frame(shortname = I(shortname),
countf = I(countf),
cell_type = cell_type,
biol_rep = biol_rep,
treatment = treatment))

#make DESeq2 dataset object
dds = DESeqDataSetFromHTSeqCount(sampleTable = coldata, directory = "/datadirectory/mydata", design = ~ cell_type + treatment + cell_type:treatment)

#collapse technical replicates
ddscollapsed <- collapseReplicates ( dds, groupby = dds$biol_rep) #relevel samples so that uninfected, BMDM are control samples ddscollapsed$cell_type <- relevel (ddscollapsed$cell_type, "BMDM") ddscollapsed$treatment <- relevel (ddscollapsed\$treatment, "uninfected")

#do differential expression analysis
ddscollapsed <- DESeq(ddscollapsed)

#final data frame:

> as.data.frame(colData(ddscollapsed))
cell_type biol_rep  treatment sizeFactor
BMDM_I_1      BMDM BMDM_I_1   infected  0.9241779
BMDM_I_2      BMDM BMDM_I_2   infected  0.6117723
BMDM_I_3      BMDM BMDM_I_3   infected  0.3563606
BMDM_I_4      BMDM BMDM_I_4   infected  0.8114125
BMDM_U_1      BMDM BMDM_U_1 uninfected  0.8375975
BMDM_U_2      BMDM BMDM_U_2 uninfected  2.2902482
BMDM_U_3      BMDM BMDM_U_3 uninfected  1.9502123
BMDM_U_4      BMDM BMDM_U_4 uninfected  2.1769115
ESDM_I_1      ESDM ESDM_I_1   infected  0.4880022
ESDM_I_2      ESDM ESDM_I_2   infected  0.8463753
ESDM_I_3      ESDM ESDM_I_3   infected  1.1217204
ESDM_I_4      ESDM ESDM_I_4   infected  0.6988828
ESDM_U_1      ESDM ESDM_U_1 uninfected  1.7999121
ESDM_U_2      ESDM ESDM_U_2 uninfected  0.7995184
ESDM_U_3      ESDM ESDM_U_3 uninfected  1.6323795
ESDM_U_4      ESDM ESDM_U_4 uninfected  1.0969576

#display what the different results names are
resultsNames(ddscollapsed)

#At this point, I get the following:

> resultsNames(ddscollapsed)
[1] "Intercept"                        "cell_type_ESDM_vs_BMDM"
[3] "treatment_infected_vs_uninfected" "cell_typeESDM.treatmentinfected" 

What I would like to be able to do is identify genes that are differentially expressed ONLY between infected ESDM cells and infected BMDM cells.

I've also tried the design formula "~ cell_type + treatment " and get:

> resultsNames(ddscollapsed)
[1] "Intercept"           "cell_typeBMDM"       "cell_typeESDM"
[4] "treatmentuninfected" "treatmentinfected"  

Here I've tried to implement the contrast list function, but I'm unsure of how to make it accept two terms for either the numerator or the denominator. I would use the matrix (numeric function) to input this data, but I'm unsure of how to indicate that "treatmentinfected" should be both +1 and -1.

Could you please let me know what to do so that I can make specific comparisons between samples that have two conditions and two groups?

Erin

deseq2 • 2.4k views
modified 5.0 years ago by Steve Lianoglou12k • written 5.0 years ago by erin.gill8140

Hi Erin,

I have read your post and since it is the same situation of my experiment, so I'm wondering if now you have  clearer idea than initial one after the Simon and Michael' s suggestion for the statistical analysis of your data with DESeq2.

I will appreciate if you will share with me the sequence of DESeq code you have used to answer your question:

"What I would like to be able to do is identify genes that are differentially expressed ONLY between infected ESDM cells and infected BMDM cells"

thanks

Annalisa

Answer: DESeq2: making specific comparisons between samples that have two conditions and
1
5.0 years ago by
Michael Love25k
United States
Michael Love25k wrote:

You should use the first design, ~ cell_type + treatment + cell_type:treatment, and then you can test for genes which are differentially expressed between infected ESDM cells and infected BMDM cells by adding the cell_type_ESDM_vs_BMDM effect (which is the difference of ESDM over BMDM in uninfected cells) with the interaction term cell_typeESDM.treatmentinfected (which is the additional difference of EDSM over BMBM in infected cells beyond the main effect). Adding these together gives the ESDM vs BMDM effect in infected cells:

results(dds, contrast=c(0,1,0,1))

Or you can also use the list() style of contrast to add these two effects:

results(dds, contrast=list(c("cell_type_ESDM_vs_BMDM","cell_typeESDM.treatmentinfected"), character()))

edit: the contrast can also be specified in the upcoming release (v1.6) by leaving off the second element of the list:

results(dds, contrast=list(c("cell_type_ESDM_vs_BMDM","cell_typeESDM.treatmentinfected")))

"What I would like to be able to do is identify genes that are differentially expressed ONLY between infected ESDM cells and infected BMDM cells."

So, you are looking for genes for which you

(i) have evidence that the expression differs between infected ESDM and infected BMDM cells, and

(ii) have evidence that the expression does not differ between non-infected ESDM and non-infected BMDM cells.

Note that (ii) is not the same as (ii') "have no evidence that the expression differs between non-infected ESDM and non-infected BMDM cells".

Mike's solution will give you only genes that fulfil (i), with no regard to (ii) or (ii').

DESeq2 can answer all three possible questions, so which one do you want?

Ah, good point! I didn't catch that the first time I read the question.

Hi Simon,

I am in exactly the same situation having exactly the same biological question with the same design (control_vs_case ,male_vs_female). Could you please give all the possible answers to these 3 questions?

I am particularly interested to find genes that are differentially expressed in cases_male_VS cases_females and not differentially expressed in control_male_VS control_female with no evidence that the expression differs between control_male VS control_female

Thanks,

Nikos