Hello,

Here is some code I am using to run topGO. Essentially starting from text
files which contain one gene list each (no headers, list name from file
nama). This may be a bit excessive, but it seems to work for me. You need a
feature data (fd) matrix / data frame which has at least the Columns Entrez
(that contains Entrez ids for your array) and Symbol (that contains gene
symbols corresponding to the Entrez ids) for your entire array.

BR, Pekka

## Reading in gene lists, forming a glist

readGeneList <- function(folder, pattern, recursive=T)
{
  gl.files <- list.files(path=folder, pattern=pattern, recursive=recursive)
  names <- gsub(".txt", "",basename(gl.files))

  glist <- lapply(gl.files, function(f){
    cat("Reading:", f, "\n")
    read.delim(f, stringsAsFactors=F, header=F)[,1]
  })
  names(glist) <- names
  glist
}

# Removing white-spaces and other problematic characters
cleanIds <- function(ids)
{
  ids <- gsub(" ", "", ids)
  ids <- gsub(",", "", ids)
  ids <- gsub(";", "", ids)
  ids <- gsub("\'", "", ids)
  ids <- gsub("\"", "", ids)
  ids <- ids[!is.na(ids)]
  ids <- ids[!ids =="0"]
  ids <- as.character(unique(ids))
}

## function assumes that using gene symbols: ##
# fd = feature data containing at least mapping between allg for Symbol,
Entrez
# GObranch = "BP", "MF", "CC"

mapGenesToTerms <- function(fd, GObranch="BP")
{
  require(org.Hs.eg.db)
  keys <- cleanIds(as.character(unique(fd$Entrez)))
  cols <- c("GO", "SYMBOL")
  mapping <- select(org.Hs.eg.db, keys, cols, keytype="ENTREZID")
  mapping.bp <- mapping[mapping$ONTOLOGY==GObranch,]
  go.list <- by(mapping.bp[!mapping.bp$SYMBOL=="",],
INDICES=mapping.bp$SYMBOL[!mapping.bp$SYMBOL==""], function(x){x$GO})
  go.list
}

## function assumes that using gene symbols: ##
# Variables:
# deg = differentially expressed genes in gene symbols
# allg = all genes that were measured or are left after independent
filtering
# go.list from mapGenesToTerms()

setUpGOdata <- function(deg, allg, go.list, GObranch="BP")
{
  require(topGO)
  allg <- cleanIds(allg)
  deg <- cleanIds(deg)
  geneList <- factor(as.integer(allg %in% deg))
  names(geneList) <- allg
  cat("Starting mapping of identifiers to GO-DAG... \n")
  GOdata <- new("topGOdata", ontology=GObranch, allGenes=geneList,
annot=annFUN.gene2GO, nodeSize = 5, gene2GO=go.list)
  cat("Done \n")
  GOdata
}

# Running the function: only the weight01 contains all the results,
# so in effect only weight01 can be used, others are only for comparison.
# Especially the ParentChild and Classic can give results not in the top
500 of weight01

runTopGO <- function(godata, nnodes=500)
{
  require(topGO)
  cat("Running Classic \n")
  go.test <- new("classicCount", testStatistic=GOFisherTest, name="Fisher
test")
  resultFis <- getSigGroups(godata, go.test)
  cat("Running weight \n")
  go.weight <- new("weightCount", testStatistic=GOFisherTest, name="Fisher
test", sigRatio="ratio")
  resultWeight <- getSigGroups(godata, go.weight)
  cat("Running parent-child \n")
  go.parentChild <- new("parentChild", testStatistic=GOFisherTest,
name="Fisher test")
  resultparentChild <- getSigGroups(godata, go.parentChild)
  cat("Running weight01 \n")
  go.weight01 <- new("weight01Count", testStatistic=GOFisherTest,
name="Fisher test")
  resultweight01 <- getSigGroups(godata, go.weight01)
  cat("Running elim \n")
  go.elimCount <- new("elimCount", testStatistic=GOFisherTest, name="Fisher
test")
  resultelimCount <- getSigGroups(godata, go.elimCount)

  allres <- GenTable(godata, classic=resultFis, p.c=resultparentChild,
weight=resultWeight, weight01=resultweight01, elim=resultelimCount,
orderBy="weight01", ranksOf="classic", topNodes=nnodes)
  allres
}

# Function cycles through a list of gene lists and assumes that all is in
SYMBOL format:#
doTopGO <- function(glist, allg, fd, GObranch="BP", nnodes=500)
{
  require(topGO)
  glnames <-names(glist)
  go.list <- mapGenesToTerms(fd=fd, GObranch=GObranch)
  resultslist <- lapply(1:length(glnames), function(n) {
    cat("############## \n")
    cat("Running:", glnames[n], "\n")
    g <- glist[[n]]
    cat("Number of genes:", length(g), "\n")
    gd <- setUpGOdata(deg=g, allg=allg, go.list=go.list, GObranch=GObranch)
    tg <- runTopGO(godata=gd, nnodes=nnodes)
    tg
  })
  names(resultslist) <- names(glist)
  resultslist
}

### Using the functions: ###
#######################

## Reading in the data: ##
glist <- readGeneList(folder=f, pattern=".txt")
allg <- fd$Symbol

## Running the topGO analysis on all the lists:#
## Gives a list of first 250 "weight01" topGO results: ##
tgl <- doTopGO(glist=glist, allg=allg, fd=fd, GObranch="BP", nnodes=250)




2013/9/13 Marc Carlson <mcarlson@fhcrc.org>

> Hi Raj,
>
> The problem seems to be that instead of using the appropriate annotation
> package you are grabbing the one used by the ALL example object. As a 1st
> step, you might want to look at our website and see if there is a different
> annotation package that better matches your array. You can see the list of
> packages on our website here:
>
> http://www.bioconductor.org/**packages/release/BiocViews.**
> html#___AnnotationData<http://www.bioconductor.org/packages/release/BiocViews.html#___AnnotationData>
>
>
> Marc
>
>
>
> On 09/11/2013 08:15 AM, kaushal Raj Chaudhary [guest] wrote:
>
>> Hi List,
>>
>> I have illumina expression data from mouse chip.  I was using topGO
>> package with probe id and P values .I am not sure if I can use hgu95av2.db
>> for annotation.  I encountered follwing error.  I appreciate your help.
>>  Thank you.
>>
>>   -- output of sessionInfo():
>>
>>   data1=read.csv ("go.csv", header=T)
>> data1
>>         ProbeID           p
>> 1   2810273 0.022920083
>> 2   5720014 0.042961586
>> 3   2600424 0.025965275
>> 4   3520176 0.032705434
>> 5   1240739 0.000668000
>> 6   1190735 0.016157765
>> …………………………**………………..
>> rn<-paste(data1[,1], sep="")
>> P_values=data1[,-1]
>> names(P_values)<-rn  ####creating named vectors
>> P_values    ### this contains probeID with corresponding p-values
>>
>> 2810273     5720014     2600424     3520176     1240739     1190735
>> 0.022920083 0.042961586 0.025965275 0.032705434 0.000668000 0.016157765
>>      4280164     2480441     5870504     5390019     5220576     4060465
>> 0.002103902 0.017717950 0.034700550 0.021078015 0.012964631 0.014807984
>>      4050653     1470619     4610131     2650730     2690035     4670288
>> 0.036829505 0.001426032 0.010399277 0.010104929 0.003046091 0.005084754
>>      4900064     6400377     5560204     2120377     4150037     1240370
>> 0.004906018 0.020934649 0.020285547 0.009254484 0.038157470 0.001123237
>>      2350068     3400246     4730630     3850167     4860048     2000674
>> …………………………**…………………………**
>> …………………………**…………..
>> affyLib <- paste(annotation (ALL), "db", sep = ".")
>> library(package = affyLib, character.only = TRUE)
>> [1] "hgu95av2.db"
>> sum(topDiffGenes(P_values))
>> 75
>>
>> sampleGOdata <- new("topGOdata",
>>   description = "Simple session", ontology = "BP",
>>   allGenes = P_values, geneSel = topDiffGenes,
>>   nodeSize = 10,
>>   annot = annFUN.db, affyLib = affyLib)
>>
>>
>> Building most specific GOs .....        ( 0 GO terms found. )
>>
>> Build GO DAG topology ..........        ( 0 GO terms and 0 relations. )
>> Error in if (is.na(index) || index < 0 || index > length(nd))
>> stop("vertex is not in graph: ",  :
>>    missing value where TRUE/FALSE needed
>>
>>
>> --
>> Sent via the guest posting facility at bioconductor.org.
>>
>> ______________________________**_________________
>> Bioconductor mailing list
>> Bioconductor@r-project.org
>> https://stat.ethz.ch/mailman/**listinfo/bioconductor<https://stat.ethz.ch/mailman/listinfo/bioconductor>
>> Search the archives: http://news.gmane.org/gmane.**
>> science.biology.informatics.**conductor<http://news.gmane.org/gmane.science.biology.informatics.conductor>
>>
>
> ______________________________**_________________
> Bioconductor mailing list
> Bioconductor@r-project.org
> https://stat.ethz.ch/mailman/**listinfo/bioconductor<https://stat.ethz.ch/mailman/listinfo/bioconductor>
> Search the archives: http://news.gmane.org/gmane.**
> science.biology.informatics.**conductor<http://news.gmane.org/gmane.science.biology.informatics.conductor>
>

	[[alternative HTML version deleted]]

