[BioC] Designing a model with blocking and other interactions
Gordon K Smyth
smyth at wehi.EDU.AU
Fri Apr 4 01:31:22 CEST 2014
It is, but that would require to a shift from edgeR to voom and limma, and
one still can't test for a main effect in the presence of an interaction.
I have refered Eleanor to a section in the limma User's Guide which
explains the use of duplicateCorrelation for this purpose.
Gordon
On Thu, 3 Apr 2014, Aaron Mackey wrote:
> Is this not a case for using duplicateCorrelation across the Family units,
> which would account for the correlation within a family, while still
> allowing for a mitoHap main effect?
>
> -Aaron
>
>
> On Thu, Apr 3, 2014 at 7:04 PM, Gordon K Smyth <smyth at wehi.edu.au> wrote:
>
>> On Thu, 3 Apr 2014, Eleanor Su wrote:
>>
>> Hi Gordon,
>>>
>>> When I enter the design that you've suggested,
>>>
>>> design1 <- model.matrix(~Family)
>>> design2 <- model.matrix(~mitoHap*Treatment)
>>> design <- cbind(design1,design2[,3:4])
>>>
>>> and test for the last coefficient, I see that I get DE for the interaction
>>> between Treatment:mitoHap (which is what I wanted to look at). As I look
>>> through the other columns in the design matrix, I see that I have data for
>>> Treatment (coef=6) but not for mitoHap.
>>>
>>
>> No, naturally you can't have a mitoHap column because that factor is
>> confounded with Family.
>>
>> If I use an equivalent formula for design2
>>>
>>> design2<-model.matrix(~mitoHap+Treatment+mitoHap:Treatment)
>>>
>>> would this allow me to see both factors (treatment and mitoHap)
>>> independently in other columns of the design matrix AND the interaction
>>> between the two in the last coefficient? I'd like to be able to look at
>>> differential expression in each factor independently and the interaction
>>> between the two.
>>>
>>
>> Well, now you are ignoring Family, which previously you felt it was
>> important to account for.
>>
>> You are also asking for things that are impossible. It isn't meaningful
>> to test for factors independently of their interaction.
>>
>> You might find it helpful to read the sections on "multi level designs" in
>> the edgeR and limma User's Guides.
>>
>> Gordon
>>
>> If so, how would the last "design" formula change?
>>>
>>> design<-cbind(design1,design2[?])
>>>
>>> Thanks for you help.
>>>
>>> Best,
>>> Eleanor
>>>
>>>
>>> On Wed, Apr 2, 2014 at 8:25 PM, Gordon K Smyth <smyth at wehi.edu.au> wrote:
>>>
>>> Dear Eleanor,
>>>>
>>>> design1 <- model.matrix(~Family)
>>>> design2 <- model.matrix(~mitoHap*Treatment)
>>>> design <- cbind(design1,design2[,3:4])
>>>>
>>>> Then test for the last coefficient.
>>>>
>>>> Best wishes
>>>> Gordon
>>>>
>>>> Date: Tue, 1 Apr 2014 11:24:52 -0700
>>>>
>>>>> From: Eleanor Su <eleanorjinsu at gmail.com>
>>>>> To: "bioconductor at stat.math.ethz.ch" <bioconductor at stat.math.ethz.ch>
>>>>> Subject: [BioC] Designing a model with blocking and other interactions
>>>>>
>>>>> Hi All,
>>>>>
>>>>> I'm trying to set up a model matrix where I can look at the interaction
>>>>> between Treatment and mitochondrial haplotypes in my paired samples.
>>>>> These
>>>>> are the preliminary commands that I've set up:
>>>>>
>>>>> rawdata<-read.delim("piRNAtotalcount<10.txt", check.names=FALSE,
>>>>>
>>>>>>
>>>>>> stringsAsFactors=FALSE)
>>>>>
>>>>> y <- DGEList(counts=rawdata[,2:11], genes=rawdata[,1])
>>>>>> Family<-factor(c(6,6,9,9,11,11,26,26,28,28))
>>>>>> Treatment<-factor(c("C","H","C","H","C","H","C","H","C","H"))
>>>>>> mitoHap<-factor(c("S","S","S","S","S","S","D","D","D","D"))
>>>>>> data.frame(Sample=colnames(y),Family,Treatment,mitoHap)
>>>>>>
>>>>>> Sample Family Treatment mitoHap
>>>>> 1 6C (S) 6 C S
>>>>> 2 6H (S) 6 H S
>>>>> 3 9C (S) 9 C S
>>>>> 4 9H (S) 9 H S
>>>>> 5 11C (S) 11 C S
>>>>> 6 11H (S) 11 H S
>>>>> 7 26C (D) 26 C D
>>>>> 8 26H (D) 26 H D
>>>>> 9 28C (D) 28 C D
>>>>> 10 28H (D) 28 H D
>>>>>
>>>>> design<-model.matrix(?)
>>>>>
>>>>>>
>>>>>>
>>>>> I have 10 sequencing samples from 5 different families (a treatment and
>>>>> control sample from each family) and two different types of
>>>>> mitochondrial
>>>>> haplotypes. How do I set up a design where I can look at the interaction
>>>>> between the Treatments and mitoHap while still accounting for Family?
>>>>>
>>>>> Any help would be greatly appreciated. Thank you for your time.
>>>>>
>>>>> Best,
>>>>> Eleanor
>>>>>
>>>>
>> ______________________________________________________________________
>> The information in this email is confidential and intend...{{dropped:4}}
>>
>> _______________________________________________
>> Bioconductor mailing list
>> Bioconductor at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/bioconductor
>> Search the archives: http://news.gmane.org/gmane.
>> science.biology.informatics.conductor
>>
>
______________________________________________________________________
The information in this email is confidential and intend...{{dropped:4}}
More information about the Bioconductor
mailing list