This is a test version of Biostars. For the public version, visit https://www.biostars.org.
heatmap for differential expressed genes

how to generate a heatmap after the limma analysis (for differential expressed genes)?

fit<-lmFit(selExpr,design)
contrast.matrix<-makeContrasts(BT20_1_day-BT20_0, levels=design)
fits<-contrasts.fit(fit,contrast.matrix)
ebfit<-eBayes(fits)
topTable(ebfit,number=10000,coef=1,adjust="fdr",p.value=0.05)
differential-expressed-genes heatmap microarray r
• 12,147 views • link
updated by Ram Ram • written by minni9234 • minni9234

Have you checked the posts on the right of the screen? Just save the data.frame produced by topTable without sorting, subset the expressionSet accordingly, and make a heatmap of that.

• 2 views • link
updated by Ram Ram • written by Devon Ryan Devon Ryan

I tried that, and it was not working

• 2 views • link
written by minni9234 • minni9234

Define "not working".

• 2 views • link
written by Devon Ryan Devon Ryan

It keeps showing me an error even I have already set it to be a numeric matrix

Error in heatmap.plus(top) : 'x' must be a numeric matrix

• 3 views • link
written by minni9234 • minni9234

It would be rather more useful if you also showed the entirety of the relevant code...

• 2 views • link
written by Devon Ryan Devon Ryan

> top<-as.numeric(rownames(topTable(ebfit,n=10000,coef=1,adjust="fdr",p.value=0.05)))

>heatmap.plus(top)

• 3 views • link
written by minni9234 • minni9234

I'm not sure how you could hope for that to ever work. You're trying to make a heatmap of probe names, which doesn't make any sense.

tt <- topTable(ebfit,n=Inf,coef=1,adjust="fdr",p.value=0.05, sort.by="none")
idx <- which(tt$fdr <0.1)
subEset <- selExpr[idx]
heatmap(assayData(subEset))

or something like that.

• 3 views • link
written by Devon Ryan Devon Ryan

as.numeric of the rownames? Stopped reading there.

• 2 views • link
written by karl.stamm karl.stamm

2 answers 

fit<-lmFit(selExpr,design)
contrast.matrix<-makeContrasts(BT20_1_day-BT20_0, levels=design)
fits<-contrasts.fit(fit,contrast.matrix)
ebfit<-eBayes(fits)
# use top 200 probes
tt = topTable(ebfit,number=200,coef=1,adjust="fdr",p.value=0.05)
library(gplots)
idx = rownames(tt)
# assume that selExpr is an ExpressionSet
heatmap.2(exprs(selExpr)[idx,],trace='none',scale='row')
# OR, if selExpr is a matrix
heatmap.2(selExpr[idx,],trace='none',scale='row')
• 3 views • link
updated by Ram Ram • written by Sean Davis Sean Davis

Hi, Sean Davis. Thanks for your useful explanation.

However, I have this same errors again.

> heatmap.2(exprs(selExpr)[index_top,],trace='none',scale='row')
Error in function (classes, fdef, mtable)  : 
  unable to find an inherited method for function "exprs", for signature "matrix"
• 0 views • link
updated by Ram Ram • written by minni9234 • minni9234

Note the "assuming that selExpr is an ExpressionSet" comment... Try the next non-comment line.

• 0 views • link
written by Devon Ryan Devon Ryan

Sorry, now I have another error.

> heatmap.2(selExpr[idx,],trace='none',scale='row')

Error in heatmap.2(selExpr[idx, ], trace = "none", scale = "row") : 
  subscript out of bounds
• 0 views • link
updated by Ram Ram • written by minni9234 • minni9234

Change the tt=... and idx=... lines to these:

tt = topTable(ebfit,n=Inf,coef=1,adjust="fdr", sort.by="none")
idx = which(tt$adj.P.Value <0.1)

The original version from Sean Davis won't work if the rownames are probe ID numbers.

• 1 views • link
written by Devon Ryan Devon Ryan

Hi Devon, this doesn't work either.

> heatmap.2(selExpr[idx,],trace='none',scale='row')
Error in heatmap.2(selExpr[idx, ], trace = "none", scale = "row") : 
  `x' must have at least 2 rows and 2 columns

Here is the result I got from topTable:

     row.names    ID    logFC    AveExpr    t    P.Value    adj.P.Val    B
1    10603    HSD3B1    -7.0675833    5.68016    -48.837013    3.100100e-22    1.341568e-17    38.539428
2    40338    LOC100505633    7.2969167    3.78948    45.883284    1.064792e-21    1.616203e-17    37.668604
3    14880    DUSP6    -6.0338333    7.95048    -45.765236    1.120418e-21    1.616203e-17    37.631914
4    11308    HCAR3    5.3886667    3.31332    38.487688    3.416048e-20    3.695737e-16    35.041866
5    9566    COL15A1    -5.3935833    5.30452    -33.774838    4.453433e-19    3.854446e-15    32.944817
• 1 views • link
updated by Ram Ram • written by minni9234 • minni9234

Well, adj.P.Value should just be adj.P.Val, then.

• 0 views • link
written by Devon Ryan Devon Ryan

Hi Devon,

I have used the same to get heatmap for differential expressed genes. But, got an error.

design <- model.matrix(~ -1+factor(c(1,1,2,2,3,3)))
colnames(design) <- c("group1", "group2", "group3")
fit <- lmFit(celfiles.rma, design)
contrast.matrix <- makeContrasts(group2-group1, group3-group2, group3-group1, levels=design)
fit2 <- contrasts.fit(fit, contrast.matrix)
fit2 <- eBayes(fit2)
tab <- topTable(fit2, coef=1, adjust="fdr", sort.by="B", number=Inf)
head(tab)
             logFC   AveExpr          t      P.Value adj.P.Val         B
17100820  3.813465 10.084444  13.990013 1.477686e-05 0.3544622 -4.142070
17100824  3.813465 10.084444  13.990013 1.477686e-05 0.3544622 -4.142070
16744689 -3.037688  1.577746 -13.252989 1.983301e-05 0.3544622 -4.143913
16842465 -6.954801  4.096664 -11.968499 3.444300e-05 0.4616826 -4.147931
16843156  3.222599  2.013231  11.375102 4.530001e-05 0.4857702 -4.150236
16706416 -2.014341  2.072348  -8.395743 2.270946e-04 1.0000000 -4.169346

idx = which(tab$P.Val < 0.9)
heatmap.2(tab[idx,],trace='none',scale='row')

Error in heatmap.2(tab[idx, ], trace = "none", scale = "row") : 
  `x' must be a numeric matrix

Can you help me in this?

• 1 views • link
updated by Ram Ram • written by Vasu Vasu

Make a heatmap of the normalized counts, not the output of topTable().

• 1 views • link
written by Devon Ryan Devon Ryan

Thanks for the reply Devon. I have used oligo for my Affymetrix Human gene 2.0 ST array.

head(tab)
               logFC   AveExpr           t    P.Value adj.P.Val         B
16650001 -0.67713091 0.7013795 -2.11246203 0.07463104 0.8414740 -3.923065
16650003  0.01775323 0.8173575  0.05752413 0.95581904 0.9922585 -5.247542
16650005  0.20727356 1.2926154  0.28654513 0.78318782 0.9636612 -5.214740
16650007  0.13198594 0.8028477  0.37395467 0.72008320 0.9502862 -5.191016
16650009 -0.31117759 0.7187522 -0.92432176 0.38765340 0.8815813 -4.917247
16650011 -0.18142907 0.7450577 -0.46483921 0.65689248 0.9382330 -5.160082

results <- tab[which(tab$logFC >= 1.5 & tab$P.Value <= 0.05),]

idx = which(tab$P.Value < 0.05 & tab$logFC > 1.5)

heatmap(exprs(celfiles.rma[idx,]),trace='none',scale='row')

But I'm getting the Unit numbers. How to add annotation to get the gene names for this?

• 1 views • link
updated by Ram Ram • written by Vasu Vasu

No clue on that part, I haven't used microarrays in a number of years.

• 1 views • link
written by Devon Ryan Devon Ryan

Ok. Thank you very much !!

• 2 views • link
written by Vasu Vasu

Do you have a question that we can answer? When you dump code on us and ask "how to do" there is no specific thing we can give you. You are asking the internet to do your job. Search "making heatmaps in R" on the google or on biostars and you'll see many tutorials. You have then to read them and try it. When you get a specific question, search that, and we may be able to help.

• 1 views • link
written by karl.stamm karl.stamm

I tried that, and it was not working

• 0 views • link
written by minni9234 • minni9234

Log in to answer this question.

More posts like this