CysPresso.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "source": [
        "###Load the CDP dataset and embeddings"
      ],
      "metadata": {
        "id": "ZnTV_ToUwSXu"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "md7AAPdaZfBG"
      },
      "outputs": [],
      "source": [
        "import numpy as np\n",
        "import pandas as pd"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "AuK5HfESI6z2"
      },
      "outputs": [],
      "source": [
        "# Load the cysteine-dense peptides dataset (CDPs.csv) from GitHub \n",
        "df = pd.read_csv('https://raw.githubusercontent.com/Zebreu/cyspresso/main/CDPs.csv')"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "14C0ue1NtGdv"
      },
      "outputs": [],
      "source": [
        "# Process the dataset (omit extraneous columns and label \"Expressibility\" as True or False) \n",
        "express = {'+': True, '+-PR': True, '-': False}\n",
        "df['Expressibility'] = df['Expressibility'].replace(express)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 424
        },
        "id": "hBD3OJT_3Jrb",
        "outputId": "6c91492c-74b7-49d5-b765-097c7042091b"
      },
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "     Uniprot                                           Sequence  \\\n",
              "0     P01030          AKRCCQDGLTRLPMARTCEQRAARVQQPACREPFLSCCQFA   \n",
              "1     P46162               PQSCRWNMGVCIPFLCRVGMRQIGTCFGPRVPCCRR   \n",
              "2     P46163              PQSCRWNMGVCIPISCPGNMRQIGTCFGPRVPCCRRW   \n",
              "3     P46167                FVTCRINRGFCVPIRCPGHRRQIGTCLAPQIKCCR   \n",
              "4     P01223        GLACGQAMSFCIPTEYMMHVERKECAYCLTINTTVCAGYCMTR   \n",
              "...      ...                                                ...   \n",
              "1244  B6UHE2  ADLCVTRSRTFKGWCHQSENCITVCKSEGNTGGFCKLGACMCTKECVRS   \n",
              "1245  P0C1Y5  GGGCGYKDVNKAPFNSMGACGNVPIFKDGLGCGSCFEIKCDKPAECSGK   \n",
              "1246  B6SJ49  ARTCQSQSHRFRGPCLRRSNCANVCRTEGFPGGRCRGFRRRCFCTTHCH   \n",
              "1247  B6SQK6  AQICYSRSKTFKGWCYHSTNCISVCITEGEISGFCQHGICMCTYECLTG   \n",
              "1248  A1DMY0                    RDHCGQVCLNKTGCGGKCPKCDMRSLTCKKA   \n",
              "\n",
              "      Expressibility Is Knottin? Uniprot  \n",
              "0              False                   N  \n",
              "1              False                   N  \n",
              "2              False                   N  \n",
              "3              False                   N  \n",
              "4              False                   N  \n",
              "...              ...                 ...  \n",
              "1244            True                   N  \n",
              "1245           False                   N  \n",
              "1246           False                   N  \n",
              "1247           False                   N  \n",
              "1248            True                   Y  \n",
              "\n",
              "[1249 rows x 4 columns]"
            ],
            "text/html": [
              "\n",
              "  <div id=\"df-f056fa76-4f0a-43a8-8ea8-7819f576f5e5\">\n",
              "    <div class=\"colab-df-container\">\n",
              "      <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>Uniprot</th>\n",
              "      <th>Sequence</th>\n",
              "      <th>Expressibility</th>\n",
              "      <th>Is Knottin? Uniprot</th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>0</th>\n",
              "      <td>P01030</td>\n",
              "      <td>AKRCCQDGLTRLPMARTCEQRAARVQQPACREPFLSCCQFA</td>\n",
              "      <td>False</td>\n",
              "      <td>N</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>1</th>\n",
              "      <td>P46162</td>\n",
              "      <td>PQSCRWNMGVCIPFLCRVGMRQIGTCFGPRVPCCRR</td>\n",
              "      <td>False</td>\n",
              "      <td>N</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2</th>\n",
              "      <td>P46163</td>\n",
              "      <td>PQSCRWNMGVCIPISCPGNMRQIGTCFGPRVPCCRRW</td>\n",
              "      <td>False</td>\n",
              "      <td>N</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>3</th>\n",
              "      <td>P46167</td>\n",
              "      <td>FVTCRINRGFCVPIRCPGHRRQIGTCLAPQIKCCR</td>\n",
              "      <td>False</td>\n",
              "      <td>N</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>4</th>\n",
              "      <td>P01223</td>\n",
              "      <td>GLACGQAMSFCIPTEYMMHVERKECAYCLTINTTVCAGYCMTR</td>\n",
              "      <td>False</td>\n",
              "      <td>N</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>...</th>\n",
              "      <td>...</td>\n",
              "      <td>...</td>\n",
              "      <td>...</td>\n",
              "      <td>...</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>1244</th>\n",
              "      <td>B6UHE2</td>\n",
              "      <td>ADLCVTRSRTFKGWCHQSENCITVCKSEGNTGGFCKLGACMCTKECVRS</td>\n",
              "      <td>True</td>\n",
              "      <td>N</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>1245</th>\n",
              "      <td>P0C1Y5</td>\n",
              "      <td>GGGCGYKDVNKAPFNSMGACGNVPIFKDGLGCGSCFEIKCDKPAECSGK</td>\n",
              "      <td>False</td>\n",
              "      <td>N</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>1246</th>\n",
              "      <td>B6SJ49</td>\n",
              "      <td>ARTCQSQSHRFRGPCLRRSNCANVCRTEGFPGGRCRGFRRRCFCTTHCH</td>\n",
              "      <td>False</td>\n",
              "      <td>N</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>1247</th>\n",
              "      <td>B6SQK6</td>\n",
              "      <td>AQICYSRSKTFKGWCYHSTNCISVCITEGEISGFCQHGICMCTYECLTG</td>\n",
              "      <td>False</td>\n",
              "      <td>N</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>1248</th>\n",
              "      <td>A1DMY0</td>\n",
              "      <td>RDHCGQVCLNKTGCGGKCPKCDMRSLTCKKA</td>\n",
              "      <td>True</td>\n",
              "      <td>Y</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "<p>1249 rows × 4 columns</p>\n",
              "</div>\n",
              "      <button class=\"colab-df-convert\" onclick=\"convertToInteractive('df-f056fa76-4f0a-43a8-8ea8-7819f576f5e5')\"\n",
              "              title=\"Convert this dataframe to an interactive table.\"\n",
              "              style=\"display:none;\">\n",
              "        \n",
              "  <svg xmlns=\"http://www.w3.org/2000/svg\" height=\"24px\"viewBox=\"0 0 24 24\"\n",
              "       width=\"24px\">\n",
              "    <path d=\"M0 0h24v24H0V0z\" fill=\"none\"/>\n",
              "    <path d=\"M18.56 5.44l.94 2.06.94-2.06 2.06-.94-2.06-.94-.94-2.06-.94 2.06-2.06.94zm-11 1L8.5 8.5l.94-2.06 2.06-.94-2.06-.94L8.5 2.5l-.94 2.06-2.06.94zm10 10l.94 2.06.94-2.06 2.06-.94-2.06-.94-.94-2.06-.94 2.06-2.06.94z\"/><path d=\"M17.41 7.96l-1.37-1.37c-.4-.4-.92-.59-1.43-.59-.52 0-1.04.2-1.43.59L10.3 9.45l-7.72 7.72c-.78.78-.78 2.05 0 2.83L4 21.41c.39.39.9.59 1.41.59.51 0 1.02-.2 1.41-.59l7.78-7.78 2.81-2.81c.8-.78.8-2.07 0-2.86zM5.41 20L4 18.59l7.72-7.72 1.47 1.35L5.41 20z\"/>\n",
              "  </svg>\n",
              "      </button>\n",
              "      \n",
              "  <style>\n",
              "    .colab-df-container {\n",
              "      display:flex;\n",
              "      flex-wrap:wrap;\n",
              "      gap: 12px;\n",
              "    }\n",
              "\n",
              "    .colab-df-convert {\n",
              "      background-color: #E8F0FE;\n",
              "      border: none;\n",
              "      border-radius: 50%;\n",
              "      cursor: pointer;\n",
              "      display: none;\n",
              "      fill: #1967D2;\n",
              "      height: 32px;\n",
              "      padding: 0 0 0 0;\n",
              "      width: 32px;\n",
              "    }\n",
              "\n",
              "    .colab-df-convert:hover {\n",
              "      background-color: #E2EBFA;\n",
              "      box-shadow: 0px 1px 2px rgba(60, 64, 67, 0.3), 0px 1px 3px 1px rgba(60, 64, 67, 0.15);\n",
              "      fill: #174EA6;\n",
              "    }\n",
              "\n",
              "    [theme=dark] .colab-df-convert {\n",
              "      background-color: #3B4455;\n",
              "      fill: #D2E3FC;\n",
              "    }\n",
              "\n",
              "    [theme=dark] .colab-df-convert:hover {\n",
              "      background-color: #434B5C;\n",
              "      box-shadow: 0px 1px 3px 1px rgba(0, 0, 0, 0.15);\n",
              "      filter: drop-shadow(0px 1px 2px rgba(0, 0, 0, 0.3));\n",
              "      fill: #FFFFFF;\n",
              "    }\n",
              "  </style>\n",
              "\n",
              "      <script>\n",
              "        const buttonEl =\n",
              "          document.querySelector('#df-f056fa76-4f0a-43a8-8ea8-7819f576f5e5 button.colab-df-convert');\n",
              "        buttonEl.style.display =\n",
              "          google.colab.kernel.accessAllowed ? 'block' : 'none';\n",
              "\n",
              "        async function convertToInteractive(key) {\n",
              "          const element = document.querySelector('#df-f056fa76-4f0a-43a8-8ea8-7819f576f5e5');\n",
              "          const dataTable =\n",
              "            await google.colab.kernel.invokeFunction('convertToInteractive',\n",
              "                                                     [key], {});\n",
              "          if (!dataTable) return;\n",
              "\n",
              "          const docLinkHtml = 'Like what you see? Visit the ' +\n",
              "            '<a target=\"_blank\" href=https://colab.research.google.com/notebooks/data_table.ipynb>data table notebook</a>'\n",
              "            + ' to learn more about interactive tables.';\n",
              "          element.innerHTML = '';\n",
              "          dataTable['output_type'] = 'display_data';\n",
              "          await google.colab.output.renderOutput(dataTable, element);\n",
              "          const docLink = document.createElement('div');\n",
              "          docLink.innerHTML = docLinkHtml;\n",
              "          element.appendChild(docLink);\n",
              "        }\n",
              "      </script>\n",
              "    </div>\n",
              "  </div>\n",
              "  "
            ]
          },
          "metadata": {},
          "execution_count": 5
        }
      ],
      "source": [
        "small = df[['Uniprot','Sequence','Expressibility','Is Knottin? Uniprot']]\n",
        "small"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "Download AlphaFold2 embeddings from Hugging Face (https://huggingface.co/datasets/TonyKYLim/CysPresso/tree/main). \n",
        "\n",
        "Alternatively, you can use your own embeddings in a '.npy' format, one representation per peptide: (number of residues, number of dimensions) e.g. (38, 256) for Alphafold2's MSA embedding"
      ],
      "metadata": {
        "id": "9sqDUpnka2ZV"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "import glob"
      ],
      "metadata": {
        "id": "FNbhSBTyejnT"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "XqINYMUnOspb"
      },
      "outputs": [],
      "source": [
        "cembs1 = sorted(glob.glob('msa_first_row_reps/*.npy')) #Path to folder AlphaFold2 containing MSA embeddings"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "TTFxpLct5NIR"
      },
      "outputs": [],
      "source": [
        "cembs2 = sorted(glob.glob('struct_mod_reps/*.npy')) #Path to folder contianing AlphaFold2 structure embeddings"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "VYMPXYISPK5i"
      },
      "outputs": [],
      "source": [
        "cembs3 = sorted(glob.glob('pair_reps/*.npy'))#Path to folder containing AlphaFold2 pair embeddings"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "7wZrZhOI5NIR"
      },
      "outputs": [],
      "source": [
        "cembs4 = sorted(glob.glob('single_embeddings/*.npy')) #Path to folder containing AlphaFold2 single embeddings"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "RHwZgK0C5NIS"
      },
      "outputs": [],
      "source": [
        "cembs = [cembs1, cembs2, cembs3, cembs4]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "i0P7o2TR5NIT"
      },
      "outputs": [],
      "source": [
        "# Concatenate the four AlphaFold2 embeddings (MSA, structure, pair, single) for each peptide to generate the combined AlphaFold2 representation (combined embeddings, 'cembs')\n",
        "cembs_array = []\n",
        "names = []\n",
        "for c,c1,c2,c3 in zip(cembs1,cembs2,cembs3, cembs4):\n",
        "    array = np.load(c)\n",
        "    array1 = np.load(c1)\n",
        "    array2 = np.load(c2).mean(axis=1)\n",
        "    array3 = np.load(c3)\n",
        "    array = np.concatenate([array, array1, array2, array3], axis=1)\n",
        "    names.append(c[-10:-4])\n",
        "    cembs_array.append(array)\n",
        "raw_array = cembs_array.copy()                     "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "FpQUk_8f5NIT"
      },
      "outputs": [],
      "source": [
        "raw_array[0].shape #Output: (number of residues, number of dimensions)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "AWuGUIfLPM1N"
      },
      "source": [
        "### Embedding preprocessing"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "FrAT-I2r5NIW"
      },
      "outputs": [],
      "source": [
        "# Rocket works best when zero-padded\n",
        "\n",
        "features = [[] for i in range(384+256+128+384)]\n",
        "new_names = []\n",
        "for name, a in zip(names, raw_array):\n",
        "  for i in range(384+256+128+384):\n",
        "    f = a[:,i]\n",
        "    d = 50 - len(f) #Pad values with 0s to ensure equal length (50, in this case)\"\n",
        "    f = np.pad(f, (0,d)) \n",
        "    features[i].append(pd.Series(f))\n",
        "  new_names.append(name)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "cVMM9be9PM1N",
        "outputId": "a8edf170-b4df-4ba8-f439-1bf20dd1b469"
      },
      "outputs": [
        {
          "data": {
            "text/plain": [
              "1249"
            ]
          },
          "execution_count": 45,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "len(new_names)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "n4fBDAuw5NIW",
        "outputId": "95d9f343-d53d-4d02-c7a7-2553755127c4"
      },
      "outputs": [
        {
          "data": {
            "text/plain": [
              "1152"
            ]
          },
          "execution_count": 16,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "len(features)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "tHVcYW2xsfk1"
      },
      "outputs": [],
      "source": [
        "cc = pd.DataFrame(new_names, columns=['Uniprot'])"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "YmVv8girslWz",
        "tags": []
      },
      "outputs": [],
      "source": [
        "for i,f in enumerate(features):\n",
        "  cc[i] = pd.Series(f)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "1eMtBIjlPhbr"
      },
      "outputs": [],
      "source": [
        "combined = pd.merge(cc, small, left_on='Uniprot', right_on='Uniprot')"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "Llmoxv8sKfBb",
        "outputId": "508c8f30-776c-4e94-ebcc-12d3f7329028"
      },
      "outputs": [
        {
          "data": {
            "text/plain": [
              "True     678\n",
              "False    549\n",
              "Name: Expressibility, dtype: int64"
            ]
          },
          "execution_count": 47,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# Omit duplicates from dataset\n",
        "combined = combined.drop_duplicates('Uniprot').sort_values('Uniprot')\n",
        "combined['Expressibility'].value_counts()"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "###Import required modules"
      ],
      "metadata": {
        "id": "6vl1pV6awIbh"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "BR9shVurNdzz"
      },
      "outputs": [],
      "source": [
        "from sklearn import preprocessing\n",
        "from sklearn.utils import shuffle"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "r-6gqeGFN-JS"
      },
      "outputs": [],
      "source": [
        "from sklearn.ensemble import RandomForestClassifier\n",
        "from sklearn.linear_model import LogisticRegression, RidgeClassifier, LogisticRegressionCV, RidgeClassifierCV\n",
        "from sklearn.svm import SVC"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "xoJ70lJ7_VdU"
      },
      "outputs": [],
      "source": [
        "#import lightgbm\n",
        "from sklearn.metrics import auc, roc_auc_score\n",
        "from sklearn import model_selection\n",
        "from sklearn import metrics"
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "!pip install sktime"
      ],
      "metadata": {
        "id": "lMEcBJhKarCO"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "VmtboR2gs05Z"
      },
      "outputs": [],
      "source": [
        "from sktime.transformations.panel.rocket import MiniRocketMultivariate, Rocket, MultiRocketMultivariate"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "TuTnCTMWPM1O"
      },
      "source": [
        "### Split the dataset between knottins and non-knottins"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "nmfHQGy3PM1O"
      },
      "outputs": [],
      "source": [
        "# Select all rows that correspond to knottin proteins based on Uniprot identification \n",
        "knot = combined[combined['Is Knottin? Uniprot'] == 'Y']\n",
        "knot = knot.sort_values('Uniprot') "
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "###ROCKET transformation and model training"
      ],
      "metadata": {
        "id": "fEutIrvYoeuN"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "WuKYSzMWPM1O"
      },
      "outputs": [],
      "source": [
        "mega = {} # Dictionary holding all results for later analysis"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "ZRxT9OGz5NIa"
      },
      "outputs": [],
      "source": [
        "# ROCKET transformation on the knottin dataset \n",
        "\n",
        "feat_cols = list(range(384+256+128+384))\n",
        "\n",
        "rocket = Rocket(num_kernels=10000, random_state=42)\n",
        "rocket.fit(knot[feat_cols])\n",
        "transformed = rocket.transform(knot[feat_cols])"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "tags": [],
        "id": "m59huiPSPM1O"
      },
      "outputs": [],
      "source": [
        "# Train a logistic regression model\n",
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    scaler = preprocessing.StandardScaler()\n",
        "    scaler.fit(X_train)\n",
        "    tr = scaler.transform(X_train)\n",
        "\n",
        "    model = LogisticRegression(C=0.0001)\n",
        "    model.fit(tr, y_train)\n",
        "\n",
        "    tt = scaler.transform(X_test)\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['knottin_rocket']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "YmxRH8MGPM1O"
      },
      "source": [
        "### Alternative methods\n",
        "#### The following cells use a similar methodology as above, but with different embeddings and/or preprocessing, and with the non-knottin partition of the dataset\n"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "#### Proteinfer embeddings"
      ],
      "metadata": {
        "id": "q0jOCmQLpBL7"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "fdEQcEChPM1O"
      },
      "outputs": [],
      "source": [
        "# Load Proteinfer embeddings\n",
        "st = pd.read_parquet('peptide-data.parquet')\n",
        "stt = pd.read_parquet('proteInfer_embeddings_peptides.parquet') # Computable by proteinfer as directed by https://github.com/google-research/proteinfer\n",
        "\n",
        "stt['Uniprot'] = st[\"('Uniprot',)\"].values\n",
        "\n",
        "knot = pd.merge(combined[['Uniprot', 'Expressibility', 'Is Knottin? Uniprot']], stt, left_on='Uniprot', right_on='Uniprot')\n",
        "knot = knot[knot['Is Knottin? Uniprot'] == 'Y']\n",
        "knot = knot.drop_duplicates('Uniprot').sort_values('Uniprot').sort_values('Uniprot')\n",
        "transformed = knot[[str(i) for i in range(1100)]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "dJZWnUWJPM1O"
      },
      "outputs": [],
      "source": [
        "# Train Random Forest classifier \n",
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    tr = X_train\n",
        "    model = RandomForestClassifier(n_estimators=300)\n",
        "    model.fit(tr, y_train)\n",
        "    \n",
        "    tt = X_test\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['knottin_proteinfer']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Te_mVbXLPM1O"
      },
      "source": [
        "#### Averaged representations from Alphafold2"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "##### All representations are averaged"
      ],
      "metadata": {
        "id": "uOvcdU7SqsW4"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "r8zgbayePM1P"
      },
      "outputs": [],
      "source": [
        "parts = []\n",
        "for a in raw_array:\n",
        "    parts.append(a.mean(axis=0))\n",
        "averaged = pd.DataFrame(parts)\n",
        "averaged['Uniprot'] = names"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "nMk85tNAPM1P"
      },
      "outputs": [],
      "source": [
        "knot = pd.merge(combined[['Uniprot', 'Expressibility', 'Is Knottin? Uniprot']], averaged, left_on='Uniprot', right_on='Uniprot')\n",
        "knot = knot[knot['Is Knottin? Uniprot'] == 'Y']\n",
        "knot = knot.drop_duplicates('Uniprot').sort_values('Uniprot').sort_values('Uniprot')\n",
        "transformed = knot[[i for i in range(1152)]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "kfEZ8dqnPM1P"
      },
      "outputs": [],
      "source": [
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    tr = X_train\n",
        "    model = RandomForestClassifier(n_estimators=300)\n",
        "    model.fit(tr, y_train)\n",
        "    \n",
        "    tt = X_test\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['knottin_averagedfour']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "##### Averaging MSA representations only"
      ],
      "metadata": {
        "id": "NRCdP9TAsGim"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "C-AybiQaPM1P"
      },
      "outputs": [],
      "source": [
        "parts = []\n",
        "for a in raw_array:\n",
        "    parts.append(a[:,0:256].mean(axis=0))\n",
        "averaged = pd.DataFrame(parts)\n",
        "averaged['Uniprot'] = names"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "rEJW_jsKPM1P"
      },
      "outputs": [],
      "source": [
        "knot = pd.merge(combined[['Uniprot', 'Expressibility', 'Is Knottin? Uniprot']], averaged, left_on='Uniprot', right_on='Uniprot')\n",
        "knot = knot[knot['Is Knottin? Uniprot'] == 'Y']\n",
        "knot = knot.drop_duplicates('Uniprot').sort_values('Uniprot').sort_values('Uniprot')\n",
        "transformed = knot[[i for i in range(256)]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "sP313IdgPM1P"
      },
      "outputs": [],
      "source": [
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    tr = X_train\n",
        "    model = RandomForestClassifier(n_estimators=300)\n",
        "    model.fit(tr, y_train)\n",
        "    \n",
        "    tt = X_test\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['knottin_averagedmsa']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "##### Averaging structure module representations only"
      ],
      "metadata": {
        "id": "3Oh6BZMasKs6"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "5tlTguZsPM1P"
      },
      "outputs": [],
      "source": [
        "parts = []\n",
        "for a in raw_array:\n",
        "    parts.append(a[:,256:256+384].mean(axis=0))\n",
        "averaged = pd.DataFrame(parts)\n",
        "averaged['Uniprot'] = names"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "VwtnOX9uPM1P"
      },
      "outputs": [],
      "source": [
        "knot = pd.merge(combined[['Uniprot', 'Expressibility', 'Is Knottin? Uniprot']], averaged, left_on='Uniprot', right_on='Uniprot')\n",
        "knot = knot[knot['Is Knottin? Uniprot'] == 'Y']\n",
        "knot = knot.drop_duplicates('Uniprot').sort_values('Uniprot').sort_values('Uniprot')\n",
        "transformed = knot[[i for i in range(384)]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "XbzD8GOLPM1P"
      },
      "outputs": [],
      "source": [
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    tr = X_train\n",
        "    model = RandomForestClassifier(n_estimators=300)\n",
        "    model.fit(tr, y_train)\n",
        "    \n",
        "    tt = X_test\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['knottin_averagedstruct']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "##### Averaging pair representations only"
      ],
      "metadata": {
        "id": "FaEpTBHTsN9t"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "D4srLcsIPM1P"
      },
      "outputs": [],
      "source": [
        "parts = []\n",
        "for a in raw_array:\n",
        "    parts.append(a[:,256+384:256+384+128].mean(axis=0))\n",
        "averaged = pd.DataFrame(parts)\n",
        "averaged['Uniprot'] = names"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "IoF-W-1APM1P"
      },
      "outputs": [],
      "source": [
        "knot = pd.merge(combined[['Uniprot', 'Expressibility', 'Is Knottin? Uniprot']], averaged, left_on='Uniprot', right_on='Uniprot')\n",
        "knot = knot[knot['Is Knottin? Uniprot'] == 'Y']\n",
        "knot = knot.drop_duplicates('Uniprot').sort_values('Uniprot').sort_values('Uniprot')\n",
        "transformed = knot[[i for i in range(128)]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "CxatMIsIPM1P"
      },
      "outputs": [],
      "source": [
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    tr = X_train\n",
        "    model = RandomForestClassifier(n_estimators=300)\n",
        "    model.fit(tr, y_train)\n",
        "    \n",
        "    tt = X_test\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['knottin_averagedpair']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "##### Averaging single representations only"
      ],
      "metadata": {
        "id": "6oazmCHnsSLi"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "s7dkwpOQPM1P"
      },
      "outputs": [],
      "source": [
        "parts = []\n",
        "for a in raw_array:\n",
        "    parts.append(a[:,256+384+128:256+384+128+384].mean(axis=0))\n",
        "averaged = pd.DataFrame(parts)\n",
        "averaged['Uniprot'] = names"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "zYm9-gLvPM1P"
      },
      "outputs": [],
      "source": [
        "knot = pd.merge(combined[['Uniprot', 'Expressibility', 'Is Knottin? Uniprot']], averaged, left_on='Uniprot', right_on='Uniprot')\n",
        "knot = knot[knot['Is Knottin? Uniprot'] == 'Y']\n",
        "knot = knot.drop_duplicates('Uniprot').sort_values('Uniprot').sort_values('Uniprot')\n",
        "transformed = knot[[i for i in range(384)]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "Fa_U00nePM1P"
      },
      "outputs": [],
      "source": [
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    tr = X_train\n",
        "    model = RandomForestClassifier(n_estimators=300)\n",
        "    model.fit(tr, y_train)\n",
        "    \n",
        "    tt = X_test\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['knottin_averagedsingle']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "VTBSZVkhPM1P"
      },
      "source": [
        "#### Non-knottins experiments"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "AyvpGvnRPM1P"
      },
      "outputs": [],
      "source": [
        "knot = combined[combined['Is Knottin? Uniprot'] == 'N'] # Selecting non-knottins from dataset\n",
        "knot = knot.sort_values('Uniprot')"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "7mwk1q6JPM1Q"
      },
      "outputs": [],
      "source": [
        "# Rocket transformation on the non-knottin dataset\n",
        "feat_cols = list(range(384+256+128+384))\n",
        "\n",
        "rocket = Rocket(num_kernels=10000, random_state=42)\n",
        "rocket.fit(knot[feat_cols])\n",
        "transformed = rocket.transform(knot[feat_cols])"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "tags": [],
        "id": "cY6_8-mcPM1Q"
      },
      "outputs": [],
      "source": [
        "# Train a logistic regression model\n",
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    scaler = preprocessing.StandardScaler()\n",
        "    scaler.fit(X_train)\n",
        "    tr = scaler.transform(X_train)\n",
        "\n",
        "    model = LogisticRegression(C=0.0001)\n",
        "    model.fit(tr, y_train)\n",
        "\n",
        "    tt = scaler.transform(X_test)\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['nonknottin_rocket']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "##### Proteinfer embeddings for non-knottins"
      ],
      "metadata": {
        "id": "g2_y9qOysb92"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "mIPbzBq4PM1Q"
      },
      "outputs": [],
      "source": [
        "# Load Proteinfer embeddings\n",
        "st = pd.read_parquet('peptide-data.parquet')\n",
        "\n",
        "stt = pd.read_parquet('proteInfer_embeddings_peptides.parquet')\n",
        "\n",
        "stt['Uniprot'] = st[\"('Uniprot',)\"].values\n",
        "\n",
        "knot = pd.merge(combined[['Uniprot', 'Expressibility', 'Is Knottin? Uniprot']], stt, left_on='Uniprot', right_on='Uniprot')\n",
        "knot = knot[knot['Is Knottin? Uniprot'] == 'N'] # Select non-knottins from dataset\n",
        "knot = knot.drop_duplicates('Uniprot').sort_values('Uniprot').sort_values('Uniprot')\n",
        "transformed = knot[[str(i) for i in range(1100)]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "fZKiClaIPM1Q"
      },
      "outputs": [],
      "source": [
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    tr = X_train\n",
        "    model = RandomForestClassifier(n_estimators=300)\n",
        "    model.fit(tr, y_train)\n",
        "    \n",
        "    tt = X_test\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['nonknottin_proteinfer']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "##### Averaged representations from AlphaFold2"
      ],
      "metadata": {
        "id": "E6AqCksasl0i"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "###### All representations are averaged"
      ],
      "metadata": {
        "id": "RvCVAGqIv2KX"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "qK1m3eJmPM1Q"
      },
      "outputs": [],
      "source": [
        "parts = []\n",
        "for a in raw_array:\n",
        "    parts.append(a.mean(axis=0))\n",
        "averaged = pd.DataFrame(parts)\n",
        "averaged['Uniprot'] = names"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "mNybhw42PM1Q"
      },
      "outputs": [],
      "source": [
        "knot = pd.merge(combined[['Uniprot', 'Expressibility', 'Is Knottin? Uniprot']], averaged, left_on='Uniprot', right_on='Uniprot')\n",
        "knot = knot[knot['Is Knottin? Uniprot'] == 'N']\n",
        "knot = knot.drop_duplicates('Uniprot').sort_values('Uniprot').sort_values('Uniprot')\n",
        "transformed = knot[[i for i in range(1152)]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "UkZw6XFgPM1Q"
      },
      "outputs": [],
      "source": [
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    tr = X_train\n",
        "    model = RandomForestClassifier(n_estimators=300)\n",
        "    model.fit(tr, y_train)\n",
        "    \n",
        "    tt = X_test\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['nonknottin_averagedfour']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "###### Averaging MSA representations only"
      ],
      "metadata": {
        "id": "Jgtl2nxUv6jE"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "DREz2YIXPM1Q"
      },
      "outputs": [],
      "source": [
        "parts = []\n",
        "for a in raw_array:\n",
        "    parts.append(a[:,0:256].mean(axis=0))\n",
        "averaged = pd.DataFrame(parts)\n",
        "averaged['Uniprot'] = names"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "CP4lMca4PM1Q"
      },
      "outputs": [],
      "source": [
        "knot = pd.merge(combined[['Uniprot', 'Expressibility', 'Is Knottin? Uniprot']], averaged, left_on='Uniprot', right_on='Uniprot')\n",
        "knot = knot[knot['Is Knottin? Uniprot'] == 'N']\n",
        "knot = knot.drop_duplicates('Uniprot').sort_values('Uniprot').sort_values('Uniprot')\n",
        "transformed = knot[[i for i in range(256)]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "gWeCQLxPPM1Q"
      },
      "outputs": [],
      "source": [
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    tr = X_train\n",
        "    model = RandomForestClassifier(n_estimators=300)\n",
        "    model.fit(tr, y_train)\n",
        "    \n",
        "    tt = X_test\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['nonknottin_averagedmsa']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "###### Averaging structure module representations only"
      ],
      "metadata": {
        "id": "AyCl2_fcv-B1"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "XiqzeBHjPM1Q"
      },
      "outputs": [],
      "source": [
        "parts = []\n",
        "for a in raw_array:\n",
        "    parts.append(a[:,256:256+384].mean(axis=0))\n",
        "averaged = pd.DataFrame(parts)\n",
        "averaged['Uniprot'] = names"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "RBFaOos0PM1Q"
      },
      "outputs": [],
      "source": [
        "knot = pd.merge(combined[['Uniprot', 'Expressibility', 'Is Knottin? Uniprot']], averaged, left_on='Uniprot', right_on='Uniprot')\n",
        "knot = knot[knot['Is Knottin? Uniprot'] == 'N']\n",
        "knot = knot.drop_duplicates('Uniprot').sort_values('Uniprot').sort_values('Uniprot')\n",
        "transformed = knot[[i for i in range(384)]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "j9sUiu_UPM1Q"
      },
      "outputs": [],
      "source": [
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    tr = X_train\n",
        "    model = RandomForestClassifier(n_estimators=300)\n",
        "    model.fit(tr, y_train)\n",
        "    \n",
        "    tt = X_test\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['nonknottin_averagedstruct']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "###### Averaging pair representations only"
      ],
      "metadata": {
        "id": "p82f6HyrwDhO"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "NY9zMbAZPM1Q"
      },
      "outputs": [],
      "source": [
        "parts = []\n",
        "for a in raw_array:\n",
        "    parts.append(a[:,256+384:256+384+128].mean(axis=0))\n",
        "averaged = pd.DataFrame(parts)\n",
        "averaged['Uniprot'] = names"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "qpRzxneDPM1Q"
      },
      "outputs": [],
      "source": [
        "knot = pd.merge(combined[['Uniprot', 'Expressibility', 'Is Knottin? Uniprot']], averaged, left_on='Uniprot', right_on='Uniprot')\n",
        "knot = knot[knot['Is Knottin? Uniprot'] == 'N']\n",
        "knot = knot.drop_duplicates('Uniprot').sort_values('Uniprot').sort_values('Uniprot')\n",
        "transformed = knot[[i for i in range(128)]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "Th_7hjmoPM1Q"
      },
      "outputs": [],
      "source": [
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    tr = X_train\n",
        "    model = RandomForestClassifier(n_estimators=300)\n",
        "    model.fit(tr, y_train)\n",
        "    \n",
        "    tt = X_test\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['nonknottin_averagedpair']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "###### Averaging single representations only"
      ],
      "metadata": {
        "id": "KgcUU9SuwJgR"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "PfmWPmKpPM1Q"
      },
      "outputs": [],
      "source": [
        "parts = []\n",
        "for a in raw_array:\n",
        "    parts.append(a[:,256+384+128:256+384+128+384].mean(axis=0))\n",
        "averaged = pd.DataFrame(parts)\n",
        "averaged['Uniprot'] = names"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "Tp7Ab1ysPM1Q"
      },
      "outputs": [],
      "source": [
        "knot = pd.merge(combined[['Uniprot', 'Expressibility', 'Is Knottin? Uniprot']], averaged, left_on='Uniprot', right_on='Uniprot')\n",
        "knot = knot[knot['Is Knottin? Uniprot'] == 'N']\n",
        "knot = knot.drop_duplicates('Uniprot').sort_values('Uniprot').sort_values('Uniprot')\n",
        "transformed = knot[[i for i in range(384)]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "10WobRY_PM1Q"
      },
      "outputs": [],
      "source": [
        "scores=[]\n",
        "skf = model_selection.ShuffleSplit(n_splits=50, test_size=0.1, random_state=10)\n",
        "\n",
        "for train_index, test_index in skf.split(transformed.values, knot['Expressibility'].values):\n",
        "    X_train, X_test = transformed.values[train_index], transformed.values[test_index]\n",
        "    y_train, y_test = knot['Expressibility'].values[train_index], knot['Expressibility'].values[test_index]\n",
        "\n",
        "    tr = X_train\n",
        "    model = RandomForestClassifier(n_estimators=300)\n",
        "    model.fit(tr, y_train)\n",
        "    \n",
        "    tt = X_test\n",
        "    truth = y_test\n",
        "    \n",
        "    scores.append(roc_auc_score(truth, model.predict_proba(tt)[:,1]))  \n",
        "mega['nonknottin_averagedsingle']=list(scores)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "851EeJdXPM1R"
      },
      "source": [
        "### Prepare the data to draw a critical difference diagram according to https://github.com/hfawaz/cd-diagram "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "pUjrgO5DPM1R"
      },
      "outputs": [],
      "source": [
        "classifiers = []\n",
        "datasets = []\n",
        "accuracies = []\n",
        "for key in mega:\n",
        "    if key.startswith('knottin'):\n",
        "        classifier = key.split('_')[1].replace('averaged','')\n",
        "        for e, value in enumerate(mega[key]):\n",
        "            accuracy = value\n",
        "            dataset = str(e)\n",
        "\n",
        "            classifiers.append(classifier.replace('four','combined'))\n",
        "            datasets.append(dataset)\n",
        "            accuracies.append(accuracy)\n",
        "    "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "YKQn4e0lPM1R"
      },
      "outputs": [],
      "source": [
        "cdddata = pd.DataFrame()\n",
        "cdddata['classifier_name'] = classifiers\n",
        "cdddata['dataset_name'] = datasets\n",
        "cdddata['accuracy'] = accuracies"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "ngUU3AXwPM1R",
        "outputId": "ec3f7841-4649-4a96-c1e2-1ddad67b511c"
      },
      "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>classifier_name</th>\n",
              "      <th>dataset_name</th>\n",
              "      <th>accuracy</th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>0</th>\n",
              "      <td>rocket</td>\n",
              "      <td>0</td>\n",
              "      <td>0.862500</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>1</th>\n",
              "      <td>rocket</td>\n",
              "      <td>1</td>\n",
              "      <td>0.745553</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2</th>\n",
              "      <td>rocket</td>\n",
              "      <td>2</td>\n",
              "      <td>0.828526</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>3</th>\n",
              "      <td>rocket</td>\n",
              "      <td>3</td>\n",
              "      <td>0.786963</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>4</th>\n",
              "      <td>rocket</td>\n",
              "      <td>4</td>\n",
              "      <td>0.801675</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>...</th>\n",
              "      <td>...</td>\n",
              "      <td>...</td>\n",
              "      <td>...</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>345</th>\n",
              "      <td>single</td>\n",
              "      <td>45</td>\n",
              "      <td>0.714052</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>346</th>\n",
              "      <td>single</td>\n",
              "      <td>46</td>\n",
              "      <td>0.754427</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>347</th>\n",
              "      <td>single</td>\n",
              "      <td>47</td>\n",
              "      <td>0.775054</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>348</th>\n",
              "      <td>single</td>\n",
              "      <td>48</td>\n",
              "      <td>0.787854</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>349</th>\n",
              "      <td>single</td>\n",
              "      <td>49</td>\n",
              "      <td>0.720031</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "<p>350 rows × 3 columns</p>\n",
              "</div>"
            ],
            "text/plain": [
              "    classifier_name dataset_name  accuracy\n",
              "0            rocket            0  0.862500\n",
              "1            rocket            1  0.745553\n",
              "2            rocket            2  0.828526\n",
              "3            rocket            3  0.786963\n",
              "4            rocket            4  0.801675\n",
              "..              ...          ...       ...\n",
              "345          single           45  0.714052\n",
              "346          single           46  0.754427\n",
              "347          single           47  0.775054\n",
              "348          single           48  0.787854\n",
              "349          single           49  0.720031\n",
              "\n",
              "[350 rows x 3 columns]"
            ]
          },
          "execution_count": 12,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "cdddata"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "IXUbxZwHPM1R"
      },
      "outputs": [],
      "source": [
        "cdddata.to_csv('peptidecdd.csv', index=False)"
      ]
    }
  ],
  "metadata": {
    "colab": {
      "provenance": [],
      "collapsed_sections": [
        "ZnTV_ToUwSXu",
        "AWuGUIfLPM1N",
        "6vl1pV6awIbh"
      ]
    },
    "kernelspec": {
      "display_name": "Python 3 (ipykernel)",
      "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.8.5"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}