{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# coproID report"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "import datetime\n",
    "with open(\"version.txt\", \"r\") as f:\n",
    "    for line in f:\n",
    "        print(f\"Version {line.rstrip()}\")\n",
    "print(f\"Report generated on {datetime.datetime.now()}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "[![](coproid_logo.png)](https://github.com/nf-core/coproID)"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introduction\n",
    "[coproID](https://github.com/nf-core/coproID) is a pipeline to identify the source of coprolites, and in general, of a metagenomic sample.\n",
    "\n",
    "If you read these lines, coproID successfully finished running and you can find your results below.  \n",
    "You can find more informations about the different result files in the coproID documentation: [coproid.readthedocs.io](https://coproid.readthedocs.io/en/latest/output.html)"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import os\n",
    "from IPython.display import display, Markdown, Latex, HTML\n",
    "import base64\n",
    "from bokeh.plotting import figure, show, output_notebook\n",
    "from bokeh.models import ColumnDataSource\n",
    "from bokeh.transform import factor_cmap\n",
    "from bokeh.palettes import Set1\n",
    "from bokeh.models.tools import HoverTool\n",
    "from bokeh.layouts import widgetbox\n",
    "from bokeh.models.widgets import Button\n",
    "from bokeh.models.widgets import DataTable, DateFormatter, TableColumn\n",
    "from bokeh.models import CustomJS\n",
    "import ipywidgets as widgets\n",
    "from plotnine import *\n",
    "import warnings\n",
    "warnings.simplefilter('ignore')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "def coproid_summary_plot(df):\n",
    "    df = pd.read_csv(df, index_col=0)\n",
    "    organisms = [i.replace(\"normalized_bp_proportion_aligned_\",\"\") for i in list(df.columns) if \"normalized_bp_proportion_aligned_\" in i]\n",
    "    organisms_clean = [i.replace(\"_\",\" \") for i in organisms]\n",
    "    if len(organisms_clean) < 3:\n",
    "        display(Markdown(\"### coproID summary plot\"))\n",
    "        species_text = pd.DataFrame()\n",
    "        species_text['x'] = [0.25, 0.75, 0.75, 0.25]\n",
    "        species_text['y'] = [0.25, 0.25, 0.75, 0.75]\n",
    "        species_text['text'] = ['Unknown', organisms_clean[0], 'Unknown', organisms_clean[1]]\n",
    "    \n",
    "        df['samp_name'] = df.index\n",
    "        df['coproID_prediction'] = ['Unknown'] * df.shape[0]\n",
    "        df['coproID_prediction'] = np.where(df[f\"coproID_proba_{organisms[0]}\"] > 0.5, organisms_clean[0], df['coproID_prediction'])\n",
    "        df['coproID_prediction'] = np.where(df[f\"coproID_proba_{organisms[1]}\"] > 0.5, organisms_clean[1], df['coproID_prediction'])\n",
    "    \n",
    "        p = ggplot(df, aes(x = f\"coproID_proba_{organisms[0]}\",y = f\"coproID_proba_{organisms[1]}\"))\n",
    "        p = p + geom_point(aes(color='coproID_prediction'), size=2)\n",
    "        p = p + geom_label(aes(label=\"samp_name\", color='coproID_prediction'), size=8, nudge_x = 0.02, ha='left', va='bottom')\n",
    "        p = p + theme_classic() + labs(x=f\"coproID proba {organisms_clean[0]}\",y = f\"coproID proba {organisms_clean[1]}\")\n",
    "        p = p + geom_text(data=species_text, mapping=aes(x='x',y='y', label='text'), alpha=0.3, color='grey')\n",
    "        p = p + geom_hline(yintercept=0.5, linetype='dashed', alpha=0.1) \n",
    "        p = p + geom_vline(xintercept=0.5, linetype='dashed', alpha=0.1)\n",
    "        p = p + scale_color_manual(name='Predicted Organism', values = {organisms_clean[0]:'#ef7576', organisms_clean[1]:'#c194c8', 'Unknown':'#a2a3a1'})\n",
    "        p = p + coord_cartesian(xlim=[0,1],ylim=[0,1])\n",
    "        p.draw()\n",
    "        return(organisms)\n",
    "    else:\n",
    "        return"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "def plot_bokeh(df, organisms = []):\n",
    "    d = pd.read_csv(df, index_col = 0).fillna('sink')\n",
    "    d = d.rename(columns={'PC1':'DIM1', 'PC2':'DIM2'})\n",
    "    \n",
    "    orga_in_endo = [i for i in list(sorted(set(d['labels']))) if i in organisms]\n",
    "    orga_not_in_endo = [i for i in list(sorted(set(d['labels']))) if i not in organisms]\n",
    "    colors = {k:v for k, v in zip(orga_in_endo, Set1[9][0:2]) if k in organisms}\n",
    "    cnt=2\n",
    "    for i in orga_not_in_endo:\n",
    "        colors[i] = Set1[9][cnt]\n",
    "        cnt +=1\n",
    "        \n",
    "    d['colors'] = [colors[i] for i in list(d['labels'])]\n",
    "\n",
    "    TOOLS=\"pan,wheel_zoom,zoom_in,zoom_out,box_zoom,reset,save\"\n",
    "    source = ColumnDataSource(d)\n",
    "#     labels = list(set(list(set(source.data['labels'])).append(organisms)))\n",
    "#     colors = [\"#e41a1c\",\"#377eb8\",\"#4daf4a\",\"#984ea3\",\"#ff7f00\",\"\"]\n",
    "#     color_map = factor_cmap(field_name='labels', palette=colors, factors=labels)\n",
    "    p = figure(tools=TOOLS)\n",
    "    p.scatter(x = 'DIM1', y='DIM2', color='colors', alpha = 0.6, size = 6, legend='labels', source=d)\n",
    "    hover = HoverTool()\n",
    "    hover.tooltips = [(\"Organism\", \"@labels\"),('Sample',\"@name\")]\n",
    "    p.add_tools(hover)\n",
    "    output_notebook(hide_banner=True)\n",
    "    show(p)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "def bokeh_table(df):\n",
    "    \n",
    "    d = pd.read_csv(df, index_col=0)\n",
    "    d = d.rename(columns={'PC1':'DIM1', 'PC2':'DIM2'})\n",
    "    d.insert(0, \"sample\", d.index)\n",
    "    source = ColumnDataSource(d)\n",
    "\n",
    "    Columns = [TableColumn(field=Ci, title=Ci) for Ci in d.columns] # bokeh columns\n",
    "    data_table = DataTable(columns=Columns, source=ColumnDataSource(d)) # bokeh table\n",
    "\n",
    "    button = Button(label=\"Download\", button_type=\"success\")\n",
    "\n",
    "    javaScript=\"\"\"\n",
    "        function table_to_csv(source) {\n",
    "        const columns = Object.keys(source.data)\n",
    "        const nrows = source.get_length()\n",
    "        const lines = [columns.join(',')]\n",
    "\n",
    "        for (let i = 0; i < nrows; i++) {\n",
    "            let row = [];\n",
    "            for (let j = 0; j < columns.length; j++) {\n",
    "                const column = columns[j]\n",
    "                row.push(source.data[column][i].toString())\n",
    "            }\n",
    "            lines.push(row.join(','))\n",
    "        }\n",
    "        return lines.join('\\\\n').concat('\\\\n')\n",
    "    }\n",
    "\n",
    "\n",
    "    const filename = 'coproID_result.csv'\n",
    "    filetext = table_to_csv(source)\n",
    "    const blob = new Blob([filetext], { type: 'text/csv;charset=utf-8;' })\n",
    "\n",
    "    //addresses IE\n",
    "    if (navigator.msSaveBlob) {\n",
    "        navigator.msSaveBlob(blob, filename)\n",
    "    } else {\n",
    "        const link = document.createElement('a')\n",
    "        link.href = URL.createObjectURL(blob)\n",
    "        link.download = filename\n",
    "        link.target = '_blank'\n",
    "        link.style.visibility = 'hidden'\n",
    "        link.dispatchEvent(new MouseEvent('click'))\n",
    "    }\n",
    "    \"\"\"\n",
    "    button.js_event_callbacks = CustomJS(args=dict(source=source),code=javaScript)\n",
    "    output_notebook(hide_banner=True)\n",
    "    show(button)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "umap = \"sourcepredict_embedding.csv\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## coproID summary\n",
    "### Table"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "d = \"coproID_result.csv\"\n",
    "df = pd.read_csv(d, index_col=0)\n",
    "df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hideCode": false,
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "def create_download_link(df=df, title = \"Download table as CSV\", filename = \"coproid_result.csv\"):  \n",
    "    csv = df.to_csv(index =False)\n",
    "    b64 = base64.b64encode(csv.encode())\n",
    "    payload = b64.decode()\n",
    "    html = '<a download=\"{filename}\" href=\"data:text/csv;base64,{payload}\" target=\"_blank\">{title}</a>'\n",
    "    html = html.format(payload=payload,title=title,filename=filename)\n",
    "    return HTML(html)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "create_download_link(df)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "orga = coproid_summary_plot(d)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Microbiome profile embedding"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true,
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "plot_bokeh(umap, organisms=orga)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Damage plots"
   ],
   "execution_count": null,
   "outputs": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "\n",
    "allplots = os.listdir(\".\")\n",
    "allplots = [i for i in allplots if \"_freq.txt\" in i]\n",
    "fwd_ext = \".5pCtoT_freq.txt\"\n",
    "rev_ext = \".3pGtoA_freq.txt\"\n",
    "samples = list(set([i.split(\"_otu_\")[0] for i in allplots]))\n",
    "species = list(set([\"_\".join(\".\".join(i.split(\".\")[0:-2]).split(\"_otu_\")[1:]) for i in allplots]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "def plot_damage(fwd, rev, sample, organism):\n",
    "    f = pd.read_csv(fwd, skiprows=3, delimiter=\"\\t\")\n",
    "    r =  pd.read_csv(rev, skiprows=3, delimiter=\"\\t\")\n",
    "    r['pos'] = list(r['pos']*-1)[::-1]\n",
    "    fig = plt.figure(figsize=(18,3))\n",
    "    fig.suptitle(f\"Sample: {sample} - Species: {organism.replace('_',' ')}\", fontsize=\"x-large\", fontweight='bold', y = 1.2)\n",
    "    display(Markdown(f\"- **Sample**: {sample} - **Species**: *{organism.replace('_',' ')}*\"))\n",
    "    plt.title('Test')\n",
    "    plt.subplot(1, 2, 1)\n",
    "    plt.plot(f['pos'],f['5pC>T'])\n",
    "    plt.title('5pC>T')\n",
    "    plt.ylabel('Frequency')\n",
    "    plt.xticks(f['pos'])\n",
    "\n",
    "    plt.subplot(1, 2, 2)\n",
    "    plt.plot(r['pos'][::-1],r['3pG>A'], color = 'red')\n",
    "    plt.title('3pG>A')\n",
    "    plt.ylabel('Frequency')\n",
    "    plt.xticks(r['pos'])\n",
    "    plt.show()\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "remove_cell"
    ]
   },
   "outputs": [],
   "source": [
    "for sa in samples:\n",
    "    for sp in species:\n",
    "        try:\n",
    "            plot_damage(f\"{sa}_otu_{sp}{fwd_ext}\", f\"{sa}_otu_{sp}{rev_ext}\", sa, sp)\n",
    "        except FileNotFoundError:\n",
    "            print(f\"{sa}_{sp}{fwd_ext}\")\n",
    "            continue"
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Tags",
  "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.7.3-final"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": false,
   "sideBar": false,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": false,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}