RidgeRegressionMultipleAlphaModels4CDK2_Ki_MVD.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "provenance": []
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "source": [
        "# **Ridge Regression Models for the CDK2 K<sub>i</sub> Dataset Generated with Molegro Virtual Docker**\n",
        "\n",
        "This Jupyter Notebook builds regression models using the Ridge method (also known as Tikhonov regularization) ([Tikhonov, 1963](https://scholar.google.com/scholar?q=A.+N.+Tikhonov,+Solution+of+incorrectly+formulated+problems+and+the+regularization+method), [Géron, 2023](https://www.isbns.net/isbn/9781098125974/)) for data from the CDK2 K<sub>i</sub> dataset ([de Azevedo et al., 2024](https://doi.org/10.1002/jcc.27449)). This dataset employs the ligands described for case study 1 on the SAnDReS 2.0 paper but with energy terms and poses determined using Molegro Virtual Docker (MVD) ([Thomsen & Christensen, 2006](https://doi.org/10.1021/jm051197e)). MVD determined poses using differential evolution [Storn & Price, 1997](https://doi.org/10.1023/A:1008202821328) as a search algorithm and MolDock Score to estimate binding affinity. This code employs the [Ridge](https://scikit-learn.org/stable/modules/generated/sklearn.linear_model.Ridge.html) class from the [Scikit-Learn](https://scikit-learn.org/stable/) library ([Pedregosa et al., 2011](https://doi.org/10.48550/arXiv.1201.0490)) to generate the regression models. We have a loop to generate regression models with different values for alpha. It evaluates the model's predictive performance using metrics recommended by [Walsh et al., 2021](https://doi.org/10.1038/s41592-021-01205-4) and [de Azevedo et al., 2024](https://doi.org/10.1002/jcc.27449).\n",
        "This code reads two CSV files to generate the Ridge models. It is necessary to upload them into the Colab folder. We may create these CSV files using [Ridge_CDK2_Ki_MVD.ipynb](https://colab.research.google.com/drive/1Cfqn-mGO7mX6rJjPp7ytALGFIrar7CYH?usp=sharing).\n",
        "<br> </br>\n",
        "<img src=\"https://drive.usercontent.google.com/download?id=1qQSdzcOcoMwszytCxIRSe5oq7KlxHM6s&export=view&authuser=0\" width=600 alt=\"PDB: 2A4L\">\n",
        "<br><i>Structure of a protein-ligand complex ([de Azevedo et al., 1997](https://doi.org/10.1111/j.1432-1033.1997.0518a.x)) with an inhibitor bound to the macromolecule (PDB access code: [2A4L](https://www.rcsb.org/structure/2A4L)).</i></br>\n",
        "<br></br>\n",
        "**References**\n",
        "<br></br>\n",
        "de Azevedo WF Jr, Quiroga R, Villarreal MA, da Silveira NJF, Bitencourt-Ferreira G, da Silva AD, Veit-Acosta M, Oliveira PR, Tutone M, Biziukova N, Poroikov V, Tarasova O, Baud S. SAnDReS 2.0: Development of machine-learning models to explore the scoring function space. J Comput Chem. 2024; 45(27): 2333-2346. PMID: 38900052. [DOI: 10.1002/jcc.27449](https://doi.org/10.1002/jcc.27449) [PubMed](https://pubmed.ncbi.nlm.nih.gov/38900052/)\n",
        "<br></br>\n",
        "Géron, A. Hands-on Machine Learning with Scikit-Learn, Keras, and TensorFlow 3e: Concepts, Tools, and Techniques to Build Intelligent Systems, 3rd ed.; O’Reilly Media: Sebastopol, CA, 2023. [ISBN: 978-1-098-12597-4](https://www.isbns.net/isbn/9781098125974/)\n",
        "<br></br>\n",
        "Pedregosa F, Varoquaux G, Gramfort A, Michel V, Thirion B, Grisel O, Blondel M, Prettenhofer P, Weiss R, Dubourg V, Verplas J, Passos A, Cournapeau D, Brucher M, Perrot M, Duchesnay E. Scikitlearn: Machine Learning in Python. J Mach Learn Res., 2011; 12:2825–2830. [DOI: 10.48550/arXiv.1201.0490](https://doi.org/10.48550/arXiv.1201.0490)\n",
        "<br></br>\n",
        "Storn R, Price K. Differential Evolution – A Simple and Efficient Heuristic for Global Optimization over Continuous Spaces. J Glob Optim. 1997; 11: 341–359. [DOI: 10.1023/A:100820282132](https://doi.org/10.1023/A:100820282132)\n",
        "<br></br>\n",
        "Thomsen R, Christensen MH. MolDock: a new technique for high-accuracy molecular docking. J Med Chem. 2006; 49(11): 3315-21. [DOI: 10.1021/jm051197e](https://doi.org/10.1021/jm051197e) [PubMed](https://pubmed.ncbi.nlm.nih.gov/16722650/)\n",
        "<br></br>\n",
        "Tikhonov, A. N. Solution of Incorrectly Formulated Problems and the Regularization Method. Doklady Akademii Nauk SSSR 1963, 151, 1035–1038. [Google Scholar](https://scholar.google.com/scholar?q=A.+N.+Tikhonov,+Solution+of+incorrectly+formulated+problems+and+the+regularization+method)\n",
        "<br></br>\n",
        "Walsh I, Fishman D, Garcia-Gasulla D, Titma T, Pollastri G; ELIXIR Machine Learning Focus Group; Harrow J, Psomopoulos FE, Tosatto SCE. DOME: recommendations for supervised machine learning validation in biology. Nat Methods. 2021; 18(10): 1122-1127. [DOI: 10.1038/s41592-021-01205-4](https://doi.org/10.1038/s41592-021-01205-4) [PubMed](https://pubmed.ncbi.nlm.nih.gov/34316068/)\n",
        "<br></br>\n",
        "\n",
        "It follows the complete Python code."
      ],
      "metadata": {
        "id": "wx6I8oUGP8kz"
      }
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 963
        },
        "id": "9UjXmYGTPykl",
        "outputId": "db177a23-a062-446a-ebc8-ba7cefc21643"
      },
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "\n",
            "Downloading the file bind_Ki.csv and splitting it into training and test sets...done!\n",
            "\n",
            "File bind_Ki_training_set.csv has 68 instances.\n",
            "File bind_Ki_test_set.csv has 29 instances.\n",
            "\n",
            "Downloading the file cdk2_Ki_mvd.csv to add SAnDReS-generated data...done!\n",
            "\n",
            "Number of instances written to file cdk2_Ki_mvd_binding_training_set.csv : 68\n",
            "\n",
            "Number of instances written to file cdk2_Ki_mvd_binding_test_set.csv : 29\n",
            "\n",
            "Reading and preprocessing file: cdk2_Ki_mvd_binding_training_set.csv...done!\n",
            "\n",
            "Reading and preprocessing file: cdk2_Ki_mvd_binding_test_set.csv...done!\n",
            "\n",
            "Generating regression models with Ridge...done!\n",
            "\n",
            "Minimum RMSE: 0.9309 for alpha = 0.0100\n",
            "\n",
            "Saving scatter plot to file: test_set_scatter_plot.pdf..."
          ]
        },
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "<Figure size 640x480 with 1 Axes>"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkQAAAHSCAYAAAANN9SZAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABWF0lEQVR4nO3deVxU5f4H8M+wowIKIoICoqIYKu4KYyy5b9005VqWW5mlmeZNrVvmvpe5lJpeQ1ssy65mphYaarKJIrik4q5YooKCyuIAz+8PfzOXcQaYDWaG83m/XrxunPU759CdT895zvPIhBACRERERBJmY+4CiIiIiMyNgYiIiIgkj4GIiIiIJI+BiIiIiCSPgYiIiIgkj4GIiIiIJI+BiIiIiCSPgYiIiIgkj4GIiIiIJI+BiIiIiCSPgYiIiIgkj4GIqIZZunQpgoKCUFpaau5SLNK6devg5+eHoqIivfbjda35DP3boJqBgYis0sWLFzF+/Hg0bdoUTk5OcHV1hVwux8qVK1FQUAAA2LRpE2QymerHyckJPj4+6NOnD1atWoX79+9rHFe5z9GjR9WW5+bmokuXLnBycsLevXsBACkpKXjzzTcRHByM2rVrw8/PD9HR0cjIyNDpM+hzLl3l5eVhyZIlmDFjBmxsbNQ+f0U/Bw4c0Os8lUlISMDs2bNx7949nfc5efIkhg4dCn9/fzg5OaFRo0bo1asXVq9ebdIaRo8ejUePHuHzzz/X+VhPXlelJ//G7Ozs0KhRI4wePRo3btxQO0bZbQ8fPqxxDiEEfH19IZPJMHDgQI31ulyfJ+t58icpKUnnz1xWUVERZsyYAR8fHzg7O6Nr166IjY3Vef9jx46hb9++cHV1hYuLC3r37o20tDSN7U6fPo1hw4ahadOmqFWrFurXr4/w8HD8/PPPRtWl6/kN+dugmsPO3AUQ6euXX37BsGHD4OjoiJEjR6J169Z49OgRDh8+jGnTpuH06dNYv369avu5c+ciICAACoUCN2/exIEDBzBlyhQsX74cO3fuRNu2bSs8X15eHnr37o0TJ05g+/bt6Nu3LwBgyZIliI+Px7Bhw9C2bVvcvHkTn376KTp06ICkpCS0bt1a789W3rl09cUXX6C4uBgvvPACAOCrr75SW//ll18iNjZWY3mrVq30rrUiCQkJmDNnDkaPHo26devqtH1UVBT8/Pwwbtw4NGzYENevX0dSUhJWrlyJSZMmmawGJycnjBo1CsuXL8ekSZMgk8kqPdaT1/VJyr+xwsJCJCUlYdOmTTh8+DBOnToFJycntW2dnJywZcsWdO/eXW35wYMHkZmZCUdHR62fRZ/ro6znSc2bN6/0s2ozevRobNu2DVOmTEFgYCA2bdqE/v37Iy4uTuNzPCk1NRXdu3eHr68vZs2ahdLSUqxZswYRERE4cuQIWrZsqdr26tWruH//PkaNGgUfHx/k5+fjxx9/xLPPPovPP/8cr732mt516XN+Q/42qAYRRFbk0qVLok6dOiIoKEj89ddfGuvPnz8vVqxYIYQQIiYmRgAQKSkpGtvt379fODs7C39/f5Gfn69a/uQ+eXl5olu3bsLBwUHs2rVL7Rjx8fGiqKhIbVlGRoZwdHQUI0aMqPSz6HMuXbVt21a89NJL5a6fOHGiqI5/7ZctWyYAiMuXL+u0ff/+/YWnp6e4e/euxrqsrCyT13D06FEBQOzfv1+nY5V3Xcv7G5sxY4YAILZu3aqx7ZAhQ0T9+vWFQqFQ22fcuHGiY8eOwt/fXwwYMEBtna7Xp6K/eUMlJycLAGLZsmWqZQUFBaJZs2YiNDS00v379+8v6tWrJ+7cuaNa9tdff4k6deqIIUOGVLp/cXGxCAkJES1btjSoLn3Pr+/fBtUcfGRGVmXp0qV48OABNm7cCG9vb431zZs3x+TJkys9zjPPPIOZM2fi6tWr+Prrr7Vu8+DBA/Tt2xepqan48ccfMWDAALX1YWFhcHBwUFsWGBiI4OBgnDlzRo9PVfm5dHH58mWcOHECPXv21HtfpRs3bmDs2LHw8vKCo6MjgoOD8cUXX2hsd//+fUyZMgVNmjSBo6MjGjRogF69eiE1NRWzZ8/GtGnTAAABAQGqxzVXrlwp97wXL15EcHCw1takBg0a6F1nZTV07NgR7u7u+Omnnyq9JoZc16efflr1uZ70wgsvIDs7W+3RzqNHj7Bt2za8+OKLWo+n7/XR1dmzZ3Ht2rUKt9m2bRtsbW3VWmecnJzwyiuvIDExEdevX69w/z/++AM9e/aEh4eHapm3tzciIiKwa9cuPHjwoML9bW1t4evrq/HoU9e69D2/Pn8bVLMwEJFV+fnnn9G0aVOEhYUZfayXX34ZAPDbb79prHv48CH69euHlJQU/PDDD1r7dGgjhEBWVhbq16+vcx2GnutJCQkJAIAOHToYtH9WVha6deuGffv24c0338TKlSvRvHlzvPLKK1ixYoXatq+//jrWrl2L559/HmvWrME777wDZ2dnnDlzBkOGDFE9Wvrkk0/w1Vdf4auvvoKnp2e55/b398exY8dw6tQpk9SpSw0dOnRAfHx8pecz5Loqg1e9evU01jVp0gShoaH49ttvVcv27NmD3NxcDB8+XOvx9Lk+wON+aHfu3FH7yc7O1tiuVatWGDlyZIXHOn78OFq0aAFXV1e15V26dAEArX1xyioqKoKzs7PG8lq1auHRo0daP9PDhw9x584dXLx4EZ988gn27NmDHj16GFSXIefX9W+DahhzN1ER6So3N1cAEP/4xz902l6Xxwdubm6iffv2Gvv4+/sLe3t7sWPHDr1q/OqrrwQAsXHjRp3rM/RcT/rggw8EAHH//v1yt6nokdkrr7wivL291R4tCCHE8OHDhZubm9qjRTc3NzFx4sRyz6PvI7PffvtN2NraCltbWxEaGiqmT58ufv31V/Ho0SOD66yshtdee004OztXWltF11V5D/ft2ydu374trl+/LrZt2yY8PT2Fo6OjuH79usa2KSkp4tNPPxUuLi6qWocNGyaioqKEEELrIzNdr4/yHNp+HB0dNeoHICIiIir8/MHBweKZZ57RWH769GkBQKxbt67C/du0aSNatGghiouLVcuKioqEn5+fACC2bdumsc/48eNVddvY2IihQ4eKnJwcg+oy5Py6/m1QzcIWIrIaeXl5AAAXFxeTHbNOnTpa3zbLysqCk5MTfH19dT7W2bNnMXHiRISGhmLUqFE676fLuYqKiuDl5aW6BtpkZ2fDzs4OderU0fncSkII/Pjjjxg0aBCEEGotC3369EFubi5SU1NV29etWxfJycn466+/9D6XNr169UJiYiKeffZZpKenY+nSpejTpw8aNWqEnTt3GlxnRerVq4eCggLk5+dXuJ0u17Vnz57w9PSEr68vhg4ditq1a2Pnzp1o3Lix1u2jo6NRUFCAXbt24f79+9i1a1e5j8sA3a+P0meffYbY2Fi1nz179mhsJ4So9A3DgoICrR29lZ3FlW91lmfChAnIyMjAK6+8gj///BOnTp3CyJEj8ffff5e7/5QpUxAbG4vNmzejX79+KCkpwaNHjwyqy5Dz6/q3QTULAxFZDWXTuLYAY6gHDx5oDViff/45HBwc0LdvX5w7d67S49y8eRMDBgyAm5ubqm+DrnQ5l6OjI7KysjQeD5jK7du3ce/ePaxfvx6enp5qP2PGjAEA3Lp1S7X90qVLcerUKfj6+qJLly6YPXs2Ll26ZFQNnTt3xn//+1/cvXsXR44cwXvvvYf79+9j6NCh+PPPPw2qsyJCCAAwyZtEygCybds29O/fH3fu3NH6Za3k6emJnj17YsuWLfjvf/+LkpISDB06tMJz6HJ9lLp06YKePXuq/URFRRn02ZydnbWOy1NYWKhaX5HXX38d//73v7FlyxYEBwejTZs2uHjxIqZPnw4AWoNmUFAQevbsiZEjR6r6+ShDsL51GXJ+U/5tkPVgICKr4erqCh8fH537UVQmMzMTubm5Wl9Ffuqpp7B7924UFBSgV69eFXYczc3NRb9+/XDv3j3s3bsXPj4+etWhz7kq4uHhgeLiYoMCo3KwwZdeekmjZUH5I5fLVdtHR0fj0qVLWL16NXx8fLBs2TIEBwdrbYXQl4ODAzp37oyFCxdi7dq1UCgU+OGHHwyqsyJ3795FrVq1Kv1C1+W6KgPI888/j507d6J169Z48cUXK+ww/OKLL2LPnj1Yt24d+vXrp9PwBEDF16cqeHt7q1pTylIu0+XvfcGCBcjKysIff/yBEydOICUlRXUvW7RoUen+Q4cORUpKitoYX/rUpe/5df3boJqFgYisysCBA3Hx4kUkJiYafSzlWDx9+vTRur5Lly7YsWMHbt26hV69euH27dsa2xQWFmLQoEHIyMjArl278NRTTxlUS2XnWrlyJV555ZUKjxEUFATg8VtR+vL09ISLiwtKSko0WhaUP0++zeTt7Y0JEyZgx44duHz5Mjw8PLBgwQIApvsv606dOgH435ecPnVWVsPly5d1Gn9J3+tqa2uLRYsW4a+//sKnn35a7naDBw+GjY0NkpKSKnxcVpEnr09VaNeuHTIyMjQe1yYnJ6vW66JevXro3r072rRpAwDYt28fGjdurLq+FVE+1srNzTW4Ln3Or+vfBtUsDERkVaZPn47atWvj1VdfRVZWlsb6ixcvYuXKlZUe5/fff8e8efMQEBCAESNGlLtdjx498O233+LChQvo27ev2v/5lpSU4J///CcSExPxww8/IDQ01LAPpcO5Tpw4UekAksrzPznytS5sbW3x/PPP48cff9TaAlc2oJWUlKh9MQGPX/328fFRPcKoXbs2AOg8UnVcXJza4xCl3bt3A4Bq8Dx96qyshtTUVJ3eVjTkukZGRqJLly5YsWKF6hHOk+rUqYO1a9di9uzZGDRoUIXH0/X66EuX1+6HDh2KkpIStcFOi4qKEBMTg65du6r1fcvPz8fZs2dx586dCo+5detWpKSkYMqUKWojf2t73KlQKPDll1/C2dlZ7T849KlL1/Mr6fq3QTULR6omq9KsWTNs2bIF//znP1WvDCtHqk5ISMAPP/yA0aNHq+2zZ88enD17FsXFxcjKysLvv/+O2NhY+Pv7Y+fOnRojCT9p8ODB2LBhA8aOHYtnn30We/fuhZOTE/71r39h586dGDRoEHJycjTGM3rppZf0/nzlnevEiROVHq9p06Zo3bo19u3bh7Fjx+p97sWLFyMuLg5du3bFuHHj8NRTTyEnJwepqanYt28fcnJyADzuw9W4cWMMHToUISEhqFOnDvbt24eUlBR8/PHHAB6P5QIA77//PoYPHw57e3sMGjRIFVKeNGnSJOTn52Pw4MEICgpS3c+tW7eiSZMmqv5B+tRZUQ3Hjh1DTk4O/vGPf1R6XQy9rtOmTcOwYcOwadMmvP7661q30bXzvT7XB/jf3/yTwsLC0LRpU9XvrVq1QkRERIUdq7t27Yphw4bhvffew61bt9C8eXNs3rwZV65cwcaNG9W2PXLkCKKiojBr1izMnj0bAHDo0CHMnTsXvXv3hoeHB5KSkhATE4O+fftqjBk2fvx45OXlITw8HI0aNcLNmzfxzTff4OzZs/j444/V+vvoWpc+5weg198G1TDmer2NyBgZGRli3LhxokmTJsLBwUG4uLgIuVwuVq9eLQoLC4UQmq8gOzg4iIYNG4pevXqJlStXiry8PI3jVvSq/kcffSQAiIEDBwqFQiEiIiLKfcVZl3+1dD1XUVGRqFWrlsZr5tosX75c1KlTR+0V+bIqG6k6KytLTJw4Ufj6+gp7e3vRsGFD0aNHD7F+/XrVNkVFRWLatGkiJCREuLi4iNq1a4uQkBCxZs0atWPNmzdPNGrUSNjY2FT6Cv6ePXvE2LFjRVBQkKhTp45wcHAQzZs3F5MmTdI6UrUudVZUw4wZM4Sfn58oLS0tt6ayyruuFd3DkpIS0axZM9GsWTNRXFys8yjS2l671/X6VPTaPQARExOjdlzo8Nq9EI9HgH7nnXdEw4YNhaOjo+jcubPYu3evxnZxcXECgJg1a5Zq2YULF0Tv3r1F/fr1haOjowgKChKLFi3SGOVdCCG+/fZb0bNnT+Hl5SXs7OxEvXr1RM+ePcVPP/1kcF36nF8I/f82qOaQCaGlHZaILEZGRgYiIyN1esU9NzcXTZs2xdKlSyvtcyRVRUVFaNKkCd59912dRjUHeF2lwpC/Dao52IeIyMLp0n9Iyc3NDdOnT8eyZctUb9GQupiYGNjb25f7GEsbXldpMORvg2oOthARWbhZs2ahoKAAS5cuNXcpREQ1FgMRERERSZ7FPzI7dOgQBg0aBB8fH8hkMuzYsUO1TqFQYMaMGWjTpg1q164NHx8fjBw50mTTCRAREZE0WHwgevjwIUJCQvDZZ59prMvPz0dqaipmzpyJ1NRU/Pe//8W5c+fw7LPPmqFSIiIislZW9chMJpNh+/bteO6558rdJiUlBV26dMHVq1fh5+dXfcURERGR1apxAzPm5uZCJpNVOC9QUVGR2qSApaWlyMnJgYeHByfzIyIishJCCNy/fx8+Pj5aRx3XR40KRIWFhZgxYwZeeOGFCmcFX7RoEebMmVONlREREVFVuX79Oho3bmzUMWrMIzOFQoHnn38emZmZOHDgQIWB6MkWotzcXPj5+SEjIwPu7u5VUTrpSKFQIC4uDlFRUbC3tzd3OZLGe2FZeD8sB++F5cjJyUGLFi1w7949uLm5GXWsGtFCpFAoEB0djatXr+L333+vMAwBgKOjIxwdHTWWu7u7w8PDo6rKJB0oFArUqlULHh4e/D8aM+O9sCy8H5aD98LymKK7i9UHImUYOn/+POLi4hhoiIiISG8WH4gePHiACxcuqH6/fPky0tLS4O7uDm9vbwwdOhSpqanYtWsXSkpKcPPmTQCPW3scHBzMVTYRERFZEYsPREePHkVUVJTq96lTpwIARo0ahdmzZ2Pnzp0AgHbt2qntFxcXh8jIyOoqk4iIiKyYxQeiyMhIVNTv24r6hBMREZGFsviRqomIiIiqGgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUmexQeiQ4cOYdCgQfDx8YFMJsOOHTvU1gsh8OGHH8Lb2xvOzs7o2bMnzp8/b55iiYiIyCpZfCB6+PAhQkJC8Nlnn2ldv3TpUqxatQrr1q1DcnIyateujT59+qCwsLCaKyUiIiJrZWfuAirTr18/9OvXT+s6IQRWrFiBDz74AP/4xz8AAF9++SW8vLywY8cODB8+vDpLJSIiIitl8YGoIpcvX8bNmzfRs2dP1TI3Nzd07doViYmJ5QaioqIiFBUVqX7Py8sDACgUCigUiqotmiqkvP68D+bHe2FZeD8sB++F5TDlPbDqQHTz5k0AgJeXl9pyLy8v1TptFi1ahDlz5mgsj4uLQ61atUxbJBkkNjbW3CXQ/+O9sCy8H5aD98L88vPzTXYsqw5EhnrvvfcwdepU1e95eXnw9fVFVFQUPDw8zFgZKRQKxMbGolevXrC3tzd3OZLGe2FZeD8sB++F5cjOzjbZsaw6EDVs2BAAkJWVBW9vb9XyrKwstGvXrtz9HB0d4ejoqLHc3t6ef9wWgvfCcvBeWBbeD8vBe2F+prz+Fv+WWUUCAgLQsGFD7N+/X7UsLy8PycnJCA0NNWNlREREZE0svoXowYMHuHDhgur3y5cvIy0tDe7u7vDz88OUKVMwf/58BAYGIiAgADNnzoSPjw+ee+458xVNREREVsXiA9HRo0cRFRWl+l3Z92fUqFHYtGkTpk+fjocPH+K1117DvXv30L17d+zduxdOTk7mKpmIiIisjMUHosjISAghyl0vk8kwd+5czJ07txqrIiIioprEqvsQEREREZkCAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJnp0xO8fFxWH//v2Ij49HZmYm7ty5g1q1asHT0xNt2rRBREQEBg4ciIYNG5qqXiIiIiKT0zsQPXz4EKtWrcKGDRtw9epVCCEAAE5OTnB3d0dBQQFOnTqFEydO4JtvvoG9vT0GDRqEt99+G3K53OQfgIiIiMhYej0yW7duHZo3b473338frq6umDdvHvbv34/c3Fzk5+cjMzMT2dnZUCgUOHv2LDZv3ozo6Gj89ttvCA8Px5AhQ3D58uWq+ixEREREBtErEE2aNAm9evXCiRMnkJaWhn//+9+IioqCi4uL2nYymQwtWrTAyy+/jK+++gpZWVn4/PPPkZ6ejq+++sqkH4CIiIjIWHo9Mjt9+jRatGih90mcnZ3x6quvYsyYMbh27Zre+xMRERFVJb1aiAwJQ2XZ2toiICDAqGMQERERmRpfuyciIiLJ0+uR2dixYw06iUwmw8aNGw3al4iIiKiq6RWINm3apHW5TCZTvX6vbTkDEREREVkyvQLRk6/Ml5aWYvLkyUhKSsLkyZPx9NNPw8vLC1lZWTh06BBWrVqF0NBQfPLJJyYtmoiIiMiU9ApE/v7+ar8vXrwYycnJSE9Ph7e3t2p5y5YtER4ejjFjxqB9+/bYtm0bpk+fbpqKiYiIiEzMqE7VGzduRHR0tFoYKqtRo0aIjo7Ghg0bjDkNERERUZUyKhBlZmbCycmpwm2cnJyQmZlpzGmIiIiIqpRRgahx48bYvn07CgsLta7Pz8/H9u3b0bhxY2NOQ0RERFSljApEr776Ki5dugS5XI6ffvoJ2dnZAIDs7Gzs2LED3bt3x5UrVzBu3DiTFEtERERUFfSe7b6sadOmISMjAzExMRgyZAgAwMbGBqWlpQAAIQTGjBmDadOmGV8pERERURUxKhDZ2Nhg48aNGDlyJDZv3owTJ04gNzcXbm5uCAkJwcsvv4zIyEgTlUpERERUNYwKREoRERGIiIgwxaGIiIiIqh3nMiMiIiLJMzoQFRcX45NPPkGXLl3g6uoKO7v/NTqlpaVhwoQJyMjIMPY0RERERFXGqEdmBQUF6N27NxISElC/fn24urri4cOHqvUBAQGIiYmBu7s75s+fb3SxRERERFXBqBaihQsXIj4+HosWLcLNmzfx6quvqq13c3NDREQEfv31V6OKJCIiIqpKRgWirVu3IioqCtOnT4dMJoNMJtPYpmnTprh27ZoxpyEiIiKqUkYFomvXrqFTp04VbuPi4oLc3FxjTkNERERUpYwKRC4uLrh161aF21y8eBGenp7GnIaIiIioShkViLp164aff/4Z9+7d07r++vXr2L17N8LDw405DREREVGVMioQTZs2DXfv3kWPHj0QHx+P4uJiAI8ndd2/fz/69OmD4uJiTJ061STFEhEREVUFo167Dw8Px6efforJkyertQK5uLgAAGxtbbFmzRp07NjRuCqJiIiIqpDRAzO+8cYbSE9Px5tvvonOnTujWbNmaN++PV5//XUcP35c41V8UyspKcHMmTMREBAAZ2dnNGvWDPPmzYMQokrPS0RERDWHSeYya9WqFVauXGmKQ+ltyZIlWLt2LTZv3ozg4GAcPXoUY8aMgZubG9566y2z1ERERETWxahAdO3aNdStWxeurq7lbnP//n3cvXsXfn5+xpyqXAkJCfjHP/6BAQMGAACaNGmCb7/9FkeOHKmS8xEREVHNY1QgCggIwOzZszFz5sxyt1m1ahU+/PBDlJSUGHOqcoWFhWH9+vXIyMhAixYtkJ6ejsOHD2P58uXl7lNUVISioiLV73l5eQAAhUIBhUJRJXWSbpTXn/fB/HgvLAvvh+XgvTCvgoICHD16FAkJCTh48KDJjmtUIBJCVNpXp6r78rz77rvIy8tDUFAQbG1tUVJSggULFmDEiBHl7rNo0SLMmTNHY3lcXBxq1apVleWSjmJjY81dAv0/3gvLwvthOXgvqsfdu3dx5swZnD17FmfPnsWlS5dUb7Wbkkn6EFUkMzNT9dZZVfj+++/xzTffYMuWLQgODkZaWhqmTJkCHx8fjBo1Sus+7733ntpQAHl5efD19UVUVBQ8PDyqrFaqnEKhQGxsLHr16gV7e3tzlyNpvBeWhffDcvBeVJ2SkhL8+eefSExMREJCAhITE3H58mWN7Ro2bIjQ0FCEhIRg9uzZJjm33oFo7ty5ar8fOHBA63YlJSW4fv06vvvuO3Tr1s2g4nQxbdo0vPvuuxg+fDgAoE2bNrh69SoWLVpUbiBydHSEo6OjxnJ7e3v+cVsI3gvLwXthWXg/LAfvhfEePHiA5ORkJCQkID4+HomJiapuLEoymQxt2rSBXC5HWFgY5HI5mjRpAplMhuzsbPMForInlslkOHDgQLmhCAB8fHywZMkSQ2rTSX5+Pmxs1EcPsLW1RWlpaZWdk4iIiPQjhMD169cRHx+vCkDp6eka39d16tRBt27dVOGnW7duFb68ZSp6B6K4uDgAjz/YM888g9GjR2ttibG1tYW7uzuCgoI0AospDRo0CAsWLICfnx+Cg4Nx/PhxLF++HGPHjq2ycxIREVHFiouLkZaWpgo/8fHxuHHjhsZ2/v7+qvAjl8vRunVr2NlVeY8eDXqfMSIiQvXPs2bNQlRUlFnnKlu9ejVmzpyJCRMm4NatW/Dx8cH48ePx4Ycfmq0mIiIiqbl3756q7098fDySk5ORn5+vto2trS3at2+v9virUaNGZqpYnVERrGfPnvjxxx/RokULNGzYUGP933//jWXLliE6OrrK+hG5uLhgxYoVWLFiRZUcn4iIiNQJIXDp0iVVy09CQgJOnz6t8WZ53bp1ERoaqmr96dy5M2rXrm2mqitmVCBavnw5Tpw4Ue6YP97e3ti1axdu3LiBrVu3GnMqIiIiMpOioiKkpqaqWn8SEhKQlZWlsV3z5s3VHn+1atWqSrvNmJJRgSglJQU9evSocJvw8HCO1UBERGRF7ty5oxZ+UlJS1AY0BgAHBwd07NhR9fgrLCwMXl5eZqrYeEYFolu3blX67K9hw4a4deuWMachIiKiKiKEwNmzZ9U6P2dkZGhsV79+fVXrT1hYGDp16gQnJyczVFw1jApEdevWxbVr1yrc5urVq6hTp44xpyEiIiITKSgoQEpKiqr1JyEhATk5ORrbtWrVSq3zc2BgIGQymRkqrh5GBaJu3bph+/btuH79Onx9fTXWX7t2DTt27MAzzzxjzGmIiIjIQDdv3lTr/JyamqoxD5uTkxO6dOmi6vvTrVs3yc3cYFQgmjp1Kn7++WfI5XLMnz8fvXr1gre3N/7++2/89ttv+OCDD1BQUIB//etfpqqXiIiIylFSUoLTp0+rDX5Y3tQXyvAjl8vRrl07ODg4mKFiy2FUIAoPD8fy5cvxr3/9C2PGjAHwePRq5Wt3NjY2WLlypVnHKSIiIqqp7t+/j+TkZFUASkpK0pj6wsbGRmPqC39//xr9+MsQRg8FOXnyZERFRWHdunVISUlBbm4u6tatiy5duuD1119H69atTVEnERGRpAkhcO3aNbXOzydOnCh36gtl60/Xrl2rZeoLa2eSsbHbtm2LNWvWmOJQREREBEChUCA9PV3t8Vd5U1+Ubf1p06YNbG1tzVCxdav+yUKIiIhIw927d5GUlKRq/Tly5EiFU18oQ5ClTH1h7RiIiIiIqpkQAhcvXlR7/HX69GmN7erWrasa9NDSp76wdnoForFjx0Imk2HhwoXw8vLSeUZ5mUyGjRs3GlQgERGRtVNOfVF27J+aNvWFtdMrEG3atAkymQwzZsyAl5cXNm3apNN+DERERCQlZae+iI+Px9GjRyuc+kL5+KtBgwZmqpj0CkTKsQyUzyu1jW1AREQkJaWlpTh37pxa5+fypr4o2/m5Y8eONWrqC2unVyDy9/ev8HciIqKaLj8/H6dPn8bJkyeRlJSExMTECqe+UIagmj71hbVjp2oiIqIK/P3332qtP6mpqSguLlbb5smpL0JDQ+Hu7m6miskQeneqNgT7EBERkTUoKSnBqVOnVOEnISFBa/eQevXqISoqCk8//TTCwsI49UUNoHenam3KTtehbTkDERERWSLl1BfKAKRt6guZTIa2bduqHn116dIFp0+fxoABA2Bvb2+mysnUDOpUrVRaWorJkycjKSkJkydPxtNPPw0vLy9kZWXh0KFDWLVqFUJDQ/HJJ5+YtGgiIqoemZnA+fNAYCDQuLG5qzGOIVNfhIWFoVu3bmpTXygUCvz555/VXT5VMaM6VS9evBjJyclIT0+Ht7e3annLli0RHh6OMWPGoH379ti2bRumT59umoqJiKhabNwIvPYaUFoK2NgA69cDr7xi7qp0V3bqC+XjL059QeUxqlP1xo0bER0drRaGymrUqBGio6OxYcMGBiIiIiuSmfm/MAQ8/t/x44E+fSy3pUiXqS/s7OzQrl07Tn1BGowKRJmZmZWOoeDk5ITMzExjTkNERNXs/Pn/hSGlkhLgwgXLCETKqS/Kvv1V1VNfKB8fNmliZPFkkYwKRI0bN8b27dsxb948rcEoPz8f27dvR2NL+LeHiKgK1KQ+NmUFBj5+TFY2FNnaAs2bm6eeslNfKEPQrVu3NLZr3ry5WuuPqaa+UH98aIc33vBD//5GH5YsiFGB6NVXX8V7770HuVyODz/8EN27d4eHhweys7Pxxx9/YO7cubhy5QoWLVpkqnqJiCyGtfexqUjjxo8/z/jxj1uGbG2Bzz+vvtB3+/ZtJCYmVjr1RadOnVThp6qmvtB8fCjD2rUh+Ne/ShAQYPLTkZkYFYimTZuGjIwMxMTEYMiQIQAAGxsbVY99IQTGjBmDadOmGV8pEVm1mtaSYo19bPT1yiuPP8+FC49bhqrqc+k79YUyAFXX1BfaHh+Wltrg4sVSBqIaxKhAZGNjg40bN2LkyJHYvHkzTpw4gdzcXLi5uSEkJAQvv/wyIiMjTVQqEVmrmtiSYul9bEylcWPTf578/HwcPXpU1fpT3tQXTz31lNrM782bNzfL1BfaHh/a2JSiWTPN8ffIeplk6o6IiAhERESY4lBEVMPU1JYUS+tjY8l0mfrC2dkZXbp0UQUgS5r6QvPxocDrr6ejcePW5i6NTIhzmRFRlaqpLSnm7mNjqZRTX5QNQFeuXNHYzsfHR+3xV7t27Sx61Oeyjw/9/Ytx4sQ1AAxENYnRgai4uBirV6/Gt99+i7NnzyI/P1+V/NPS0rB+/XpMmTIFLVq0MLpYIrI++rakWFNfo+rqY2PJ7t+/j6SkJLWpL+7fv6+2jY2NDdq0aaMKQHK5HH5+flY387vy8aFCAZw4Ye5qyNSMCkQFBQXo3bs3EhISUL9+fbi6uuLhw4eq9QEBAYiJiYG7uzvmz59vdLFEZH30aUmxxr5GVdHHxlIpp74o2/pT3tQXoaGhqtafrl27qk19QWSJjApECxcuRHx8PBYvXoxp06Zhzpw5mDdvnmq9m5sbIiIi8OuvvzIQEUmYLi0pNbWvkTXTd+oL5U/r1q059QVZHaMC0datWxEVFaWalkNb82fTpk1x/PhxY05DRDVAZS0pNbWvkTW5e/cuEhMTVa0/5U190b59e1XnZ059QTWFUYHo2rVrGDx4cIXbuLi4IDc315jTEJEE8K2t6lV26gtl609FU18ow0+XLl1Qq1YtM1RMVLWMCkQuLi5ah04v6+LFi/D09DTmNEQkAXxrq2rpOvVFYGCg2szvQUFBJpn6whjW1NGerJdRgahbt274+eefce/ePdStW1dj/fXr17F79+5KW5GIqGbS94uMb22Zzu3bt5GQkKB6/GXOqS+MYY0d7ck6GT11R1RUFHr06IFVq1apXrfPz89HYmIiJk2ahOLiYkydOtUkxRKR9TD0i0xKb22ZStmpL5Q/58+f19jO09NTbeTnjh07wtHR0QwV64Yd7ak6GRWIwsPD8emnn2Ly5MkIDw9XLXdxcQEA2NraYs2aNejYsaNxVRKRVeEXWdXKz89HSkqK6tFXQkIC7t69q7Fdq1at1N7+MtfUF4ZiR3uqTkYPzPjGG28gMjIS69atQ3JyMnJycuDq6oquXbtiwoQJCA4ONkWdRGRFquuLTCp9S/766y8cOnQIW7ZswYIFC3D8+HGtU1907dpV9fgrNDQU9erVM1PFpsGO9lSdTDJ1R6tWrbBy5UpTHIqIaoDq+CKzlL4lpg5lJSUlOHnypKrvT3x8PK5evaqxXdmpL+RyOUJCQix66gtDsKM9VSejApGtrS2GDx+Ob775xlT1EJEW1tYSUtVfZJbySM4UoUyfqS8aNWqEf/7zn4iIiLDKqS8MwY72VF2MCkSurq7w9fU1VS1EpIWltIToqyq/yCyhb4khoazs1BfKn5MnT2pMfeHi4oJu3bqpTX3h7OyM3bt3o3///jWuJagy7GhP1cGoQNSlSxekp6ebqhYieoKltIQYqqq+yCyhb4kuoUyhUCAtLU3t8ddff/2lcawmTZqojf2jbeoLhUJRVR+FiGBkIJo9ezYiIyPx5ZdfYuTIkaaqiYj+nyW0hFgiS+hboi2U2djcxbVriXj//cdvf+ky9YVcLoePj0/1FU5EWhkViGJjYxEZGYkxY8Zg9erV6Ny5M7y8vDSea8tkMsycOdOoQomkyBJaQiyVufuWNGokMG/eBXzwQQKEiAcQj9LSPzFqlPp29erVU838LpfL0blzZ059QWSBjG4hUjp27BiOHTumdTsGIiLDWEJLiCWrzr4lRUVFOHbsmNrYP9Yy9QURVc6oQBQXF2eqOoioHOZuCakp9H1TTzn1hTIAVTb1hVwuR2hoqMVNfUFEujEqEEVERJiqDiKqAN+y0aRPwKnsTb3S0lKcPXtWrfNzeVNflG39sfSpL4hIdyYZmJGIqDrpMxSBtjf1XnstH66uKTh/vuKpL4KDg1XhJywszOqmviAi3RkUiBITE/H+++8jJSUFMpkMXbt2xYIFC9ClSxdT10dEVs7Ug0rqOhSB8ry3bwOlpTcAJACIB5CA0tLjiI5+cuqLWmjZsgueflqOvn1rxtQXRKQ7vQPRyZMn0aNHDxQWFqqW7d+/X/WKKecuI7Ic5h7h2pBBJSurubKhCEpKSjB37knMmxcPIZQhSHPqi4YNGyE8/HHLz507cixYEIK0NHucOAGEhADMQkTSoverD4sXL0ZhYSHef/993Lx5Ezdv3sTMmTNRUFCAJUuWVEWNRGSAjRsBf3/gmWce/+/GjdV7/vJacjIzy99Hl5qVQxH8z33Y2MRi587Z6NWrF9zc6mLu3PYQ4k0AW/A4DNkAaAdgImSyLViy5Ar++us6tm7diuefn4yFCztBCHud6ySimkfvFqI//vgD3bt3x7x581TL5syZgwMHDuDgwYMmLY6IDGMJI1zrO6ikLjULIVBcfBWvvJKA//wn/v/H/3k89cUnn5Q9mguAbgDk///TFd9/7wJPT8039Tj4JREBBgSirKwsDB8+XGN5165dkZycbJKiiMg4lvAlr++gktprVmDXrjQUFsar3gCraOoLuVyOZs3C0Ldvawjxv6kvbG2B0FDtn52DXxIRYEAgUigUqFOnjsby2rVrc64dIgthCV/y+g4qGRgIyGR3IUQilJ2fgSN44w3tU18oA1BYWJjG1BcbNqif9+23TVcnEdVMfO2eqAaylC/5igaVFELgwoULqoEP4+PjIcSfGseoV68ewsLCVK+/6zL1hfK8K1cCH38MfPQRsHx5+Z26OfglERkUiL7++mskJSWpLbtw4QIAoH///hrby2Qy/PLLL4aciogMZClf8spBJYuKipCQcEw18GFCQgJu376tsX1AQCBatpTjmWfkGDAgzKipL5YvB4R4/M+V9aPi4JdE0mZQILpw4YIqAD1p7969Gss4kBmReZjzS/7WrVtITExUBaCjR4/i0aNHats4ODigc+fOaoMfenp6muT8ltCPioish96B6PLly1VRBxFZsdLSUpw5c0bt8Ze2/2hq0KCB2uOvqpz6whL6URGR9dA7EPn7+1dFHUSkJ3MOupifn48jR46oWn8SExNx7949je2Cg4PVOj83a9ZM1WKcmQkkJFRd/ZbSj4qIrAM7VRNZIUNGgDbGjRs31Fp/0tLSUFysPvVFrVq10KVLF1UA6tatW7lTX1RX/ZbSj4qILB8DEZGVqepBF0tKSnDy5EkcOnQI27Ztw+TJk3H1qubUF40aNVJr/QkJCYG9vb3Z638SO0sTkS70CkR9+/bFvHnz0LlzZ71P9PDhQ6xevRouLi6YOHGi3vsT0WOm7iycl5eH5ORk1eOvpKQkPHjwQG0bGxsbhISEqPr+yOVy+Pr6GvTCBDs7E5El0isQ3b59G926dUN4eDhGjhyJIUOGwM3NrcJ9kpKS8PXXX+O7775DQUEBNm/ebFTBRFJnTGdhIQSuXr2qevQVHx+PkycfT31RlouLC7p27Yr69etj1KhRkMvlcHFxMXv9RERVRa9AdOzYMWzevBlz5szBK6+8gnHjxqFly5bo2LEjvLy8ULduXRQWFiInJwfnzp3D0aNHcf/+fdja2mL48OGYP38+/Pz8TP4hbty4gRkzZmDPnj3Iz89H8+bNERMTg06dOpn8XETmpk9nYYVCgbS0NLWxf7RNfREQEKD26nvr1q1RWlqK3bt3o0ePHjo9CquK+omIqovefYhGjRqFkSNHYvfu3YiJicGBAwfw9ddfa2xnY2ODtm3bYvDgwXj11Vfh7e1tkoKfdPfuXcjlckRFRWHPnj3w9PTE+fPny+3MSVQTlNdZOCcnB0lJSaoAdOTIERQUFKjta2dnhw4dOqjCj7apLwBotBpVR/1EROZiUKdqmUyGAQMGYMCAAQCAM2fOIDMzE9nZ2XB2doanpyeCg4MrfZxmCkuWLIGvry9iYmJUywICAqr8vETm1qiRQEHBBezb97/Wnz//NM3UF9WBnZ2JyJKY5C2zVq1aoVWrVqY4lN527tyJPn36YNiwYTh48CAaNWqECRMmYNy4ceXuU1RUhKKiItXveXl5AB4/XuAEtealvP68D5oKCwuRmpqKxMREJCQkICkpSevUF4GBgQgNDUVYWBhCQ0PRsmVLjakvdLm+vBeWhffDcvBeWA5T3gOZEMqZfqyTk5MTAGDq1KkYNmwYUlJSMHnyZKxbtw6jRo3Sus/s2bMxZ84cjeVbtmyxiP9yJgKAe/fu4ezZs6qfCxcuaIz9Y29vj+bNmyMoKEj1Ux0ts0REliA/Px8vvvgicnNz4erqatSxrD4QOTg4oFOnTkhISFAte+utt5CSkoLExESt+2hrIfL19cXff/8NDw+PKq+ZyqdQKBAbG4tevXqZtCOvpVNOfVG29ae8qS9CQ0NVLUDt27fXeeqLzEzgwgUZmjcXOj2qkuq9sFS8H5aD98JyZGdnw9vb2ySByOoHZvT29sZTTz2ltqxVq1b48ccfy93H0dFR65eIvb09/7gtRE2/Fw8fPsSRI0dUr79rm/pCJpOppr5Q9v9p2rSpQWP/GDMytL29PbKy7M02TQipq+n/blgT3gvzM+X1t/pAJJfLce7cObVlGRkZnHONLMqTU18cP34cJSUlatvUqlULXbt2VQWg0NBQ1K1b1+hzGzsydEyMDG+8UX3ThBARmYPVB6K3334bYWFhWLhwIaKjo3HkyBGsX78e69evN3dpJFHKqS+Ur77Hx8fj2rVrGtsZOvWFvowZGfrOHSe88YZttU2zQURkLlYfiDp37ozt27fjvffew9y5cxEQEIAVK1ZgxIgR5i6NJCIvLw9JSUmq1p+Kpr4o+/irKgYp1caYkaH//rsOSkvVH9FZ2zQbmZng4z4iqpTJAtGkSZOwbNkyODk5ITc3t1rfdBk4cCAGDhxYbecj6VJOfVF25Ofypr4IDQ1VtQB16dLFZFNf6MuYkaG9vR/AxkaohSJrmmbDmL5TRCQtJgtEubm5KCwshJOTE+rVqwc/Pz+EhISgbdu2qv8NDAw0qEMokbkoFAocP35cbe6vv//+W2O7slNfyOVyBAcHw9bW1gwVa2foyND16xdi7doSTJhgZ3XTbBjbd4qIpMVkgejLL79U/fPZs2eRnp6OEydOID09HV999RWuXbuGWrVqITg4GMnJyaY6LZFJ5eTkIDExUdX6o8vUF3K5vMqmpjElQ0eGHjNGoH9/65tmw5i+U0QkPVXSh6hFixZo0aIFhg0bplqWm5urCklElkAIgQsXLqg9/qpo6gtl64+3dydkZtaSVJ8Ua5xmw5i+U0QkPXoHovPnzyMwMFCvfaZOnYrly5cjPDwc4eHh+p6SyCQKCwtx7NgxVfhJSEjQOvVFixYt1Fp/yk59sXEjEBVlvj4p7CCsO2P6ThGR9OgdiHr16oWEhASts2Nr884772DlypVYvny53sURGSMrK0vV9ychIQHHjh3Do0eP1LZxdHRE586dVeEnNDQUnp6eWo9n7j4p7CCsP0P7ThGR9OgdiK5fv47evXvj0KFDcHd3r3Dbd999F8uXL0eLFi0MLpBIF6Wlpfjzzz/VOj9fvHhRY7sGDRqojf3ToUMHnae+MGefFF3DGFuQNFnj4z4iqn56B6Lly5fj7bffRv/+/bF//37Url1b63bvv/8+li5diubNmyMuLs7oQonKUk59oWz9qeqpLwDz9knRJYwZ2oJkjhDF4EZElkbvQDR58mRkZ2dj/vz5GDx4MH755ReN0XU//PBDLFq0CE2bNkVcXJxVvIFDlk059YXyJy0trcKpL+RyObp162aSqS+UzNknpbIwZujjPHM8huOjPyKyRAa9ZTZ37lzk5ORgzZo1ePHFF/H999+r/qt7zpw5mD9/Ppo0aYK4uDg0atTIpAVTzVdSUoLjx48jJSVF56kv5HI5QkJCYGdXtYOvm6tPSmVhzJDHeeboE2XuflhEROUx+Nvj008/xd27d/Htt99i/PjxWL9+PebNm4c5c+bA398fcXFx8PX1NWWtVEMpp74o2wJUWFioto05p754krn6pFQUxgx5nGeOPlEcG4iILJVR/zm9efNm3Lt3Dxs3blR1aG3cuDF+//13zjZPWgkhcOXKFbXOzydPnoQQQm07V1dXdOvWzSKmvrAk5YUxQx7nmaNPFMcGIiJLZVQgsrOzw48//ojevXvj8OHDaNSoEQ4cOICAgABT1UdWTjn1hbLzc0VTXyj7/RQXF2P8+PFwcnIyQ8XWS9/HeeboE8WxgYjIUukdiPr376+xrLS0FDKZDHXr1sWbb76psV4mk+GXX34xrEKyKvpOfaF8BKbseK9QKLB7926LmgfMmihbkDIzgbi4yt/iMkefKI4NRESWSO9AtHfv3nLXnT59GqdPn9ZYzgldayYhBM6fP6/2+OvMmTMa27m7uyMsLEzV96dz585wdnY2Q8XSoO9bXOboE8WxgYjI0ugdiC5fvlwVdZAVKCwsxNGjR9Wmvrhz547Gdi1btlQLQGWnvqCqxbe4iIgMo3cgYmdp6Sg79UV8fDyOHTsGhUKhto2Tk5PG1Bf169c3U8XEt7iIiAxjskFbsrOzkZ6ejtzcXLi5uSEkJAQeHh6mOjxVMeXUF2U7P2ub+sLLy0tt5vcOHTrAwcHBDBWTNnyLi4jIMEYHoitXrmDy5Mn45Zdf1F6dlslkGDhwIFasWIEmTZoYexoysbJTX8THxyMxMRG5ublq25Sd+kLZ+dmYqS+o6vEtLiIiwxgViC5evAi5XI5bt24hMDAQcrkcXl5eqkctO3fuRFJSEhISEtC0aVNT1UwGyMzMVHv8pW3qi9q1a6Nr166qFiBTT31B1YNvcRER6c+oQDRjxgzcvn0b69atw7hx49RaDoQQWL9+PSZMmIAZM2bghx9+MLpY0k1xcTFOnjyp9vhL29QXjRs3Vmv9qY6pL6h68C0uIiL9GPXtt3//fjz77LN47bXXNNbJZDKMHz8eu3fvxr59+4w5DVUiNzcXycnJqtaf5ORkPHjwQG0bGxsbtGvXTtX6ExYWZrapL4iIiCyNUYGopKQEwcHBFW7TunVrxMXFGXMaKkOfqS9CQ0NVAahr166oU6eOmaomIiKybEYFog4dOmgdiLGs06dPo1OnTsacRtIePXqE48ePqwJQQkJChVNfKANQcHAwR3smIiLSkVGBaMGCBejRowf+85//4NVXX9VYv379evz666/Yv3+/MaeRlJycHLXwc+TIEY2Z3+3t7dGhQwe1x1/KqS+IiIhIf0b3IYqKisL48ePx8ccfq71lFh8fj4yMDPTp0wf79u1T60ckk8kwc+ZMo4u3dkIIZGRkqAUgbVNfeHh4IDQ0VNUBulOnTpz6goiIyISMCkSzZ89W/fO5c+dw7tw5jW327t2rMf+ZVAORPlNflH381bJlS479Q0REVIWMCkTsLF0xXaa+cHR0ROfOnVWtP5z6goiIqPoZFYgiIiJMVYfV02fqi7Jj/3DqCyIiIvPjKHwG0nXqi9atW6t1fubUF0RERJaHgUhHmZmZaq0/FU19oQw/nPqCiIjIOjAQaVF26gtlCKps6gu5XI62bdty6gsiIiIrxG/vMhYtWoS0tDSdpr6Qy+Xw9fU1U6VERERkSgxEZXz88ceqf1ZOfaF8/MWpL4iIiGouBqIyhg0bhmeeeQZyuRxPPfUUp74gIiKSCAaiMtauXQsPDw9zl0FERETVzMbcBRARERGZGwMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJHgMRERERSR4DEREREUkeAxERERFJXo0LRIsXL4ZMJsOUKVPMXQoRERFZiRoViFJSUvD555+jbdu25i6FiIiIrEiNCUQPHjzAiBEjsGHDBtSrV8/c5RAREZEVsTN3AaYyceJEDBgwAD179sT8+fMr3LaoqAhFRUWq3/Py8gAACoUCCoWiSuukiimvP++D+fFeWBbeD8vBe2E5THkPakQg+u6775CamoqUlBSdtl+0aBHmzJmjsTwuLg61atUydXlkgNjYWHOXQP+P98Ky8H5YDt4L88vPzzfZsWRCCGGyo5nB9evX0alTJ8TGxqr6DkVGRqJdu3ZYsWKF1n20tRD5+vri77//hoeHR3WUTeVQKBSIjY1Fr169YG9vb+5yJI33wrLwflgO3gvLkZ2dDW9vb+Tm5sLV1dWoY1l9C9GxY8dw69YtdOjQQbWspKQEhw4dwqeffoqioiLY2tqq7ePo6AhHR0eNY9nb2/OP20LwXlgO3gvLwvthOXgvzM+U19/qA1GPHj1w8uRJtWVjxoxBUFAQZsyYoRGGiIiIiJ5k9YHIxcUFrVu3VltWu3ZteHh4aCwnIiIi0qbGvHZPREREZCirbyHS5sCBA+YugYiIiKwIW4iIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPIYiIiIiEjyGIiIiIhI8hiIiIiISPKsPhAtWrQInTt3houLCxo0aIDnnnsO586dM3dZREREZEWsPhAdPHgQEydORFJSEmJjY6FQKNC7d288fPjQ3KURERGRlbAzdwHG2rt3r9rvmzZtQoMGDXDs2DGEh4ebqSoiIiKyJlbfQvSk3NxcAIC7u7uZKyEiIiJrYfUtRGWVlpZiypQpkMvlaN26dbnbFRUVoaioSPW7MkTl5ORUeY1UMYVCgfz8fGRnZ8Pe3t7c5Uga74Vl4f2wHLwXlkP5vS2EMPpYNSoQTZw4EadOncLhw4cr3G7RokWYM2eOxvIWLVpUVWlERERURbKzs+Hm5mbUMWTCFLHKArz55pv46aefcOjQIQQEBFS47ZMtRPfu3YO/vz+uXbtm9AUl4+Tl5cHX1xfXr1+Hq6urucuRNN4Ly8L7YTl4LyxHbm4u/Pz8cPfuXdStW9eoY1l9C5EQApMmTcL27dtx4MCBSsMQADg6OsLR0VFjuZubG/+4LYSrqyvvhYXgvbAsvB+Wg/fCctjYGN8l2uoD0cSJE7Flyxb89NNPcHFxwc2bNwE8DjfOzs5mro6IiIisgdW/ZbZ27Vrk5uYiMjIS3t7eqp+tW7eauzQiIiKyElbfQmSKLlCOjo6YNWuW1sdoVL14LywH74Vl4f2wHLwXlsOU96LGdKomIiIiMpTVPzIjIiIiMhYDEREREUkeAxERERFJHgMRERERSZ6kA9GiRYvQuXNnuLi4oEGDBnjuuedw7tw5c5cleYsXL4ZMJsOUKVPMXYpk3bhxAy+99BI8PDzg7OyMNm3a4OjRo+YuS3JKSkowc+ZMBAQEwNnZGc2aNcO8efNM8nYtVezQoUMYNGgQfHx8IJPJsGPHDrX1Qgh8+OGH8Pb2hrOzM3r27Inz58+bp1gJqOh+KBQKzJgxA23atEHt2rXh4+ODkSNH4q+//tLrHJIORAcPHsTEiRORlJSE2NhYKBQK9O7dGw8fPjR3aZKVkpKCzz//HG3btjV3KZJ19+5dyOVy2NvbY8+ePfjzzz/x8ccfo169euYuTXKWLFmCtWvX4tNPP8WZM2ewZMkSLF26FKtXrzZ3aTXew4cPERISgs8++0zr+qVLl2LVqlVYt24dkpOTUbt2bfTp0weFhYXVXKk0VHQ/8vPzkZqaipkzZyI1NRX//e9/ce7cOTz77LP6nUSQyq1btwQAcfDgQXOXIkn3798XgYGBIjY2VkRERIjJkyebuyRJmjFjhujevbu5yyAhxIABA8TYsWPVlg0ZMkSMGDHCTBVJEwCxfft21e+lpaWiYcOGYtmyZapl9+7dE46OjuLbb781Q4XS8uT90ObIkSMCgLh69arOx5V0C9GTcnNzAQDu7u5mrkSaJk6ciAEDBqBnz57mLkXSdu7ciU6dOmHYsGFo0KAB2rdvjw0bNpi7LEkKCwvD/v37kZGRAQBIT0/H4cOH0a9fPzNXJm2XL1/GzZs31f6/ys3NDV27dkViYqIZKyOl3NxcyGQyvSZ8tfqRqk2ltLQUU6ZMgVwuR+vWrc1djuR89913SE1NRUpKirlLkbxLly5h7dq1mDp1Kv79738jJSUFb731FhwcHDBq1Chzlycp7777LvLy8hAUFARbW1uUlJRgwYIFGDFihLlLkzTlnJleXl5qy728vFTryHwKCwsxY8YMvPDCC3pNvstA9P8mTpyIU6dO4fDhw+YuRXKuX7+OyZMnIzY2Fk5OTuYuR/JKS0vRqVMnLFy4EADQvn17nDp1CuvWrWMgqmbff/89vvnmG2zZsgXBwcFIS0vDlClT4OPjw3tBpIVCoUB0dDSEEFi7dq1e+/KRGYA333wTu3btQlxcHBo3bmzuciTn2LFjuHXrFjp06AA7OzvY2dnh4MGDWLVqFezs7FBSUmLuEiXF29sbTz31lNqyVq1a4dq1a2aqSLqmTZuGd999F8OHD0ebNm3w8ssv4+2338aiRYvMXZqkNWzYEACQlZWltjwrK0u1jqqfMgxdvXoVsbGxerUOARIPREIIvPnmm9i+fTt+//13BAQEmLskSerRowdOnjyJtLQ01U+nTp0wYsQIpKWlwdbW1twlSopcLtcYfiIjIwP+/v5mqki68vPzYWOj/n/Ttra2KC0tNVNFBAABAQFo2LAh9u/fr1qWl5eH5ORkhIaGmrEy6VKGofPnz2Pfvn3w8PDQ+xiSfmQ2ceJEbNmyBT/99BNcXFxUz37d3Nzg7Oxs5uqkw8XFRaPfVu3ateHh4cH+XGbw9ttvIywsDAsXLkR0dDSOHDmC9evXY/369eYuTXIGDRqEBQsWwM/PD8HBwTh+/DiWL1+OsWPHmru0Gu/Bgwe4cOGC6vfLly8jLS0N7u7u8PPzw5QpUzB//nwEBgYiICAAM2fOhI+PD5577jnzFV2DVXQ/vL29MXToUKSmpmLXrl0oKSlRfZ+7u7vDwcFBt5MY9e6blQOg9ScmJsbcpUkeX7s3r59//lm0bt1aODo6iqCgILF+/XpzlyRJeXl5YvLkycLPz084OTmJpk2bivfff18UFRWZu7QaLy4uTuv3w6hRo4QQj1+9nzlzpvDy8hKOjo6iR48e4ty5c+Ytugar6H5cvny53O/zuLg4nc8hE4JDnhIREZG0SboPERERERHAQERERETEQERERETEQERERESSx0BEREREksdARERERJLHQERERESSx0BEREREksdARERERJLHQEREGoQQ6NixI3r37m3uUqrVlStXIJPJMHr0aHOXYjG0XZNz587Bzs4Oa9asMV9hRCbGQERkZZRfUH379tW6fuXKlbCxsYGfn5/GrPW6+vLLL5Gamoq5c+dqPXdFP02aNDHonGQaBw4cgEwmw+zZs6vsHC1btsQLL7yAOXPm4P79+1V2HqLqJOnZ7olqmg8//BDz5s1DUFAQfvvtN/j6+up9jNLSUsyePRtPP/00unXrpnWbZs2a4aWXXtK6rm7dunqf01I0atQIZ86cgZubm7lLsXjTp0/H119/jVWrVuH99983dzlERmMgIqoBhBCYNGkSPvvsM3Tq1Al79uxB/fr1DTrWnj17cOXKlQq/5Jo3b16lLRDmYm9vj6CgIHOXYRXatGmDtm3bYsOGDXjvvfdgY8MHDmTd+BdMZGZlH3EcPnwYkZGRcHFxQd26dfH888/jwoULFe6vUCjw0ksv4bPPPsMzzzyD33//3eAwBAAxMTGQyWR4/vnnDT4G8Dik9e/fHzKZDFu3btVY169fP411hl6LQ4cOYdCgQahfvz4cHR0RGBiIDz74APn5+WrblT1+QkICevfujbp160ImkwEovw/Rk/tFRUXBxcUFnp6emDBhAgoKCgAAv/zyC0JDQ1G7dm14eXlh+vTpKC4uNlnNR48eRa9eveDi4gI3NzcMHjwYV65cUW07e/ZsREVFAQDmzJmj9ihTud2jR4+wevVq9OnTB76+vnB0dESDBg0wZMgQHD9+XGut5YmOjsbVq1cRFxen135EloiBiMhCJCUloUePHnBzc8OkSZMQERGB7du3IywsDJcuXdK6T0FBAZ577jls2bIFgwcPxu7du+Hi4mJwDUIIxMXFoWXLlqhXr57BxwEAmUyGmJgYNGjQAOPHj8fVq1dV61asWIG9e/di9OjR+Oc//6mxrz7XYu3atYiMjER8fDwGDBiAt956C40bN8aCBQvQq1cvPHr0SOP4CQkJiIyMhEwmw2uvvaa1Bm2Sk5NVdY0fPx5+fn5Yu3Ytxo0bh61bt2Lo0KHw9/fH+PHjUbduXSxbtgwLFy7UOI4hNaekpCA8PBwODg4YP348OnXqhB07dqBnz54oLCwEAERGRmLUqFEAgIiICMyaNUv1o3yUmZOTgylTpqCoqAj9+/fH22+/jcjISOzevRthYWFISUnR6VoAQGhoKABg//79Ou9DZLEEEZlVXFycACAAiHXr1qmtW7dunQAgBg4cqFp2+fJlAUCEhoaK7t27CwBi7Nixori42OhaTp8+LQCIESNGaF2vPHezZs3ErFmztP7s2bNHbZ89e/YImUwmwsLCRHFxsTh+/LhwcHAQgYGB4v79+0Zdi9OnTws7OzsREhIi7ty5o7b9okWLBADx0UcfaT3+F198Ue7nGzVqVLl17dixQ7X80aNHom3btkImk4n69euLI0eOqNbl5eWJBg0aCHd3d/Ho0SOT1Pzdd9+pbf/yyy8LAOLbb7/V2H7WrFkan08IIQoLC0VmZqbG8lOnTok6deqInj176nRNhBAiNzdXABDh4eFaz0VkTRiIiMxM+QXWokULUVJSoraupKREBAYGCplMJm7duiWE+N8XlPInNDRUp/O89NJLYsKECRVu8+uvvwoAYurUqVrXP3lubT+TJ0/W2G/KlCkCgHjnnXdEUFCQsLe3F0ePHjX6Wrz11lsCgDh06JDGsUpKSoSnp6fo2LGjxvE7dOhQ4ecrLxBFRUVp7DN37lwBQIwZM0Zj3dixYwUAcenSJdUyQ2vWFjqU68rer8oCUUUGDRokHBwc1AJcRYFICCGcnJxE06ZN9T4XkaVhp2oiCyGXyzU6ptrY2EAul+P8+fNIT09Hz549Veueeuop3Lt3D4mJiZg7dy4+/PDDCo+/cuVKODk5VbhNdnY2gMrfFOvTpw/27t1b4TZlLV68GAcOHMBHH30EAFiyZAk6duxY7va6XoukpCQAwK+//qr1sY29vT3Onj2rsbxz5846115Wu3btNJZ5e3tXuu6vv/5CQEAAABhcs7br1bhxYwDAvXv3dKpfKS0tDUuXLsXhw4dx8+ZNKBQKtfV37txR1V4Zd3d33LlzR6/zE1kiBiIiC+Hl5VXh8tzcXLXlvr6++OmnnxAVFYVZs2ahpKQEc+bMKff47u7uldbg7OwMAKo+Kabi6OiIfv36IS0tDU5OTnj11Vcr3F7Xa5GTkwMAWLBggV71lHf8yri6umoss7Ozq3Rd2cBhaM0VHb+kpETn4yQkJOCZZ54BAPTu3RuBgYGoU6cOZDIZduzYgfT0dBQVFel8vIKCAtSqVUvn7YksFQMRkYXIysqqcLm2sXGaN2+OgwcPIioqCnPnzkVJSQnmz5+vsd25c+cQFBSEnJycCjtLe3p6Avjfl7apJCcnY9myZfDw8EB2djbeeOMNjTfPytL1WihDQl5enl6dyZVvlZmDoTWbyoIFC1BUVIQ//vgD3bt3V1uXlJSE9PR0nY9VWlqK3NxcBAcHm7pMomrHt8yILER8fDxKS0vVlpWWliIhIQEymQwhISFa92vatCkOHDgAf39/LFiwAO+9957GNunp6fD19a30zbHg4GDY2NgYPMK1Nvfv38eLL74IOzs7HDhwAM8//zy+//57fPHFF+Xuo+u16Nq1K4D/PYayBlVds62tLYDyW40uXrwId3d3jTCUn5+P1NRUvc51/vx5lJaWok2bNoYVS2RBGIiILERGRgY2bNigtmzDhg3IyMjAgAEDVK032gQEBODgwYMICAjA4sWLMX36dLX1J06cKDdQlVW3bl20bdsWR48e1QgkhpowYQIuXbqEjz76CK1bt8aGDRvg6+uLt956CxkZGVr30fVaTJgwAXZ2dpg0aRKuXbumcZx79+7pPbZOVavqmpWPRq9fv651vb+/P+7evYvTp0+rlpWUlOCdd97B7du39TpXcnIygMev+BNZOz4yI7IQffr0wVtvvYXdu3cjODgYp0+fxs8//4z69etj5cqVle7v7++veny2bNkylJSU4OOPPwbwuIWobdu2OtUxePBgzJo1C0lJSQgLC9O6zYULFyocqfrdd9+Fk5MTvv76a3z99dcYNGgQJk6cCACoV68evv76a0RFReHFF19EYmIi7O3tDboWrVu3xpo1a/DGG2+gZcuW6N+/P5o1a4b79+/j0qVLOHjwIEaPHo1169bp9NmrQ1XXHBQUBB8fH3z33XdwdHRE48aNIZPJMGnSJNW4Tr/99hu6d++O6OhoODk54cCBA7hx4wYiIyNx4MABnc8VGxsLOzs7DBw40KBaiSyKuV9zI5K6sq9J//HHHyIiIkLUrl1buLq6isGDB4vz58+rba98DbpPnz5aj5eZmSkCAwPVXoH38/MTW7du1ameGzduCDs7O/HGG29orNPltXsA4u7du+LSpUvC1dVVeHt7i9u3b2sc6/3331e9im/otVA6cuSIGD58uPDx8RH29vaifv36okOHDuLdd98VZ86c0Xp8bSp77V7bfjExMQKAiImJ0Vg3a9YsAUDExcVVSc3l1ZuUlCQiIiKEi4uL6p5cvnxZtX7btm2iQ4cOolatWqJ+/foiOjpaXLx4UYwaNUpj2/LO8fDhQ1GnTh3x3HPPadRFZI1kQghRbemLiDQcOHBA9aZYVcwPdu/ePdSrVw9nzpzReZ6ul19+Gb/88guuXr1arR1/q/pakOn85z//wbhx43Dw4EGEh4ebuxwio7EPEVENd+LECTg7OyMwMFDnfebPn4+CggKsXr26Cisja1VcXIyFCxfi2WefZRiiGoN9iIhquPT0dAQHB6vePtKFv78/Nm/eXO7r7yRt165dw8iRI/Hyyy+buxQik2EgIqrhJk2ahEmTJum9X3R0dBVUQzVB06ZN+UiTahz2ISIiIiLJYx8iIiIikjwGIiIiIpI8BiIiIiKSPAYiIiIikjwGIiIiIpI8BiIiIiKSPAYiIiIikjwGIiIiIpI8BiIiIiKSPAYiIiIikjwGIiIiIpK8/wOMWrppYK9nIwAAAABJRU5ErkJggg==\n"
          },
          "metadata": {}
        },
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "done!\n",
            "\n",
            "Number of generated regression models:  10\n",
            "\n",
            "Metrics saved to file: regression_metrics.csv\n"
          ]
        }
      ],
      "source": [
        "#!/usr/bin/env python3\n",
        "#\n",
        "################################################################################\n",
        "# Dr. Walter F. de Azevedo, Jr.                                                #\n",
        "# [Scopus](https://www.scopus.com/authid/detail.uri?authorId=7006435557)       #\n",
        "# [GitHub](https://github.com/azevedolab)                                      #\n",
        "# July 20, 2024                                                                #\n",
        "################################################################################\n",
        "#\n",
        "################################################################################\n",
        "# Import section                                                               #\n",
        "################################################################################\n",
        "import sys\n",
        "import csv\n",
        "import numpy as np\n",
        "from urllib.request import urlretrieve\n",
        "import pandas as pd\n",
        "from sklearn import preprocessing\n",
        "from sklearn.linear_model import Ridge\n",
        "import matplotlib.pyplot as plt\n",
        "\n",
        "################################################################################\n",
        "# Define merge() function                                                      #\n",
        "################################################################################\n",
        "def merge(san_in,mvd_in):\n",
        "    \"\"\"Function to merge a SAnDReS-generated file with a MVD-generated file.\"\"\"\n",
        "\n",
        "    # Import section\n",
        "    import csv\n",
        "    import sys\n",
        "\n",
        "    # Define mvd_out string\n",
        "    if \"training\" in san_in:\n",
        "        mvd_out = mvd_in.replace(\".csv\",\"_binding_training_set.csv\")\n",
        "    elif \"test\" in san_in:\n",
        "        mvd_out = mvd_in.replace(\".csv\",\"_binding_test_set.csv\")\n",
        "    else:\n",
        "        mvd_out = mvd_in.replace(\".csv\",\"_binding.csv\")\n",
        "\n",
        "    # Define columns of the SAnDReS CSV file\n",
        "    ref_col_san1 = 9\n",
        "    ref_col_san2 = 11\n",
        "    ref_col_san3 = 13\n",
        "    ref_col_san4 = 16\n",
        "    ref_col_san5 = 17\n",
        "    ref_col_san6 = 18\n",
        "    ref_col_san7 = 24\n",
        "\n",
        "    # Try to open san_in\n",
        "    try:\n",
        "        fo_san = open(san_in,\"r\")\n",
        "        csv_san = csv.reader(fo_san)\n",
        "    except IOError:\n",
        "        msg_out = \"\\nIOError! I can't find file \"+san_in\n",
        "        sys.exit(msg_out)\n",
        "\n",
        "    # Looping through csv_san\n",
        "    for line1 in csv_san:\n",
        "        header1 = str(line1[:ref_col_san1])+\",\"\n",
        "        header1 += str(line1[ref_col_san2:ref_col_san3])+\",\"\n",
        "        header1 += str(line1[ref_col_san4:ref_col_san5])+\",\"\n",
        "        header1 += str(line1[ref_col_san6:ref_col_san7])\n",
        "        header1 = header1.replace(\"[\",\"\").replace(\"]\",\"\").replace(\"'\",\"\").\\\n",
        "                    replace(\" \",\"\").replace(\"AverageQ\",\"Average Q\")\n",
        "\n",
        "        break\n",
        "\n",
        "    # Try to open mvd_in\n",
        "    try:\n",
        "        fo_mvd = open(mvd_in,\"r\")\n",
        "        csv_mvd = csv.reader(fo_mvd)\n",
        "    except IOError:\n",
        "        msg_out = \"\\nIOError! I can't find file \"+mvd_in\n",
        "        sys.exit(msg_out)\n",
        "\n",
        "    # Looping through csv_mvd\n",
        "    for line2 in csv_mvd:\n",
        "        i = 0\n",
        "        for aux in line2:\n",
        "            if aux.strip() == \"SMILES\":\n",
        "                i_SMILES = i\n",
        "            elif aux.strip() == \"Complex\":\n",
        "                i_Complex = i\n",
        "            elif aux.strip() == \"Filename\":\n",
        "                i_Filename = i\n",
        "            elif aux.strip() == \"Ligand\":\n",
        "                i_Ligand = i\n",
        "            elif aux.strip() == \"Path\":\n",
        "                i_Path = i\n",
        "            elif aux.strip() == \"RMSD\":\n",
        "                i_RMSD = i\n",
        "            elif aux.strip() == \"SimilarityScore\":\n",
        "                i_S = i\n",
        "\n",
        "            i += 1\n",
        "\n",
        "        # Clean header2\n",
        "        # Warning!\n",
        "        # This code keeps the labels used in MVD except for\n",
        "        # \"E-Intra (tors, ligand atoms)\".\n",
        "        # We replaced it for \"E-Intra(tors-ligand atoms)\".\n",
        "        header2 = str(line2[i_SMILES+1:i_Complex])+\",\"\n",
        "        header2 += str(line2[i_Complex+1:i_Filename])+\",\"\n",
        "        header2 += str(line2[i_Filename+1:i_Ligand])+\",\"\n",
        "        header2 += str(line2[i_Ligand+1:i_Path])+\",\"\n",
        "        header2 += str(line2[i_Path+1:i_RMSD])+\",\"\n",
        "        header2 += str(line2[i_RMSD+1:i_S])+\",\"\n",
        "        header2 += str(line2[i_S+1:])\n",
        "        header2 = header2.replace(\"[\",\"\").replace(\" \",\"\").\\\n",
        "        replace(\"]\",\"\").replace(\"'\",\"\").\\\n",
        "        replace(\"Cofactor(VdW)\",\"Cofactor (VdW)\").\\\n",
        "        replace(\"Cofactor(elec)\",\"Cofactor (elec)\").\\\n",
        "        replace(\"Cofactor(hbond)\",\"Cofactor (hbond)\").\\\n",
        "        replace(\"E-Inter(cofactor-ligand)\",\"E-Inter (cofactor - ligand)\").\\\n",
        "        replace(\"E-Inter(protein-ligand)\",\"E-Inter (protein - ligand)\").\\\n",
        "        replace(\"E-Inter(water-ligand)\",\"E-Inter (water - ligand)\").\\\n",
        "        replace(\"E-Intertotal\",\"E-Inter total\").\\\n",
        "        replace(\"E-Intra(clash)\",\"E-Intra (clash)\").\\\n",
        "        replace(\"E-Intra(elec)\",\"E-Intra (elec)\").\\\n",
        "        replace(\"E-Intra(hbond)\",\"E-Intra (hbond)\").\\\n",
        "        replace(\"E-Intra(sp2-sp2)\",\"E-Intra (sp2-sp2)\").\\\n",
        "        replace(\"E-Intra(steric)\",\"E-Intra (steric)\").\\\n",
        "        replace(\"E-Intra(tors)\",\"E-Intra (tors)\").\\\n",
        "        replace(\"E-Intra(tors,ligandatoms)\",\"E-Intra (tors-ligand atoms)\").\\\n",
        "        replace(\"E-Intra(vdw)\",\"E-Intra (vdw)\").\\\n",
        "        replace(\"E-SoftConstraintPenalty\",\"E-Soft Constraint Penalty\").\\\n",
        "        replace(\"VdW(LJ12-6)\",\"VdW (LJ12-6)\")\n",
        "\n",
        "        break\n",
        "\n",
        "    # Close fo_mvd to re-open it into the next loop\n",
        "    fo_mvd.close()\n",
        "\n",
        "    # New header\n",
        "    data_out = header1+\",\"+header2+\"\\n\"\n",
        "\n",
        "    # Looping through csv_san\n",
        "    count_instances = 0\n",
        "    for line1 in csv_san:\n",
        "        line_out1 = str(line1[:ref_col_san1])+\",\"\n",
        "        line_out1 += str(line1[ref_col_san2:ref_col_san3])+\",\"\n",
        "        line_out1 += str(line1[ref_col_san4:ref_col_san5])+\",\"\n",
        "        line_out1 += str(line1[ref_col_san6:ref_col_san7])\n",
        "        line_out1 = line_out1.replace(\"[\",\"\").replace(\"]\",\"\").\\\n",
        "                                        replace(\"'\",\"\").replace(\" \",\"\")\n",
        "\n",
        "        # Open mvd_in\n",
        "        fo_mvd = open(mvd_in,\"r\")\n",
        "        csv_mvd = csv.reader(fo_mvd)\n",
        "\n",
        "        count_instances += 1\n",
        "\n",
        "        # Looping through csv_mvd\n",
        "        for line2 in csv_mvd:\n",
        "            if line1[2].strip() in str(line2):\n",
        "\n",
        "                # Clean line_out2\n",
        "                line_out2 = str(line2[i_SMILES+1:i_Complex])+\",\"\n",
        "                line_out2 += str(line2[i_Complex+1:i_Filename])+\",\"\n",
        "                line_out2 += str(line2[i_Filename+1:i_Ligand])+\",\"\n",
        "                line_out2 += str(line2[i_Ligand+1:i_Path])+\",\"\n",
        "                line_out2 += str(line2[i_Path+1:i_RMSD])+\",\"\n",
        "                line_out2 += str(line2[i_RMSD+1:i_S])+\",\"\n",
        "                line_out2 += str(line2[i_S+1:])\n",
        "                line_out2 = line_out2.\\\n",
        "                 replace(\"[\",\"\").replace(\"]\",\"\").replace(\"'\",\"\").replace(\" \",\"\")\n",
        "                data_out += line_out1+\",\"+line_out2+\"\\n\"\n",
        "                fo_mvd.close()\n",
        "                break\n",
        "\n",
        "    # Close fo_san\n",
        "    fo_san.close()\n",
        "\n",
        "    # Open a new file and write data_out\n",
        "    fo_new = open(mvd_out,\"w\")\n",
        "    fo_new.write(data_out)\n",
        "\n",
        "    # Close fo_new\n",
        "    fo_new.close()\n",
        "\n",
        "    # Show message\n",
        "    msg_out = \"\\nNumber of instances written to file \"+mvd_out+\" : \"\n",
        "    msg_out += str(count_instances)\n",
        "    print(msg_out)\n",
        "\n",
        "\n",
        "################################################################################\n",
        "# Define cross_validation() function                                           #\n",
        "################################################################################\n",
        "# Function to set up k-fold class. Kfold class to build a n fold\n",
        "# cross-validation loop and test the generalization ability of regression. With\n",
        "# cross-validation, we generally obtain a more conservative estimate(that is,\n",
        "# the error is larger). The cross-validation estimate is a better estimate of\n",
        "# how well we could generalize to predict on unseen data.\n",
        "#\n",
        "################################################################################\n",
        "# Reference                                                                    #\n",
        "################################################################################\n",
        "# Coelho LP, Richert W. (2015) Building Machine Learning Systems with\n",
        "# Python. 2nd ed. Packt Publishing Ltd. Birmingham UK. 301 pp. See page 162\n",
        "# (Cross-validation for regression)\n",
        "def cross_validation(model,X,y,n_splits,random_state,verbose):\n",
        "\n",
        "    # Import section\n",
        "    from sklearn.model_selection import KFold, cross_val_score\n",
        "    from warnings import simplefilter\n",
        "    import numpy as np\n",
        "\n",
        "    # Set up k-fold class\n",
        "    kf = KFold(n_splits=n_splits,shuffle=True, random_state=random_state)\n",
        "\n",
        "    # Ignore all future warnings\n",
        "    simplefilter(action='ignore', category=DeprecationWarning)\n",
        "\n",
        "    # Looping through kf.split()\n",
        "    for train,test in kf.split(X):\n",
        "\n",
        "        # Generate regression model\n",
        "        model.fit(X[train],y[train])\n",
        "\n",
        "    # Show Walsh metrics if requestes\n",
        "    if verbose:\n",
        "        # Show average coefficient of determination using n-fold crossvalidation\n",
        "        scores = cross_val_score(model,X,y,cv=kf)\n",
        "        msg_out = \"Average coefficient of determination using n-fold \"\n",
        "        msg_out += \"cross-validation\"\n",
        "        print(\"\\n\"+msg_out+\":\",np.mean(scores))\n",
        "\n",
        "    # Return model\n",
        "    return model\n",
        "\n",
        "################################################################################\n",
        "# Define Walsh_and_de_Azevedo_regression_metrics() function                    #\n",
        "################################################################################\n",
        "def Walsh_and_de_Azevedo_regression_metrics(y,y_pred):\n",
        "    \"\"\"Function to calculate regression metrics as recommended by\n",
        "    [Walsh et 2021](https://pubmed.ncbi.nlm.nih.gov/34316068/)\n",
        "    and\n",
        "    [de Azevedo Jr et al., 2024](https://pubmed.ncbi.nlm.nih.gov/38900052/))\n",
        "\n",
        "    Parameters:\n",
        "    y: array with experimental data (dimensions: (n,) )\n",
        "    y_pred: array with predicted values (dimensions: (n,) )\n",
        "\n",
        "    where n is the number of instances\n",
        "    \"\"\"\n",
        "\n",
        "    # Import section for this function\n",
        "    from sklearn.metrics import median_absolute_error\n",
        "    from scipy import stats\n",
        "    import numpy as np\n",
        "\n",
        "    # Determine metrics\n",
        "    n = len(y)\n",
        "    y_bar = y.mean()\n",
        "    rss = np.sum( np.square(y_pred - y) )\n",
        "    rmse =  np.sqrt( (np.sum( (y - y_pred)**2 ))/n )\n",
        "    ss_res = np.sum( (y - y_pred)**2 )\n",
        "    ss_tot = np.sum( (y - y_bar)**2 )\n",
        "    R2 = 1 - ss_res/ss_tot\n",
        "    mae = median_absolute_error(y,y_pred)\n",
        "    r_pearson,p_pearson = stats.pearsonr(y,y_pred)\n",
        "    rho,p_spearman = stats.spearmanr(y,y_pred)\n",
        "    core_dome = rmse**2 + mae**2 + (R2 -1)**2\n",
        "    dome =    np.sqrt( core_dome )\n",
        "    edomer2 = np.sqrt( core_dome + (r_pearson**2 - 1)**2 )\n",
        "    edomerho =np.sqrt( core_dome + (rho**2 - 1)**2 )\n",
        "    edome =   np.sqrt( core_dome + (r_pearson**2 - 1)**2 + (rho**2 - 1)**2 )\n",
        "\n",
        "    # Return metrics\n",
        "    return r_pearson,p_pearson,rho,p_spearman,rss,rmse,mae,R2,dome,edomer2,\\\n",
        "    edomerho,edome\n",
        "\n",
        "################################################################################\n",
        "# Define scatter_plot() function                                               #\n",
        "################################################################################\n",
        "def scatter_plot(y,y_pred,data_in,color_in,color_reg,plt_out,x_l,y_l,title,\n",
        "                 xlim,ylim,dpi):\n",
        "    \"\"\"Function to generate a basic scatter plot\"\"\"\n",
        "\n",
        "    # Import section\n",
        "    import numpy as np\n",
        "    import matplotlib.pyplot as plt\n",
        "\n",
        "    # Create basic scatterplot\n",
        "    plt.plot(y, y_pred,data_in,color=color_in)\n",
        "\n",
        "    # Obtain m (slope) and b(intercept) of Ridge regression line\n",
        "    m, b = np.polyfit(y,y_pred, 1)\n",
        "\n",
        "    y_in = np.linspace(xlim[0], xlim[1], 100)\n",
        "\n",
        "    # Use color_reg as color for regression line\n",
        "    plt.plot(y_in, m*y_in+b, color=color_reg)\n",
        "    plt.xlabel(x_l,fontsize = 14)\n",
        "    plt.ylabel(y_l,fontsize = 14)\n",
        "    plt.title(title,fontsize = 12)\n",
        "    plt.xlim( xlim )\n",
        "    plt.ylim( ylim )\n",
        "    plt.grid()\n",
        "    plt.savefig(plt_out,dpi=dpi)\n",
        "    plt.show()\n",
        "\n",
        "################################################################################\n",
        "# main function                                                                #\n",
        "################################################################################\n",
        "def main():\n",
        "\n",
        "    ############################################################################\n",
        "    # Split dataset\n",
        "    ############################################################################\n",
        "    sandres_bind_in = \"bind_Ki.csv\"              # Update here for a new dataset\n",
        "    url4sandres_bind_in = \\\n",
        "              \"https://bit.ly/cdk2_Ki_sandres\"   # Update here for a new dataset\n",
        "    msg_out = \"\\nDownloading the file \"+sandres_bind_in+\" and splitting it \"\n",
        "    msg_out += \"into training and test sets\"\n",
        "    print(msg_out,end=\"...\")\n",
        "    urlretrieve(url4sandres_bind_in, sandres_bind_in)\n",
        "\n",
        "    # Define PDB access codes for the test set\n",
        "    test_set = [\"EL4C\",\"38QX\",\"ARNY\",\"EB8B\",\"MHNA\",\"0UTY\",\"5GFN\",\"M50X\",\"BN07\",\n",
        "    \"7WSJ\",\"RBRI\",\"38R7\",\"XBUL\",\"1XN5\",\"SKL3\",\"ICRX\",\"VWKX\",\"INK5\",\"45B2\",\n",
        "    \"19M2\",\"3SM7\",\"06QM\",\"5ZEJ\",\"GKUV\",\"8AT0\",\"44WP\",\"8YJR\",\"RUXV\",\n",
        "                \"7PXQ\"]                          # Update here for a new dataset\n",
        "\n",
        "    # Try to open a csv file\n",
        "    try:\n",
        "        fo_bind = open(sandres_bind_in,\"r\")\n",
        "        csv_bind = csv.reader(fo_bind)\n",
        "    except IOError:\n",
        "        msg_out = \"\\nIOError! I can't find file \"+sandres_bind_in+\"!\"\n",
        "        sys.exit(msg_out)\n",
        "\n",
        "    # Get header\n",
        "    for line in csv_bind:\n",
        "        training_out = \"\"\n",
        "\n",
        "        # Looping through column labels\n",
        "        for line1 in line:\n",
        "            training_out += line1.strip()+\",\"\n",
        "\n",
        "        training_out = training_out[:len(training_out)-1]+\"\\n\"\n",
        "        test_out = training_out\n",
        "        break\n",
        "\n",
        "    # Assign zero to count_train\n",
        "    count_train = 0\n",
        "\n",
        "    # Looping through csv_bind\n",
        "    for line in csv_bind:\n",
        "        if line[0].strip() in test_set:\n",
        "            test_out += str(line).replace(\"[\",\"\").replace(\"]\",\"\").\\\n",
        "                        replace(\"'\",\"\").replace(\" \",\"\")+\"\\n\"\n",
        "        else:\n",
        "            training_out += str(line).replace(\"[\",\"\").replace(\"]\",\"\").\\\n",
        "                        replace(\"'\",\"\").replace(\" \",\"\")+\"\\n\"\n",
        "            count_train += 1\n",
        "\n",
        "    # Open new files\n",
        "    training_set_file = sandres_bind_in.replace(\".csv\",\"_training_set.csv\")\n",
        "    test_set_file = sandres_bind_in.replace(\".csv\",\"_test_set.csv\")\n",
        "    fo_training = open(training_set_file,\"w\")\n",
        "    fo_test = open(test_set_file,\"w\")\n",
        "\n",
        "    # Write data\n",
        "    fo_training.write(training_out)\n",
        "    fo_test.write(test_out)\n",
        "\n",
        "    # Close files\n",
        "    fo_bind.close()\n",
        "    fo_training.close()\n",
        "    fo_test.close()\n",
        "\n",
        "    # Show message\n",
        "    print(\"done!\")\n",
        "    msg_out = \"\\nFile \"+training_set_file+\" has \"+str(count_train)+\" instances.\"\n",
        "    msg_out += \"\\nFile \"+test_set_file+\" has \"+str(len(test_set))+\" instances.\"\n",
        "    print(msg_out)\n",
        "\n",
        "    ############################################################################\n",
        "    # Merge dataset\n",
        "    ############################################################################\n",
        "    # Docking results generated with MVD\n",
        "    mvd_in = \"cdk2_Ki_mvd.csv\"                   # Update here for a new dataset\n",
        "    url_mvd_in =  \"https://bit.ly/cdk2_Ki_mvd\"   # Update here for a new dataset\n",
        "    msg_out = \"\\nDownloading the file \"+mvd_in+\" to add SAnDReS-generated data\"\n",
        "    print(msg_out,end=\"...\")\n",
        "    urlretrieve(url_mvd_in, mvd_in)\n",
        "    print(\"done!\")\n",
        "\n",
        "    # Define SAnDReS-generated files for the training and test sets\n",
        "    san_training_in = \"bind_Ki_training_set.csv\"  # File generated using SAnDReS\n",
        "    san_test_in = \"bind_Ki_test_set.csv\"          # File generated using SAnDReS\n",
        "\n",
        "    # Call merge() function for the training set\n",
        "    merge(san_training_in,mvd_in)\n",
        "\n",
        "    # Call merge() function for the test set\n",
        "    merge(san_test_in,mvd_in)\n",
        "\n",
        "    ############################################################################\n",
        "    # Training set\n",
        "    ############################################################################\n",
        "    training_set_file_in = \\\n",
        "    \"cdk2_Ki_mvd_binding_training_set.csv\"       # Update here for a new dataset\n",
        "    print(\"\\nReading and preprocessing file: \"+training_set_file_in,end=\"...\")\n",
        "    feature_list = [\"Average Q\",\"C\",\"Csp2\",\"HeavyAtoms\",\n",
        "                    \"Steric\"]                    # Update here for a new dataset\n",
        "    training_set_data = pd.read_csv(training_set_file_in,delimiter=\",\")\n",
        "    y_train = training_set_data.iloc[:,8]        # Update here for a new dataset\n",
        "    print(\"done!\")\n",
        "\n",
        "    ############################################################################\n",
        "    # Test set\n",
        "    ###########################################################################\n",
        "    test_set_file_in = \\\n",
        "    \"cdk2_Ki_mvd_binding_test_set.csv\"           # Update here for a new dataset\n",
        "    print(\"\\nReading and preprocessing file: \"+test_set_file_in,end=\"...\")\n",
        "    test_set_data = pd.read_csv(test_set_file_in,delimiter=\",\")\n",
        "    y_test = test_set_data.iloc[:,8]             # Update here for a new dataset\n",
        "    print(\"done!\")\n",
        "\n",
        "    ############################################################################\n",
        "    # Define additional inputs\n",
        "    ############################################################################\n",
        "    plt_scatter_X_label_str,plt_scatter_y_label_str = \"pK$_i$ (Experimental)\",\\\n",
        "                        \"pK$_i$ (Predicted)\"            # Define labels for axes\n",
        "    plt_dpi = 1500                                      # Plot dpi\n",
        "\n",
        "    ############################################################################\n",
        "    # Build a regression model with Ridge\n",
        "    ############################################################################\n",
        "    msg_out = \"\\nGenerating regression models with Ridge\"\n",
        "    print(msg_out,end=\"...\")\n",
        "\n",
        "    # For training set\n",
        "    X_train = training_set_data[feature_list]\n",
        "    scaler_train = preprocessing.StandardScaler().fit(X_train)\n",
        "    X_train = scaler_train.transform(X_train)\n",
        "\n",
        "    # For test set\n",
        "    X_test = test_set_data[feature_list]\n",
        "    scaler_test = preprocessing.StandardScaler().fit(X_test)\n",
        "    X_test = scaler_test.transform(X_test)\n",
        "\n",
        "    # Set up string, list, and array for regression models\n",
        "    alpha_lst = [0.01,0.02,0.03,0.04,0.05,0.06,0.07,0.08,0.09,0.1]\n",
        "    rmse_array = np.zeros(len(alpha_lst))\n",
        "    data_out=\"Features,Data,n,r,p-value(r),r2,rho,p-value(rho),RSS,RMSE,MAE,R2,\"\n",
        "    data_out += \"DOME,EDOMEr2,EDOMErho,EDOME\\n\"\n",
        "\n",
        "    # Looping through alpha_lst\n",
        "    for i,a in enumerate(alpha_lst):\n",
        "\n",
        "        # Instantiate an object of Ridge() class\n",
        "        lin_reg = Ridge(alpha = a,\n",
        "                    solver = \"cholesky\",\n",
        "                    random_state = 46).fit(X_train, y_train)\n",
        "\n",
        "        # Call cross_validation() function\n",
        "        lin_reg = cross_validation(lin_reg,X_train,y_train,5,46,False)\n",
        "        y_train_pred_reg = lin_reg.predict(X_train)\n",
        "        y_test_pred_reg = lin_reg.predict(X_test)\n",
        "\n",
        "        # Call Walsh_and_de_Azevedo_regression_metrics() function\n",
        "        r_pearson,p_pearson,rho,p_spearman,rss,rmse,mae,R2,dome,\\\n",
        "        edomer2,edomerho,edome = \\\n",
        "        Walsh_and_de_Azevedo_regression_metrics(y_train,y_train_pred_reg)\n",
        "\n",
        "        # Set up output line\n",
        "        line_out = str(feature_list).replace(\"'\",\"\").replace(\"[\",\"\").\\\n",
        "                                                replace(\"]\",\"\").replace(\",\",\" \")\n",
        "        line_out += \",Training Set,\"+str(len(y_train))+\",\"\n",
        "        line_out += \"{:.4f}\".format(r_pearson)+\",\"\n",
        "        line_out += \"{:.4e}\".format(p_pearson)+\",\"\n",
        "        line_out += \"{:.4f}\".format(r_pearson**2)+\",\"\n",
        "        line_out += \"{:.4f}\".format(rho)+\",\"\n",
        "        line_out += \"{:.4e}\".format(p_spearman)+\",\"\n",
        "        line_out += \"{:.4e}\".format(rss)+\",\"\n",
        "        line_out += \"{:.4f}\".format(rmse)+\",\"+\"{:.4f}\".format(mae)+\",\"\n",
        "        line_out += \"{:.4f}\".format(R2)+\",\"+\"{:.4f}\".format(dome)+\",\"\n",
        "        line_out += \"{:.4f}\".format(edomer2)+\",\"\n",
        "        line_out += \"{:.4f}\".format(edomerho)+\",\"\n",
        "        line_out += \"{:.4f}\".format(edome)\n",
        "\n",
        "        # Update data_out\n",
        "        data_out += line_out+\"\\n\"\n",
        "\n",
        "        # Call Walsh_and_de_Azevedo_regression_metrics() function\n",
        "        r_pearson,p_pearson,rho,p_spearman,rss,rmse,mae,R2,dome,\\\n",
        "        edomer2,edomerho,edome = \\\n",
        "        Walsh_and_de_Azevedo_regression_metrics(y_test,y_test_pred_reg)\n",
        "\n",
        "        # Set up output line\n",
        "        line_out = str(feature_list).replace(\"'\",\"\").replace(\"[\",\"\").\\\n",
        "                                                replace(\"]\",\"\").replace(\",\",\" \")\n",
        "        line_out += \",Test Set,\"+str(len(y_test))+\",\"\n",
        "        line_out += \"{:.4f}\".format(r_pearson)+\",\"\n",
        "        line_out += \"{:.4e}\".format(p_pearson)+\",\"\n",
        "        line_out += \"{:.4f}\".format(r_pearson**2)+\",\"\n",
        "        line_out += \"{:.4f}\".format(rho)+\",\"\n",
        "        line_out += \"{:.4e}\".format(p_spearman)+\",\"\n",
        "        line_out += \"{:.4e}\".format(rss)+\",\"\n",
        "        line_out += \"{:.4f}\".format(rmse)+\",\"+\"{:.4f}\".format(mae)+\",\"\n",
        "        line_out += \"{:.4f}\".format(R2)+\",\"+\"{:.4f}\".format(dome)+\",\"\n",
        "        line_out += \"{:.4f}\".format(edomer2)+\",\"\n",
        "        line_out += \"{:.4f}\".format(edomerho)+\",\"\n",
        "        line_out += \"{:.4f}\".format(edome)\n",
        "\n",
        "        # Update data_out\n",
        "        data_out += line_out+\"\\n\"\n",
        "\n",
        "        # Add rmse to rmse_array\n",
        "        rmse_array[i] = rmse\n",
        "\n",
        "    print(\"done!\")\n",
        "\n",
        "    # Find minimum\n",
        "    min_i = np.argmin(rmse_array)\n",
        "    best_alpha = alpha_lst[min_i]\n",
        "    best_rmse = rmse_array[min_i]\n",
        "    msg_out = \"\\nMinimum RMSE: {:.4f}\".format(best_rmse)\n",
        "    msg_out += \" for alpha = {:.4f}\".format(best_alpha)\n",
        "    print(msg_out)\n",
        "\n",
        "    # Regression for best alpha\n",
        "    # Instantiate an object of Ridge() class\n",
        "    lin_reg = Ridge(alpha = best_alpha,\n",
        "                    solver = \"cholesky\",\n",
        "                    random_state = 46).fit(X_train, y_train)\n",
        "\n",
        "    # Call cross_validation() function\n",
        "    lin_reg = cross_validation(lin_reg,X_train,y_train,5,46,False)\n",
        "    y_train_pred_reg = lin_reg.predict(X_train)\n",
        "    y_test_pred_reg = lin_reg.predict(X_test)\n",
        "\n",
        "    # Call scatter_plot() function\n",
        "    print(\"\\nSaving scatter plot to file: test_set_scatter_plot.pdf\",end=\"...\")\n",
        "    plt_scatter_title_str = \"CDK2 K$_i$ (Test Set)\"\n",
        "                                                 # Update here for a new dataset\n",
        "    plt_scatter_title_str += \" (RMSE: {:.4f}\".format(best_rmse)+\")\"\n",
        "    plot_out = \"test_set_scatter_plot.pdf\"\n",
        "    xlim = (2.0, 12.0)                           # Update here for a new dataset\n",
        "    ylim = (2.0, 12.0)                           # Update here for a new dataset\n",
        "    scatter_plot(y_test,y_test_pred_reg,\".\",\"blue\",\"black\",plot_out,\n",
        "    plt_scatter_X_label_str,plt_scatter_y_label_str,plt_scatter_title_str,xlim,\n",
        "             ylim,plt_dpi)\n",
        "    print(\"done!\")\n",
        "\n",
        "    # Open a new file and write data_out\n",
        "    fo_out = open(\"regression_metrics.csv\",\"w\")\n",
        "    fo_out.write(data_out)\n",
        "    fo_out.close()\n",
        "    print(\"\\nNumber of generated regression models: \",len(alpha_lst))\n",
        "    print(\"\\nMetrics saved to file: regression_metrics.csv\")\n",
        "\n",
        "\n",
        "# Call main() function\n",
        "main()\n",
        "################################################################################"
      ]
    }
  ]
}