Peptide_Extraction_Workflow_1C.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# SARS-ARENA: Structure-based identification of SARS-derived peptides with potential to induce broad protective immunity\n",
    "\n",
    "## *Workflow 1C* - Sequence Alignment and Peptide Selection"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Welcome to the Sequence Alignment and Peptide Selection **Workflow 1C**. This notebook will allow you to recover beta-coronaviruses protein sequences already aligned and search for conserved regions according to your especifications. In the first part you will recover the consensus sequence from SARS-CoV-2 from a pre-computed alignment. The protein sequences are aligned in a weekly basis, so you have information that is updated. After that, you will be able to alingn and recover conserved peptides using other beta-coronaviruses proteins in your input. Just like Workflow 1A and 1B, the peptide's list output can be used in the subsequent workflow. (**Workflow 2: Peptide-HLA Prediction for Conserved SARS-CoV-2 Peptides**) \n",
    "\n",
    "This workflow consists of four steps: \n",
    "    1. Fetch Pre-computed Multiple Sequence Alignment (MSA) from SARS-CoV-2,\n",
    "    2. Multiple Sequence Alignment,\n",
    "    3. Computing conservation score,\n",
    "    4. Computing conserved peptides.\n",
    "    \n",
    "**In order to run a cell, first click on the cell, then press shift-enter. The code inside the cell will then be executed. Note that the content of the cell can be executed as Code or Markdown. Also, inside the cell you may find comments to explain a specific command. These comments are marked with \"#\"**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Step 1) Fetch Pre-computed Multiple Sequence Alignment (MSA) from SARS-CoV-2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1.1. Necessary imports:\n",
    "Run this cell to make the necessary imports. This cell should be run only one time, unless you close this session and open it again."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "a4a3d83a23684dee9cc0158cc6b98425",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "_ColormakerRegistry()"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# System-based imports\n",
    "import os\n",
    "import glob\n",
    "\n",
    "# Data processing\n",
    "import pandas as pd\n",
    "\n",
    "# Subprocess module\n",
    "from subprocess import PIPE, run\n",
    "import multiprocessing\n",
    "\n",
    "#For visualization and interaction purposes\n",
    "from tqdm.notebook import tqdm\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "from ipywidgets import *\n",
    "sns.set(rc={'figure.figsize':(18, 8)}) # Use seaborn style defaults and set the default figure size\n",
    "\n",
    "# For utility functions used within the code\n",
    "from SARS_Arena import *"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1.2. Setting a working directory:\n",
    "Choose an appropriate directory for storing all the files or use the default (*Peptide_Extraction_Workflow*)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "dir_of_workflow = \"./Peptide_Extraction_Workflow_1C\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.makedirs(dir_of_workflow, exist_ok=True)\n",
    "os.chdir(dir_of_workflow)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1.3. Setting month/year to recover sequences:\n",
    "\n",
    "Set up desired month and year that you want to draw sequences from, until today. Use 3-letter abbreviation for the month (Jan, Feb, ...)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "month = \"Jan\"\n",
    "year = \"2021\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Fetch the sequences and print the consensus sequence:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'aligned.faa'"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "fetch_precomputed_sequences(year, month)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Consensus sequence:\n",
      "MSDNGPQNQRNAPRITFGGPSDSTGSNQNGERSGARSKQRRPQGLPNNTASWFTALTQHGKEDLKFPRGQGVPINTNSSPDDQIGYYRRATRRIRGGDGKMKDLSPRWYFYYLGTGPEAGLPYGANKDGIIWVATEGALNTPKDHIGTRNPANNAAIVLQLPQGTTLPKGFY-AEGSRGGSQASSRSSSRSRNSSRNSTPG-SSKRTSPARMAGNGGDAALALLLLDRLNQLESKMS-------GK------G-QQQQ---GQTVTK---KSAAEASKKPRQKRTATKAYN--VTQAFGRRGPEQTQGNFGDQELIRQGTDYKHWPQIAQFAPSASAFFGMSRIGMEVTPSGTWLTYTGAIKLDDKDPNFKDQVILLNKHIDAYKTFPPT---EPKKD--------KKKKADETQALPQRQKKQQ----TVT-LLPAADLDDFSKQLQQSMSSADSTQA\n"
     ]
    }
   ],
   "source": [
    "print(\"Consensus sequence:\")\n",
    "with open('Consensus_sequence.txt', 'r') as f:\n",
    "    cons = f.read().replace('\\n', '')\n",
    "    print(cons)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Include the SARS-CoV-2 consensus sequence with other beta-coronaviruses sequences:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "filtered_sequence_file = \"protein_refs.faa\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.system(\"cp ../utils/Coronaviruses_sequences.fasta \" + filtered_sequence_file)\n",
    "\n",
    "cons = '\\n'.join(''.join(filter(lambda n: n is not None, chunk)) for chunk in helper_grouper(str(cons), 70))\n",
    "\n",
    "with open(filtered_sequence_file, \"a\") as file_object:\n",
    "    file_object.write(\"\\n>gi||cons|| N gene product [SARS-CoV-2]\\n\")\n",
    "    file_object.write(cons)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Step 2) Multiple Sequence Alignment"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this step you will perform the Multiple Sequence Alignment based on the sequences you have chosen. This is performed with the software [MAFFT](https://mafft.cbrc.jp/alignment/software/). "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2.1 Run in parallel (optional)\n",
    "In case your machine has multiple cores, you can select a specific number of cores to run the alignment. If your machine has a single core, you can move to step 2.2."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Run the cell below in case you don't know how many cores your machine has."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of cores : 12\n"
     ]
    }
   ],
   "source": [
    "print(\"Number of cores :\", multiprocessing.cpu_count())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, define the number of cores to be used (`nthreads`)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "nthreads = 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2.2 Run the Multiple Sequence Alignment \n",
    "After run this cell, you will be able to see a consensus sequence for this alignment. Choose a threshold for calculating the consensus sequence (frequencies below this threshold will have an unknown amino acid)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "threshold = 0.5\n",
    "consensus_sequence = run_msa(filtered_sequence_file, nthreads, threshold)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Step 3) Computing conservation score"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using the file with the aligned sequences, you can now score each position in terms of conservation. We offer four scoring method options:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- *Jensen-Shannon divergence score* (used as 'js_divergence') (**Recommended**)\n",
    "- *Shannon Entropy* (used as 'shannon_entropy')\n",
    "- *Property entropy* (used as 'property_entropy')\n",
    "- *Von Neumann entropy* (used as 'vn_entropy')\n",
    "\n",
    "For details on the scoring method options, please consult [Capra & Singh (2007)](https://academic.oup.com/bioinformatics/article/23/15/1875/203579)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For scoring matrices you can choose one of the BLOSUM options:\n",
    "- BLOSUM62 (**Recommended**)\n",
    "- BLOSUM35\n",
    "- BLOSUM40\n",
    "- BLOSUM45\n",
    "- BLOSUM50\n",
    "- BLOSUM80\n",
    "- BLOSUM100"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "scoring_method = 'js_divergence'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "scoring_matrix = 'blosum62' #This only applies to methods that actually use a scoring matrix for calculating conservation, like JS-divergence, else, it is ignored (e.g. Shannon Entropy)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now that the arguments have been defined, run conservation analysis, and store the conservation results in the `conservation_file` variable:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Scoring Completed!\n",
      "Results written to conservation.csv\n",
      " \n"
     ]
    }
   ],
   "source": [
    "conservation_file = conservation_analysis(scoring_method, scoring_matrix)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Step 4) Computing conserved peptides"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the final step of this workflow you will be able to compute the conservation of peptides based on residue conservation."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 4.1 Retrieve information on conservation residues:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "conservation_df = pd.read_csv(filepath_or_buffer = conservation_file,\n",
    "                              header = 0,\n",
    "                              names = ['Position', 'Score', 'Alignment'],\n",
    "                              converters={'Score': lambda x: max(float(x)*100, 0)})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "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>Position</th>\n",
       "      <th>Score</th>\n",
       "      <th>Alignment</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0</td>\n",
       "      <td>91.8892</td>\n",
       "      <td>MMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMM...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>-------SS----LS-S---SS------S---------S--SS---...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>-------YF----EF-----FF----------------H--FF---...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>-------TT----VT-----VT----------------T--VV---...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>4</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>-------PP----EP-----PP----------------P--PP---...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>711</th>\n",
       "      <td>711</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>TSSS-STSTT-ST-TS-N--DT----SS--N-T--SESSSSDDS-T...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>712</th>\n",
       "      <td>712</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>QEEE--QVSQ-FQ-SE-E--SS----VE--E-G--AASSAASSE-F...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>713</th>\n",
       "      <td>713</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>AAAA--AAEA-AN-EI-L--NE----AI--L-N--A-AEAANNI-A...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>714</th>\n",
       "      <td>714</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>--------I-----I-----VI----------------I--VV---...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>715</th>\n",
       "      <td>715</td>\n",
       "      <td>0.0000</td>\n",
       "      <td>----------------------------------------------...</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>716 rows × 3 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "     Position    Score                                          Alignment\n",
       "0           0  91.8892  MMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMMM...\n",
       "1           1   0.0000  -------SS----LS-S---SS------S---------S--SS---...\n",
       "2           2   0.0000  -------YF----EF-----FF----------------H--FF---...\n",
       "3           3   0.0000  -------TT----VT-----VT----------------T--VV---...\n",
       "4           4   0.0000  -------PP----EP-----PP----------------P--PP---...\n",
       "..        ...      ...                                                ...\n",
       "711       711   0.0000  TSSS-STSTT-ST-TS-N--DT----SS--N-T--SESSSSDDS-T...\n",
       "712       712   0.0000  QEEE--QVSQ-FQ-SE-E--SS----VE--E-G--AASSAASSE-F...\n",
       "713       713   0.0000  AAAA--AAEA-AN-EI-L--NE----AI--L-N--A-AEAANNI-A...\n",
       "714       714   0.0000  --------I-----I-----VI----------------I--VV---...\n",
       "715       715   0.0000  ----------------------------------------------...\n",
       "\n",
       "[716 rows x 3 columns]"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "conservation_df #Show the conservation by residue"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now fetch the aligned sequences, where we extract the peptides from:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "aligned_sequences_df = pd.read_csv(filepath_or_buffer = \"aligned.csv\",\n",
    "                                   header = 0,\n",
    "                                   names = ['Sequence_ID', 'Aligned_Sequences'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "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>Sequence_ID</th>\n",
       "      <th>Aligned_Sequences</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1</td>\n",
       "      <td>M--------TDNG-QSNSRN----APRITF---GVSDTSD----NN...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2</td>\n",
       "      <td>M--------ATPAPP----------RAVVF-----ANDNE----TP...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3</td>\n",
       "      <td>M--------ATPAAP----------RTISF-----ADNND----NQ...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>4</td>\n",
       "      <td>M--------ASTSGKGKNPAD----KSVKF----------------...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>5</td>\n",
       "      <td>M--------A----------------TVNW-----GDAVE------...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>60</th>\n",
       "      <td>61</td>\n",
       "      <td>M--------A------------------------------------...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>61</th>\n",
       "      <td>62</td>\n",
       "      <td>M--------T------------------------------------...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>62</th>\n",
       "      <td>63</td>\n",
       "      <td>M--------A------------------------------------...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>63</th>\n",
       "      <td>64</td>\n",
       "      <td>M--------A------------------------------------...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>64</th>\n",
       "      <td>65</td>\n",
       "      <td>M--------SDNGPQ-NQRN----APRITF--GGPSDSTG----SN...</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>65 rows × 2 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "    Sequence_ID                                  Aligned_Sequences\n",
       "0             1  M--------TDNG-QSNSRN----APRITF---GVSDTSD----NN...\n",
       "1             2  M--------ATPAPP----------RAVVF-----ANDNE----TP...\n",
       "2             3  M--------ATPAAP----------RTISF-----ADNND----NQ...\n",
       "3             4  M--------ASTSGKGKNPAD----KSVKF----------------...\n",
       "4             5  M--------A----------------TVNW-----GDAVE------...\n",
       "..          ...                                                ...\n",
       "60           61  M--------A------------------------------------...\n",
       "61           62  M--------T------------------------------------...\n",
       "62           63  M--------A------------------------------------...\n",
       "63           64  M--------A------------------------------------...\n",
       "64           65  M--------SDNGPQ-NQRN----APRITF--GGPSDSTG----SN...\n",
       "\n",
       "[65 rows x 2 columns]"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "aligned_sequences_df #Show the aligned sequences"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, before using the interactive plot to filter the peptides by conservation, we will pre-compute all the peptides in the sequences. For that, define the peptide length boundaries you want to analyze and extract the peptides:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "max_len = 10 #Maximum length of the peptide\n",
    "min_len = 8 #Minimum length of the peptide"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Extracting all peptides from sequences\n"
     ]
    },
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "013aea49b8a949d4a07c1909a530bb55",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "HBox(children=(FloatProgress(value=0.0, max=65.0), HTML(value='')))"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Post-processing for all peptide lengths\n"
     ]
    },
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "215ca8bf0292488989b91d58bd4c2c8c",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "HBox(children=(FloatProgress(value=0.0, max=2.0), HTML(value='')))"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n"
     ]
    }
   ],
   "source": [
    "extracted_peptides_from_sequences = extract_peptides(min_len, max_len, aligned_sequences_df)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 4.2 Choose the peptides based on conservation values:\n",
    "Use the sliders below the cell (after run) to set the following parameters:\n",
    "\n",
    "- *Conservation threshold (CV_cutoff)*: Conservation degree of the peptides.\n",
    "- *Rolling Median Window length (RMW_cutoff)*: As conservation values are different and not homogeneous for each position, the regions can be smoothed based on this filter. Alternatively, you can set to 1 to take conservation as it is. \n",
    "- *Peptide Length (Pep_length)*: Fetch peptides of desired length for post-processing."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "e2916514c1fd4bac91070e9d22815fd3",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "interactive(children=(FloatSlider(value=21.65286494413408, continuous_update=False, description='CV_cutoff', m…"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "interactive_plot_selection(conservation_df, extracted_peptides_from_sequences, min_len, max_len)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 4.3 Print the peptides:\n",
    "\n",
    "Print the peptides sequence based on the threshold set above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "['HFYYHGTGP', 'YFYYTGTGR', 'YFYFLGTGP', 'YFYYTGTGP', 'YFYYLGTGP', 'YFYPLGHGP', 'HFYYLGTGP', 'EFPLLSEYP', 'AFYYTGTGP', 'HFYFLGTGP', 'FFYFLGTGP', 'SFYYTGTGP', 'FFYYTGTGP', 'FFYYLGTGP']\n"
     ]
    }
   ],
   "source": [
    "peptide_file = open(\"peptides.list\", \"r\")\n",
    "peptide_list = peptide_file.readlines()\n",
    "peptide_list = [peptide.strip() for peptide in peptide_list] \n",
    "print(peptide_list)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<font size=\"+2\"><center><b>This is the end of Workflow 1C</font></center></b>\n",
    "\n",
    "\n",
    "You will find a file named *peptides.list* in your folder that can be used as input for the [Workflow 2](Peptide-HLA_Binding_Prediction_Workflow.ipynb)."
   ]
  }
 ],
 "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.6.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}