Limma package for gene expression

1.1k views Asked by At

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"))

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 lmfit(design, t(M)) : 0 (non-NA) cases

1

There are 1 answers

0
Djork On

Should just be:

Group <- factor(c("p1", "p1", "p2", "p2", "p3", "p3", "p3", "p4", "p4"))

The group category order should correspond to the order of samples in colnames(modified)

In your original code, you assigned different unique levels to each factor in Group so essentially you have no grouping category, no replicates and cannot estimate coefficients.

Create your design matrix:

design <- model.matrix(~0 + Group)

You then want to assign name of samples to rownames of design matrix, not colnames as columns represent the grouping.

rownames(design) <- c("GSM362180","GSM362181","GSM362188","GSM362189","GSM362192","GSM362193","GSM362194","GSM362197","GSM362198")
# or simply
rownames(design) <- colnames(modified)

Your design matrix should look like this, with a 1 where a sample belongs to specific group category.

design
          Groupp1 Groupp2 Groupp3 Groupp4
GSM362180       1       0       0       0
GSM362181       1       0       0       0
GSM362188       0       1       0       0
GSM362189       0       1       0       0
GSM362192       0       0       1       0
GSM362193       0       0       1       0
GSM362194       0       0       1       0
GSM362197       0       0       0       1
GSM362198       0       0       0       1
attr(,"assign")
[1] 1 1 1 1
attr(,"contrasts")
attr(,"contrasts")$Group
[1] "contr.treatment"