0_0-Yao2023_Brain_signatures.ipynb

{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "id": "ca0f8968-07f1-4764-8634-2a567e7dcaf4",
   "metadata": {},
   "source": [
    "Download data from BrainAtlas_links.txt\n",
    "\n",
    "DE genes from [Allen whole brain atlas](https://www.nature.com/articles/s41586-023-06812-z) and [here](https://alleninstitute.github.io/abc_atlas_access/descriptions/WMB-taxonomy.html)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "20f92b78-c30e-468e-9626-a5ce52b9ecd1",
   "metadata": {},
   "outputs": [],
   "source": [
    "import urllib.request\n",
    "from pathlib import Path\n",
    "\n",
    "import anndata as ad\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import scanpy as sc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "d4b475cb-b41e-4de6-b200-7e08dc545547",
   "metadata": {},
   "outputs": [],
   "source": [
    "data_path = Path(\"data\") / \"ABC_atlas\"\n",
    "\n",
    "de_genes_file = data_path / \"ABC_atlas\" / \"DE_genes.xlsx\"\n",
    "de_genes_atlas_file = data_path / \"ABC_atlas\" / \"DE_genes_atlas.xlsx\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "419705c4-1838-448d-9f34-f51f9e93cfcb",
   "metadata": {},
   "outputs": [],
   "source": [
    "pd.read_excel(\n",
    "    de_genes_file, sheet_name=\"DE_gene_list\", header=None, names=[\"gene\"]\n",
    ").loc[:, \"gene\"].to_csv(Path(\".\") / \"8k_genes.txt\", header=False, index=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "488393cb-a159-48b6-b14d-2b96ea144fd8",
   "metadata": {},
   "outputs": [],
   "source": [
    "de_genes_file.parent.mkdir(exist_ok=True, parents=True)\n",
    "\n",
    "url = \"https://static-content.springer.com/esm/art%3A10.1038%2Fs41586-023-06812-z/MediaObjects/41586_2023_6812_MOESM6_ESM.xlsx\"\n",
    "_ = urllib.request.urlretrieve(url, de_genes_file)\n",
    "\n",
    "url = \"https://allen-brain-cell-atlas.s3-us-west-2.amazonaws.com/metadata/WMB-taxonomy/20231215/cl.df_CCN202307220.xlsx\"\n",
    "_ = urllib.request.urlretrieve(url, de_genes_atlas_file)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "0aa1c741-b2c9-4a0b-ab17-a3803b8d55fc",
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "region_annotation = pd.read_csv(data_path / \"region_of_interest_metadata.csv\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "e1d36d49-2e11-4b9a-8cf1-62cb0265a995",
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "genes = pd.read_csv(data_path / \"gene.csv\", index_col=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "57bcb66e-02aa-4114-b200-7d0f44242637",
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "cell_metadata = pd.read_csv(\n",
    "    data_path / \"cell_metadata_with_cluster_annotation.csv\",\n",
    "    usecols=[\"cell_label\", \"class\", \"subclass\"],\n",
    "    index_col=\"cell_label\",\n",
    "    dtype={\"class\": \"category\", \"subclass\": \"category\"},\n",
    "    low_memory=False,\n",
    ").loc[:, lambda df: ~df.columns.str.contains(\"color\")]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "e5451e18-aa70-4934-8912-33e5a4330d26",
   "metadata": {},
   "outputs": [],
   "source": [
    "adata = ad.concat([ad.read_h5ad(f) for f in (data_path / \"h5ad_files\").glob(\"*.h5ad\")])\n",
    "\n",
    "adata.obs = adata.obs.join(cell_metadata[[\"class\", \"subclass\"]])\n",
    "\n",
    "adata = adata[adata.obs[\"class\"].notna()]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "a538411a-7caa-485b-ab05-5789b4b518ed",
   "metadata": {},
   "outputs": [],
   "source": [
    "# subsample to 5,000 cells per class\n",
    "target_cells = 5_000\n",
    "\n",
    "subsampled = []\n",
    "for ct in adata.obs[\"class\"].cat.categories:\n",
    "    subsampled.append(adata[adata.obs[\"class\"] == ct])\n",
    "    if subsampled[-1].n_obs > target_cells:\n",
    "        sc.pp.subsample(subsampled[-1], n_obs=target_cells, random_state=1)\n",
    "\n",
    "adata = ad.concat(subsampled)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "190cbba4-b474-4662-bcda-5ec25d13586c",
   "metadata": {},
   "outputs": [],
   "source": [
    "adata.write_h5ad(data_path / \"ABC_subsampled.h5ad\")"
   ]
  },
  {
   "cell_type": "raw",
   "id": "0fa0b7b3-73b4-48e9-bfc8-5e6b48287071",
   "metadata": {},
   "source": [
    "adata = ad.read_h5ad(data_path / \"ABC_subsampled.h5ad\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "e3f7fc58-45fe-4df2-a06e-cdded805e6cd",
   "metadata": {},
   "outputs": [],
   "source": [
    "# drop empty genes\n",
    "adata = adata[:, adata.X.sum(axis=0) != 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "664aa206-d7ea-4b61-98d4-e597a01ea867",
   "metadata": {},
   "outputs": [],
   "source": [
    "anno = adata.obs[\"class\"].astype(str)\n",
    "anno.loc[adata.obs[\"subclass\"] == \"149 PVT-PT Ntrk1 Glut\"] = \"149 PVT-PT Ntrk1 Glut\"\n",
    "anno = anno.astype(\"category\")\n",
    "anno.name = \"manual_anno_pvt\"\n",
    "\n",
    "adata.obs = adata.obs.join(anno)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "9ab5a77e-e77f-4232-abda-8fff123d3bbf",
   "metadata": {},
   "outputs": [],
   "source": [
    "adata.var = adata.var.join(genes[\"gene_symbol\"]).reset_index().set_index(\"gene_symbol\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "319859fe-93fa-474c-abff-f5fe2511069d",
   "metadata": {},
   "source": [
    "# Gene signatures"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "75686b98-8d30-4a63-b8f5-4e31a89e4cee",
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pp.normalize_total(adata, target_sum=1e4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "1bc7899e-747f-4af2-8b23-0e6e2f7deaf4",
   "metadata": {},
   "outputs": [],
   "source": [
    "from utils import celltype_signatures"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "6a125eef-5d5b-403a-b54a-4782f63db81c",
   "metadata": {},
   "outputs": [],
   "source": [
    "celltype_signatures(adata, celltype_col=\"class\").to_csv(\n",
    "    \"abc_brain_signatures.tsv\", sep=\"\\t\"\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "c47402b0-185a-4d62-b1a2-d35dfddd0260",
   "metadata": {},
   "outputs": [],
   "source": [
    "celltype_signatures(adata, celltype_col=\"manual_anno_pvt\").to_csv(\n",
    "    \"abc_brain_signatures_pvt.tsv\", sep=\"\\t\"\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "70456a3d-75b1-4c00-8890-2fc5b96ea202",
   "metadata": {},
   "outputs": [],
   "source": [
    "sc.pp.log1p(adata)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "f7762c54-2b3a-4da9-8079-bdf5b9b4a0e7",
   "metadata": {},
   "outputs": [],
   "source": [
    "celltype_signatures(adata, celltype_col=\"class\").to_csv(\n",
    "    \"abc_brain_signatures_log.tsv\", sep=\"\\t\"\n",
    ")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [conda env:muellni-abc_atlas]",
   "language": "python",
   "name": "conda-env-muellni-abc_atlas-py"
  },
  "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.12.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}