annFUN                 package:topGO                 R Documentation

_f_u_n_c_t_i_o_n_s _t_o _m_a_p _g_e_n_e _I_D_s _t_o _G_O _t_e_r_m_s

_D_e_s_c_r_i_p_t_i_o_n:

     These functions are used to compile a list of GO terms and their
     mappings to gene identifiers.

_U_s_a_g_e:

     annFUN.db(whichOnto, feasibleGenes = NULL, affyLib)
     annFUN(whichOnto, feasibleGenes = NULL, affyLib)
     annFUN.gene2GO(whichOnto, feasibleGenes = NULL, gene2GO)
     annFUN.GO2genes(whichOnto, feasibleGenes = NULL, GO2genes)

_A_r_g_u_m_e_n_t_s:

whichOnto: character string specifying one of the three GO ontologies:
          '"BP"', '"MF"', '"CC"'

feasibleGenes: character vector containing a subset of gene
          identifiers. Only these genes will be used to annotate GO
          terms. Default value is 'NULL' which means all gene
          identifiers will be used.

 affyLib: character string containing the name of the Affymetrix chip.

 gene2GO: named list of character vectors. The list names are genes
          identifiers. For each gene the character vector contains the
          GO terms IDs it maps to. Only the most specific annotations
          are required.

GO2genes: named list of character vectors. The list names are GO terms
          IDs. For each GO the character vector contains the genes
          identifiers which are mapped to it. Only the most specific
          annotations are required.

_D_e_t_a_i_l_s:

     The function 'annFUN.db' uses the mappings provided in the
     Bioconductor annotation data packages. For example, if the
     Affymetrix hgu133a chip it is used, then the user should set
     'affyLib = "hgu133a.db"'.

     The functions 'annFUN.gene2GO' and 'annFUN.GO2genes' are used when
     the user provide his own annotations.

     All these function restrict the GO terms to the ones belonging to
     the specified ontology.

_V_a_l_u_e:

     A named(GO terms IDs) list of character vectors.

_A_u_t_h_o_r(_s):

     Adrian Alexa

_S_e_e _A_l_s_o:

     'topGOdata-class'

_E_x_a_m_p_l_e_s:

     library(hgu133a.db)
     set.seed(111)

     ## generate a gene list and the GO annotations
     numGenes <- 50
     selGenes <- sample(ls(hgu133aGO), numGenes)
     gene2GO <- lapply(mget(selGenes, envir = hgu133aGO), names)
     gene2GO[sapply(gene2GO, is.null)] <- NA

     ## the annotation for the first three genes
     gene2GO[1:3]

     ## inverting the annotations
     go2genes <- annFUN.gene2GO(whichOnto = "CC", gene2GO = gene2GO)


     ## generate a GO list with the genes annotations
     numGO <- 30
     selGO <- sample(ls(hgu133aGO2PROBE), numGO)
     GO2gene <- lapply(mget(selGO, envir = hgu133aGO2PROBE), as.character)

     GO2gene[1:3]

     ## select only the GO terms for a specific ontology
     go2gene <- annFUN.GO2genes(whichOnto = "CC", GO2gene = GO2gene)

