08_Viral_replication_and_targets_Influenza.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import cobra\n",
    "import pandas as pd\n",
    "from cobra import Model, Reaction, Metabolite\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "\n",
    "import os\n",
    "import copy\n",
    "import time\n",
    "\n",
    "from cobra.flux_analysis import (\n",
    "    single_gene_deletion, single_reaction_deletion, double_gene_deletion,\n",
    "    double_reaction_deletion)\n",
    "\n",
    "from libsbml import *\n",
    "reader = SBMLReader()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "5227\n",
      "0.2925755033557048\n"
     ]
    }
   ],
   "source": [
    "new_model = cobra.io.read_sbml_model('Models/Recon2_Influenza_consistent.xml')\n",
    "print(len(new_model.reactions))\n",
    "print(new_model.optimize().objective_value)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def analyse_specific_model(directory, filename): \n",
    "    from cobra import Model, Reaction, Metabolite\n",
    "    '''\n",
    "    filename = e.g. lung_15.csv !!with ending!!\n",
    "    '''\n",
    "    \n",
    "    #Reading in the files from matlab\n",
    "    try:\n",
    "        \n",
    "        tissue_reacts = pd.read_csv(directory + filename, header=None)\n",
    "        \n",
    "        list_of_tissue_reacts = list(tissue_reacts[0])\n",
    "        list_of_tissue_reacts = [x - 1 for x in list_of_tissue_reacts]\n",
    "        \n",
    "        new_tissue_model = Model('Consistent_Recon2')\n",
    "\n",
    "        for i in list_of_tissue_reacts:\n",
    "            new_tissue_model.add_reactions([new_model.reactions[i]])\n",
    "\n",
    "        NumReacts = len(new_tissue_model.reactions)\n",
    "        \n",
    "        new_tissue_model.objective = 'biomass_reaction'\n",
    "        BOF = new_tissue_model.optimize().objective_value\n",
    "        \n",
    "        #reactions\n",
    "        #react_ko = single_reaction_deletion(new_tissue_model)\n",
    "        #genes\n",
    "        react_ko = single_gene_deletion(new_tissue_model)\n",
    "        react_ko = react_ko.rename(columns={'growth': 'growth_host', 'status': 'status_host'})\n",
    "\n",
    "            \n",
    "        with new_tissue_model:\n",
    "            new_tissue_model.objective='VBOF'\n",
    "            VBOF = new_tissue_model.optimize().objective_value\n",
    "            \n",
    "            #reactions\n",
    "            #react_ko_vbof = single_reaction_deletion(new_tissue_model)\n",
    "            #genes\n",
    "            react_ko_vbof = single_gene_deletion(new_tissue_model)\n",
    "            react_ko_vbof = react_ko_vbof.rename(columns={'growth': 'growth_virus', 'status': 'status_virus'})\n",
    "            \n",
    "        #Analysis-result\n",
    "        merged = react_ko.merge(react_ko_vbof, left_on='ids', right_on='ids')\n",
    "        \n",
    "        new_indices = [list(i)[0] for i in merged.index]\n",
    "        merged = merged.set_index([pd.Index(new_indices)])\n",
    "        \n",
    "        percent_bof = [merged.loc[i]['growth_host']/BOF*100 for i in merged.index]\n",
    "        merged.insert(loc = 0, column='%growth_host', value=percent_bof)\n",
    "        \n",
    "        percent_vbof = [merged.loc[i]['growth_virus']/VBOF*100 for i in merged.index]\n",
    "        merged.insert(loc = 0, column='%growth_virus', value=percent_vbof)\n",
    "        \n",
    "        possible_targets = []\n",
    "        relevant = []\n",
    "        \n",
    "        for i in merged.index:\n",
    "            #if merged.loc[i]['%growth_virus'] < 50 and merged.loc[i]['%growth_virus'] < merged.loc[i]['%growth_host']:\n",
    "            if merged.loc[i]['%growth_virus'] < 50 and merged.loc[i]['%growth_host'] > 80 and i != 'VBOF':\n",
    "                possible_targets.append([i, merged.loc[i]['%growth_virus'], merged.loc[i]['%growth_host']])\n",
    "                #print(i, merged.loc[i]['%growth_virus'], merged.loc[i]['%growth_host'])\n",
    "            if merged.loc[i]['%growth_virus'] <= 0.05:\n",
    "                relevant.append(i)\n",
    "                \n",
    "        if len(possible_targets) == 0:\n",
    "            possible_targets = 0\n",
    "        if len(relevant) == 0:\n",
    "            relevant = 0\n",
    "            \n",
    "    except:\n",
    "        BOF = np.nan\n",
    "        VBOF = np.nan\n",
    "        NumReacts = np.nan\n",
    "        possible_targets = np.nan\n",
    "        relevant = np.nan\n",
    "\n",
    "    \n",
    "    return([BOF, VBOF, NumReacts, possible_targets, relevant])\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Time:  2.97479510307312\n",
      "0.004483041295214671 0.01892189288013699 1574 0 ['HGNC:758', 'HGNC:2292', 'HGNC:2279', 'HGNC:24380', 'HGNC:12582', 'HGNC:16232', 'HGNC:2289', 'HGNC:30863', 'HGNC:2269', 'HGNC:30862', 'HGNC:2579', 'HGNC:2285', 'HGNC:10682', 'HGNC:4433', 'HGNC:18170', 'HGNC:24381', 'HGNC:4693', 'HGNC:2291', 'HGNC:12590', 'HGNC:11067', 'HGNC:2265', 'HGNC:12585', 'HGNC:7892', 'HGNC:2267', 'HGNC:11063', 'HGNC:3700', 'HGNC:10683', 'HGNC:29594', 'HGNC:2288', 'HGNC:12587', 'HGNC:2287', 'HGNC:2280', 'HGNC:16484', 'HGNC:7419', 'HGNC:24382', 'HGNC:2277', 'HGNC:11064', 'HGNC:746', 'HGNC:12586', 'HGNC:7427', 'HGNC:9020', 'HGNC:2294', 'HGNC:7421', 'HGNC:10971', 'HGNC:12563', 'HGNC:10680', 'HGNC:2501', 'HGNC:10681', 'HGNC:7422']\n"
     ]
    }
   ],
   "source": [
    "#testing\n",
    "directory = 'CoreModelsInfluenza/'\n",
    "start_time = time.time()\n",
    "bof, vbof, reacts, targets, relevant = analyse_specific_model(directory, 'GSM5740432_S2_AAACCCACAATCAAGA-1.csv')#('A10_B000232_B008361_S10_mm10_plus_7_0')#('K22_MAA000891_3_10_M_1_1_1_1')\n",
    "print('Time: ', time.time() - start_time)\n",
    "print(bof, vbof, reacts, targets, relevant)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0\n",
      "200\n",
      "400\n",
      "600\n",
      "800\n",
      "1000\n",
      "1200\n",
      "1400\n",
      "1600\n",
      "1800\n",
      "2000\n",
      "2200\n",
      "2400\n",
      "2600\n",
      "2800\n",
      "3000\n",
      "3200\n",
      "3400\n"
     ]
    }
   ],
   "source": [
    "### All models\n",
    "directory = 'CoreModelsInfluenza/'\n",
    "bof_dict = {}\n",
    "vbof_dict = {}\n",
    "react_dict = {}\n",
    "target_dict = {}\n",
    "relevant_dict = {}\n",
    "counter = 0\n",
    "\n",
    "for filename in os.listdir(directory): #infected_epithelial: \n",
    "    if counter %200 == 0:\n",
    "        print(counter)\n",
    "    bof, vbof, reacts, target, relevant = analyse_specific_model(directory, filename)\n",
    "    bof_dict[filename] = bof\n",
    "    vbof_dict[filename] = vbof\n",
    "    react_dict[filename]= reacts\n",
    "    target_dict[filename] = target\n",
    "    relevant_dict[filename] = relevant\n",
    "    counter += 1\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3503\n",
      "3066\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>BOF</th>\n",
       "      <th>VBOF</th>\n",
       "      <th>Number_Reactions</th>\n",
       "      <th>Targets (Gene,%VBOF, %BOF)</th>\n",
       "      <th>Virus relevant</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>GSM5740432_S2_GGGTTTACACAAGTTC-1.csv</th>\n",
       "      <td>0.00654094</td>\n",
       "      <td>0.0858842</td>\n",
       "      <td>1818</td>\n",
       "      <td>[[HGNC:8016, 35.890905400818994, 100.0], [HGNC...</td>\n",
       "      <td>[HGNC:10979, HGNC:4175, HGNC:18170, HGNC:8091,...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>GSM5740432_S2_ACGTTCCAGCTGAGCA-1.csv</th>\n",
       "      <td>0.00358027</td>\n",
       "      <td>0.00946095</td>\n",
       "      <td>1518</td>\n",
       "      <td>0</td>\n",
       "      <td>[HGNC:10452, HGNC:12586, HGNC:7427, HGNC:2277,...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>GSM5740432_S2_ACTACGAAGATGACCG-1.csv</th>\n",
       "      <td>0.00448304</td>\n",
       "      <td>0.0189219</td>\n",
       "      <td>1711</td>\n",
       "      <td>0</td>\n",
       "      <td>[HGNC:2265, HGNC:24380, HGNC:2289, HGNC:2269, ...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>GSM5740432_S2_AGTACCAAGGGCGAGA-1.csv</th>\n",
       "      <td>0.00448304</td>\n",
       "      <td>0.0103348</td>\n",
       "      <td>817</td>\n",
       "      <td>[[HGNC:2995, 35.55704366170918, 100.0], [HGNC:...</td>\n",
       "      <td>[HGNC:4693, HGNC:2291, HGNC:4378, HGNC:848, HG...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>GSM5740432_S2_ACTTCCGAGTATCTGC-1.csv</th>\n",
       "      <td>0.00448304</td>\n",
       "      <td>0.0189219</td>\n",
       "      <td>1124</td>\n",
       "      <td>0</td>\n",
       "      <td>[HGNC:8907, HGNC:10293, HGNC:746, HGNC:2205, H...</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                                             BOF        VBOF Number_Reactions  \\\n",
       "GSM5740432_S2_GGGTTTACACAAGTTC-1.csv  0.00654094   0.0858842             1818   \n",
       "GSM5740432_S2_ACGTTCCAGCTGAGCA-1.csv  0.00358027  0.00946095             1518   \n",
       "GSM5740432_S2_ACTACGAAGATGACCG-1.csv  0.00448304   0.0189219             1711   \n",
       "GSM5740432_S2_AGTACCAAGGGCGAGA-1.csv  0.00448304   0.0103348              817   \n",
       "GSM5740432_S2_ACTTCCGAGTATCTGC-1.csv  0.00448304   0.0189219             1124   \n",
       "\n",
       "                                                             Targets (Gene,%VBOF, %BOF)  \\\n",
       "GSM5740432_S2_GGGTTTACACAAGTTC-1.csv  [[HGNC:8016, 35.890905400818994, 100.0], [HGNC...   \n",
       "GSM5740432_S2_ACGTTCCAGCTGAGCA-1.csv                                                  0   \n",
       "GSM5740432_S2_ACTACGAAGATGACCG-1.csv                                                  0   \n",
       "GSM5740432_S2_AGTACCAAGGGCGAGA-1.csv  [[HGNC:2995, 35.55704366170918, 100.0], [HGNC:...   \n",
       "GSM5740432_S2_ACTTCCGAGTATCTGC-1.csv                                                  0   \n",
       "\n",
       "                                                                         Virus relevant  \n",
       "GSM5740432_S2_GGGTTTACACAAGTTC-1.csv  [HGNC:10979, HGNC:4175, HGNC:18170, HGNC:8091,...  \n",
       "GSM5740432_S2_ACGTTCCAGCTGAGCA-1.csv  [HGNC:10452, HGNC:12586, HGNC:7427, HGNC:2277,...  \n",
       "GSM5740432_S2_ACTACGAAGATGACCG-1.csv  [HGNC:2265, HGNC:24380, HGNC:2289, HGNC:2269, ...  \n",
       "GSM5740432_S2_AGTACCAAGGGCGAGA-1.csv  [HGNC:4693, HGNC:2291, HGNC:4378, HGNC:848, HG...  \n",
       "GSM5740432_S2_ACTTCCGAGTATCTGC-1.csv  [HGNC:8907, HGNC:10293, HGNC:746, HGNC:2205, H...  "
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "new_dict = {}\n",
    "for key in bof_dict.keys():\n",
    "    if key != '.DS_Store':\n",
    "        new_dict[key] = [bof_dict[key], vbof_dict[key], react_dict[key], target_dict[key], relevant_dict[key]]\n",
    "evaluation = pd.DataFrame.from_dict(new_dict).T\n",
    "evaluation = evaluation.rename(columns={0: \"BOF\", 1: \"VBOF\", 2: \"Number_Reactions\", 3: 'Targets (Gene,%VBOF, %BOF)', \n",
    "                                                   4: 'Virus relevant'})\n",
    "print(len(evaluation))\n",
    "evaluation = evaluation.dropna()\n",
    "print(len(evaluation))\n",
    "evaluation.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "evaluation.to_csv('csvs/Influenza_evaluation.csv')\n",
    "\n"
   ]
  }
 ],
 "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.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}