{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Gene Ontology Enrichment Analysis (GOEA)\n",
    "This is the same GOEA as in goea_nbt3102.ipynb, but the GOEA results can be obtained by calling a single function.\n",
    "\n",
    "We use data from a 2014 Nature paper:    \n",
    "[Computational analysis of cell-to-cell heterogeneity\n",
    "in single-cell RNA-sequencing data reveals hidden \n",
    "subpopulations of cells\n",
    "](http://www.nature.com/nbt/journal/v33/n2/full/nbt.3102.html#methods)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_goeaobj_nbt3102(method='fdr_bh'):\n",
    "    \"\"\"Return GOEA Object ready to run Nature data.\"\"\"\n",
    "    from goatools.obo_parser import GODag\n",
    "    from genes_ncbi_10090_proteincoding import GENEID2NT as GeneID2nt_mus\n",
    "    from goatools.base import download_go_basic_obo, download_ncbi_associations\n",
    "    from goatools.anno.genetogo_reader import Gene2GoReader\n",
    "    from goatools.goea.go_enrichment_ns import GOEnrichmentStudyNS\n",
    "    # Load Ontologies\n",
    "    obo_fname = download_go_basic_obo()\n",
    "    obodag = GODag(\"go-basic.obo\")\n",
    "    # Load Associations\n",
    "    download_ncbi_associations() # Get ftp://ftp.ncbi.nlm.nih.gov/gene/DATA/gene2go.gz\n",
    "    # Read NCBI's gene2go. Store annotations in a list of namedtuples\n",
    "    objanno = Gene2GoReader(\"gene2go\", taxids=[10090])\n",
    "    # Get associations for each branch of the GO DAG (BP, MF, CC)\n",
    "    ns2assoc = objanno.get_ns2assc()\n",
    "    # GOE Object holds Ontologies, Associations, and Background gene set\n",
    "    return GOEnrichmentStudyNS(\n",
    "        GeneID2nt_mus.keys(), # Background gene set: mouse protein-coding genes\n",
    "        ns2assoc, # geneid/GO Associations for BP, MF, anc CC GODAG branches\n",
    "        obodag, # Ontologies\n",
    "        propagate_counts = False,\n",
    "        alpha = 0.05, # default significance cut-off\n",
    "        methods = [method]) # defult multipletest correction method"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def read_data_nbt3102():\n",
    "    \"\"\"Read data from Nature paper.\"\"\"\n",
    "    import os\n",
    "    # Data will be stored in this variable\n",
    "    geneid2symbol = {}\n",
    "    # Get xlsx filename where data is stored\n",
    "    ROOT = os.path.dirname(os.getcwd()) # go up 1 level from current working directory\n",
    "    din_xlsx = os.path.join(ROOT, \"goatools/test_data/nbt_3102/nbt.3102-S4_GeneIDs.xlsx\")\n",
    "    # Read data\n",
    "    if os.path.isfile(din_xlsx):  \n",
    "        import xlrd\n",
    "        book = xlrd.open_workbook(din_xlsx)\n",
    "        pg = book.sheet_by_index(0)\n",
    "        for r in range(pg.nrows):\n",
    "            symbol, geneid, pval = [pg.cell_value(r, c) for c in range(pg.ncols)]\n",
    "            if geneid:\n",
    "                geneid2symbol[int(geneid)] = symbol\n",
    "    return geneid2symbol"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Copyright (C) 2016-present, DV Klopfenstein, H Tang. All rights reserved."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}