{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Annotations in NCBI's gene2go file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1) Download the NCBI annotations\n",
    "The NCBI annotations are stored in the file named \"gene2go\"."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  EXISTS: gene2go\n"
     ]
    }
   ],
   "source": [
    "from goatools.base import download_ncbi_associations\n",
    "# fin -> Filename of input file (file to be read)\n",
    "fin_gene2go = download_ncbi_associations()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2) Read NCBI annotation file, \"gene2go\"\n",
    "### 2a) Read one taxid: human"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "HMS:0:00:05.892614 323,107 annotations, 19,649 genes, 18,246 GOs, 1 taxids READ: gene2go \n"
     ]
    }
   ],
   "source": [
    "from goatools.anno.genetogo_reader import Gene2GoReader\n",
    "objanno_hsa = Gene2GoReader(fin_gene2go, taxids=[9606])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2b) Read all taxids"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "HMS:0:00:20.773931 2,057,323 annotations, 205,158 genes, 26,896 GOs, 46 taxids READ: gene2go \n"
     ]
    }
   ],
   "source": [
    "objanno_all = Gene2GoReader(fin_gene2go, taxids=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3) Get associations, split by namespace (Only human annotations loaded)\n",
    "Namespaces include *biological_process* (BP), *molecular_function* (MF), *cellular_component* (CC).\n",
    "\n",
    "The *taxid* argument in function, *get_ns2assc*, is ignored if there only one taxid is loaded."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "ns2assc_hsa1 = objanno_hsa.get_ns2assc()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "from itertools import chain\n",
    "\n",
    "def prt_assc_counts(ns2assc):\n",
    "    \"\"\"Print the number of genes and GO IDs in an association\"\"\"\n",
    "    for nspc, gene2goids in sorted(ns2assc.items()):\n",
    "        print(\"{NS} {N:6,} genes, {GOs:6,} GOs\".format(\n",
    "            NS=nspc, N=len(gene2goids), GOs=len(set.union(*gene2goids.values()))))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "BP 17,541 genes, 12,285 GOs\n",
      "CC 18,648 genes,  1,737 GOs\n",
      "MF 17,384 genes,  4,170 GOs\n"
     ]
    }
   ],
   "source": [
    "prt_assc_counts(ns2assc_hsa1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4) Get associations, split by namespace (Many taxids loaded)\n",
    "### 4a) Get associations for one species (human)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "BP 17,541 genes, 12,285 GOs\n",
      "CC 18,648 genes,  1,737 GOs\n",
      "MF 17,384 genes,  4,170 GOs\n"
     ]
    }
   ],
   "source": [
    "ns2assc_hsa2 = objanno_all.get_ns2assc(9606)\n",
    "\n",
    "prt_assc_counts(ns2assc_hsa2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4b) Get associations for one species (mouse)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "BP 17,859 genes, 12,282 GOs\n",
      "CC 18,824 genes,  1,726 GOs\n",
      "MF 16,721 genes,  4,156 GOs\n"
     ]
    }
   ],
   "source": [
    "ns2assc_mmu = objanno_all.get_ns2assc(10090)\n",
    "\n",
    "prt_assc_counts(ns2assc_mmu)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4c) Combine associations for multiple species (human and mouse)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "BP 35,400 genes, 13,017 GOs\n",
      "CC 37,472 genes,  1,798 GOs\n",
      "MF 34,105 genes,  4,373 GOs\n"
     ]
    }
   ],
   "source": [
    "ns2assc_two = objanno_all.get_ns2assc({9606, 10090})\n",
    "\n",
    "prt_assc_counts(ns2assc_two)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4d) Combine all associations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "BP 147,644 genes, 17,777 GOs\n",
      "CC 155,564 genes,  2,546 GOs\n",
      "MF 140,246 genes,  6,514 GOs\n"
     ]
    }
   ],
   "source": [
    "ns2assc_all = objanno_all.get_ns2assc(True)\n",
    "prt_assc_counts(ns2assc_all)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4e) Try getting unspecified taxids\n",
    "If annotations have been loaded for multiple or all taxids, the taxid argument is required.\n",
    "\n",
    "If the taxid arg is unused, an error will be printed and an empty dict will be returned."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "**ERROR: ARG taxid MUST BE AN int, list of ints, OR True\n"
     ]
    }
   ],
   "source": [
    "ns2assc_all = objanno_all.get_ns2assc()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{}\n"
     ]
    }
   ],
   "source": [
    "print(ns2assc_all)"
   ]
  }
 ],
 "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": 2
}