Thanks a lot Naomi for your help. There are two more things that I'm not
sure about regarding the batch effect. Is it enough to include the batch in
the design matrix for it to be corrected, or do I need to correct it using
"duplicateCorrelation"  (fit <- lmFit(..., block=...,
correlation=dupcor$consensus)?
 Also in few posts in the mailing lists there is a recommendation to include
in the contrast matrix only those comparisons one is really interested with.
What is the effect of including more/unnecessary contrasts? would it change
something to the results I get from decideTests?
Thanks again
Mali

On Thu, Mar 10, 2011 at 3:04 AM, Naomi Altman <naomi@stat.psu.edu> wrote:

> Dear Mali,
> I have a bit more time now.  I think I answered your first question in my
> last email.  Limma does not directly do a "main effects" and interactions
> type of F-test, although you can can force the issue by using the F.p.value
> and several different sets of contrasts.
>
> The second approach is not valid.  The t-test denominator will be incorrect
> because you did not account for batch and treatment.
>
> Best wishes,
> Naomi
>
>
> Regarding
>
>
> At 07:05 AM 3/7/2011, you wrote:
>
>> Dear Naomi
>> First I would like to apologize for writing you directly.
>> My name is Mali Salmon, I am a researcher in Tel Hashomer hospital in
>> Israel.
>> I am currently analyzing microarray data using limma package, and I would
>> like to have some advice from a statistician who has an experience with
>> microarray analysis using Limma and bioconductor.
>> I was searching the mailing list, and saw your name as one who assist many
>> people with their analysis.
>> I have a simple factorial design which I build the design and contrast
>> matrix for, and I just want to be sure that what I did is correct, I would
>> be very grateful is you could give me your opinion.
>>
>> I have two factors, one is strain with three levels (C,D,S), and a
>> treatment factor with two levels: treated (DSR) and untreated (NO), I also
>> have batch effect which I included in the matrix
>>
>> > targets
>>
>>   FileName Treatment Strain batch
>> 1     NO2_C        NO      C     1
>> 2    DSR2_C       DSR      C     1
>> 3     NO2_D        NO      D     1
>> 4    DSR2_D       DSR      D     1
>> 5     NO2_S        NO      S     1
>> 6    DSR2_S       DSR      S     1
>> 7     NO3_C        NO      C     2
>> 8    DSR3_C       DSR      C     2
>> 9     NO3_D        NO      D     2
>> 10   DSR3_D       DSR      D     2
>> 11    NO3_S        NO      S     2
>> 12   DSR3_S       DSR      S     2
>>
>> I am trying to follow example 8.7 in the user guide (factorial design)
>> In order to build the design matrix this is what I did:
>>
>> TS <- paste(targets$Treatment, targets$Strain, sep=".")
>> TS <- factor(TS, unique(TS))
>> design<-model.matrix(~0+TS)
>> colnames(design) <- levels(TS)
>> batch<-c(0,0,0,0,0,0,1,1,1,1,1,1)
>> design<-cbind(design,batch)
>>
>> > design
>>   NO.C DSR.C NO.D DSR.D NO.S DSR.S batch
>> 1     1     0    0     0    0     0     0
>> 2     0     1    0     0    0     0     0
>> 3     0     0    1     0    0     0     0
>> 4     0     0    0     1    0     0     0
>> 5     0     0    0     0    1     0     0
>> 6     0     0    0     0    0     1     0
>> 7     1     0    0     0    0     0     1
>> 8     0     1    0     0    0     0     1
>> 9     0     0    1     0    0     0     1
>> 10    0     0    0     1    0     0     1
>> 11    0     0    0     0    1     0     1
>> 12    0     0    0     0    0     1     1
>>
>>
>> Now the questions I would like to answer are:
>> 1. what is the main effect of the strain (here I want to ignore the
>> treatment, just to find the difference between different strains)
>> 2. what is the main effect of the treatment (genes that change their
>> expression because of treatment, again I want to ignore the strain here).
>> 3. genes that are respond to treatment in each strain + interaction
>>
>> Below is how I created the contrast.matrix, the first three contrasts are
>> to answer the first question, contrast number 4 is to answer question 2, and
>> the last 6 contrasts to answer question 3
>>
>> cont.matrix <- makeContrasts(
>>     DvsC=(DSR.D+NO.D)-(DSR.C+NO.C),
>>
>>     SvsC=(DSR.S+NO.S)-(DSR.C+NO.C),
>>     DvsS=(DSR.D+NO.D)-(DSR.S+NO.S),
>>     DSRvNO=(DSR.C+DSR.D+DSR.S)-(NO.C+NO.D+NO.S),
>>     C.DSRvsNO=DSR.C-NO.C,
>>     D.DSRvsNO=DSR.D-NO.D,
>>     S.DSRvsNO=DSR.S-NO.S,
>>     DiffDvsC=(DSR.D-NO.D)-(DSR.C-NO.C),
>>     DiffSvsC=(DSR.S-NO.S)-(DSR.C-NO.C),
>>     DiffSvsD=(DSR.S-NO.S)-(DSR.D-NO.D),
>>     levels=design)
>> fit2 <- contrasts.fit(fit, cont.matrix)
>> fit2 <- eBayes(fit2)
>>
>> Is this matrix correct?
>> The batch is included in the design matrix, is that mean that the batch
>> will be removed?
>>
>> One of the comparisons I am interested with is the different between
>> strains (first question). Would you suggest to do a separate "several group"
>> analysis (example 8.6 in limma user guide) instead of the one indicated
>> above?
>> The analysis for answer the first question would be:
>>
>> design.Strain <- model.matrix(~ 0+factor(c(1,1,2,2,3,3,1,1,2,2,3,3)))
>> colnames(design.Strain) <- c("C", "D", "S")
>> #correct for batch effect
>> batch<-c(0,0,0,0,0,0,1,1,1,1,1,1)
>> design.Strain<-cbind(design.Strain,batch)
>> fit.Strain <- lmFit(selDataMatrix, design.Strain)
>> contrast.matrix.Strain <- makeContrasts(D-C, S-D, S-C,
>> levels=design.Strain)
>> fit2.Strain <- contrasts.fit(fit.Strain, contrast.matrix.Strain)
>> fit2.Strain <- eBayes(fit2.Strain)
>>
>> I have tried the two approaches, and I got different list of differential
>> expressed genes, so which one is preferred?
>>
>> Looking forward to your reply, and apologising for bothering you
>> Thank you
>> Mali
>>
>>
>>
>
>
>

	[[alternative HTML version deleted]]

