Multiple samples on the same graph
1
0
Entering edit mode
Luiz • 0
@luiz-13209
Last seen 3.3 years ago
Brazil

Hi, good morning!

I hope you have a nice weekend.

I'm working on pathview today and i can't put different matrixes on the graph to compare. Only one is working and my version is 1.14.0. Already run the demo code and everything is fine accord with vignette.

I'm trying to do this : https://bioconductor.org/packages/devel/bioc/vignettes/pathview/inst/doc/pathview.pdf pg14

Workflow

setwd("~/R")
require(dplyr)
require(DESeq2)
require(ggplot2)
require(ggrepel)
require(pheatmap)
library(pathview)
library(gage)
library(gageData)
data(kegg.sets.hs)
data(sigmet.idx.hs)
data(gse16873.d)

kegg.sets.hs = kegg.sets.hs[sigmet.idx.hs]
head(kegg.sets.hs, 3)
bckCountTable <- read.table("Matrix.txt", header=TRUE, row.names=1)
samples <- data.frame(row.names=c("R1", "R2", "R3","IR1","IR2", "IR3"), condition=as.factor(c(rep("R",3),rep("IR",3))))
samples
bckCDS <- DESeqDataSetFromMatrix(countData = bckCountTable, colData=samples, design=~condition)
bckCDS <- bckCDS[ rowSums(counts(bckCDS)) > 1, ]
bckCDS_1 <- DESeq(bckCDS)
res <- results(bckCDS_1, alpha=0.05, betaPrior = F, contrast=c("condition", "IR", "R"), pAdjustMethod = "BH")
res$symbol = mapIds(org.Hs.eg.db,
                    keys=row.names(res), 
                    column="SYMBOL",
                    keytype="REFSEQ",
                    multiVals="first")
res$entrez = mapIds(org.Hs.eg.db,
                    keys=row.names(res), 
                    column="ENTREZID",
                    keytype="REFSEQ",
                    multiVals="first")
res$name = mapIds(org.Hs.eg.db,
                    keys=row.names(res), 
                    column="GENENAME",
                    keytype="REFSEQ",
                    multiVals="first")
foldchanges = res$log2FoldChange
names(foldchanges) = res$entrez
head(foldchanges)
head(res,10)
foldchanges = res$log2FoldChange
names(foldchanges) = res$entrez
head(foldchanges)
keggres = gage(foldchanges, gsets=kegg.sets.hs, same.dir=TRUE)
lapply(keggres, head)

bckCountTable2 <- read.table("Matrix2.txt", header=TRUE, row.names=1)
samples2 <- data.frame(row.names=c("CTL1", "CTL2", "CTL3","HT1","HT2", "HT3"), condition=as.factor(c(rep("CTL",3),rep("HT",3))))
samples2
bckCDS2 <- DESeqDataSetFromMatrix(countData = bckCountTable2, colData=samples2, design=~condition)
bckCDS2 <- bckCDS2[ rowSums(counts(bckCDS2)) > 1, ]
bckCDS_12 <- DESeq(bckCDS2)
res2 <- results(bckCDS_12, alpha=0.05, betaPrior = F, contrast=c("condition", "HT", "CTL"), pAdjustMethod = "BH")
res2$symbol = mapIds(org.Hs.eg.db,
                    keys=row.names(res2), 
                    column="SYMBOL",
                    keytype="REFSEQ",
                    multiVals="first")
res2$entrez = mapIds(org.Hs.eg.db,
                    keys=row.names(res2), 
                    column="ENTREZID",
                    keytype="REFSEQ",
                    multiVals="first")
res2$name = mapIds(org.Hs.eg.db,
                    keys=row.names(res2), 
                    column="GENENAME",
                    keytype="REFSEQ",
                    multiVals="first")
foldchanges2 = res2$log2FoldChange
names(foldchanges2) = res2$entrez
head(foldchanges2)
head(res2,10)
foldchanges2 = res2$log2FoldChange
names(foldchanges2) = res2$entrez
head(foldchanges2)
keggres2 = gage(foldchanges2, gsets=kegg.sets.hs, same.dir=TRUE)
lapply(keggres2, head)

Then i proceed to:

pv.out <- pathview(gene.data = foldchanges, cpd.data = foldchanges2, pathway.id = "any", + species = "hsa", out.suffix = "any", keys.align = "y", + kegg.native = T, match.data = F, multi.state = T, same.layer = T) 

But i get a graph with only the gene.data, like cpd.data don't exist.

Thank you for your attention!

 

Pathview probe • 834 views
ADD COMMENT
0
Entering edit mode
Luo Weijun ★ 1.5k
@luo-weijun-1783
Last seen 4 weeks ago
United States

You are not calling pathview the right way. In your analysis, you have two fold change vectors from two analyses, i.e. foldchanges and foldchanges2. You should combined them into a single matrix with two columns, say fc.mat, then pass that into argument gene.data, something like:

pv.out <- pathview(gene.data = fc.mat, pathway.id = "any", species = "hsa", out.suffix = "any")

 

the row names of fc.mat should be gene ID (Entrez in your case).

note cpd.data argument is for compound or metabolite data, not for addition gene data (as you did in your code). For details of pathview, do:

?pathview

 

You may always use pathview web server as an alternative:

https://pathview.uncc.edu/

Pathview Web server provides a intuitive GUI and API to Pathview, a full pathway analysis workflow supporting multiple omics data and integrated analysis, and much more.

 

ADD COMMENT

Login before adding your answer.

Traffic: 302 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