SEQanswers

Go Back   SEQanswers > Applications Forums > RNA Sequencing



Similar Threads
Thread Thread Starter Forum Replies Last Post
DESEQ2: Which of the following 2 contrast methods does what I want ? EVELE RNA Sequencing 0 08-24-2015 06:46 AM
Advanced interaction contrast in DESeq2 Skiaphrene Bioinformatics 2 04-28-2015 03:25 PM
DESeq2 Multiple Contrast Correction jsage8 Bioinformatics 1 04-18-2015 07:12 AM
Paired design versus unpaired design in DESeq2 KristenC RNA Sequencing 1 05-29-2014 11:05 AM
DESeq2 Contrast Ruhi Bioinformatics 6 04-17-2014 06:12 AM

Reply
 
Thread Tools
Old 01-27-2017, 04:05 AM   #1
Michel Meyer
Junior Member
 
Location: France

Join Date: Jan 2017
Posts: 3
Default Extract contrast from DESeq2 advanced design

Dear all,

I wish to analyze RNAseq experiment data with complex matrix design. The matrix I wish to analyze is the same as displayed in the edgeR documentation section 3.5:
Code:
   Disease  Patient Treatment
1  Healthy  1       None
2  Healthy  1       Hormone
3  Healthy  2       None
4  Healthy  2       Hormone
5  Healthy  3       None
6  Healthy  3       Hormone
7  Disease1 4       None
8  Disease1 4       Hormone
9  Disease1 5       None
10 Disease1 5       Hormone
11 Disease1 6       None
12 Disease1 6       Hormone
13 Disease2 7       None
14 Disease2 7       Hormone
15 Disease2 8       None
16 Disease2 8       Hormone
17 Disease2 9       None
18 Disease2 9       Hormone
So, to analyse this kind of matrix design with DESeq2, I have followed the instructions in the DESeq2 vignette section 3.12.1. Creating a "nested" Patient variable (called PN) like that:
Code:
> dm
   Disease  PN Treatment
1  Healthy  P1 None
2  Healthy  P1 Hormone
3  Healthy  P2 None
4  Healthy  P2 Hormone
5  Healthy  P3 None
6  Healthy  P3 Hormone
7  Disease1 P1 None
8  Disease1 P1 Hormone
9  Disease1 P2 None
10 Disease1 P2 Hormone
11 Disease1 P3 None
12 Disease1 P3 Hormone
13 Disease2 P1 None
14 Disease2 P1 Hormone
15 Disease2 P2 None
16 Disease2 P2 Hormone
17 Disease2 P3 None
18 Disease2 P3 Hormone
and using the following design formulae to analyse with DESeq2:
Code:
> df
~ Disease + Disease:PN + Disease:Treatment
Using this design matrix and formulae, I got the following output from the resultsNames function:
Code:
> dds <- DESeqDataSetFromMatrix(countData=rawData, colData=dm, design=df) 
> dds <- DESeq(dds, quiet=TRUE)
> resultsNames(dds)
[1] "Intercept"             "Disease_Disease2_vs_Disease1"  "Disease_Healthy_vs_Disease1"   "DiseaseDisease1.PNP2"         
[5] "DiseaseDisease2.PNP2"  "DiseaseHealthy.PNP2"           "DiseaseDisease1.PNP3"          "DiseaseDisease2.PNP3"         
[9] "DiseaseHealthy.PNP3"   "DiseaseDisease1.TraitmentNone" "DiseaseDisease2.TraitmentNone" "DiseaseHealthy.TraitmentNone"
So, now I can extract comparisons DiseaseDisease1.TraitmentNone vs DiseaseHealthy.TraitmentNone and DiseaseDisease2.TraitmentNone vs DiseaseHealthy.TraitmentNone, using contrast option of results function.

However, I'm also interested in comparisons DiseaseDisease1.TraitmentHormone vs DiseaseHealthy.TraitmentHormone and DiseaseDisease2.TraitmentHormone vs DiseaseHealthy.TraitmentHormone. To get these comparisons, I have rerun a DESeq2 analyses after inverting the factor levels from the Treatment condition. So, now I can access to the other ones:
Code:
> dm$Traitment
[1] None Hormone None Hormone None Hormone None Hormone None Hormone None Hormone None Hormone None Hormone None Hormone
Levels: Hormone None
> dm$Traitment <- relevel(dm$Traitment,'None')
> dm$Traitment
[1] None Hormone None Hormone None Hormone None Hormone None Hormone None Hormone None Hormone None Hormone None Hormone
Levels: None Hormone
> dds <- DESeqDataSetFromMatrix(countData=rawdata, colData=dm, design=df) 
> dds <- DESeq(dds, quiet=TRUE)
> resultsNames(dds)
[1] "Intercept"             "Disease_Disease2_vs_Disease1"     "Disease_Healthy_vs_Disease1"      "DiseaseDisease1.PNP2"            
[5] "DiseaseDisease2.PNP2"  "DiseaseHealthy.PNP2"              "DiseaseDisease1.PNP3"             "DiseaseDisease2.PNP3"            
[9] "DiseaseHealthy.PNP3"   "DiseaseDisease1.TraitmentHormone" "DiseaseDisease2.TraitmentHormone" "DiseaseHealthy.TraitmentHormone"
I wish to know if it is the right way to get the requested comparisons?
Because I have seen that the results between (1) DiseaseDisease1.TraitmentNone vs DiseaseHealthy.TraitmentNone and (2) DiseaseDisease2.TraitmentHormone vs DiseaseHealthy.TraitmentHormone are nearly the same:

Code:
 (1)        (2)
 0,268438  -0,268438
 0,141639  -0,141639
 0,000000   0,000000
-0,696627   0,696634
 0,000000  -0,973195
 0,034611  -0,034611
 0,530830  -0,530829
-0,118112   0,118113
-0,016033   0,016033
 0,270009  -0,270009
 3,021271  -3,021392
 0,000000   0,000000
 0,000000   0,000000
 1,617347  -1,617340
-0,491919   0,491938
-0,380177   0,380001
 0,528211  -0,528204
-0,052266   0,052266
 0,204210  -0,204210
 0,267654  -0,267591
 0,012433  -0,012432
-0,046289   0,046290
 ...
Thanks in advance for any help.
For information, I use R version 3.3.0 (2016-05-03) and DESeq2_1.12.4.
Michel Meyer is offline   Reply With Quote
Old 01-31-2017, 07:34 AM   #2
shunyip
Member
 
Location: San Francisco South Bay Area

Join Date: Oct 2013
Posts: 20
Default

Hello Michel,

While your way of doing it can be correct, maybe you can use a much simpler way to perform this analysis.

For example, you can extract DiseaseDisease1.TreatmentNone and DiseaseHealthy.TreatmentNone's samples into one matrix. Then, create a much less complicated design matrix to compare them. (1 condition vs 1 condition) This way, you will definitely know that you are doing the right thing.

Best,
shunyip is offline   Reply With Quote
Old 02-06-2017, 07:55 AM   #3
Michel Meyer
Junior Member
 
Location: France

Join Date: Jan 2017
Posts: 3
Default

Thank you for your answer.
Michel Meyer is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 10:49 AM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO