[BioC] How to create a heatmap of a gene family (using subset of exprs(vsd) in DESeq)

[Wolfgang Huber]

Dear Theresa

you can use the function 'grepl' to create a logical vector that indicates whether or not you would like to include a gene in a particular plot, e.g.

  isOR = grepl( "^OR",  xx)

where xx is, for instance, fData(cds)$SYMBOL, or whereever you have stored the gene names. Then

  select = order( isOR, rowMeans(counts(cds)), decreasing=TRUE)[1:100]

will select you the top 100 most highly expressed genes among those for which isOR is TRUE.
(More precisely, it will first take the most highly expressed genes among those for which isOR is TRUE, until these are exhausted, then continue with the ones for which it is FALSE.) Another option is thus 

  select = order(rowMeans(counts(cds)), decreasing=TRUE)[1:100]
  select = intersect(select, which(isOR))

Hope this helps
	Wolfgang 


On Jul 2, 2013, at 11:43 pm, Theresa Hodges [guest] <guest at bioconductor.org> wrote:

> 
> I have been using DESeq to look at differential gene expression and have run into a problem creating heatmap plots of my data.   I am able to create a heatmap of my top 100 expressed genes using the following code:   
> 
> select=order(rowMeans(counts(cds)), decreasing=TRUE)[1:100]
> hmcol=colorRampPalette(brewer.pale(9, "GnBu"))(100)
> heatmap.2(exprs(vsd)[select, ], col=hmcol, trace="none", margin=c(10,6))
> 
> I would also like to create a separate heatmap for 4 gene families, such as the olfactory receptors, and have already replaced the gene id with the gene name for my families of interest.  For example, all olfactory receptors' ids in my dataset of 13,500 genes now start with OR (OR1, OR50, etc). Can you tell me how I can modify the above code to create a subset of the exprs(vsd) containing just the ORs and a heatmap to display the ORs?  I am new to R so I'm not sure if this is possible.
> 
> Thank you in advance.  
> 
> Theresa
> 
> 
> -- output of sessionInfo(): 
> 
> select=order(rowMeans(counts(cds)), decreasing=TRUE)[1:100]
> hmcol=colorRampPalette(brewer.pale(9, "GnBu"))(100)
> heatmap.2(exprs(vsd)[select, ], col=hmcol, trace="none", margin=c(10,6))
> 
> --
> 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