SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
DESeq data analysis kasutubh Bioinformatics 2 08-07-2012 11:00 AM
negative p-value in DEseq analysis Marianna85 Bioinformatics 7 06-07-2012 04:05 AM
error in DESeq analysis stephenhart General 5 11-08-2011 02:55 AM
DESeq-statistical analysis without replicate lynn012 RNA Sequencing 0 10-27-2011 02:47 AM
RNAseq analysis using DESeq katussa10 Bioinformatics 9 08-29-2011 06:32 AM

Reply
 
Thread Tools
Old 12-14-2015, 09:12 AM   #1
andrewelamb
Junior Member
 
Location: Boston

Join Date: Dec 2015
Posts: 4
Default eEtting up DESeq 2 analysis

I have some RNA-seq data I'm trying to do DESeq2 analysis on, and it's more complicated than what I've done before.

Patients came in two different times and had their blood drawn and sequenced. There is some missing data, so that some patients only have data for visit 1 or 2. In between visits patients were given either a placebo or one of two treatments.

The answer we're trying to answer is what expression levels changed in between visits in a different way between the treatment groups. Based on the vignettes and other answered questions it seems the design would be:

~ patient + visit + treatment + visit:treatment

the interaction visit:treatment would be the 'difference in differences' of each treatment over the time between visits. Is this correct?
andrewelamb is offline   Reply With Quote
Old 12-14-2015, 10:11 AM   #2
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

You can simply exclude the patients for whom you only have a single sample. They'll get ignored in the analysis anyway. Anyway, yes your model is correct and you do indeed care most about the "visit:treatment" term.
dpryan is offline   Reply With Quote
Old 12-15-2015, 10:57 AM   #3
andrewelamb
Junior Member
 
Location: Boston

Join Date: Dec 2015
Posts: 4
Default

Thanks for the answer!

I did get this error however:

Error in checkFullRank(modelMatrix) :
the model matrix is not full rank, so the model cannot be fit as specified.
One or more variables or interaction terms in the design formula are linear
combinations of the others and must be removed.

my pheno file looks like:

sampleName visit condition patient
1 V2 control 1
2 V5 control 1
3 V2 treatment 2
4 V5 treatment 2
5 V5 treatment 3
6 V2 treatment 3
7 V5 treatment 4
8 V2 treatment 4
9 V2 control 5
10 V5 control 5

Removing patients from the experimental design worked. Is there any way, or value, to preserve the patient data?
andrewelamb is offline   Reply With Quote
Old 12-15-2015, 11:42 PM   #4
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

Indeed, I should have foreseen that :P

If you were to instead use "~patient+condition:visit+visit" and got rid of the "conditiontreatment:visitV2" column in the model matrix then the result would work. The original problem was that each condition is comprised of a set of patients, so you can't have patient coefficients and a "condition" coefficient (which is just the average of the patient coefficients!).

Sorry that that's so confusing.

Last edited by dpryan; 12-16-2015 at 05:34 AM.
dpryan is offline   Reply With Quote
Old 12-16-2015, 05:30 AM   #5
andrewelamb
Junior Member
 
Location: Boston

Join Date: Dec 2015
Posts: 4
Default

Thank you for the help!

I apologize, I'm not entirely clear on how to set up my model matrix based on your answer. It seems I would still need every column if I were to use "~patient+condition:visit+visit".
andrewelamb is offline   Reply With Quote
Old 12-16-2015, 05:35 AM   #6
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

I had a typo in my reply, I meant to remove the "conditiontreatment:visitV2" from the model matrix. That'll make it full rank,
dpryan is offline   Reply With Quote
Old 12-16-2015, 06:31 AM   #7
andrewelamb
Junior Member
 
Location: Boston

Join Date: Dec 2015
Posts: 4
Default

Ahh I see, I'm getting my sample table and the model matrix confused.

So is this the correct way to use my own model matrix?

design_string <- "~patient+condition:visit+visit"
sample_table <- read.table(input_file, row.names = NULL, header = T, sep = ",")
deseq_object <- DESeqDataSetFromHTSeqCount(sampleTable = sample_table,
design = ~condition, #have to have something here
directory = count_folder)
mm <- model.matrix(as.formula(design_string), sample_table)
mm <- mm[,-19] # gets rid of conditiontreatment:visitV2
deseq_object <- DESeq(deseq_object, full=mm, betaPrior=FALSE)
andrewelamb is offline   Reply With Quote
Old 12-16-2015, 10:51 AM   #8
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

Something along those lines at least.
dpryan 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 05:27 AM.


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