Hi! 

Sorry if this is a very basic question, but how do I add additional information from a data frame (i.e. gene symbols, chromosomal coordinates, etc) to a DESeqDataSet which I obtained from a count table with gene names as row names?

I.e. my count table looks like this: 

 head(counttable)

I'm using the following code to create the object:

groups <- c(rep("A", 3), rep("B", 3))
dds <- DESeqDataSetFromMatrix(countData = counttable, colData = as.data.frame(groups) ,design = ~ groups)

I've used biomaRt to query ensembl and obtain an annotation data frame which has MGI gene symbols and other information I'm interested in:

head(annotation)

     ensembl_gene_id mgi_symbol   gene_biotype
1 ENSMUSG00000064336      mt-Tf        Mt_tRNA
2 ENSMUSG00000064337    mt-Rnr1        Mt_rRNA
3 ENSMUSG00000064338      mt-Tv        Mt_tRNA
4 ENSMUSG00000064339    mt-Rnr2        Mt_rRNA
5 ENSMUSG00000064340     mt-Tl1        Mt_tRNA
6 ENSMUSG00000064341     mt-Nd1 protein_coding

What I am trying to do (in the end) is to add this information to the DESeqDataSet so that (for example)

- when I plot the heatmap my gene names are MGI symbols, not ensembl IDs

- when I use ReportingTools to generate the report, I have human - readable gene names and biotype classes in the output table as well

Thanks in advance !!!

R version 3.1.1 (2014-07-10)
Platform: x86_64-apple-darwin13.1.0 (64-bit)

locale:
[1] en_AU.UTF-8/en_AU.UTF-8/en_AU.UTF-8/C/en_AU.UTF-8/en_AU.UTF-8

attached base packages:
[1] parallel  grid      stats     graphics  grDevices datasets  utils     methods  
[9] base     

other attached packages:
 [1] gplots_2.14.2             RColorBrewer_1.0-5        Biobase_2.24.0           
 [4] pasilla_0.4.0             BiocInstaller_1.14.3      DESeq2_1.4.5             
 [7] RcppArmadillo_0.4.450.1.0 Rcpp_0.11.3               GenomicRanges_1.16.4     
[10] GenomeInfoDb_1.0.2        IRanges_1.22.10           BiocGenerics_0.10.0      
[13] gridExtra_0.9.1           ggplot2_1.0.0             edgeR_3.6.8              
[16] limma_3.20.9              biomaRt_2.20.0           

loaded via a namespace (and not attached):
 [1] annotate_1.42.1      AnnotationDbi_1.26.1 bitops_1.0-6         caTools_1.17.1      
 [5] colorspace_1.2-4     DBI_0.3.1            DESeq_1.16.0         digest_0.6.4        
 [9] gdata_2.13.3         genefilter_1.46.1    geneplotter_1.42.0   gtable_0.1.2        
[13] gtools_3.4.1         KernSmooth_2.23-13   lattice_0.20-29      locfit_1.5-9.1      
[17] MASS_7.3-35          munsell_0.4.2        plyr_1.8.1           proto_0.3-10        
[21] RCurl_1.95-4.3       reshape2_1.4         RSQLite_1.0.0        scales_0.2.4        
[25] splines_3.1.1        stats4_3.1.1         stringr_0.6.2        survival_2.37-7     
[29] tools_3.1.1          XML_3.98-1.1         xtable_1.7-4         XVector_0.4.0       

hi Darya,

DESeqDataSet is a subclass of SummarizedExperiment, and you can add to the metadata columns of the rowData with code like

mcols(dds) <- cbind(mcols(dds), newMcols)

However, you have to make sure that these new columns line up row for row with the DESeqDataSet:

all(rownames(dds) == annotation$ensembl_gene_id)

If this is not the case you need to reorder the annotation data.frame using the match() function, so that the Ensembl IDs match up, then check again that everything is lined up:

annotation <- annotation[match(rownames(dds), annotation$ensembl_gene_id),]

all(rownames(dds) == annotation$ensembl_gene_id)

The heatmap.2 function looks for the rownames of the matrix you give it. So you could assign the mgi_symbols as the rownames of the DESeqDataSet, but you need to make sure you have unique mgi_symbols for every ensembl_gene_id.

sum(is.na(annotation$mgi_symbol))
sum(duplicated(annotation$mgi_symbol))

If these sums above are both 0, then:

rownames(dds) <- mcols(dds)$mgi_symbol

should work.

If these sums are not 0, you can manually supply the row labels to the heatmap.2() function using the "labRow" argument. See ?heatmap.2

You can add this mgi symbol column to the results table like so:

res <- results(dds)
res$mgi_symbol <- mcols(dds)$mgi_symbol

I don't know the details on how to tell ReportingTools which columns to output into the final report. The ReportingTools vignette might have some more information.

I tried the above steps but i cant get it to work to add my features. Here is what I tried:

dds <- DESeqDataSetFromMatrix(countData = countData, colData = colData, design = ~ condition)
dds$condition <- relevel(dds$condition, "WT")
dds <- DESeq(dds)

#Add Feature data
featureData <-read.delim(file="FeatureData.txt", sep="\t", header=T)
all(rownames(dds) ==featureData$GeneID) #to confirm that they are in the same order
-->TRUE
mcols(dds) <- cbind(mcols(dds), featureData)
names(mcols(dds)) #I see the new columns added at the end

res <- results(dos) #When i do this, i get the following error:

Error in results(dds) : 
  cannot find results columns in object, first call DESeq, nbinomWaldTest, or nbinomLRT

If I do NOT add the feature, i dont get an error. hence, how can I run this with the added features, and how can i output a table with the added features include in it?

Further, how can I output a matrix with the normalized counts of each sample along with the log2fold change and p-value results, + added features? 

.