{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "0",
   "metadata": {},
   "source": [
    "# Visualizing tools’ results\n",
    "\n",
    "graph stuff:\n",
    "- https://holoviews.org/gallery/demos/bokeh/network_graph.html\n",
    "- https://holoviews.org/reference/elements/bokeh/Graph.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1",
   "metadata": {},
   "outputs": [],
   "source": [
    "%load_ext autoreload\n",
    "%autoreload 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2",
   "metadata": {},
   "outputs": [],
   "source": [
    "import holoviews as hv\n",
    "\n",
    "from hv_anndata import A, register\n",
    "from hv_anndata import scanpy as hv_sc\n",
    "\n",
    "register()\n",
    "\n",
    "hv.extension(\"bokeh\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import scanpy as sc\n",
    "from anndata import AnnData"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4",
   "metadata": {},
   "outputs": [],
   "source": [
    "adata = sc.datasets.pbmc68k_reduced()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5",
   "metadata": {},
   "source": [
    "## Embeddings"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6",
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pp.pca(adata)\n",
    "sc.pp.neighbors(adata)\n",
    "sc.tl.umap(adata)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7",
   "metadata": {},
   "source": [
    "{func}`scanpy.pl.embedding` (i.e. {func}`scanpy.pl.pca`, {func}`scanpy.pl.umap`, {func}`scanpy.pl.diffmap`, {func}`scanpy.pl.tsne`, {func}`scanpy.pl.draw_graph`)\n",
    "\n",
    "missing:\n",
    "- `add_outline`\n",
    "- `legend_position` doesn’t work for graphs\n",
    "\n",
    "missing convenience:\n",
    "- `groups` to restrict to a subset easily"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8",
   "metadata": {},
   "outputs": [],
   "source": [
    "# without edges, see `scatter` in Basic notebook\n",
    "\n",
    "hv_sc.draw_graph(adata, A.obsm[\"X_umap\"], \"distances\", [A.obs[\"bulk_labels\"]]).opts(\n",
    "    node_color=A.obs[\"bulk_labels\"],\n",
    "    node_cmap=\"tab10\",\n",
    "    aspect=\"square\",\n",
    "    show_legend=True,\n",
    "    legend_position=\"right\",\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9",
   "metadata": {},
   "source": [
    "{func}`scanpy.pl.ranking` ({func}`scanpy.pl.pca_loadings`, {func}`scanpy.pl.pca_variance_ratio`)\n",
    "\n",
    "missing:\n",
    "- `xoffset`/`yoffset` taken from dimension or `text_xoffset`/`text_yoffset` available for `Labels`: https://github.com/holoviz/holoviews/issues/3884 "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "10",
   "metadata": {},
   "outputs": [],
   "source": [
    "hv.Layout([\n",
    "    hv_sc.ranking(adata, A.varm[\"PCs\"][0]).opts(aspect=1.2),\n",
    "    hv_sc.ranking(adata, A.varm[\"PCs\"][0], include_lowest=False).opts(aspect=0.6),\n",
    "]).opts(shared_axes=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "11",
   "metadata": {},
   "source": [
    "{func}`scanpy.pl.pca_overview` is just a layout of multiple `pca_` plots"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "12",
   "metadata": {},
   "source": [
    "{func}`scanpy.pl.embedding_density`\n",
    "\n",
    "TODO: maybe better approach: https://holoviews.org/gallery/demos/bokeh/iris_density_grid.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "13",
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.tl.embedding_density(adata, basis=\"umap\", groupby=\"phase\")\n",
    "# sc.pl.embedding_density(adata, \"umap\", groupby=\"phase\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "14",
   "metadata": {},
   "outputs": [],
   "source": [
    "hv_sc.embedding_density(adata, A.obsm[\"X_umap\"], groupby=\"phase\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "15",
   "metadata": {},
   "source": [
    "## Pseudotime and clustering"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "16",
   "metadata": {},
   "outputs": [],
   "source": [
    "adata.uns[\"iroot\"] = np.flatnonzero(adata.obs[\"bulk_labels\"] == \"CD56+ NK\")[0]\n",
    "sc.tl.diffmap(adata)\n",
    "sc.tl.dpt(adata, n_branchings=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "17",
   "metadata": {},
   "source": [
    "{func}`scanpy.pl.dpt_timeseries` (`as_heatmap=True`)\n",
    "\n",
    "missing:\n",
    "- figure out why the weird extra vdim is needed"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "18",
   "metadata": {},
   "outputs": [],
   "source": [
    "markers: list[str] = [\"C1QA\", \"PSAP\", \"CD79A\", \"CD79B\", \"CST3\", \"LYZ\"]\n",
    "display(\n",
    "    hv.HeatMap(adata[:, markers], [A.obs.index, A.var.index], [A.X[:, :]])\n",
    "    * hv.VLines(adata.uns[\"dpt_changepoints\"])\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "19",
   "metadata": {},
   "outputs": [],
   "source": [
    "markers: list[str] = [\"C1QA\", \"PSAP\", \"CD79A\", \"CD79B\", \"CST3\", \"LYZ\"]\n",
    "hv.HeatMap(\n",
    "    adata[adata.obs[\"dpt_order_indices\"], markers],\n",
    "    [A.obs.index, A.var.index],\n",
    "    [A.X[:, :]],\n",
    ").opts(xticks=0, colorbar=True, width=400, height=200) * hv.VLines(\n",
    "    (adata.uns[\"dpt_changepoints\"], [0] * len(adata.uns[\"dpt_changepoints\"])),\n",
    "    vdims=[\"?\"],\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "20",
   "metadata": {},
   "source": [
    "{func}`scanpy.pl.dpt_timeseries` (`as_heatmap=False`)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "21",
   "metadata": {},
   "outputs": [],
   "source": [
    "def dpt_timeseries(adata: AnnData):\n",
    "    adata = adata[adata.obs[\"dpt_order_indices\"]].copy()\n",
    "    return (\n",
    "        hv.Overlay([\n",
    "            hv.Points(adata, [A.obs.index, A.X[:, gene]], label=gene)\n",
    "            for gene in adata.var_names\n",
    "        ])\n",
    "        * hv.VLines(\n",
    "            (\n",
    "                adata.uns[\"dpt_changepoints\"],\n",
    "                [0] * len(adata.uns[\"dpt_changepoints\"]),\n",
    "            ),\n",
    "            vdims=[\"?\"],\n",
    "        )\n",
    "    ).opts(xticks=0)\n",
    "\n",
    "\n",
    "markers = [\"C1QA\", \"PSAP\", \"CD79A\", \"CD79B\", \"CST3\", \"LYZ\"]\n",
    "dpt_timeseries(adata[:, markers]).opts(width=800, legend_position=\"right\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "22",
   "metadata": {},
   "source": [
    "{func}`scanpy.pl.dpt_groups_pseudotime`\n",
    "\n",
    "Skip for now, since it’s broken in scanpy: https://github.com/scverse/scanpy/issues/3086"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "23",
   "metadata": {},
   "outputs": [],
   "source": [
    "try:\n",
    "    sc.pl.dpt_groups_pseudotime(adata)\n",
    "except Exception:  # noqa: BLE001\n",
    "    import traceback\n",
    "\n",
    "    traceback.print_exc()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "24",
   "metadata": {},
   "source": [
    "{func}`scanpy.pl.correlation_matrix`\n",
    "\n",
    "missing:\n",
    "- again: dendrogram on heatmaps\n",
    "- invalid warning:\n",
    "  > “WARNING:param.RasterPlot02026: aspect value was ignored because absolute width and height values were provided. Either supply explicit frame_width and frame_height to achieve desired aspect OR supply a combination of width or height and an aspect value.”"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "25",
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.tl.dendrogram(adata, \"bulk_labels\")\n",
    "# sc.pl.correlation_matrix(adata, groupby=\"bulk_labels\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "26",
   "metadata": {},
   "outputs": [],
   "source": [
    "def correlation_matrix(adata: AnnData, groupby: str):\n",
    "    # from scipy.sparse import coo_array\n",
    "\n",
    "    dendrogram_key = f\"dendrogram_{groupby}\"  # _get_dendrogram_key(...)\n",
    "    mat = adata.uns[dendrogram_key][\"correlation_matrix\"]\n",
    "    index = adata.uns[dendrogram_key][\"categories_idx_ordered\"]\n",
    "    mat = mat[index[::-1], :][:, index]  # TODO: why do we have to flip?  # noqa: TD003\n",
    "\n",
    "    labels = list(adata.obs[groupby].cat.categories)\n",
    "    labels = np.array(labels).astype(\"str\")[index]\n",
    "    ticks = list(enumerate(labels))\n",
    "\n",
    "    return hv.Image(\n",
    "        mat, bounds=(-0.5, -0.5, len(labels) - 0.5, len(labels) - 0.5)\n",
    "    ).opts(\n",
    "        xticks=ticks,\n",
    "        yticks=ticks,\n",
    "        xrotation=45,\n",
    "        frame_width=400,\n",
    "        frame_height=400,\n",
    "        aspect=1,\n",
    "    )\n",
    "\n",
    "\n",
    "correlation_matrix(adata, \"bulk_labels\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "27",
   "metadata": {},
   "source": [
    "## Marker genes"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "28",
   "metadata": {},
   "source": [
    "{func}`scanpy.pl.rank_genes_groups`, {func}`scanpy.pl.rank_genes_groups_dotplot`, {func}`scanpy.pl.rank_genes_groups_heatmap`, {func}`scanpy.pl.rank_genes_groups_matrixplot`, {func}`scanpy.pl.rank_genes_groups_tracksplot`, {func}`scanpy.pl.rank_genes_groups_stacked_violin`, {func}`scanpy.pl.rank_genes_groups_violin`\n",
    "\n",
    "are all just pre-parametrized versions of other plots; `rank_genes_groups` is `ranking`, the others are just the suffix."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "hv-anndata",
   "language": "python",
   "name": "hv-anndata"
  },
  "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.13.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}