# Load required libraries
source("https://bioconductor.org/biocLite.R")
biocLite("BgeeDB")
library(BgeeDB)
biocLite("biomaRt")
library(biomaRt)

# Connection to Biomart for mouse data in Ensembl release 87
mouse_mart <- useMart("ENSEMBL_MART_ENSEMBL", dataset="mmusculus_gene_ensembl",
    host="Dec2016.archive.ensembl.org");

# Retrieve mouse genes annotated to GO term GO:0060537 "muscle tissue development",
# along with cattle orthologous gene IDs
mouse_muscle_genes <- getBM(attributes=c("ensembl_gene_id", "btaurus_homolog_ensembl_gene"),
    filters="go_parent_term", values=c("GO:0060537"), mart=mouse_mart);
# Retrieve all mouse genes with any GO annotation, along with cattle orthologous gene IDs, 
# for proper definition of gene universe
mouse_GO_genes <- getBM(attributes=c("ensembl_gene_id", "btaurus_homolog_ensembl_gene"),
    filters="with_go_id", values=c(T), mart=mouse_mart);
# Get cattle orthologous gene IDs of mouse genes mapped to GO term GO:0060537 "muscle tissue development"
cow_muscle_gene_IDs <- unique(mouse_muscle_genes$btaurus_homolog_ensembl_gene);
# Get cattle orthologous gene IDs of mouse genes mapped to any GO term
cow_GO_gene_IDs <- unique(mouse_GO_genes$btaurus_homolog_ensembl_gene);

# Prepare the gene list vector
gene_list <- factor(as.integer(cow_GO_gene_IDs %in% cow_muscle_gene_IDs));
names(gene_list) <- cow_GO_gene_IDs;

# Load TopAnat object with cattle expression data
bgee_topanat <- Bgee$new(species="Bos_taurus");
topanat_data <- loadTopAnatData(bgee_topanat);
topanat_obj  <- topAnat(topanat_data, gene_list);
# Run test
cow_results <- runTest(topanat_obj, algorithm='weight', statistic='fisher');
cow_table_over <- makeTable(topanat_data, topanat_obj, cow_results, cutoff=0.2);