notebooks/tutorial_scRNA-seq.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "42af5d72-c713-4374-bb64-c70494a284cf",
   "metadata": {
    "tags": []
   },
   "source": [
    "# 1. Import modules"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "4a2f04dd-c40e-4a9d-b1e1-3dfcbc4dee80",
   "metadata": {},
   "outputs": [],
   "source": [
    "import scanpy as sc\n",
    "import scGeneClust as gc"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "800fad03-e86f-4e72-b41d-245be13e470b",
   "metadata": {},
   "source": [
    "# 2. Load the example PBMC3k dataset\n",
    "GeneClust expects raw counts as input."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "40df6b66-9197-4a69-9fab-d94eac5c525e",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "AnnData object with n_obs × n_vars = 2700 × 32738\n",
       "    var: 'gene_ids'"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pbmc3k = sc.read_h5ad(\"../data/pbmc3k_raw.h5ad\")\n",
    "pbmc3k"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2fb47b40-ec19-437a-9bb2-3ea4e3a461ab",
   "metadata": {},
   "source": [
    "# 3. Quality control"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "c5e1f3bc-22e4-4401-88e4-4aac4a285961",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "AnnData object with n_obs × n_vars = 2685 × 11135\n",
       "    obs: 'n_genes'\n",
       "    var: 'gene_ids', 'n_cells'"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sc.pp.filter_cells(pbmc3k, min_genes=300)\n",
    "sc.pp.filter_genes(pbmc3k, min_cells=10)\n",
    "pbmc3k.X = pbmc3k.X.toarray()\n",
    "pbmc3k"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "49e242cc-f6de-4fe8-89b6-a55df4aaa31e",
   "metadata": {},
   "source": [
    "# 4. Feature selection\n",
    "## Use GeneClust-fast to select features"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "ec0d31d0-cb4e-4305-881c-aa337ffbfa2c",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[32m2023-02-12 15:14:18.107\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mPerforming \u001b[35mGeneClust-fast\u001b[0m\u001b[1m on \u001b[35mscRNA-seq\u001b[0m\u001b[1m data, with \u001b[33m95\u001b[0m\u001b[1m workers.\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:18.156\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mPreprocessing data...\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:19.207\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mData preprocessing done.\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:19.208\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mClustering genes...\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:19.657\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mGene clustering done!\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:20.182\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mSelected \u001b[33m396\u001b[0m\u001b[1m genes.\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:20.183\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mGeneClust-fast finished.\n",
      "\u001b[0m"
     ]
    }
   ],
   "source": [
    "info, selected_genes_fast = gc.scGeneClust(pbmc3k, n_var_clusters=200, version='fast', return_info=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "0db2cc0b-3678-4033-a076-05b90d18168e",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(['TSFM', 'MTF1', 'LONRF1', 'WDR13', 'RPL39', 'SEPT8', 'S100A11',\n",
       "       'TTC7A', 'KDELC2', 'GNPTG'], dtype=object)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "selected_genes_fast[:10]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "feaa6566-b0fb-4239-933f-31cf3a9bb581",
   "metadata": {},
   "source": [
    "## Use GeneClust-ps to select features"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "d7608308-3518-4886-b440-c64de9eeffe9",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[32m2023-02-12 15:14:20.563\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mPerforming \u001b[35mGeneClust-ps\u001b[0m\u001b[1m on \u001b[35mscRNA-seq\u001b[0m\u001b[1m data, with \u001b[33m95\u001b[0m\u001b[1m workers.\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:20.613\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mPreprocessing data...\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:21.913\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mData preprocessing done.\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:21.914\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mClustering genes...\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:21.915\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mFinding high-confidence cells...\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:24.594\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mFound \u001b[33m1641\u001b[0m\u001b[1m (\u001b[33m61.0%\u001b[0m\u001b[1m) high-confidence cells.\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:24.609\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mFinding relevant genes...\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:28.836\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1m\u001b[33m2227\u001b[0m\u001b[1m (\u001b[33m20%\u001b[0m\u001b[1m) genes are marked as relevant genes.\n",
      "\u001b[0m\u001b[32m2023-02-12 15:14:28.839\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mComputing gene redundancy...\n",
      "\u001b[0m\u001b[32m2023-02-12 15:15:52.055\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mGene redundancy computed.\n",
      "\u001b[0m\u001b[32m2023-02-12 15:15:52.058\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mComputing gene complementarity...\n",
      "\u001b[0m\u001b[32m2023-02-12 15:15:56.123\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mGene complementarity computed.\n",
      "\u001b[0m\u001b[32m2023-02-12 15:15:56.188\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mGene clustering done!\n",
      "\u001b[0m\u001b[32m2023-02-12 15:15:56.375\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mSelected \u001b[33m251\u001b[0m\u001b[1m genes.\n",
      "\u001b[0m\u001b[32m2023-02-12 15:15:56.376\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[1mGeneClust-ps finished.\n",
      "\u001b[0m"
     ]
    }
   ],
   "source": [
    "info, selected_genes_ps = gc.scGeneClust(pbmc3k, n_obs_clusters=7, version='ps', return_info=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "c15528ba-5729-471d-b860-97bd6ea133b0",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(['TNFRSF4', 'MIB2', 'PIK3CD', 'RBP7', 'EXOSC10', 'EFHD2', 'CDA',\n",
       "       'RCAN3', 'LDLRAP1', 'FGR'], dtype=object)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "selected_genes_ps[:10]"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "bioinfo",
   "language": "python",
   "name": "bioinfo"
  },
  "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.9.15"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}