Pathview : why are some genes not colored?
0
1
Entering edit mode
okna0215 ▴ 10
@okna0215-23735
Last seen 4.4 years ago

Hello all. I am using bioconductor pathview package, and have a question. My goal is to input 'up-regulated genes' and get output of enriched pathway maps. I generated one pathview map as an example compared it to my Differentially Expressed Genes statistics, just to be sure that the map is drawn from my input DEGs.

Then I faced some doubts : - Some genes are represented in the pathway map, but not mapped (all.mapped = NA). How can this be possible? - Mapped genes are colored (ex - mol.data = 6.85822838471096, mol.col = #FF0000). Red colors seem to be based on foldchanges, but 1 gene data does not match between pathway map datat and DEG statistics. (ex - mol.data = 48.21641641, log2FoldChange = 5.477477252)

Is this kind of a bug?

Here is my script :

## Let's annotate entrez gene id, for KEGG pathway is based on entrez id
library("AnnotationDbi")
library("org.Hs.eg.db")
columns(org.Hs.eg.db)

res$symbol = mapIds(org.Hs.eg.db,
                    keys=row.names(res), 
                    column="SYMBOL",
                    keytype="ENSEMBL",
                    multiVals="first")

res$entrez = mapIds(org.Hs.eg.db,
                    keys=row.names(res), 
                    column="ENTREZID",
                    keytype="ENSEMBL",
                    multiVals="first")

res$name =   mapIds(org.Hs.eg.db,
                    keys=row.names(res), 
                    column="GENENAME",
                    keytype="ENSEMBL",
                    multiVals="first")

head(res, 10)

## setup the KEGG data-sets
library(pathview)
library(gage)
library(gageData)

data(kegg.sets.hs)
data(sigmet.idx.hs)

kegg.sets.hs = kegg.sets.hs[sigmet.idx.hs]
head(kegg.sets.hs, 3)

##### Input up DEG pathview #####################################################

## get significant res : logFoldchange >=2, p.val <= 0.05
## subset removes NAs.

res_sig = subset(res, (log2FoldChange >= 2) & pvalue <= 0.05)

foldchanges = res_sig$log2FoldChange
names(foldchanges) = res_sig$entrez
head(foldchanges)

##### Run GAGE ################################################################
keggres = gage(foldchanges, gsets=kegg.sets.hs, same.dir=TRUE)

##### Draw an example pathway map ################################################
## 1) Get data that are used to draw a pathway map
pv.out <- pathview(gene.data = foldchanges, pathway.id = '00040', speices = 'hsa', same.layer = F, limit = c(-10,10))
list.files(pattern = 'hsa04110', full.names = T)
write.csv(pv.out$plot.data.cpd, 'hsa00040_cpd_data.csv')
write.csv(pv.out$plot.data.gene, 'hsa00040_gene_data.csv')

## 2) Draw example pathway map
pathview(gene.data = foldchanges, pathway.id = '00040', speices = 'hsa', same.layer = F, limit = c(-10,10))
pathview rna-seq • 1.1k views
ADD COMMENT

Login before adding your answer.

Traffic: 702 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6