-
Notifications
You must be signed in to change notification settings - Fork 0
/
revision_ORA_COSMOS.R
72 lines (53 loc) · 2.38 KB
/
revision_ORA_COSMOS.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
library(readr)
library(biomaRt)
library(piano)
setwd("~/Dropbox/COSMOS/")
source("scripts/revision_COSMOS_functions.R")
source("scripts/revision_piano_function.R")
#Load inputs
meta_network <- as.data.frame(
read_csv("support/meta_network_carnival_ready_exch_solved_fullomni_metfiltered_expfiltered.csv"))
prots <- unique(c(meta_network$source,meta_network$target))
prots <- prots[!grepl("XMetab",prots)]
prots <- gsub("^X","",prots)
prots <- gsub("Gene[0-9]+__","",prots)
prots <- gsub("_reverse","",prots)
prots <- gsub("EXCHANGE.*","",prots)
prots <- unique(prots)
prots <- unlist(sapply(prots, function(x){unlist(strsplit(x,split = "_"))}))
gene_mapping <- "ensembl"
if(gene_mapping == "ensembl")
{
ensembl = useEnsembl(biomart="ensembl", dataset="hsapiens_gene_ensembl")
G_list <- getBM(filters = "entrezgene_id",
attributes = c('hgnc_symbol','entrezgene_id', "description"),
values = prots, mart = ensembl)
} else
{
G_list <- as.data.frame(read_csv(gene_mapping))
}
gene_mapping <- G_list[,1]
names(gene_mapping) <- G_list[,2]
background <- sapply(prots, function(x, translation_dictionary){
return(translation_dictionary[x])
},translation_dictionary = gene_mapping)
names(background) <- 1:length(background)
background <- unique(background)
hallmarks <- gmt_to_df("support/h.all.v7.1.symbols.gmt_20200709.gmt")
cosmos_att <- as.data.frame(
read_csv("results_revisions/COSMOS_cluster_runs/run_fullomni/met_exp_filtered_continuous_input_sigfilter/subnet_combined_att_full_newDoro_long.csv"))
cosmos_att <- cosmos_att[cosmos_att$Nodes %in% background,]
cosmos_att <- unique(cosmos_att$Nodes)
hallmarks_ora <- runGSAhyper(genes = cosmos_att, universe = background, gsc = loadGSC(hallmarks))
hallmarks_ora <- as.data.frame(hallmarks_ora$resTab)
hallmarks_ora$pathway <- row.names(hallmarks_ora)
cgp <- gmt_to_df("support/c2.cgp.v7.1.symbols.gmt")
cgp <- cgp[cgp$gene %in% background,]
cgp_ora <- runGSAhyper(genes = cosmos_att, universe = background, gsc = loadGSC(cgp))
cgp_ora <- as.data.frame(cgp_ora$resTab)
cgp_ora$pathway <- row.names(cgp_ora)
cp <- gmt_to_df("support/c2.cp.v7.1.symbols.gmt")
cp <- cp[cp$gene %in% background,]
cp_ora <- runGSAhyper(genes = cosmos_att, universe = background, gsc = loadGSC(cp))
cp_ora <- as.data.frame(cp_ora$resTab)
cp_ora$pathway <- row.names(cp_ora)