noriakis / ggkegg

Analyzing and visualizing KEGG information using the grammar of graphics
https://noriakis.github.io/software/ggkegg
MIT License
201 stars 16 forks source link

Does the "rawMap" function support "enricher" object? #13

Closed jijitoutou closed 6 months ago

jijitoutou commented 7 months ago

I collected KEGG gene sets (KEGG_gene_entrezid.txt) from the KEGG database, and I conducted a ORA analysis using the "enricher" function from the "clusterProfiler" package. My codes are as below: `data <- read.delim("KEGG_gene_entrezid.txt", header=FALSE) data msigdb <- data.frame( term = sub(".:", "",data$V2), gene = sub(".:", "", data$V1) ) %>% filter( term %in% filtered_results$ID) head(msigdb) cp<- enricher(c("1029","4171"),

universe=gene_all,

                    pAdjustMethod = "fdr",
                    TERM2GENE = msigdb)

rawMap(list(cp,cp), pid="hsa04110")`

However, no gene was highlighted in the picture.

But this works well: `if (require("clusterProfiler")) { cp <- enrichKEGG(c("1029","4171"))

Multiple class object can be passed by list

rawMap(list(cp,cp), pid="hsa04110") }`

Does the "rawMap" function support "enricher" object? If not, could you support this object in the later version of the "ggKEGG" package? Thanks a lot.

noriakis commented 7 months ago

Hi @jijitoutou, thanks for raising the issue. enricher function outputs enrichResult class, which is supported in rawMap. The code below runs without errors with the genes highlighted.

library(ggkegg)
library(clusterProfiler)
library(gson)
kegg <- read.gson("KEGG_human.gson") ## https://github.com/YuLab-SMU/gson-files
data(geneList, package="DOSE")
de <- names(geneList)[abs(geneList) > 2]
k1 <- enricher(de, gson=kegg)
rawMap(k1)

Could you please provide the KEGG_gene_entrezid.txt used in your example?

jijitoutou commented 6 months ago

KEGG_gene_entrezid.txt Sure. Here is my file. Thank you for your help.

This table was downloaded using the following Linux code: `wget -c "http://rest.kegg.jp/list/organism" -O species

d=$(awk '{print $2}' species | head -n 1)

echo $d mkdir $d wget "http://rest.kegg.jp/list/pathway/$d" -O ./$d/map_pathway wget "http://rest.kegg.jp/list/$d" -O ./$d/gene_symbol wget "http://rest.kegg.jp/link/pathway/$d" -O ./$d/gene_map wget "http://rest.kegg.jp/link/pathway/$d" -O ./$d/gene_entrezid wget "https://www.genome.jp/kegg-bin/download_htext?htext=${d}00001.keg&format=htext&filedir=" -O ./$d/00001.keg | perl -e 'while(<>){chomp;if(/^D\s+(\w+)\s.\t(\w+)\s+./){$g=$1;$ `

noriakis commented 6 months ago

Thank you for sharing the file.

The rawMap function (and append_cp function) uses organism ID to make gene IDs compatible with those in the pathway. enricher results in UNKNOWN organism ID, so we should manually insert the prefix (in this case hsa:). Also, I added infer option in rawMap, so that even when the organism ID is UNKNOWN, the function tries to infer the organism from pathway ID and append the prefix automatically (default to FALSE). This option is in devel branch.

So in your case, these two codes would work.

data <- read.delim("KEGG_gene_entrezid.txt", header=FALSE)
msigdb <- data.frame( term = sub(".*:", "",data$V2), gene =  sub(".*:", "", data$V1)) 
cp<- enricher(c("1029","4171"), pAdjustMethod = "fdr", TERM2GENE = msigdb)
rawMap(cp, infer=TRUE) ## In devel branch
data <- read.delim("KEGG_gene_entrezid.txt", header=FALSE)
## No renaming in gene column
msigdb <- data.frame( term = sub(".*:", "",data$V2), gene =  data$V1) 
cp <- enricher(c("hsa:1029","hsa:4171"), pAdjustMethod = "fdr", TERM2GENE = msigdb)
rawMap(cp)
jijitoutou commented 6 months ago

Thank you for your tutoring. I see how the codes work. Write the entreziID as the format "organism abbreviation:number".

load("kidney_target.Rdata") load("human_KEGG.Rdata")

gene_kidney = paste0("hsa:",bitr(kidney_target$Symbol, #dataset fromType="SYMBOL", #SYMBOL toType="ENTREZID", #ENTERZID OrgDb="org.Hs.eg.db")$ENTREZID) #hsa R.utils::setOption( "clusterProfiler.download.method",'auto' )

data <- read.delim("KEGG_gene_entrezid.txt", header=FALSE) msigdb <- data.frame( term = sub(".*:", "",data$V2), gene = data$V1)%>% filter( term %in% filtered_results$ID)

enrich_kidney <- enricher(gene_kidney,

universe=gene_all,

                      pAdjustMethod = "fdr",
                      TERM2GENE = msigdb)

library(ggkegg)

hsa04530 TJ pathway---

g1 <- pathway("hsa04530") |> mutate(fibrosis=append_cp(enrich_kidney, how="all"), converted_name=convert_id("hsa")) gg <- ggraph(g1, layout="manual", x=x, y=y)+ ggfx::with_outer_glow( geom_node_rect(aes(filter=fibrosis), color="gold", fill="transparent"), colour="gold", expand=5, sigma=10)+ geom_node_rect(aes(fill=fibrosis, filter=type=="gene"))+

geom_node_rect(aes(fill=rptec, xmin=x, filter=type=="gene")) +

overlay_raw_map("hsa04530", transparent_colors = c("#cccccc","#FFFFFF","#BFBFFF","#BFFFBF"))+ scale_fill_manual(values=c("deepskyblue2","tomato"), name="Renal fibrosis")+ theme_void()+ theme(legend.title = element_text(color="black",face="bold"))

gg ggsave("KEGG_TJ.svg", width=7,height = 5)

image

noriakis commented 6 months ago

Thank you for reporting back, and I'm glad that it worked. I'm closing the issue now, but feel free to reopen it if you encounter other problems.