-
Notifications
You must be signed in to change notification settings - Fork 4
Example Script
Ewoud Ewing edited this page Jul 17, 2020
·
2 revisions
require(devtools)
install_github("https://github.com/TranslationalBioinformaticsUnit/GeneSetCluster/")
require(GeneSetCluster)
Canonical.files <- c(system.file("extdata", "MM10.IPA.KO.uGvsMac.Canonical_pathways.xls", package = "GeneSetCluster"),
system.file("extdata", "MM10.IPA.WT.uGvsMac.Canonical_pathways.xls", package = "GeneSetCluster"))
require(readxl)
MM10.IPA.KO.uGvsMac.Canonical <- read_excel(path = Canonical.files[1], skip = 1, sheet = 1)
MM10.IPA.WT.uGvsMac.Canonical <- read_excel(path = Canonical.files[2], skip = 1, sheet = 1)
MM10.IPA.KO.uGvsMac.Canonical <- as.data.frame(MM10.IPA.KO.uGvsMac.Canonical)
MM10.IPA.WT.uGvsMac.Canonical <- as.data.frame(MM10.IPA.WT.uGvsMac.Canonical)
require(limma)
MM10.IPA.KO.uGvsMac.Canonical$MoleculesCount <- NA
for(can.i in 1:nrow(MM10.IPA.KO.uGvsMac.Canonical))
{
mol <- as.vector(strsplit2(MM10.IPA.KO.uGvsMac.Canonical[can.i,"Molecules"], split=","))
MM10.IPA.KO.uGvsMac.Canonical[can.i,"MoleculesCount"] <- length(mol)
}
head(MM10.IPA.KO.uGvsMac.Canonical)
MM10.IPA.WT.uGvsMac.Canonical$MoleculesCount <- NA
for(can.i in 1:nrow(MM10.IPA.WT.uGvsMac.Canonical))
{
mol <- as.vector(strsplit2(MM10.IPA.WT.uGvsMac.Canonical[can.i,"Molecules"], split=","))
MM10.IPA.WT.uGvsMac.Canonical[can.i,"MoleculesCount"] <- length(mol)
}
####################
MM10.IPA.KO.uGvsMac.Canonical.sig <- MM10.IPA.KO.uGvsMac.Canonical[MM10.IPA.KO.uGvsMac.Canonical$`-log(p-value)` > 1.31 & MM10.IPA.KO.uGvsMac.Canonical$MoleculesCount > 5,]
MM10.IPA.WT.uGvsMac.Canonical.sig <- MM10.IPA.WT.uGvsMac.Canonical[MM10.IPA.WT.uGvsMac.Canonical$`-log(p-value)` > 1.31 & MM10.IPA.WT.uGvsMac.Canonical$MoleculesCount > 5,]
####################
IPA.KOvsWT <- ObjectCreator(Pathways = c(MM10.IPA.KO.uGvsMac.Canonical.sig$`Ingenuity Canonical Pathways`,
MM10.IPA.WT.uGvsMac.Canonical.sig$`Ingenuity Canonical Pathways`),
Molecules = c(MM10.IPA.KO.uGvsMac.Canonical.sig$Molecules,
MM10.IPA.WT.uGvsMac.Canonical.sig$Molecules),
Groups = c(rep("KO", times=nrow(MM10.IPA.KO.uGvsMac.Canonical.sig)),
rep("WT", times=nrow(MM10.IPA.WT.uGvsMac.Canonical.sig))),
Source = "IPA", Type = "Canonical", structure = "SYMBOL", sep = ",", organism = "org.Mm.eg.db")
IPA.KOvsWT <- CombineGeneSets(IPA.KOvsWT)
require(ggplot2)
OptimalGeneSets(object = IPA.KOvsWT, method = "gap", max_cluster = 15, cluster_method = "kmeans", main = "IPA.KOvsWT")
IPA.KOvsWT <- ClusterGeneSets(Object = IPA.KOvsWT, clusters = 6,method = "kmeans", order = "cluster" )
PlotGeneNetworks(IPA.KOvsWT, labels = F, RRmin = 5)
PlotGeneSets(Object = IPA.KOvsWT, annotation.mol = F, main = "IPA.KOvsWT with 10 kmeans clusters")
PlotGeneSets(Object = IPA.KOvsWT, annotation.mol = F, main = "IPA.KOvsWT with 10 kmeans clusters", RR.max = 90, cluster.order = c(1,8,3,4,5,6,7,9,10,2))
########################
require(WebGestaltR)
IPA.KOvsWT.ORA.all <- ORAGeneSets(Object = IPA.KOvsWT, ORA.returned = 10, unique.per.cluster = F)
View(IPA.KOvsWT.ORA.all)
########################
IPA.KOvsWT.String.all <- GetSTRINGdbPerGeneSets(Object = IPA.KOvsWT, unique.per.cluster = F)
par(mfrow=c(2,3))
plotSTRINGdbPerGeneSets(StringObject = IPA.KOvsWT.String.all, plot.cluster = 1)
plotSTRINGdbPerGeneSets(StringObject = IPA.KOvsWT.String.all, plot.cluster = 2)
plotSTRINGdbPerGeneSets(StringObject = IPA.KOvsWT.String.all, plot.cluster = 3)
plotSTRINGdbPerGeneSets(StringObject = IPA.KOvsWT.String.all, plot.cluster = 4)
plotSTRINGdbPerGeneSets(StringObject = IPA.KOvsWT.String.all, plot.cluster = 5)
plotSTRINGdbPerGeneSets(StringObject = IPA.KOvsWT.String.all, plot.cluster = 6)
Example Script: Example
Step 1A: Loading the data
Step 1B: Creating an Object
Step 2: Combine and Cluster
Step 2B: User supplied distance function
Step 2C: Highlighting-Genes
Step 3: Exporting Data
Step 4: Functional Investigation
Video: Step-by-step user guide