indexed_gset.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import json\n",
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "fasta_file_path = './data/genome/gencode.v38.transcripts.fa'\n",
    "json_geneset_fname = './data/c2.cp.kegg_medicus.v2023.2.Hs.json'\n",
    "\n",
    "target_gene_sets = [\"M47495\", \"M47605\", \"M47715\"]\n",
    "\n",
    "\n",
    "out_indexed_gset = f'./output/gene_sets_with_indices__{\"_\".join(target_gene_sets)}.tsv'\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# M47495\n",
    "\"geneSymbols\": [\n",
    "            \"ARAF\",\n",
    "            \"BRAF\",\n",
    "            \"EGFR\",\n",
    "            \"EREG\",\n",
    "            \"GRB2\",\n",
    "            \"HRAS\",\n",
    "            \"KRAS\",\n",
    "            \"MAP2K1\",\n",
    "            \"MAP2K2\",\n",
    "            \"MAPK1\",\n",
    "            \"MAPK3\",\n",
    "            \"NRAS\",\n",
    "            \"RAF1\",\n",
    "            \"SOS1\",\n",
    "            \"SOS2\"\n",
    "        ],\n",
    "\n",
    "\n",
    "# M47605\n",
    "\"geneSymbols\": [\n",
    "            \"ALAD\",\n",
    "            \"ALAS1\",\n",
    "            \"ALAS2\",\n",
    "            \"CPOX\",\n",
    "            \"FECH\",\n",
    "            \"HMBS\",\n",
    "            \"PPOX\",\n",
    "            \"UROD\",\n",
    "            \"UROS\"\n",
    "        ],\n",
    "\n",
    "\n",
    "# M47715\n",
    "\"geneSymbols\": [\n",
    "            \"ABI1\",\n",
    "            \"ACTB\",\n",
    "            \"ACTG1\",\n",
    "            \"ACTR2\",\n",
    "            \"ACTR3\",\n",
    "            \"ARF1\",\n",
    "            \"ARF6\",\n",
    "            \"ARPC1A\",\n",
    "            \"ARPC1B\",\n",
    "            \"ARPC2\",\n",
    "            \"ARPC3\",\n",
    "            \"ARPC4\",\n",
    "            \"ARPC5\",\n",
    "            \"ARPC5L\",\n",
    "            \"BRK1\",\n",
    "            \"CYFIP1\",\n",
    "            \"CYFIP2\",\n",
    "            \"CYTH1\",\n",
    "            \"CYTH2\",\n",
    "            \"CYTH3\",\n",
    "            \"CYTH4\",\n",
    "            \"NCKAP1\",\n",
    "            \"WASF1\",\n",
    "            \"WASF2\",\n",
    "            \"WASF3\"\n",
    "        ],        "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Extract Target genes with matching index from Fasta file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "All genes were found in the reference Fasta file.\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>gene_set</th>\n",
       "      <th>target_gene</th>\n",
       "      <th>index</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>M47495</td>\n",
       "      <td>ARAF</td>\n",
       "      <td>230520</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>M47495</td>\n",
       "      <td>BRAF</td>\n",
       "      <td>94630</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>M47495</td>\n",
       "      <td>EGFR</td>\n",
       "      <td>88559</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>M47495</td>\n",
       "      <td>EREG</td>\n",
       "      <td>56700</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>M47495</td>\n",
       "      <td>GRB2</td>\n",
       "      <td>194138</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>M47495</td>\n",
       "      <td>HRAS</td>\n",
       "      <td>122636</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>M47495</td>\n",
       "      <td>KRAS</td>\n",
       "      <td>139341</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>M47495</td>\n",
       "      <td>MAP2K1</td>\n",
       "      <td>167339</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>M47495</td>\n",
       "      <td>MAP2K2</td>\n",
       "      <td>202275</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>M47495</td>\n",
       "      <td>MAPK1</td>\n",
       "      <td>224591</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>M47495</td>\n",
       "      <td>MAPK3</td>\n",
       "      <td>176266</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>M47495</td>\n",
       "      <td>NRAS</td>\n",
       "      <td>10928</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>M47495</td>\n",
       "      <td>RAF1</td>\n",
       "      <td>39279</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13</th>\n",
       "      <td>M47495</td>\n",
       "      <td>SOS1</td>\n",
       "      <td>23554</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14</th>\n",
       "      <td>M47495</td>\n",
       "      <td>SOS2</td>\n",
       "      <td>156967</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>15</th>\n",
       "      <td>M47605</td>\n",
       "      <td>ALAD</td>\n",
       "      <td>111073</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>16</th>\n",
       "      <td>M47605</td>\n",
       "      <td>ALAS1</td>\n",
       "      <td>43155</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17</th>\n",
       "      <td>M47605</td>\n",
       "      <td>ALAS2</td>\n",
       "      <td>231417</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>18</th>\n",
       "      <td>M47605</td>\n",
       "      <td>CPOX</td>\n",
       "      <td>45042</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>19</th>\n",
       "      <td>M47605</td>\n",
       "      <td>FECH</td>\n",
       "      <td>199869</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>20</th>\n",
       "      <td>M47605</td>\n",
       "      <td>HMBS</td>\n",
       "      <td>135054</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>21</th>\n",
       "      <td>M47605</td>\n",
       "      <td>PPOX</td>\n",
       "      <td>14557</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>22</th>\n",
       "      <td>M47605</td>\n",
       "      <td>UROD</td>\n",
       "      <td>6022</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>23</th>\n",
       "      <td>M47605</td>\n",
       "      <td>UROS</td>\n",
       "      <td>121910</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>24</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ABI1</td>\n",
       "      <td>115668</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ACTB</td>\n",
       "      <td>85522</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>26</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ACTG1</td>\n",
       "      <td>195589</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>27</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ACTR2</td>\n",
       "      <td>25280</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>28</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ACTR3</td>\n",
       "      <td>29321</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>29</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ARF1</td>\n",
       "      <td>19145</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>30</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ARF6</td>\n",
       "      <td>156935</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>31</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ARPC1A</td>\n",
       "      <td>91267</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>32</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ARPC1B</td>\n",
       "      <td>91272</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>33</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ARPC2</td>\n",
       "      <td>35485</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>34</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ARPC3</td>\n",
       "      <td>146965</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>35</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ARPC4</td>\n",
       "      <td>38707</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>36</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ARPC5</td>\n",
       "      <td>16215</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>37</th>\n",
       "      <td>M47715</td>\n",
       "      <td>ARPC5L</td>\n",
       "      <td>111649</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>38</th>\n",
       "      <td>M47715</td>\n",
       "      <td>BRK1</td>\n",
       "      <td>38958</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>39</th>\n",
       "      <td>M47715</td>\n",
       "      <td>CYFIP1</td>\n",
       "      <td>162942</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>40</th>\n",
       "      <td>M47715</td>\n",
       "      <td>CYFIP2</td>\n",
       "      <td>71930</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>41</th>\n",
       "      <td>M47715</td>\n",
       "      <td>CYTH1</td>\n",
       "      <td>195088</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>42</th>\n",
       "      <td>M47715</td>\n",
       "      <td>CYTH2</td>\n",
       "      <td>211695</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>43</th>\n",
       "      <td>M47715</td>\n",
       "      <td>CYTH3</td>\n",
       "      <td>85629</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>44</th>\n",
       "      <td>M47715</td>\n",
       "      <td>CYTH4</td>\n",
       "      <td>226710</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>45</th>\n",
       "      <td>M47715</td>\n",
       "      <td>NCKAP1</td>\n",
       "      <td>33304</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>46</th>\n",
       "      <td>M47715</td>\n",
       "      <td>WASF1</td>\n",
       "      <td>81492</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>47</th>\n",
       "      <td>M47715</td>\n",
       "      <td>WASF2</td>\n",
       "      <td>3556</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>48</th>\n",
       "      <td>M47715</td>\n",
       "      <td>WASF3</td>\n",
       "      <td>150187</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   gene_set target_gene   index\n",
       "0    M47495        ARAF  230520\n",
       "1    M47495        BRAF   94630\n",
       "2    M47495        EGFR   88559\n",
       "3    M47495        EREG   56700\n",
       "4    M47495        GRB2  194138\n",
       "5    M47495        HRAS  122636\n",
       "6    M47495        KRAS  139341\n",
       "7    M47495      MAP2K1  167339\n",
       "8    M47495      MAP2K2  202275\n",
       "9    M47495       MAPK1  224591\n",
       "10   M47495       MAPK3  176266\n",
       "11   M47495        NRAS   10928\n",
       "12   M47495        RAF1   39279\n",
       "13   M47495        SOS1   23554\n",
       "14   M47495        SOS2  156967\n",
       "15   M47605        ALAD  111073\n",
       "16   M47605       ALAS1   43155\n",
       "17   M47605       ALAS2  231417\n",
       "18   M47605        CPOX   45042\n",
       "19   M47605        FECH  199869\n",
       "20   M47605        HMBS  135054\n",
       "21   M47605        PPOX   14557\n",
       "22   M47605        UROD    6022\n",
       "23   M47605        UROS  121910\n",
       "24   M47715        ABI1  115668\n",
       "25   M47715        ACTB   85522\n",
       "26   M47715       ACTG1  195589\n",
       "27   M47715       ACTR2   25280\n",
       "28   M47715       ACTR3   29321\n",
       "29   M47715        ARF1   19145\n",
       "30   M47715        ARF6  156935\n",
       "31   M47715      ARPC1A   91267\n",
       "32   M47715      ARPC1B   91272\n",
       "33   M47715       ARPC2   35485\n",
       "34   M47715       ARPC3  146965\n",
       "35   M47715       ARPC4   38707\n",
       "36   M47715       ARPC5   16215\n",
       "37   M47715      ARPC5L  111649\n",
       "38   M47715        BRK1   38958\n",
       "39   M47715      CYFIP1  162942\n",
       "40   M47715      CYFIP2   71930\n",
       "41   M47715       CYTH1  195088\n",
       "42   M47715       CYTH2  211695\n",
       "43   M47715       CYTH3   85629\n",
       "44   M47715       CYTH4  226710\n",
       "45   M47715      NCKAP1   33304\n",
       "46   M47715       WASF1   81492\n",
       "47   M47715       WASF2    3556\n",
       "48   M47715       WASF3  150187"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "with open(json_geneset_fname, 'r') as json_file:\n",
    "    gene_sets = json.load(json_file)\n",
    "\n",
    "# Extract the gene names from the JSON file - gene sets\n",
    "gene_set_data = []\n",
    "for gene_set_name, gene_set_info in gene_sets.items():\n",
    "    if gene_set_info[\"systematicName\"] in target_gene_sets:\n",
    "        for gene in gene_set_info[\"geneSymbols\"]:\n",
    "            gene_set_data.append((gene_set_info[\"systematicName\"], gene))\n",
    "\n",
    "# Map gene names to their first occurrence index\n",
    "gene_index_map = {}\n",
    "with open(fasta_file_path, 'r') as fasta:\n",
    "    index = 1\n",
    "    for line in fasta:\n",
    "        if line.startswith('>'):\n",
    "            gene_name = line.strip().split('|')[5]\n",
    "            if gene_name not in gene_index_map:\n",
    "                gene_index_map[gene_name] = index\n",
    "            index += 1\n",
    "\n",
    "# Compare the gene names from the JSON file to the FASTA file\n",
    "genes_not_found = []\n",
    "output_data = []\n",
    "for gene_set_name, gene in gene_set_data:\n",
    "    if gene in gene_index_map:\n",
    "        output_data.append((gene_set_name, gene, gene_index_map[gene]))\n",
    "    else:\n",
    "        genes_not_found.append(gene)\n",
    "\n",
    "# Names of genes not found in the FASTA file\n",
    "if len(genes_not_found) > 0:\n",
    "    _go = \", \".join(genes_not_found)\n",
    "    print(f\"Genes not found in FASTA file:{_go}\")\n",
    "else:\n",
    "    print('All genes were found in the reference Fasta file.')\n",
    "\n",
    "# Write output table: gene set - target gene - index in Fasta\n",
    "output_df = pd.DataFrame(output_data, columns=['gene_set', 'target_gene', 'index'])\n",
    "output_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "min      3556\n",
       "max    231417\n",
       "Name: index, dtype: int64"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Verify index\n",
    "output_df['index'].agg(['min','max'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Filter out shared genes\n",
    "\n",
    "handle shared genes between gene sets"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(49, 3)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "output_df.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(49, 3)"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# keep first occurance\n",
    "df_filtered = output_df.drop_duplicates(subset='target_gene', keep='first')\n",
    "df_filtered.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Save resulted table"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "df_filtered.to_csv(out_indexed_gset, sep='\\t', index=False)"
   ]
  }
 ],
 "metadata": {
  "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.9.16"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}