## ---- message = FALSE, warning = FALSE-----------------------------------
library(BgeeDB)

## ------------------------------------------------------------------------
listBgeeSpecies()

## ------------------------------------------------------------------------
bgee <- Bgee$new(species = "Mus_musculus", datatype = "rna_seq")

## ---- eval=FALSE---------------------------------------------------------
## # check the path where your folder with annotation will be saved. The folder is named after your chosen species.
## getwd()

## ------------------------------------------------------------------------
annotation_bgee_mouse <- bgee$get_annotation()
# head the experiments and libraries
lapply(annotation_bgee_mouse, head)

## ------------------------------------------------------------------------
# download all RPKMs and counts for Mus musculus
data_bgee_mouse <- bgee$get_data()

## ------------------------------------------------------------------------
# the number of experiments downloaded from Bgee
length(data_bgee_mouse)

## ------------------------------------------------------------------------
# see your first experiment
data_bgee_experiment1 <- data_bgee_mouse[[1]]
head(data_bgee_experiment1)

## ------------------------------------------------------------------------
# download RPKMs and counts only for GSE30617 for Mus musculus
data_bgee_mouse_gse30617 <- bgee$get_data(experiment.id = "GSE30617")

## ------------------------------------------------------------------------
# only present calls and rpkm values for GSE30617
gene.expression.mouse.rpkm <- bgee$format_data(data_bgee_mouse_gse30617, calltype = "expressed", stats = "rpkm")
names(gene.expression.mouse.rpkm)

## ------------------------------------------------------------------------
head(gene.expression.mouse.rpkm$"Ammon's horn")

## ------------------------------------------------------------------------
# Loading calls of expression based on RNA-seq data only
myTopAnatData <- loadTopAnatData(species=10090, datatype="rna_seq")
# Loading calls observed in embryonic stages only
myTopAnatData <- loadTopAnatData(species=10090, datatype="rna_seq", stage="UBERON:0000068")

# Look at the data
lapply(myTopAnatData, head)

## ------------------------------------------------------------------------
source("https://bioconductor.org/biocLite.R")
biocLite("biomaRt")
library(biomaRt)
ensembl <- useMart("ensembl")
ensembl <- useDataset("mmusculus_gene_ensembl", mart=ensembl)

# Foreground genes are those with a GO annotation "spermatogenesis"
myGenes <- getBM(attributes= "ensembl_gene_id", filters=c("go_id"), values=list(c("GO:0007283")), mart=ensembl)

# Background are all genes with a GO annotation
universe <- getBM(attributes= "ensembl_gene_id", filters=c("with_go_go"), values=list(c(TRUE)), mart=ensembl)

# Prepares the gene list vector 
geneList <- factor(as.integer(universe[,1] %in% myGenes[,1]))
names(geneList) <- universe[,1]
head(geneList)
summary(geneList == 1)

# Prepares the topGO object
myTopAnatObject <-  topAnat(myTopAnatData, geneList)

## ------------------------------------------------------------------------
results <- runTest(myTopAnatObject, algorithm = 'classic', statistic = 'fisher')
# You can also use the topGO decorrelation methods, for example the "weight" method to get less redundant results
results <- runTest(myTopAnatObject, algorithm = 'weight', statistic = 'fisher')

## ------------------------------------------------------------------------
# Display results sigificant at a 1% FDR threshold
tableOver <- makeTable(myTopAnatData, myTopAnatObject, results, 0.01)
# Display all results
tableOver <- makeTable(myTopAnatData, myTopAnatObject, results, 1)