From 32cd6b70056b0ebe78b03b8fb3e60a8144b0f8c0 Mon Sep 17 00:00:00 2001
From: Jayoung Ryu <jayoung_ryu@g.harvard.edu>
Date: Sat, 25 Nov 2023 18:48:12 -0500
Subject: [PATCH] add per-variant stat

---
 notebooks/sample_quality_report.ipynb | 78 +++++++++++++++++++++++++--
 1 file changed, 75 insertions(+), 3 deletions(-)

diff --git a/notebooks/sample_quality_report.ipynb b/notebooks/sample_quality_report.ipynb
index 70b7c44..4d056ea 100644
--- a/notebooks/sample_quality_report.ipynb
+++ b/notebooks/sample_quality_report.ipynb
@@ -20,11 +20,13 @@
    "metadata": {},
    "outputs": [],
    "source": [
+    "import numpy as np\n",
+    "import matplotlib.pyplot as plt\n",
     "import perturb_tools as pt\n",
     "import bean as be\n",
     "from bean.qc.utils import fill_in_missing_samples\n",
-    "import matplotlib.pyplot as plt\n",
-    "plt.style.use('default')"
+    "\n",
+    "plt.style.use(\"default\")"
    ]
   },
   {
@@ -282,6 +284,46 @@
     "        editable_base_end=edit_quantification_end_pos\n",
     "    )\n",
     "    be.qc.plot_sample_edit_rates(bdata)"
+<<<<<<< Updated upstream
+=======
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### 5. Variant coverage"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "if not tiling:\n",
+    "    n_guides = bdata.guides.groupby(\"target\").size()\n",
+    "    int_bins = np.arange(min(0.5, n_guides.min() - 0.5), n_guides.max() + 0.5, 1)\n",
+    "    plt.hist(n_guides, bins=int_bins)\n",
+    "    plt.xticks(np.arange(n_guides.min() - 1, n_guides.max() + 1, 1))\n",
+    "    plt.title(\"# Guides per target\")"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "if \"target_base_change\" not in bdata.uns or not base_edit_data:\n",
+    "    print(\n",
+    "        \"Not a base editing data or target base change not provided. Passing editing-related QC\"\n",
+    "    )\n",
+    "elif tiling:\n",
+    "    total_edits = bdata.guides.groupby(\"target\")[\"edit_rate\"].sum()\n",
+    "    plt.hist(n_guides)\n",
+    "    plt.title(\"Total edit rates per target\")"
+>>>>>>> Stashed changes
    ]
   },
   {
@@ -348,12 +390,42 @@
     "bdata_filtered.samples.style.background_gradient(cmap=\"coolwarm_r\")"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### 5. Variant coverage"
+   ]
+  },
   {
    "cell_type": "code",
    "execution_count": null,
    "metadata": {},
    "outputs": [],
-   "source": []
+   "source": [
+    "if not tiling:\n",
+    "    n_guides = bdata.guides.groupby(\"target\").size()\n",
+    "    int_bins = np.arange(min(0.5, n_guides.min() - 0.5), n_guides.max() + 0.5, 1)\n",
+    "    plt.hist(n_guides, bins=int_bins)\n",
+    "    plt.xticks(np.arange(n_guides.min() - 1, n_guides.max() + 1, 1))\n",
+    "    plt.title(\"# Guides per target\")"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "if \"target_base_change\" not in bdata.uns or not base_edit_data:\n",
+    "    print(\n",
+    "        \"Not a base editing data or target base change not provided. Passing editing-related QC\"\n",
+    "    )\n",
+    "elif not tiling:\n",
+    "    total_edits = bdata.guides.groupby(\"target\")[\"edit_rate\"].sum()\n",
+    "    plt.hist(n_guides)\n",
+    "    plt.title(\"Total edit rates per target\")"
+   ]
   },
   {
    "cell_type": "markdown",