Bio2Byte_Biophysical_features_applied_on_proteins_3D_structure.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "private_outputs": true,
      "provenance": [],
      "collapsed_sections": [],
      "authorship_tag": "ABX9TyMh9M7eraK5WIaBgdChHqL4",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/Bio2Byte/public_notebooks/blob/main/Bio2Byte_Biophysical_features_applied_on_proteins_3D_structure.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "## The bio2Byte Toolkit to analyze AlphaFold results\n",
        "\n",
        "This Jupyter Notebook on Google Colab aims to simplify the analysis of the predicted three-dimensional structure after running AlphaFold with your sequence of interest. \n",
        "\n",
        "### About bio2Byte\n",
        "\n",
        "<img src=\"https://pbs.twimg.com/profile_images/1247824923546079232/B9b_Yg7n.jpg\" alt=\"b2b\" width=\"48\"/>\n",
        "<img src=\"https://upload.wikimedia.org/wikipedia/commons/thumb/3/31/Vrije_Universiteit_Brussel_logo.svg/1200px-Vrije_Universiteit_Brussel_logo.svg.png\" alt=\"vub\" width=\"120\"/>\n",
        "\n",
        "Proteins are the molecular machines that make cells work. They perform a wide variety of functions through interactions with each other and many additional molecules. Traditionally, proteins are described in a single static state (a picture). It is now increasingly recognised that many proteins can adopt multiple states and move between these conformational states dynamically (a movie).\n",
        "\n",
        "We investigate how the dynamics, conformational states and available experimental data of proteins relates to their amino acid sequence. Underlying physical and chemical principles are computationally unravelled through data integration, analysis and machine learning, so connecting them to biological events and improving our understanding of the way proteins work.\n",
        "\n",
        "Online predictors are available on https://bio2byte.be/b2btools/\n",
        "\n",
        "#### Tools and predictors\n",
        "\n",
        "* **DynaMine** is a fast predictor of protein backbone dynamics using only sequence information as input. The version here also predicts side-chain dynamics and secondary structure predictors using the same principle. This tool returns the following predictions: \n",
        "  * DynaMine backbone dynamics\n",
        "  * DynaMine side chain dynamics\n",
        "  * DynaMine conformational propensities (polyproline II)\n",
        "  * DynaMine conformational propensities (coil)\n",
        "  * DynaMine conformational propensities (sheet)\n",
        "  * DynaMine conformational propensities (helix)\n",
        "\n",
        "* **EFoldMine** is a method that predicts from the primary amino acid sequence of a protein, which amino acids are likely involved in early folding events. This tool returns the following predictions:\n",
        "  * EFoldMine earlyFolding propensity\n",
        "\n",
        "* **DisoMine** method predicts protein disorder with recurrent neural networks not directly from the amino acid sequence, but instead from more generic predictions of key biophysical properties, here protein dynamics, secondary structure and early folding. This tool returns the following predictions:\n",
        "  * Disomine disorder\n",
        "\n",
        "* **AgMata** is a single-sequence based predictor of protein regions that are likely to cause beta-aggregation. This tool returns the following predictions:\n",
        "  * Agmata aggregation propensity"
      ],
      "metadata": {
        "id": "R-QOv8W52z-x"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Environment setup\n",
        "\n",
        "**Notes:** Please be patiente, this section might take several minutes to download and install the dependencies in the Jupyter context..."
      ],
      "metadata": {
        "id": "UNbGbdvUjS0i"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "cellView": "form",
        "id": "sCy393Obg3dB"
      },
      "outputs": [],
      "source": [
        "#@markdown ## Install extra Python packages\n",
        "%%capture\n",
        "\n",
        "!pip install b2bTools==3.0.5 py3Dmol biopython"
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "#@markdown ## Import required libraries and define constants\n",
        "\n",
        "%%capture\n",
        "from google.colab import files\n",
        "import py3Dmol\n",
        "import matplotlib.pyplot as plt\n",
        "import pickle\n",
        "import numpy as np\n",
        "import pandas as pd\n",
        "from math import floor\n",
        "\n",
        "from IPython.display import display\n",
        "import ipywidgets as widgets\n",
        "from ipywidgets import GridspecLayout\n",
        "from ipywidgets import Output\n",
        "\n",
        "from b2bTools import SingleSeq\n",
        "from Bio import SeqIO \n",
        "from Bio.PDB import *\n",
        "\n",
        "from IPython.utils import io\n",
        "import tqdm.notebook\n",
        "import subprocess\n",
        "\n",
        "TQDM_BAR_FORMAT = '{l_bar}{bar}| {n_fmt}/{total_fmt} [elapsed: {elapsed} remaining: {remaining}]'\n",
        "\n",
        "VALUES_CONSTANTS = {\n",
        "    # \"pLDDT\": {\n",
        "    #     \"colors\": [\n",
        "    #         (0, 50, \"#FF7D45\"),\n",
        "    #         (50, 70, \"#FFDB13\"),\n",
        "    #         (70, 90, \"#65CBF3\"),\n",
        "    #         (90, 100, \"#0053D6\"),\n",
        "    #     ],\n",
        "    #     \"legends\": [\n",
        "    #         \"Very low (pLDDT < 50)\",\n",
        "    #         \"Low (50 < pLDDT < 70)\",\n",
        "    #         \"Confident (70 < pLDDT < 90)\",\n",
        "    #         \"Very high (pLDDT > 90)\",\n",
        "    #     ],\n",
        "    #     \"title\": \"AlphaFold's Model Confidence (pLDDT)\",\n",
        "    # },\n",
        "    \"backbone\": {\n",
        "        \"colors\": [\n",
        "            (0, 0.69, \"#FF7D45\"),\n",
        "            (0.69, 0.80, \"#FFDB13\"),\n",
        "            (0.80, 1.00, \"#65CBF3\"),\n",
        "            (1.00, 2.00, \"#0053D6\"),\n",
        "        ],\n",
        "        \"legends\": [\n",
        "            \"Flexible region (value < 0.69)\",\n",
        "            \"Context dependent region (0.69 < value < 0.80)\",\n",
        "            \"Rigid region (0.80 < value < 1.00)\",\n",
        "            \"Membrane spanning region (value > 1.00)\",\n",
        "        ],\n",
        "        \"title\": \"DynaMine backbone dynamics\",\n",
        "    },\n",
        "    \"sidechain\": {\n",
        "        \"colors\": [\n",
        "            (0, -0.6, \"#FF7D45\"),\n",
        "            (-0.6, 0.0, \"#FFDB13\"),\n",
        "            (0.0, 0.6, \"#65CBF3\"),\n",
        "            (0.6, 1.2, \"#0053D6\"),\n",
        "        ],\n",
        "        \"legends\": [\n",
        "            \"Lower values\",\n",
        "            \"Low values\",\n",
        "            \"High values\",\n",
        "            \"Higher values\",\n",
        "        ],\n",
        "        \"title\": \"DynaMine sidechain dynamics\",\n",
        "    },\n",
        "    \"ppII\": {\n",
        "        \"colors\": [   \n",
        "            (0.00, 0.25, \"#ff0100\"),\n",
        "            (0.25, 0.50, \"#ffff00\"),\n",
        "            (0.50, 0.75, \"#7fff00\"),\n",
        "            (0.75, 1.00, \"#01ff00\"),\n",
        "        ],\n",
        "        \"legends\": [\n",
        "            \"Lowest propensity (value = 0.00)\",\n",
        "            \"Low propensity (0.00 < value < 0.25)\",\n",
        "            \"Mid propensity (0.25 < value < 0.50)\",\n",
        "            \"High propensity (0.50 < value < 0.75)\",\n",
        "            \"Very high propensity(0.75 < value < 1.00)\",\n",
        "        ],\n",
        "        \"title\": \"DynaMine conformational propensities (polyproline II)\",\n",
        "    },\n",
        "    \"coil\": {\n",
        "        \"colors\": [   \n",
        "            (0.00, 0.25, \"#ff0100\"),\n",
        "            (0.25, 0.50, \"#ffff00\"),\n",
        "            (0.50, 0.75, \"#7fff00\"),\n",
        "            (0.75, 1.00, \"#01ff00\"),\n",
        "        ],\n",
        "        \"legends\": [\n",
        "            \"Lowest propensity (value = 0.00)\",\n",
        "            \"Low propensity (0.00 < value < 0.25)\",\n",
        "            \"Mid propensity (0.25 < value < 0.50)\",\n",
        "            \"High propensity (0.50 < value < 0.75)\",\n",
        "            \"Very high propensity(0.75 < value < 1.00)\",\n",
        "        ],\n",
        "        \"title\": \"DynaMine conformational propensities (coil)\",\n",
        "    },\n",
        "    \"sheet\": {\n",
        "        \"colors\": [   \n",
        "            (0.00, 0.25, \"#ff0100\"),\n",
        "            (0.25, 0.50, \"#ffff00\"),\n",
        "            (0.50, 0.75, \"#7fff00\"),\n",
        "            (0.75, 1.00, \"#01ff00\"),\n",
        "        ],\n",
        "        \"legends\": [\n",
        "            \"Lowest propensity (value = 0.00)\",\n",
        "            \"Low propensity (0.00 < value < 0.25)\",\n",
        "            \"Mid propensity (0.25 < value < 0.50)\",\n",
        "            \"High propensity (0.50 < value < 0.75)\",\n",
        "            \"Very high propensity(0.75 < value < 1.00)\",\n",
        "        ],\n",
        "        \"title\": \"DynaMine conformational propensities (sheet)\",\n",
        "    },\n",
        "    \"helix\": {\n",
        "        \"colors\": [   \n",
        "            (0.00, 0.25, \"#ff0100\"),\n",
        "            (0.25, 0.50, \"#ffff00\"),\n",
        "            (0.50, 0.75, \"#7fff00\"),\n",
        "            (0.75, 1.00, \"#01ff00\"),\n",
        "        ],\n",
        "        \"legends\": [\n",
        "            \"Lowest propensity (value = 0.00)\",\n",
        "            \"Low propensity (0.00 < value < 0.25)\",\n",
        "            \"Mid propensity (0.25 < value < 0.50)\",\n",
        "            \"High propensity (0.50 < value < 0.75)\",\n",
        "            \"Very high propensity(0.75 < value < 1.00)\",\n",
        "        ],\n",
        "        \"title\": \"DynaMine conformational propensities (helix)\",\n",
        "        \"subtitle\": \"DynaMine is a fast predictor of protein backbone dynamics using only sequence information as input. The version here also predicts side-chain dynamics and secondary structure predictors using the same principle.\",\n",
        "    },\n",
        "    \"earlyFolding\": {\n",
        "        \"colors\": [\n",
        "            (0, 0.169, \"#FF7D45\"),\n",
        "            (0.169 ,2, \"#0053D6\"),\n",
        "        ],\n",
        "        \"legends\": [\n",
        "            \"Not likely to start the protein folding process (value < 0.169)\",\n",
        "            \"Likely to start the protein folding process (value > 0.169)\",\n",
        "        ],\n",
        "        \"title\": \"EFoldMine earlyFolding propensity\",\n",
        "    },\n",
        "    \"disoMine\": {\n",
        "        \"colors\": [\n",
        "            (0, 0.50, \"#FF7D45\"),\n",
        "            (0.50, 1.00, \"#0053D6\"),\n",
        "        ],\n",
        "        \"legends\": [\n",
        "            \"Not likely a disordered residue (value < 0.50)\",\n",
        "            \"Likely a disordered residue (value > 0.50)\",\n",
        "        ],\n",
        "        \"title\": \"Disomine disorder\",\n",
        "    },\n",
        "    \"agmata\": {\n",
        "        \"colors\": [\n",
        "            (0.00, 0.25, \"#ff0100\"),\n",
        "            (0.50, 0.75, \"#ffff00\"),\n",
        "            (0.75, 1.00, \"#01ff00\"),\n",
        "        ],\n",
        "        \"legends\": [\n",
        "            \"Not likely to be involved in beta-sheet aggregation (Trough)\",\n",
        "            \"Likely to be involved in beta-sheet aggregation (Mid Peak)\",\n",
        "            \"Most likely to be involved in beta-sheet aggregation (Higher Peak)\",\n",
        "        ],\n",
        "        \"title\": \"Agmata aggregation propensity\",\n",
        "    },\n",
        "}\n",
        "\n",
        "B2B_TOOLS = [\"dynamine\", \"disomine\", \"efoldmine\", \"agmata\"]\n",
        "\n",
        "B2B_TOOLS_RESULTS = np.array(\n",
        "    [\n",
        "        [\"backbone\", \"sidechain\", \"ppII\"],\n",
        "        [\"coil\", \"sheet\", \"helix\"],\n",
        "        [\"earlyFolding\", \"disoMine\", \"agmata\"],\n",
        "    ]\n",
        ")\n",
        "\n",
        "# PLDDT_COLOR_MAP = {i: bands[2] for i, bands in enumerate(VALUES_CONSTANTS[\"pLDDT\"][\"colors\"])}\n",
        "# PLDDT_DEFAULT_ATOM_STYLE = {\"cartoon\": {\"colorscheme\": {\"prop\": \"b\", \"map\": PLDDT_COLOR_MAP}}}"
      ],
      "metadata": {
        "cellView": "form",
        "id": "GaAsdwGciI0k"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "#@markdown ## Internal source code\n",
        "#@markdown Please execute this cell by pressing the _Play_ button \n",
        "\n",
        "class Atom(dict):\n",
        "    def __init__(self, line):\n",
        "        self[\"type\"] = line[0:6].strip()\n",
        "        self[\"idx\"] = line[6:11].strip()\n",
        "        self[\"name\"] = line[12:16].strip()\n",
        "        self[\"resname\"] = line[17:20].strip()\n",
        "        self[\"resid\"] = int(int(line[22:26]))\n",
        "        self[\"x\"] = float(line[30:38])\n",
        "        self[\"y\"] = float(line[38:46])\n",
        "        self[\"z\"] = float(line[46:54])\n",
        "        self[\"sym\"] = line[76:78].strip()\n",
        "\n",
        "    def __str__(self):\n",
        "        line = list(\" \" * 80)\n",
        "\n",
        "        line[0:6] = self[\"type\"].ljust(6)\n",
        "        line[6:11] = self[\"idx\"].ljust(5)\n",
        "        line[12:16] = self[\"name\"].ljust(4)\n",
        "        line[17:20] = self[\"resname\"].ljust(3)\n",
        "        line[22:26] = str(self[\"resid\"]).ljust(4)\n",
        "        line[30:38] = str(self[\"x\"]).rjust(8)\n",
        "        line[38:46] = str(self[\"y\"]).rjust(8)\n",
        "        line[46:54] = str(self[\"z\"]).rjust(8)\n",
        "        line[76:78] = self[\"sym\"].rjust(2)\n",
        "        return \"\".join(line) + \"\\n\"\n",
        "\n",
        "\n",
        "class Molecule(list):\n",
        "    def __init__(self, file):\n",
        "        for line in file:\n",
        "            if \"ATOM\" in line or \"HETATM\" in line:\n",
        "                self.append(Atom(line))\n",
        "\n",
        "    def __str__(self):\n",
        "        outstr = \"\"\n",
        "        for at in self:\n",
        "            outstr += str(at)\n",
        "\n",
        "        return outstr\n",
        "\n",
        "\n",
        "def check_min_max(sequence_df, former_min, former_max):\n",
        "    seq_max = max(sequence_df)\n",
        "    seq_min = min(sequence_df)\n",
        "\n",
        "    if seq_max + \\\n",
        "            0.1 > former_max and not np.isnan(seq_max) \\\n",
        "            and not np.isinf(seq_max):\n",
        "        former_max = seq_max + 0.1\n",
        "    if seq_min - \\\n",
        "            0.1 < former_min and not np.isnan(seq_min) \\\n",
        "            and not np.isinf(seq_min):\n",
        "        former_min = seq_min - 0.1\n",
        "    return former_min, former_max\n",
        "\n",
        "\n",
        "def plot_prediction(prediction_name, highlighting_regions, sequence_df):\n",
        "    thresholds_dict = {'backbone': {'membrane spanning': [1., 1.5],\n",
        "                                    'rigid': [0.8, 1.],\n",
        "                                    'context-dependent': [0.69, 0.8],\n",
        "                                    'flexible': [-1.0, 0.69]},\n",
        "                       'earlyFolding': {'early folds': [0.169, 2.],\n",
        "                                        'late folds': [-1., 0.169]},\n",
        "                       'disoMine': {'ordered': [-1., 0.5],\n",
        "                                    'disordered': [0.5, 2.]},\n",
        "                       }\n",
        "    ordered_regions_dict = {'backbone': ['flexible',\n",
        "                                         'context-dependent',\n",
        "                                         'rigid',\n",
        "                                         'membrane spanning'],\n",
        "                            'earlyFolding': ['late folds', 'early folds'],\n",
        "                            'disoMine': ['ordered', 'disordered'],\n",
        "                            }\n",
        "    colors = ['yellow', 'orange', 'pink', 'red']\n",
        "    ranges_dict = {\n",
        "        'backbone': [-0.2, 1.2],\n",
        "        'sidechain': [-0.2, 1.2],\n",
        "        'ppII': [-0.2, 1.2],\n",
        "        'earlyFolding': [-0.2, 1.2],\n",
        "        'disoMine': [-0.2, 1.2],\n",
        "        'agmata': [-0.2, 1.2],\n",
        "        'helix': [-1., 1.],\n",
        "        'sheet': [-1., 1.],\n",
        "        'coil': [-1., 1.],\n",
        "    }\n",
        "    fig, ax = plt.subplots(1, 1)\n",
        "    fig.set_figwidth(10)\n",
        "    fig.set_figheight(5)\n",
        "    ax.set_title(prediction_name + ' ' + 'prediction')\n",
        "    min_value, max_value = ranges_dict[prediction_name]\n",
        "\n",
        "    predictions = sequence_df[prediction_name]\n",
        "    min_value, max_value = check_min_max(predictions, min_value, max_value)\n",
        "    \n",
        "    ax.plot(range(len(predictions)), predictions, label=selected_sequence)\n",
        "    ax.set_xlim([0, len(predictions) - 1])\n",
        "    \n",
        "    legend_lines = plt.legend(\n",
        "        bbox_to_anchor=(\n",
        "            1.04,\n",
        "            1),\n",
        "        loc=\"upper left\",\n",
        "        fancybox=True,\n",
        "        shadow=True)\n",
        "    ax.add_artist(legend_lines)\n",
        "    \n",
        "    # Define regions\n",
        "    if highlighting_regions:\n",
        "        if prediction_name in ordered_regions_dict.keys():\n",
        "            for i, prediction in enumerate(\n",
        "                    ordered_regions_dict[prediction_name]):\n",
        "                lower = thresholds_dict[prediction_name][prediction][0]\n",
        "                upper = thresholds_dict[prediction_name][prediction][1]\n",
        "                color = colors[i]\n",
        "                ax.axhspan(\n",
        "                    lower,\n",
        "                    upper,\n",
        "                    alpha=0.3,\n",
        "                    color=color,\n",
        "                    label=prediction)\n",
        "            # to sort it \"from up to low\"\n",
        "            included_in_regions_legend = list(reversed(\n",
        "                [r_pred for r_pred in ordered_regions_dict[prediction_name]]))\n",
        "            # Get handles and labels\n",
        "            handles, labels = plt.gca().get_legend_handles_labels()\n",
        "            handles_dict = {label: handles[idx]\n",
        "                            for idx, label in enumerate(labels)}\n",
        "            # Add legend for regions, if available\n",
        "            lgnd_labels = [handles_dict[r] for r in included_in_regions_legend]\n",
        "            lgnd_regions = [region for region in included_in_regions_legend]\n",
        "            region_legend = ax.legend(lgnd_labels,\n",
        "                                      lgnd_regions,\n",
        "                                      fancybox=True,\n",
        "                                      shadow=True,\n",
        "                                      loc='lower left',\n",
        "                                      bbox_to_anchor=(1.04, 0))\n",
        "            ax.add_artist(region_legend)\n",
        "    \n",
        "    ax.set_ylim([min_value, max_value])\n",
        "    ax.set_xlabel('residue index')\n",
        "    ax.set_ylabel('prediction values')\n",
        "    ax.grid(axis='y')\n",
        "    plt.show()\n",
        "\n",
        "def get_value_color(result_key, value, min_value=0, max_value=1):\n",
        "    if (result_key == \"ppII\" \n",
        "        or result_key == \"coil\" \n",
        "        or result_key == \"sheet\" \n",
        "        or result_key == \"helix\" \n",
        "        or result_key == \"agmata\"):\n",
        "        return get_spectrum_color(value, min_value, max_value)\n",
        "    else:\n",
        "        for band in VALUES_CONSTANTS[result_key][\"colors\"]:\n",
        "            band_floor, band_ceiling, band_color = band\n",
        "            if value > band_floor and value <= band_ceiling:\n",
        "                return band_color\n",
        "\n",
        "def get_spectrum_color(real, min_value, max_value):\n",
        "    color_r = 0\n",
        "    color_g = 0\n",
        "    color_b = 0\n",
        "\n",
        "    value = real\n",
        "    base = max_value - min_value\n",
        "    if (value < min_value):\n",
        "        value = min_value \n",
        "    elif value > max_value:\n",
        "        value = max_value\n",
        "\n",
        "    if base == 0:\n",
        "        value = 100\n",
        "    else:\n",
        "        value = (value - min_value) / base * 100\n",
        "    \n",
        "    if value < 50:\n",
        "        color_r = 255\n",
        "        off_g = 5.1 * value\n",
        "        color_g = max(1, off_g)\n",
        "    else:\n",
        "        color_g = 255\n",
        "        off_r = 510 - 5.1 * value\n",
        "        color_r = max(1, off_r)\n",
        "  \n",
        "    return '#%02x%02x%02x' % (int(color_r), int(color_g), int(color_b))\n",
        "\n",
        "def create_py3dmol_view(molecule, color_key):\n",
        "    view = py3Dmol.view(width=400, height=300)\n",
        "    view.setBackgroundColor(\"white\")\n",
        "    view.addModelsAsFrames(str(molecule))\n",
        "\n",
        "    default_color_style = {\"cartoon\": {\"color\": \"black\"}}\n",
        "    for atom_index, atom in enumerate(molecule):\n",
        "        view.setStyle(\n",
        "            {\"model\": -1, \"serial\": atom_index + 1}, \n",
        "            atom.get(color_key, default_color_style)\n",
        "        )\n",
        "\n",
        "    view.zoomTo()\n",
        "    return view\n",
        "\n",
        "def plot_result_legend(result_key):\n",
        "    \"\"\"Plots the legend for a predicted result.\"\"\"\n",
        "\n",
        "    value_constants = VALUES_CONSTANTS[result_key]\n",
        "    colors = [x[2] for x in value_constants[\"colors\"]]\n",
        "\n",
        "    plt.figure(figsize=(2, 2))\n",
        "    for c in colors:\n",
        "        plt.bar(0, 0, color=c)\n",
        "\n",
        "    plt.legend(value_constants[\"legends\"], frameon=False, loc=\"center\", fontsize=14)\n",
        "    plt.xticks([])\n",
        "    plt.yticks([])\n",
        "    ax = plt.gca()\n",
        "    ax.spines[\"right\"].set_visible(False)\n",
        "    ax.spines[\"top\"].set_visible(False)\n",
        "    ax.spines[\"left\"].set_visible(False)\n",
        "    ax.spines[\"bottom\"].set_visible(False)\n",
        "    # plt.suptitle(value_constants[\"title\"], fontsize=18)\n",
        "    plt.title(value_constants[\"title\"], fontsize=16, pad=16)\n",
        "\n",
        "    return plt\n",
        "\n",
        "def plot_structure(pdb_file, b2b_predictions_results):\n",
        "    print(pdb_file_path)\n",
        "    with tqdm.notebook.tqdm(total=100, bar_format=TQDM_BAR_FORMAT) as pbar:\n",
        "        b2b_agmata_result = b2b_predictions_results['agmata']\n",
        "        agmata_min_value = min(b2b_agmata_result)\n",
        "        agmata_max_value = max(b2b_agmata_result)\n",
        "\n",
        "        pbar.update(2)\n",
        "\n",
        "        with open(pdb_file, \"r\") as ifile:\n",
        "            molecule = Molecule(ifile)\n",
        "    \n",
        "        for atom in molecule:\n",
        "            res_id = int(atom[\"resid\"]) - 1\n",
        "            for result_key in VALUES_CONSTANTS.keys():\n",
        "                if result_key != \"pLDDT\":\n",
        "                    atom[result_key] = b2b_predictions_results[result_key][res_id]\n",
        "\n",
        "                if result_key == \"agmata\":\n",
        "                    atom_color = get_value_color(\n",
        "                        result_key, \n",
        "                        atom[result_key], \n",
        "                        min_value=agmata_min_value, \n",
        "                        max_value=agmata_max_value\n",
        "                    )\n",
        "                else:\n",
        "                    atom_color = get_value_color(result_key, atom[result_key])\n",
        "\n",
        "                atom[\"{0}_color\".format(result_key)] = {\n",
        "                    \"cartoon\": {\n",
        "                        \"color\": atom_color\n",
        "                      }\n",
        "                }\n",
        "        \n",
        "        pbar.update(3)\n",
        "        grid = GridspecLayout(len(np.matrix.flatten(B2B_TOOLS_RESULTS)) + 1, 2)\n",
        "\n",
        "        # output_plot = Output()\n",
        "        # with output_plot:\n",
        "        #     create_py3dmol_view(molecule, \"pLDDT_color\").show()\n",
        "        # grid[0, 0] = output_plot\n",
        "\n",
        "        # output_legend = Output()\n",
        "        # with output_legend:\n",
        "        #     plot_result_legend(\"pLDDT\").show()\n",
        "        # grid[0, 1] = output_legend\n",
        "        pbar.update(5)\n",
        "\n",
        "        predictors_count = B2B_TOOLS_RESULTS.shape[0] * B2B_TOOLS_RESULTS.shape[1]\n",
        "        for index, result_key in enumerate(np.matrix.flatten(B2B_TOOLS_RESULTS)):\n",
        "            output_plot = Output()\n",
        "            with output_plot:\n",
        "                create_py3dmol_view(\n",
        "                    molecule, \n",
        "                    \"{0}_color\".format(result_key)\n",
        "                ).show()\n",
        "\n",
        "            grid[index + 1, 0] = output_plot\n",
        "\n",
        "            output_legend = Output()\n",
        "            with output_legend:\n",
        "                plot_result_legend(result_key).show()\n",
        "\n",
        "            grid[index + 1, 1] = output_legend\n",
        "            pbar.update(int(90.0 / predictors_count))\n",
        "\n",
        "        display(grid)\n",
        "\n",
        "def build_pdb_by_prediction(pdb_file_path, predictions, predictor):\n",
        "    p = PDBParser()\n",
        "    structure = p.get_structure(\"A\", pdb_file_path)\n",
        "\n",
        "    for model in structure:\n",
        "        for chain in model:\n",
        "            for residue_position, residue_object in enumerate(chain):\n",
        "                for atom in residue_object:\n",
        "                    atom.set_bfactor(round(predictions[predictor][residue_position], 3))\n",
        "\n",
        "    io = PDBIO()\n",
        "    io.set_structure(structure)\n",
        "    io.save(path.join(\"/content/\", f\"{predictor}.pdb\"))\n",
        "    return path.join(\"/content/\", f\"{predictor}.pdb\")\n",
        "\n",
        "\n",
        "def download_pdb_by_prediction(pdb_file_path, b2b_predictions_results, predictor):\n",
        "    filename = build_pdb_by_prediction(pdb_file_path, b2b_predictions_results, predictor)\n",
        "    files.download(filename)"
      ],
      "metadata": {
        "cellView": "form",
        "id": "GFyoVkkbr7z2"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Submit your files\n",
        "Open this section to see the input cells where you will be able to submit the sequence of interest in both FASTA and PDB format."
      ],
      "metadata": {
        "id": "d9gNUlssjNZ2"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "#@title Upload Sequence file in FASTA format { display-mode: \"form\" }\n",
        "#@markdown Please execute this cell by pressing the _Play_ button and upload just one FASTA file\n",
        "\n",
        "from google.colab import files\n",
        "from os import path\n",
        "from Bio import SeqIO \n",
        "import ipywidgets as widgets\n",
        "\n",
        "fasta_uploaded = files.upload()\n",
        "print(fasta_uploaded.keys())\n",
        "\n",
        "fasta_file_path = \"\"\n",
        "\n",
        "for fasta_file_name in fasta_uploaded.keys():\n",
        "  fasta_file_path = path.join('/content/', fasta_file_name)\n",
        "\n",
        "  print('User fasta_uploaded file \"{name}\" with length {length} bytes'.format(\n",
        "      name=fasta_file_path, length=len(fasta_uploaded[fasta_file_name])))\n",
        "\n",
        "sequence_keys = []\n",
        "fasta = SeqIO.parse(open(fasta_file_path), \"fasta\")\n",
        "formatted_fasta_content = \"\"\n",
        "for sequence in fasta:\n",
        "  sequence_id = sequence.id.replace(\"|\", \"_\")\n",
        "  sequence_keys.append(sequence_id)\n",
        "  formatted_fasta_content += f\">{sequence_id}\\n{sequence.seq}\"\n",
        "\n",
        "formatted_fasta_filepath = path.join(\"/content\", \"formatted_input.fasta\")\n",
        "with open(formatted_fasta_filepath, \"w\") as fasta_writter:\n",
        "    fasta_writter.write(formatted_fasta_content)\n",
        "\n",
        "sequence_picker = widgets.Dropdown(options=sequence_keys, description=\"Sequence:\")\n",
        "display(sequence_picker)"
      ],
      "metadata": {
        "id": "Kd8Qk7COhVam"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "#@markdown Run this cell to confirm the sequence selected in the previous cell\n",
        "selected_sequence = sequence_picker.value\n",
        "print(\"Selected sequence:\", selected_sequence)"
      ],
      "metadata": {
        "cellView": "form",
        "id": "qPXdGKLslbah"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "#@title Upload your structure file in PDB format { display-mode: \"form\" }\n",
        "#@markdown Please execute this cell by pressing the _Play_ button and upload just one PDB file\n",
        "\n",
        "from google.colab import files\n",
        "from os import path\n",
        "\n",
        "pdb_uploaded = files.upload()\n",
        "\n",
        "for pdb_file_name in pdb_uploaded.keys():\n",
        "  pdb_file_path = path.join(\"/content/\", pdb_file_name)\n",
        "\n",
        "  print('User pdb_uploaded file \"{name}\" with length {length} bytes'.format(\n",
        "      name=pdb_file_path, length=len(pdb_uploaded[pdb_file_name])))"
      ],
      "metadata": {
        "id": "sLVvE_lmhngw"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Processing of the input files"
      ],
      "metadata": {
        "id": "xF7QhKh8kGzS"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "#@markdown ## Step: Predict biophysical features\n",
        "\n",
        "with tqdm.notebook.tqdm(total=100, bar_format=TQDM_BAR_FORMAT) as pbar:\n",
        "    with io.capture_output() as captured:\n",
        "        try:\n",
        "            single_seq = SingleSeq(formatted_fasta_filepath)\n",
        "            pbar.update(25)\n",
        "            single_seq.predict(tools=B2B_TOOLS)\n",
        "            pbar.update(50)\n",
        "            print(single_seq.get_all_predictions().keys())\n",
        "            b2b_predictions_results = single_seq.get_all_predictions()[selected_sequence]\n",
        "            pbar.update(25)            \n",
        "        except subprocess.CalledProcessError:\n",
        "            print(captured)\n",
        "            raise"
      ],
      "metadata": {
        "cellView": "form",
        "id": "yZblWVyJkGU1"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Analysis"
      ],
      "metadata": {
        "id": "iwV-CnO-p-Yn"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "#@markdown ## Plotting the biophysical predictions\n",
        "#@markdown Please execute this cell by pressing the _Play_ button. \n",
        "\n",
        "#@markdown ### DynaMine backbone dynamics\n",
        "#@markdown Values above _0.8_ indicate **rigid conformations**, values above _1.0_ **membrane spanning regions**, values below _0.69_ **flexible regions**. \n",
        "#@markdown Values between _0.69-0.80_ are **'context' dependent** and capable of being either rigid or flexible.\n",
        "\n",
        "#@markdown ### DynaMine sidechain dynamics\n",
        "#@markdown Higher values mean more likely rigid. These values are highly dependent on the amino acid type (i.e. a Trp will be rigid, an Asp flexible).\n",
        "#@markdown ### DynaMine conformational propensities (sheet, helix, coil, ppII (polyproline II))\n",
        "#@markdown Higher values indicate higher propensities.\n",
        "#@markdown ### EFoldMine earlyFolding propensity\n",
        "#@markdown Values above _0.169_ indicate residues that are **likely to start the protein folding process**, based on only local interactions with other amino acids.\n",
        "#@markdown ### Disomine disorder\n",
        "#@markdown Values above _0.5_ indicate that this is **likely a disordered residue**.\n",
        "#@markdown ### Agmata aggregation propensity\n",
        "#@markdown These values are divided by a factor of 20 from the original. Peaks indicate residues likely to be involved in beta-sheet aggregation.\n",
        "\n",
        "residues = b2b_predictions_results['seq']\n",
        "residues_count = len(residues)\n",
        "sequence_df = pd.DataFrame(columns=b2b_predictions_results.keys(), index=range(residues_count))\n",
        "sequence_df.index.name = 'residue_index'\n",
        "for predictor in b2b_predictions_results.keys():\n",
        "    sequence_df[predictor] = b2b_predictions_results[predictor]\n",
        "sequence_df = sequence_df.rename(columns={\"seq\": \"residue\"})\n",
        "sequence_df = sequence_df.round(decimals=3)\n",
        "\n",
        "for predictor in b2b_predictions_results.keys():\n",
        "    if predictor != 'seq':\n",
        "        plot_prediction(prediction_name=predictor,\n",
        "                        highlighting_regions=True,\n",
        "                        sequence_df=sequence_df)"
      ],
      "metadata": {
        "cellView": "form",
        "id": "36dBXmx_qAOc"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "#@markdown ## Render three-dimensional structures { display-mode: \"form\" }\n",
        "#@markdown Please execute this cell by pressing the _Play_ button\n",
        "\n",
        "#@markdown ### DynaMine backbone dynamics\n",
        "#@markdown Values above _0.8_ indicate **rigid conformations**, values above _1.0_ **membrane spanning regions**, values below _0.69_ **flexible regions**. \n",
        "#@markdown Values between _0.69-0.80_ are **'context' dependent** and capable of being either rigid or flexible.\n",
        "\n",
        "#@markdown ### DynaMine sidechain dynamics\n",
        "#@markdown Higher values mean more likely rigid. These values are highly dependent on the amino acid type (i.e. a Trp will be rigid, an Asp flexible).\n",
        "#@markdown ### DynaMine conformational propensities (sheet, helix, coil, ppII (polyproline II))\n",
        "#@markdown Higher values indicate higher propensities.\n",
        "#@markdown ### EFoldMine earlyFolding propensity\n",
        "#@markdown Values above _0.169_ indicate residues that are **likely to start the protein folding process**, based on only local interactions with other amino acids.\n",
        "#@markdown ### Disomine disorder\n",
        "#@markdown Values above _0.5_ indicate that this is **likely a disordered residue**.\n",
        "#@markdown ### Agmata aggregation propensity\n",
        "#@markdown These values are divided by a factor of 20 from the original. Peaks indicate residues likely to be involved in beta-sheet aggregation.\n",
        "\n",
        "plot_structure(pdb_file_path, b2b_predictions_results)"
      ],
      "metadata": {
        "cellView": "form",
        "id": "xS33BuKqqYQM"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "#@markdown ## PDB download\n",
        "#@markdown Select a biophysical feature and press on the _Play_ button to download\n",
        "#@markdown a new PDB file containing the prediction values on the B-FACTOR column\n",
        "\n",
        "download_pdb_by_feature = \"backbone\" #@param ['backbone', 'sidechain', 'ppII', 'coil', 'sheet', 'helix', 'earlyFolding', 'agmata', 'disoMine']\n",
        "download_pdb_by_prediction(pdb_file_path, b2b_predictions_results, download_pdb_by_feature)"
      ],
      "metadata": {
        "cellView": "form",
        "id": "qMG-86w3El3J"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Appendix"
      ],
      "metadata": {
        "id": "dHXoi550qLDW"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "## Conclusion and notes\n",
        "\n",
        "**Important note:** If you need to run these analysis again, please restart the Runtime from the menu \"Runtime > Restart Runtime\"\n",
        "\n",
        "Download the Bio2Byte's tools package from our PyPi repository: https://pypi.org/project/b2bTools/.\n",
        "\n",
        "### Citations\n",
        "\n",
        "DynaMine:\n",
        "> Elisa Cilia, Rita Pancsa, Peter Tompa, Tom Lenaerts, and Wim Vranken From protein sequence to dynamics and disorder with DynaMine Nature Communications 4:2741(2013)\n",
        "doi: 10.1038/ncomms3741\n",
        "\n",
        "EFoldMine:\n",
        "> Raimondi, D., Orlando, G., Pancsa, R. et al. Exploring the Sequence-based Prediction of Folding Initiation Sites in Proteins. Sci Rep 7, 8826 (2017).\n",
        "doi: https://doi.org/10.1038/s41598-017-08366-3\n",
        "\n",
        "Disomine:\n",
        "> Gabriele Orlando, Daniele Raimondi, Francesco Codice, Francesco Tabaro, Wim Vranken bioRxiv 2020.05.25.115253;\n",
        "doi: https://doi.org/10.1101/2020.05.25.115253\n",
        "\n",
        "AgMata:\n",
        "> O.G, S.A, MR.S, R.D, V.W.; Accurate prediction of protein beta-aggregation with generalized statistical potentials. doi: 10.1093/bioinformatics/btz912.\n",
        "\n",
        "Regarding AlphaFold:\n",
        "```\n",
        "@Article{AlphaFold2021,\n",
        "  author  = {Jumper, John and Evans, Richard and Pritzel, Alexander and Green, Tim and Figurnov, Michael and Ronneberger, Olaf and Tunyasuvunakool, Kathryn and Bates, Russ and {\\v{Z}}{\\'\\i}dek, Augustin and Potapenko, Anna and Bridgland, Alex and Meyer, Clemens and Kohl, Simon A A and Ballard, Andrew J and Cowie, Andrew and Romera-Paredes, Bernardino and Nikolov, Stanislav and Jain, Rishub and Adler, Jonas and Back, Trevor and Petersen, Stig and Reiman, David and Clancy, Ellen and Zielinski, Michal and Steinegger, Martin and Pacholska, Michalina and Berghammer, Tamas and Bodenstein, Sebastian and Silver, David and Vinyals, Oriol and Senior, Andrew W and Kavukcuoglu, Koray and Kohli, Pushmeet and Hassabis, Demis},\n",
        "  journal = {Nature},\n",
        "  title   = {Highly accurate protein structure prediction with {AlphaFold}},\n",
        "  year    = {2021},\n",
        "  doi     = {10.1038/s41586-021-03819-2},\n",
        "  note    = {(Accelerated article preview)},\n",
        "}\n",
        "```"
      ],
      "metadata": {
        "id": "IBW8bp-PqNJr"
      }
    }
  ]
}