Differential expression analysis using Limma (Not on cel files)

Question

[Bioinformatist Newbie]

Hi,

Instead of using limma on cel files downloaded from GEO, I now have gene expression data for differential expression of genes under various condition (8 replicates) in txt file in following format:

Gene Name      Condition1    Condition 2 ..... Condition 8
Gene A             11.235          10.4657 ......   9.34567
Gene B             15.324           11.334  ......  16.3245

I was following this tutorial but as it uses cel files as input so I was encountering problem in fitting the linear model on my expression dataset.

Can someone guide me how I can further proceed with limma if I have data in such format?

Thanks for your suggestions.

[Manvendra Singh]

are you looking for differential expression between two conditions or relative expression of genes in each condition?

[Bioinformatist Newbie]

looking for differential expression between two conditions and for each condition I have 4 replicates.

Answer 1

[GouthamAtla]

Here is Limma tutorial for RNA-Seq or for a count matrix. Chapter 15. Here is a case study.

But what are those numbers in your data? Already normalised?

[Bioinformatist Newbie]

Yes, log2 normalized expression values.

[Bioinformatist Newbie]

Geek_y you didn't tell me what to do in case my data is already log2 normalized expression values...!!

Answer 2

[Bioinformatist Newbie]

So, here is the solution:

library(simpleaffy)
library(limma)
setwd("path to your working directory where series_matrix_file is located")

xx <-as.matrix(read.table(file.path("GSE5007_series_matrix.txt"), skip=52, header=TRUE, sep="\t",row.names=1))

minimalSet <- ExpressionSet(assayData=xx)
minimalSet
ex <- exprs(minimalSet)
head(ex)
qx <- as.numeric(quantile(ex, c(0., 0.25, 0.5, 0.75, 0.99, 1.0), na.rm=T))
LogC <- (qx[5] > 100) ||
           (qx[6]-qx[1] > 50 && qx[2] > 0) ||
          (qx[2] > 0 && qx[2] < 1 && qx[4] > 1 && qx[4] < 2)
if (LogC) { ex[which(ex <= 0)] <- NaN
 exprs(minimalSet) <- log2(ex) }
groups<-as.factor(c(rep("Control",3),rep("Treated",3)))
groups
design<-model.matrix(~0+groups)
colnames(design)=levels(groups)
design
fit<-lmFit(minimalSet, design)
cont.matrix<-makeContrasts(Treated-Control, levels=design)
fit2<-contrasts.fit(fit, cont.matrix)
ebfit<-eBayes(fit2)
topTable(ebfit, coef=1)
nrow(topTable(ebfit, coef=1, number=25000, lfc=1)) 
#nrow(topTable(ebfit, coef=1, adjust="fdr", number=25000, lfc=1))</p>