[BioC] limma for finding differentialy expressed genes of several groups

[James W. MacDonald]

Hi Priya,

On 10/17/2012 10:20 AM, priya [guest] wrote:
> I would like to find the differentially expressed genes for several variables using the limma package for several groups.
>
>   -- output of sessionInfo():
>
> I have the rma normalized matrix in the following format :
>
>
> ID_REF GSM362180    GSM362181  GSM362188    GSM362189  GSM362192
> 244901 5.094871713 4.626623079 4.554272515 4.748604391 4.759221647
> 244902 5.194528083 4.985930299 4.817426064 5.151654407 4.838741605
> 244903 5.412329253 5.352970877 5.06250609  5.305709079 8.365082403
> 244904 5.529220594 5.28134657  5.467445095 5.62968933  5.458388909
> 244905 5.024052699 4.714631878 4.792865831 4.843975286 4.657188246
> 244906 5.786557533 5.242403911 5.060605782 5.458148567 5.890061836
>
> where the different columns correspond to four different types of promoters and each of the four promoters has a biological replicate so totally there are 8 columns.
>
> I tried using the Limma package to find the differentially expressed genes across several promoters ( with replicates) and I always get an error as Iam new to r and unable to understand it fully .
>
> This is the code that I used:
>
> Group<- factor(c("p1", "p1", "p2", "p2", "p3","p3","p3","p4","p4"), levels = c("GSM362180","GSM362181","GSM362188","GSM362189","GSM362192","GSM362193","GSM362194","GSM362197","GSM362198"))

This ^^^^^^^^^^^^^ doesn't make sense. If you look at the results from 
that, you get this:

 > Group
[1] <NA> <NA> <NA> <NA> <NA> <NA> <NA> <NA> <NA>
9 Levels: GSM362180 GSM362181 GSM362188 GSM362189 GSM362192 ... GSM362198

Plus you say you have 8 samples, yet you list 9.

You most likely want something like

 > Group <- factor(c("p1", "p1", "p2", "p2", "p3","p3","p3","p4","p4"))
 > design <- model.matrix(~0+Group)
 > colnames(design) <- gsub("Group","", colnames(design))
 > design
   p1 p2 p3 p4
1  1  0  0  0
2  1  0  0  0
3  0  1  0  0
4  0  1  0  0
5  0  0  1  0
6  0  0  1  0
7  0  0  1  0
8  0  0  0  1
9  0  0  0  1
attr(,"assign")
[1] 1 1 1 1
attr(,"contrasts")
attr(,"contrasts")$Group
[1] "contr.treatment"

Then you can make your contrasts matrix and go from there.

Best,

Jim


> design<- model.matrix(~0 + Group)
>
> colnames(design)<- c("GSM362180","GSM362181","GSM362188","GSM362189","GSM362192","GSM362193","GSM362194","GSM362197","GSM362198")
> fit<- lmFit(modified, design)
>
> where modified is the rma normalized data matrix as inputted in the above format.
> I get the following error:
>
> Coefficients not estimable: GSM362180 GSM362181 GSM362188 GSM362189 GSM362192 GSM362193 GSM362194 GSM362197 GSM362198
> Error in lm.fit(design, t(M)) : 0 (non-NA) cases
>
>
> --
> Sent via the guest posting facility at bioconductor.org.
>
> _______________________________________________
> 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

-- 
James W. MacDonald, M.S.
Biostatistician
University of Washington
Environmental and Occupational Health Sciences
4225 Roosevelt Way NE, # 100
Seattle WA 98105-6099