[BioC] Limma - gls.series - design matrix?

[Gordon Smyth]

Your problem is with your code

design <- model %*% contrast.matrix

which is entirely unnecessary and produces a singular design matrix. The 
correct design matrix is computed by modelMatrix -- why don't you trust it? 
The design matrix should have only 4 columns because you have 4 distinct 
RNA sources apart from the reference. The formation of contrasts is done 
after the fitting of the linear model as per the limma documentation and 
examples.

In your example the last three constrasts are always NA because any set of 
3 contrasts span the contrast space - the last 3, whichever they are, are 
redundant.

Gordon

At 04:11 PM 11/08/2004, Jakob Hedegaard wrote:
>Hi
>
>We are studying the impact of different traits on the expression in pigs 
>using cDNA microarray. We have some problems using Limma when taking 
>within-array replicated spots into acount. Analysing the effects of a 
>single trait (like sick-healthy) works fine, but combing traits (like 
>far6syg-far2syg, far=father, syg=sick) results in the presence of NA's in 
>the output from gls.series. We are using R ver 1.9.0 and Bioconductor 1.4.
>
>What we do:
>
> > targets <- readTargets("targets_farXsyg.txt")
> > targets
>    SlideNumber Cy3      Cy5
>1     12755473 mix far6rask
>2     12755474 mix  far6syg
>3     12755475 mix far2rask
>4     12755476 mix  far2syg
>5     12755477 mix far6rask
>.......
>25    12759971 mix  far2syg
>26    12760017 mix  far6syg
>27    12760018 mix far6rask
>28    12760019 mix  far6syg
> >
> > model <- modelMatrix(targets, ref="mix")
>Found unique target names:
>  far2rask far2syg far6rask far6syg mix
> > model
>    far2rask far2syg far6rask far6syg
>1         0       0        1       0
>2         0       0        0       1
>3         1       0        0       0
>4         0       1        0       0
>5         0       0        1       0
>...........
>25        0       1        0       0
>26        0       0        0       1
>27        0       0        1       0
>28        0       0        0       1
> >
> > contrast.matrix <- makeContrasts(far6rask-far6syg, far2rask-far2syg, 
> far6rask-far2syg, far6syg-far2rask, far6rask-far2rask, far6syg-far2syg, 
> levels=model)
> > contrast.matrix
>          far6rask - far6syg far2rask - far2syg far6rask - far2syg
>far2rask                  0                  1                  0
>far2syg                   0                 -1                 -1
>far6rask                  1                  0                  1
>far6syg                  -1                  0                  0
>          far6syg - far2rask far6rask - far2rask far6syg - far2syg
>far2rask                 -1                  -1                 0
>far2syg                   0                   0                -1
>far6rask                  0                   1                 0
>
>
> > contrast.matrix
>          far6rask - far6syg far2rask - far2syg far6rask - far2syg
>far2rask                  0                  1                  0
>far2syg                   0                 -1                 -1
>far6rask                  1                  0                  1
>far6syg                  -1                  0                  0
>          far6syg - far2rask far6rask - far2rask far6syg - far2syg
>far2rask                 -1                  -1                 0
>far2syg                   0                   0                -1
>far6rask                  0                   1                 0
>far6syg                   1                   0                 1
> >
> > design <- model %*% contrast.matrix
> > design
>
>      far6rask - far6syg far2rask - far2syg far6rask - far2syg
>   1                   1                  0                  1
>   2                  -1                  0                  0
>   3                   0                  1                  0
>   4                   0                 -1                 -1
>   5                   1                  0                  1
>  ............
>   25                  0                 -1                 -1
>   26                 -1                  0                  0
>   27                  1                  0                  1
>   28                 -1                  0                  0
>
>      far6syg - far2rask far6rask - far2rask far6syg - far2syg
>   1                   0                   1                 0
>   2                   1                   0                 1
>   3                  -1                  -1                 0
>   4                   0                   0                -1
>   5                   0                   1                 0
>...........
>   25                  0                   0                -1
>   26                  1                   0                 1
>   27                  0                   1                 0
>   28                  1                   0                 1
> > cor <- duplicateCorrelation(MArep$M,design,ndups=4)
> > fitcor <- gls.series(MArep$M, design,ndups=4,correlation=cor$cor)
> > fitcor$coefficients[1:5,]
>      far6rask - far6syg far2rask - far2syg far6rask - far2syg
>[1,]        0.019897478       -0.058434294        -0.06046312
>[2,]       -0.006362523        0.061185620        -0.11697534
>[3,]       -0.003859953       -0.050241979        -0.02122334
>[4,]        0.013675174       -0.003196691        -0.02909511
>[5,]        0.005154980        0.043414444        -0.04125395
>      far6syg - far2rask far6rask - far2rask far6syg - far2syg
>[1,]                 NA                  NA                NA
>[2,]                 NA                  NA                NA
>[3,]                 NA                  NA                NA
>[4,]                 NA                  NA                NA
>[5,]                 NA                  NA                NA
> >
>
>For some reason (?) it is always the final three contrasts that results in 
>NA´s - chancing the design matrix by flipping the final three contrast 
>with the first three contrasts, still results in NA's in the final three 
>contrast (which were the the first three before flipping......). Some 
>technical problem? The design matrix ("design") are of dim 28x6 as it it 
>should be......
>Any suggestions?
>
>------------------------------------------------------
>Jakob Hedegaard
>
>Danish Institute of Agricultural Sciences
>Department of Animal Breeding and Genetics
>Research Centre Foulum
>P.O. Box 50
>DK-8830 Tjele, Denmark
>
>Tel: (+45) 8999 1363
>Fax: (+45) 8999 1300
>
>_______________________________________________
>Bioconductor mailing list
>Bioconductor at stat.math.ethz.ch
>https://stat.ethz.ch/mailman/listinfo/bioconductor