Question: DEseq2 - design help
0
8 months ago by
udp3f0
udp3f0 wrote:

Hi,

I have the following design where each sample is pooled from many samples and the control and treated(S) samples are not paired. I want genes that are differentially expressed under the effect of the 2 drugs for treatment S. I created 'sample.nested' to account for the pairing of control and treated samples. I used a design

design(dds) = ~ drug + drug:sample.nested + drug:treatment


• Firstly, is my design right for the question I mention?
• Secondly, resultsNames(dds) gives a term "drugisovs_pro". Will this tell me the differentially expressed genes for drug 'iso' compared to 'pro'? Does this account for sample differences or treatment effects while talking about this difference? I use relevel in the design for dds$treatment and dds$drug, so that the reference levels are 'control' and 'pro' respectively.
• Thirdly, one of the questions I am interested to know is the drug specific effects. I guess in this above design I am controlling for the drug specific effects with ~drug. How can I do this. Please advice
sample  treatment   drug    sample.nested
U46801  S       iso 1
U46813  control iso 1
U46802  S       iso 2
U46814  control iso 2
U46803  S       iso 3
U46815  control iso 3
U46804  S       iso 4
U46816  control iso 4
U46917  S       pro 1
U46929  control pro 1
U46918  S       pro 2
U46930  control pro 2
U46919  S       pro 3
U46931  control pro 3
U46920  S       pro 4
U46932  control pro 4

deseq2 • 174 views
modified 8 months ago by Michael Love25k • written 8 months ago by udp3f0
0
8 months ago by
Michael Love25k
United States
Michael Love25k wrote:

You say "the control and treated samples are not paired", but then you also have the numbers 1-4 in your table above. I don't understand the discrepancy.

Take the first two rows. Do these two samples have anything in common?

Hi Michael,

Sorry if that was confusing. Considering the first 2 rows, the 2 samples do not have anything in common. For the metadata, I paired the control and test groups randomly to add an extra column "sample.nested" and later accounted for them, following a post as in https://support.bioconductor.org/p/62357/#62368.

Each sample here for e.g. U46801 is liver tissue coming from 3 other samples/animals (pooled). (So for each type of drug administered, actually we are looking at 3x8 samples).

1

You shouldn't use sample.nested here if those two samples have nothing in common. That is for matched samples.

If you want to get treatment effects in each drug and then contrast the drugs, use: ~drug + drug:treatment. You can use results() with name to pull out the treatment effects, and results() with contrast=list(... , ...) to compare the two treatment effects.

Hi Michael,

Thank you for your reply. That was extremely helpful. I seemed to have overlooked the example which was for a paired design. One other question is, is there no need to account for the sample specific effects?

Also, with what you suggested resultsNames(dds) gives

resultsNames(dds)
[1] "Intercept"          "drug_pro_vs_iso"    "drugiso.treatmentS"
[4] "drugpro.treatmentS"


How are names= "drugprovs_iso" & contrast=list("drugpro.treatmentS","drugiso.treatmentS") different?

-Uma

There are no repeated measures per sample, so no. You can’t control for sample specific effects and have replication at the same level.

Ok thanks! Can you please answer how the below results for the design you suggested are different

names= "drug_pro_vs_iso" & contrast=list("drugpro.treatmentS","drugiso.treatmentS")


I see that I did not get any differentially expressed genes for any of these below

results(dds, name="drugpro.treatmentS")
results(dds, name="drugiso.treatmentS")
results(dds, contrast=list("drugiso.treatmentS", "drugpro.treatmentS"))


But when I run this analysis separately for each drug with the metadata

sample drug treatment
U46917  pro S
U46918  pro S
U46919  pro S
U46920  pro S
U46929  pro C
U46930  pro C
U46931  pro C
U46932  pro C


I have significant results.

Those should be similar results, although not identical because the estimation of dispersion uses all the samples in the dataset. Rather than just saying that there were “some” or “none” can you say if the genes are similar by rank?

results(dds, name="drugpro.treatmentS")
results(dds, name="drugiso.treatmentS")
results(dds, contrast=list("drugiso.treatmentS", "drugpro.treatmentS"))


but, have results for (I am interested to know what question this output answers?

results(dds, name="drug_pro_vs_iso")



But I get differentially expressed genes (at Padj<0.05) when I run the analysis separately. So I have no genes to compare here.