source("http://bioconductor.org/biocLite.R")
## Bioconductor version 3.5 (BiocInstaller 1.26.0), ?biocLite for help
biocLite("topGO")
## BioC_mirror: https://bioconductor.org
## Using Bioconductor 3.5 (BiocInstaller 1.26.0), R 3.4.0 (2017-04-21).
## Installing package(s) 'topGO'
## installation path not writeable, unable to update packages: foreign
## Old packages: 'gdata'
#Be sure to install these three packages############################
source("https://bioconductor.org/biocLite.R")
## Bioconductor version 3.5 (BiocInstaller 1.26.0), ?biocLite for help
biocLite("ALL")
## BioC_mirror: https://bioconductor.org
## Using Bioconductor 3.5 (BiocInstaller 1.26.0), R 3.4.0 (2017-04-21).
## Installing package(s) 'ALL'
## installation path not writeable, unable to update packages: foreign
## Old packages: 'gdata'
source("http://bioconductor.org/biocLite.R")
## Bioconductor version 3.5 (BiocInstaller 1.26.0), ?biocLite for help
biocLite("hgu95av2.db")
## BioC_mirror: https://bioconductor.org
## Using Bioconductor 3.5 (BiocInstaller 1.26.0), R 3.4.0 (2017-04-21).
## Installing package(s) 'hgu95av2.db'
## installation path not writeable, unable to update packages: foreign
## Old packages: 'gdata'
source("http://bioconductor.org/biocLite.R")
## Bioconductor version 3.5 (BiocInstaller 1.26.0), ?biocLite for help
biocLite("org.Hs.eg.db")
## BioC_mirror: https://bioconductor.org
## Using Bioconductor 3.5 (BiocInstaller 1.26.0), R 3.4.0 (2017-04-21).
## Installing package(s) 'org.Hs.eg.db'
## installation path not writeable, unable to update packages: foreign
## Old packages: 'gdata'
#####################################################################
library(topGO)
## Loading required package: BiocGenerics
## Loading required package: parallel
##
## Attaching package: 'BiocGenerics'
## The following objects are masked from 'package:parallel':
##
## clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
## clusterExport, clusterMap, parApply, parCapply, parLapply,
## parLapplyLB, parRapply, parSapply, parSapplyLB
## The following objects are masked from 'package:stats':
##
## IQR, mad, sd, var, xtabs
## The following objects are masked from 'package:base':
##
## anyDuplicated, append, as.data.frame, cbind, colMeans,
## colnames, colSums, do.call, duplicated, eval, evalq, Filter,
## Find, get, grep, grepl, intersect, is.unsorted, lapply,
## lengths, Map, mapply, match, mget, order, paste, pmax,
## pmax.int, pmin, pmin.int, Position, rank, rbind, Reduce,
## rowMeans, rownames, rowSums, sapply, setdiff, sort, table,
## tapply, union, unique, unsplit, which, which.max, which.min
## Loading required package: graph
## Loading required package: Biobase
## Welcome to Bioconductor
##
## Vignettes contain introductory material; view with
## 'browseVignettes()'. To cite Bioconductor, see
## 'citation("Biobase")', and for packages 'citation("pkgname")'.
## Loading required package: GO.db
## Loading required package: AnnotationDbi
## Loading required package: stats4
## Loading required package: IRanges
## Loading required package: S4Vectors
##
## Attaching package: 'S4Vectors'
## The following object is masked from 'package:base':
##
## expand.grid
##
## Loading required package: SparseM
##
## Attaching package: 'SparseM'
## The following object is masked from 'package:base':
##
## backsolve
##
## groupGOTerms: GOBPTerm, GOMFTerm, GOCCTerm environments built.
##
## Attaching package: 'topGO'
## The following object is masked from 'package:IRanges':
##
## members
library(ALL)
data(ALL)
data(geneList)
affyLib <- paste(annotation(ALL), "db", sep = ".")
library(package = affyLib, character.only = TRUE)
## Loading required package: org.Hs.eg.db
##
##
sum(topDiffGenes(geneList))
## [1] 50
#Changed nnot = annFUN.db to annotationFun = annFUN.db
#in response to "Error in annotationFun(ontology, .Object@allGenes, ...) :
#argument "annotationFun" is missing, with no default"
sampleGOdata <- new("topGOdata", description = "Simple session", ontology = "BP", allGenes = geneList, geneSel = topDiffGenes, nodeSize = 10, annotationFun = annFUN.db, affyLib = affyLib)
##
## Building most specific GOs .....
## ( 1555 GO terms found. )
##
## Build GO DAG topology ..........
## ( 4396 GO terms and 10270 relations. )
##
## Annotating nodes ...............
## ( 310 genes annotated to the GO terms. )
sampleGOdata
##
## ------------------------- topGOdata object -------------------------
##
## Description:
## - Simple session
##
## Ontology:
## - BP
##
## 323 available genes (all genes from the array):
## - symbol: 1095_s_at 1130_at 1196_at 1329_s_at 1340_s_at ...
## - score : 1 1 0.62238 0.541224 1 ...
## - 50 significant genes.
##
## 310 feasible genes (genes that can be used in the analysis):
## - symbol: 1095_s_at 1130_at 1196_at 1329_s_at 1340_s_at ...
## - score : 1 1 0.62238 0.541224 1 ...
## - 46 significant genes.
##
## GO graph (nodes with at least 10 genes):
## - a graph with directed edges
## - number of nodes = 1064
## - number of edges = 2388
##
## ------------------------- topGOdata object -------------------------
resultFisher <- runTest(sampleGOdata, algorithm = "classic", statistic = "fisher")
##
## -- Classic Algorithm --
##
## the algorithm is scoring 950 nontrivial nodes
## parameters:
## test statistic: fisher
resultFisher
##
## Description: Simple session
## Ontology: BP
## 'classic' algorithm with the 'fisher' test
## 1064 GO terms scored: 17 terms with p < 0.01
## Annotation data:
## Annotated genes: 310
## Significant genes: 46
## Min. no. of genes annotated to a GO: 10
## Nontrivial nodes: 950
resultKS <- runTest(sampleGOdata, algorithm = "classic", statistic = "ks")
##
## -- Classic Algorithm --
##
## the algorithm is scoring 1064 nontrivial nodes
## parameters:
## test statistic: ks
## score order: increasing
resultKS.elim <- runTest(sampleGOdata, algorithm = "elim", statistic = "ks")
##
## -- Elim Algorithm --
##
## the algorithm is scoring 1064 nontrivial nodes
## parameters:
## test statistic: ks
## cutOff: 0.01
## score order: increasing
##
## Level 15: 3 nodes to be scored (0 eliminated genes)
##
## Level 14: 8 nodes to be scored (0 eliminated genes)
##
## Level 13: 14 nodes to be scored (0 eliminated genes)
##
## Level 12: 28 nodes to be scored (0 eliminated genes)
##
## Level 11: 40 nodes to be scored (21 eliminated genes)
##
## Level 10: 57 nodes to be scored (31 eliminated genes)
##
## Level 9: 101 nodes to be scored (31 eliminated genes)
##
## Level 8: 131 nodes to be scored (43 eliminated genes)
##
## Level 7: 163 nodes to be scored (84 eliminated genes)
##
## Level 6: 177 nodes to be scored (200 eliminated genes)
##
## Level 5: 160 nodes to be scored (200 eliminated genes)
##
## Level 4: 120 nodes to be scored (214 eliminated genes)
##
## Level 3: 43 nodes to be scored (216 eliminated genes)
##
## Level 2: 18 nodes to be scored (216 eliminated genes)
##
## Level 1: 1 nodes to be scored (216 eliminated genes)
allRes <- GenTable(sampleGOdata, classicFisher = resultFisher, classicKS = resultKS, elimKS = resultKS.elim, orderBy = "elimKS", ranksOf = "classicFisher", topNodes = 10)
pValue.classic <- score(resultKS)
pValue.elim <- score(resultKS.elim)[names(pValue.classic)]
gstat <- termStat(sampleGOdata, names(pValue.classic))
gSize <- gstat$Annotated / max(gstat$Annotated) * 4
#Defined colMap, ref. https://github.com/Bioconductor-mirror/topGO/blob/master/vignettes/topGO.Rnw
colMap <- function(x) {
.col <- rep(rev(heat.colors(length(unique(x)))), time = table(x))
return(.col[match(1:length(x), order(x))])
}
gCol <- colMap(gstat$Significant)
plot(pValue.classic, pValue.elim, xlab = "p-value classic", ylab = "p-value elim", pch = 19, cex = gSize, col = gCol)
