hyperGTest not working on custom annotation package
2
0
Entering edit mode
@lucasmichel-20031
Last seen 5.8 years ago

I have created a custom annotation package using makeOrgPackage(). After succesfully creating it, installing it and loading it I have am trying to make a functional enrichment using it's annotations (which is why I created it) using GOstats. I can load it into the params object with no errors but when i call hyperGTest() i get the following error:

*Error in (function (classes, fdef, mtable)  : 
  unable to find an inherited method for function ‘columns’ for signature ‘"function"’*

I am guessing it might have something to do with the name of the custom package (if i call columns("org.PfalciparumPlasmoDB.eg") i get exactly the same error) but I don't know how to fix it.

Thanks

## CODE ##

library("GOstats")
library("GO.db")
library("org.PfalciparumPlasmoDB.eg.db") ## My custom annotation package

## Create universe; only genes with GO BP annotations ##

all_tbl <-  AnnotationDbi::select(org.PfalciparumPlasmoDB.eg.db, keys = featureNames(xgene_impute), columns = c("GO"), keytype = "SYMBOL")
ids = keys(org.PfalciparumPlasmoDB.eg.db, "GO")
df = AnnotationDbi::select(GO.db, keys(GO.db), "ONTOLOGY")

bp_ids = ids[ids %in% df$GOID[df$ONTOLOGY == "BP"]]
universeBP <- unique(all_tbl[all_tbl$GO %in% bp_ids, "SYMBOL"])


## Remove genes that have go terms that do not appear in GO.db ##

noGO <- all_tbl[!all_tbl$GO %in% keys(GO.db),]$SYMBOL
universeBP <- universeBP[!universeBP %in% noGO]


## Run FE (hyperGTest)  ##

  hgCutoff <- 0.05
  paramsBP <- new("GOHyperGParams",
               geneIds=cl_bp,
               universeGeneIds=universeBP,
               annotation="org.PfalciparumPlasmoDB.eg.db",
                 ontology="BP",
               pvalueCutoff=hgCutoff,
               conditional=FALSE,
               testDirection="over")

  hgOverBP <- hyperGTest(paramsBP)

Error in (function (classes, fdef, mtable)  : 
  unable to find an inherited method for function ‘columns’ for signature ‘"function"’

> sessionInfo()
R version 3.5.2 (2018-12-20)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 16.04.6 LTS

Matrix products: default
BLAS: /usr/lib/openblas-base/libblas.so.3
LAPACK: /usr/lib/libopenblasp-r0.2.18.so

locale:
 [1] LC_CTYPE=en_US.UTF-8          LC_NUMERIC=C                  LC_TIME=es_ES.UTF-8           LC_COLLATE=en_US.UTF-8       
 [5] LC_MONETARY=es_ES.UTF-8       LC_MESSAGES=en_US.UTF-8       LC_PAPER=es_ES.UTF-8          LC_NAME=es_ES.UTF-8          
 [9] LC_ADDRESS=es_ES.UTF-8        LC_TELEPHONE=es_ES.UTF-8      LC_MEASUREMENT=es_ES.UTF-8    LC_IDENTIFICATION=es_ES.UTF-8

attached base packages:
 [1] stats4    grid      parallel  stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] org.Pf.plasmo.db_3.7.0              org.PfalciparumPlasmoDB.eg.db_1.0.0 forcats_0.4.0                       stringr_1.4.0                      
 [5] purrr_0.3.0                         readr_1.3.1                         tibble_2.0.1                        tidyverse_1.2.1                    
 [9] Rgraphviz_2.26.0                    xtable_1.8-3                        RColorBrewer_1.1-2                  GO.db_3.7.0                        
[13] GOstats_2.48.0                      graph_1.60.0                        Category_2.48.0                     Matrix_1.2-15                      
[17] AnnotationForge_1.24.0              AnnotationDbi_1.44.0                IRanges_2.16.0                      S4Vectors_0.20.1                   
[21] kableExtra_1.0.1                    ggdendro_0.1-20                     venneuler_1.1-0                     rJava_0.9-10                       
[25] pca3d_0.10                          gridExtra_2.3                       tidyr_0.8.2                         dplyr_0.8.0.1                      
[29] VennDiagram_1.6.20                  futile.logger_1.4.3                 ggfortify_0.4.5                     gplots_3.0.1.1                     
[33] ggplot2_3.1.0                       reshape2_1.4.3                      Biobase_2.42.0                      BiocGenerics_0.28.0                

loaded via a namespace (and not attached):
 [1] colorspace_1.4-0        rstudioapi_0.9.0        bit64_0.9-7             lubridate_1.7.4         xml2_1.2.0              splines_3.5.2          
 [7] knitr_1.21              jsonlite_1.6            broom_0.5.1             annotate_1.60.0         shiny_1.2.0             compiler_3.5.2         
[13] httr_1.4.0              backports_1.1.3         assertthat_0.2.0        lazyeval_0.2.1          cli_1.0.1               later_0.8.0            
[19] formatR_1.5             htmltools_0.3.6         tools_3.5.2             gtable_0.2.0            glue_1.3.0              Rcpp_1.0.0             
[25] cellranger_1.1.0        gdata_2.18.0            nlme_3.1-137            crosstalk_1.0.0         xfun_0.5                rvest_0.3.2            
[31] mime_0.6                miniUI_0.1.1.1          gtools_3.8.1            XML_3.98-1.17           MASS_7.3-51.1           scales_1.0.0           
[37] hms_0.4.2               promises_1.0.1          RBGL_1.58.1             lambda.r_1.2.3          yaml_2.2.0              memoise_1.1.0          
[43] stringi_1.3.1           RSQLite_2.1.1           genefilter_1.64.0       caTools_1.17.1.1        manipulateWidget_0.10.0 rlang_0.3.1            
[49] pkgconfig_2.0.2         bitops_1.0-6            rgl_0.99.16             evaluate_0.13           lattice_0.20-38         htmlwidgets_1.3        
[55] bit_1.1-14              tidyselect_0.2.5        GSEABase_1.44.0         plyr_1.8.4              magrittr_1.5            R6_2.4.0               
[61] generics_0.0.2          DBI_1.0.0               pillar_1.3.1            haven_2.1.0             withr_2.1.2             survival_2.43-3        
[67] RCurl_1.95-4.11         modelr_0.1.4            crayon_1.3.4            futile.options_1.0.1    KernSmooth_2.23-15      ellipse_0.4.1          
[73] rmarkdown_1.11          readxl_1.3.0            blob_1.1.1              digest_0.6.18           webshot_0.5.1           httpuv_1.4.5.1         
[79] munsell_0.5.0           viridisLite_0.3.0    
GOstats AnnotationDbi GO.db • 1.7k views
ADD COMMENT
0
Entering edit mode

The input IDs for this function should be the central ID for your OrgDb package. Usually that's the Gene ID, but it's not the GO ID because the functions in GOstats will already be trying to map the ID to GO IDs, which obviously won't work. In general the Universe should be the set of genes that you measured (which is probably the set of genes on your microarray), again, using the central ID of your annotation package, which you can get by doing head(keys(org.PfalciparumPlasmoDB.eg.db)), and which are usually the NCBI Gene ID, but may be different for you.

The error is mysterious to me. The only place in either GOstats or Category that calls columns on anything is in a function for Yeast arrays. Which should only happen if you get YEASTCHIPDB if you do dbmeta(dbconn(org.PfalciparumPlasmoDB.eg.db), "DBSCHEMA"), but you should get NOSCHEMADB for this OrgDb.

Anyway, try fixing the IDs you use and if you get the same error, run traceback() right after the error and post the output from that.

ADD REPLY
0
Entering edit mode

Hi James, thank you very much for your reply.

Indeed my schema is NOSCHEMADB. On the other hand, the central identifier of my annotation package is ENTREZGID (it is mandatory for makeOrgPackage()) but the IDs I was passing to the function and the gene universe are gene symbols. I have changed both the gene universe and the passed IDs to ENTREZGIDs and the problem is still there though.

I attach the result of running traceback().

> hgOverBP <- hyperGTest(paramsBP)
Error in (function (classes, fdef, mtable)  : 
  unable to find an inherited method for function ‘columns’ for signature ‘"function"’
> traceback()
15: stop(gettextf("unable to find an inherited method for function %s for signature %s", 
        sQuote(fdef@generic), sQuote(cnames)), domain = NA)
14: (function (classes, fdef, mtable) 
    {
        methods <- .findInheritedMethods(classes, fdef, mtable)
        if (length(methods) == 1L) 
            return(methods[[1L]])
        else if (length(methods) == 0L) {
            cnames <- paste0("\"", vapply(classes, as.character, 
                ""), "\"", collapse = ", ")
            stop(gettextf("unable to find an inherited method for function %s for signature %s", 
                sQuote(fdef@generic), sQuote(cnames)), domain = NA)
        }
        else stop("Internal error in finding inherited methods; didn't return a unique method", 
            domain = NA)
    })(list("function"), new("nonstandardGenericFunction", .Data = function (x) 
    {
        value <- standardGeneric("columns")
        sort(value)
    }, generic = "columns", package = "AnnotationDbi", group = list(), 
        valueClass = character(0), signature = "x", default = NULL, 
        skeleton = (function (x) 
        stop("invalid call in method dispatch to 'columns' (no default method)", 
            domain = NA))(x)), <environment>)
13: columns(db)
12: map %in% columns(db)
11: getAnnMap(col, datpkg@name)
10: annMapChooser("GO", p)
9: ID2GO(datPkg)
8: ID2GO(datPkg)
7: eapply(ID2GO(datPkg), function(goids) {
       if (length(goids) == 0 || (length(goids) == 1 && is.na(goids))) 
           return(FALSE)
       if (!is.character(goids[[1]])) 
           goids <- subListExtract(goids, "GOID", simplify = TRUE)
       if (any(goids %in% ontIds)) 
           return(TRUE)
       FALSE
   })
6: getUniverseViaGo(p)
5: universeBuilder(p)
4: universeBuilder(p)
3: .hyperGTestInternal(p, "GOHyperGResult", getGoToChildGraph)
2: hyperGTest(paramsBP)
1: hyperGTest(paramsBP)

ADD REPLY
0
Entering edit mode
@james-w-macdonald-5106
Last seen 17 hours ago
United States

There are two problems here. First, there is a bug in the annotate package that is simple to fix and I will do in the next few minutes. But there is another problem that will take more time.

The annotation packages contain SQLite databases, and various functions in GOstats make SQL queries on these databases. The SQL queries rely on a consistent column naming convention for the underlying tables, so as an example, a query might be

SELECT DISTINCT go_id
FROM go_all INNER JOIN genes USING (_id)
WHERE gene_id IN ('812888','810731','814083') AND ontology='BP'

and the expectation is that the go_all table has columns called go_id and ontology. HOWEVA, the org.Pfalciparum.3d7.eg.db package I made the other day has this:

dbGetQuery(dbconn(org.Pfalciparum.3d7.eg.db), "select * from go_all limit 3;")
   _id      GOALL EVIDENCEALL ONTOLOGYALL
1    1 GO:0020013         TAS          BP
2    1 GO:0020033         TAS          BP
3    1 GO:0020035         TAS          BP

Where the column names are obviously not consistent with all the other OrgDb packages, and is a problem with makeOrgDbFromNCBI and may take a bit longer to fix. I'll post back here when I have fixed both bugs.

ADD COMMENT
0
Entering edit mode

Thank you very much for the prompt reply. Let me know if I can help in any way.

ADD REPLY
0
Entering edit mode
@james-w-macdonald-5106
Last seen 17 hours ago
United States

OK, it's fixed now:

> library(GOstats)
> library(org.Pfalciparum.3d7.eg.db)
> universeGeneIds <- keys(org.Pfalciparum.3d7.eg.db)
> geneIds <- sample(universeGeneIds, 300)
> p <- new("GOHyperGParams", geneIds = geneIds, universeGeneIds = universeGeneIds, annotation = "org.Pfalciparum.3d7.eg.db", ontology = "BP", conditional = TRUE)
> z <- hyperGTest(p)
'select()' returned 1:many mapping between keys and columns
> summary(z)
      GOBPID      Pvalue OddsRatio  ExpCount Count Size
1 GO:0071840 0.001587302  2.811628 6.6714184    15  144
2 GO:0006364 0.002190388  9.683060 0.6022808     4   13
3 GO:0046907 0.005611950  3.147527 3.4283678     9   74
4 GO:0009451 0.006528955 10.741935 0.4169636     3    9
5 GO:0022613 0.006682985  4.984848 1.2508909     5   27
                                           Term
1 cellular component organization or biogenesis
2                               rRNA processing
3                       intracellular transport
4                              RNA modification
5          ribonucleoprotein complex biogenesis

I just checked in the bugfixes to both release and devel. It takes a day or two to propagate - you will need to get annotate version 1.60.1 and Category 2.48.1 when they are available.

ADD COMMENT

Login before adding your answer.

Traffic: 542 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6