I recently posted a question regarding identifying GO terms that were mapped to GOslims and received a great answer on how to do this.

However, in the grand scheme of things, I really want to find out which genes/transcripts are mapped to a set of GOslims. Since GSEABase is able to map GO terms to GOslims, it seems like identifying corresponding gene/transcript IDs should be feasible. Unfortunately, I have no idea how to proceed.

If I have an input file with GO terms and corresponding gene/transcript IDs (note: actual format/layout is unimportant - I'm confident I can change things around to match needed input format for use in GSEABase):

GO:0016874  TRINITY_DN8634_c0_g1
GO:0016874  TRINITY_DN9386_c0_g1
GO:0019693  TRINITY_DN10297_c1_g1
GO:0019693  TRINITY_DN12835_c0_g1
GO:0019693  TRINITY_DN1421_c0_g1

I know how to map the GO terms to GOslims, as well as identify which GO terms map to GOslims.

Now, I am hoping there's a "simple" solution built-in to GSEABase that will allow me to subset gene/transcript IDs from any given GOslim.

Would anyone have any suggestions as to how to approach this within GSEABase? It seems like creating a GeneSetCollection (or, a GOCollection?) would be the first approach, but I'm not familiar enough with R and GSEABase to figure out where to go from there. Heck, I might not even be able to figure out how to properly construct a GeneSetCollection properly with the information provided above...

I can probably wrangle some command line stuff to accomplish this task, but it would be nicer (and probably easier?) to just use GSEABase to accomplish this.

Oh, one last thing to add, all of our work is with non-model organisms with limited/non-existent genomic resources. So, I do have UniProt IDs associated with these gene/transcript IDs, but won't be able to tap into the standard ENTREZ/PMID/etc databases.

hi,

one possible way to follow using the GSEABase infrastructure is to build a GeneSetCollection object with the GO annotations. Let's say we have the GO annotations you show as an example in a text file called goannot.txt, where the first columns has the GO identifier and the second column, separated by blank spaces from the first, has the gene identifier. We can build this GeneSetCollection object of these annotations as follows:

library(GSEABase)
library(GO.db)

## read GO annotations from a text file with one row per annotation
goannot <- read.table("goannot.txt", sep="", quote="", colClasses="character")

## organize GO annotations by GO term
genesbygo <- split(goannot[, 2], goannot[, 1])
genesbygo
$`GO:0016874`
[1] "TRINITY_DN8634_c0_g1" "TRINITY_DN9386_c0_g1"

$`GO:0019693`
[1] "TRINITY_DN10297_c1_g1" "TRINITY_DN12835_c0_g1" "TRINITY_DN1421_c0_g1" 

## fetch GO terms into a named vector
goterms <- select(GO.db, keys=names(genesbygo), columns="TERM")
goterms <- setNames(goterms$TERM, goterms$GOID)

## build GeneSetCollection object
gsc <- do.call("GeneSetCollection",
               mapply(function(geneSetID, geneIDs, goterms) {
                        GeneSet(NullIdentifier(),
                                geneIds=geneIDs,
                                setName=geneSetID,
                                shortDescription=goterms[geneSetID])
                     }, names(genesbygo), genesbygo,
                     MoreArgs=list(goterms=goterms), USE.NAMES=FALSE))
gsc
GeneSetCollection
  names: GO:0016874, GO:0019693 (2 total)
  unique identifiers: TRINITY_DN8634_c0_g1, TRINITY_DN9386_c0_g1, ..., TRINITY_DN1421_c0_g1 (5 total)
  types in collection:
    geneIdType: NullIdentifier (1 total)
    collectionType: NullCollection (1 total)

If you have UnitProt IDs in the GO annotations, then you could replace the constructor NullIdentifier() by UniprotIdentifier(), which will help to map with other identifiers, whenever those mappings exist. In principle, you should be able to use the solution given in your previous question to do the mapping you're looking for.

cheers,

robert.