creating_RDKIT_dataset.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Initialisations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-30T09:35:23.292935Z",
     "start_time": "2020-01-30T09:35:14.178639Z"
    },
    "code_folding": []
   },
   "outputs": [],
   "source": [
    "import pandas as pd ; pd.set_option('display.max_columns', None, 'display.max_rows', None)\n",
    "pd.set_option('expand_frame_repr', True)\n",
    "import numpy as np\n",
    "import pubchempy as pcp\n",
    "from rdkit import Chem\n",
    "import cirpy #for changing between identifiers\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "\n",
    "#for web integration\n",
    "from googlesearch import search\n",
    "import webbrowser\n",
    "\n",
    "#for molecular drawings and descriptor calculation\n",
    "from rdkit import Chem\n",
    "from rdkit.Chem import Draw\n",
    "from rdkit.Chem import rdMolDescriptors as descriptor\n",
    "from rdkit.Chem import Crippen\n",
    "from rdkit.Chem import AllChem\n",
    "#size = (120, 120)  # Smaller figures than the default.\n",
    "\n",
    "from sklearn import linear_model\n",
    "import statsmodels.api as sm\n",
    "\n",
    "#for progress bars\n",
    "from ipywidgets import widgets, Layout\n",
    "from IPython.display import display\n",
    "from ipywidgets import IntProgress \n",
    "#Run this part outside of the loop - create progress bar\n",
    "\n",
    "#max_count = len(cal_df.index.values)\n",
    "#f = IntProgress(min=0, max=max_count, description='CIDs') # instantiate the bar\n",
    "#display(f) # display the bar\n",
    "#count = 0\n",
    "    \n",
    "    \n",
    "#Place this part inside the loop\n",
    "\n",
    "    #f.value += 1 # signal to increment the progress bar\n",
    "    #count += 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Each constant is marked with a symbol denoting how the HLC was obtained\n",
    "- Some values are calculated (T for thermo calc, Q for QSAR calc)\n",
    "- Some values are estimates (E)\n",
    "- Some values are probable incorrect (W)\n",
    "- As stated above, some data was extrapolated to a different temperature or concentration range (R)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here I will remove estimated/calculated values from the dataset"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Getting data from the supplement"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-30T09:35:26.087314Z",
     "start_time": "2020-01-30T09:35:25.864932Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "noVPAS 4632 4632 4632 \n",
      " <class 'list'> <class 'list'> <class 'list'> \n",
      " <class 'str'> <class 'str'> <class 'list'>\n",
      "VPAS 4632 4632 4632 \n",
      " <class 'list'> <class 'list'> <class 'list'> \n",
      " <class 'str'> <class 'str'> <class 'list'>\n"
     ]
    }
   ],
   "source": [
    "file = open(\"henry_law_dataset/output/HcpSI.f90\", \"r\")\n",
    "species_f90 = (file.read().split(sep=\"! species:\")) # create list where each element contains one species\n",
    "species_list = [] \n",
    "constant_list_noVPAS = []\n",
    "constant_list_VPAS = []\n",
    "cas_list = []\n",
    "\n",
    "for i in species_f90[1:]: #first line is ! alpha, skip this line with [1:]\n",
    "    cas_line = i.split(sep='\\n')[3] # fourth line contains CAS number\n",
    "    cas_list.append(cas_line.split(sep='! casrn:   ')[1]) # isolate CAS number and place it in the list\n",
    "    \n",
    "    species = i.split(sep='\\n')[0] #first line contains species name\n",
    "    species_list.append(species[1:]) # isolate species name and place it in list\n",
    "    hcpsi_lines = (i.split(sep='\\n')[5:]) #lines past the 6th contain constants\n",
    "    # each species has a different number of constants\n",
    "    # loop through them and pull the constant from each line, place them into a list\n",
    "    const_per_line_VPAS = []\n",
    "    const_per_line_noVPAS = [] #list for each species' constants\n",
    "    for j in range(0, len(hcpsi_lines)):\n",
    "        if 'HcpSI' in hcpsi_lines[j]: #check if line contains constant\n",
    "            const_type = hcpsi_lines[j].split(sep='type: ')[1].split(sep=',')[0] # isolate 'type' code\n",
    "            if (const_type != 'W') and (const_type != 'E') and (const_type != 'T') and (const_type != 'Q') and (const_type != 'V'):\n",
    "                # remove calculated / incorrect constants\n",
    "                const_per_line_noVPAS.append(hcpsi_lines[j].split(sep='=  ')[1][0:8]) # add remaining constants to list\n",
    "            if (const_type != 'W') and (const_type != 'E') and (const_type != 'T') and (const_type != 'Q'):\n",
    "                const_per_line_VPAS.append(hcpsi_lines[j].split(sep='=  ')[1][0:8]) # add remaining constants to list\n",
    "        else:\n",
    "            pass\n",
    "    \n",
    "    if 'HcpSI' in i:\n",
    "        new_value_list_noVPAS = []\n",
    "        for const in const_per_line_noVPAS:\n",
    "            # convert constants from strings to floats\n",
    "            if type(const) == str:\n",
    "                number = float(const.split(sep='E')[0])\n",
    "                power = float((const.split(sep='E')[1]))\n",
    "                new_value = (number*np.power(10, power)) # create float\n",
    "                log_value = np.log(1/new_value) # inverse of H gives Kh, then take log.\n",
    "            new_value_list_noVPAS.append(log_value) # put new value into list per species\n",
    "        constant_list_noVPAS.append(new_value_list_noVPAS) # put new lists into final list of constant lists\n",
    "    else:\n",
    "        constant_list_noVPAS.append(np.nan)\n",
    "\n",
    "    if 'HcpSI' in i:\n",
    "        new_value_list_VPAS = []\n",
    "        for const in const_per_line_VPAS:\n",
    "            # convert constants from strings to floats\n",
    "            if type(const) == str:\n",
    "                number = float(const.split(sep='E')[0])\n",
    "                power = float((const.split(sep='E')[1]))\n",
    "                new_value = (number*np.power(10, power)) # create float\n",
    "                log_value = np.log(1/new_value) # inverse of H gives Kh, then take log.\n",
    "            new_value_list_VPAS.append(log_value) # put new value into list per species\n",
    "        constant_list_VPAS.append(new_value_list_VPAS) # put new lists into final list of constant lists\n",
    "    else:\n",
    "        constant_list_VPAS.append(np.nan)\n",
    "        \n",
    "print(\"noVPAS\", len(cas_list), len(species_list), len(constant_list_noVPAS), '\\n',\n",
    "     type(cas_list), type(species_list), type(constant_list_noVPAS), '\\n',\n",
    "     type(cas_list[0]), type(species_list[0]), type(constant_list_noVPAS[0])) #checking lists are full and correct type\n",
    "print(\"VPAS\", len(cas_list), len(species_list), len(constant_list_VPAS), '\\n',\n",
    "     type(cas_list), type(species_list), type(constant_list_VPAS), '\\n',\n",
    "     type(cas_list[0]), type(species_list[0]), type(constant_list_VPAS[0])) #checking lists are full and correct type"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Conversion of CAS numbers to SMILES strings"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T17:25:55.908037Z",
     "start_time": "2020-01-26T16:47:34.488124Z"
    }
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "82c3a4f42b5a459a866e4e127629756f",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "IntProgress(value=0, description='CAS->SMILES', max=4632)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "smiles_list = []\n",
    "\n",
    "max_count = len(cas_list)\n",
    "f = IntProgress(min=0, max=max_count, description='CAS->SMILES') # instantiate the bar\n",
    "display(f) # display the bar\n",
    "count = 0\n",
    "    \n",
    "for cas in cas_list:\n",
    "    smiles_list.append(cirpy.resolve(cas, 'smiles')) #convert CAS to SMILES\n",
    "    f.value += 1 # signal to increment the progress bar\n",
    "    count += 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T17:40:25.664608Z",
     "start_time": "2020-01-26T17:40:25.035938Z"
    }
   },
   "outputs": [],
   "source": [
    "for dataset in ['VPAS', 'noVPAS']:\n",
    "    constant_list = eval('constant_list_%s'%dataset)\n",
    "    min_list = []; max_list = []; mean_list = []; median_list = []; num_constants_list = []\n",
    "    for clist in constant_list:\n",
    "        if (type(clist) == list):\n",
    "            if (len(clist) != 0):\n",
    "                #get various stats for each species' constants\n",
    "                min_list.append(np.min(clist))\n",
    "                max_list.append(np.max(clist))\n",
    "                mean_list.append(np.mean(clist))\n",
    "                median_list.append(np.median(clist))\n",
    "            else:\n",
    "                min_list.append(0)\n",
    "                max_list.append(0)\n",
    "                mean_list.append(0)\n",
    "                median_list.append(0)\n",
    "        else:\n",
    "            min_list.append(0)\n",
    "            max_list.append(0)\n",
    "            mean_list.append(0)\n",
    "            median_list.append(0)\n",
    "        if type(clist) == list:\n",
    "            num_constants_list.append(len(clist))\n",
    "        else:\n",
    "            num_constants_list.append(0)\n",
    "    const_list_df = pd.DataFrame(species_list)\n",
    "    const_list_df.rename({0:'Species'}, axis=1, inplace=True)\n",
    "    const_list_df['Num'] = num_constants_list\n",
    "    const_list_df['min'] = min_list    \n",
    "    const_list_df['max'] = max_list    \n",
    "    const_list_df['mean'] = mean_list    \n",
    "    const_list_df['median'] = median_list\n",
    "    const_list_df['max-min'] = np.array(max_list) - np.array(min_list)\n",
    "    const_list_df.set_index('Species', inplace=True)\n",
    "    exec('const_list_df_%s = const_list_df'%dataset)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Visualise data distribution, look at varying constants"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T17:41:19.581962Z",
     "start_time": "2020-01-26T17:41:19.192929Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEcCAYAAAAoSqjDAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAdk0lEQVR4nO3de7wU5Z3n8c9XUPEWlYAGFT1qjPEOBB0TkxmjxiiOEjNr1HEioybEV3ATN5m8RM0Y1407uuslcTKrwdENuJqId2bEUXRiHDPxggS5iEZUFIQA8YZ3A/ntH/V0WTTd5/ThnOruc/i+X69+naqnqp769UNTv66nqp9SRGBmZgawUasDMDOz9uGkYGZmOScFMzPLOSmYmVnOScHMzHJOCmZmlnNSMDOznJOCmZnlnBSsT5A0X9KhPdh+saQRvRDHP0g6uxfqeUzSPj2tp6/o6b+fNY+TguUkLZL0gaQhVeWzJYWkjlbtMyL2iYgH13Mf2wIfAxb0MNahwKnAT4t1pzh3SfOSdFF6Xwd0Ut1lwEU9iacv6cm/nzWXk4JVewE4uTIjaT9gsz6+z/2A30XE+z2s52+B6RHxbqFsBPBaRLwoaQvgVuAw4KCIeLKTuqYBn5c0rIcxdUnSwLL3Yf2Hk4JVu4Hs23DFOGBKcQVJEyU9J+lNSU9JOj6V7y7pVUmj0vwOkv7QQLdBI/tcJOmIwvTfSZoj6Q1JN0sa1En9+wNz07aSdI6kFyW9LmmqpK3Tso0kXZC6mpZKOjadxWyb6jka+FVV3SOA2ZJ2Bh4G3gAOi4gVnb3hiHgPeAI4st466X2em9r4NUn/t/I+U9veJmmlpBckfavGtudImgO8XSsxpOUvp3/HZyQd3gv7HS7p9rT8FUk/KdR5RGG9uvXUi8uaJCL88ouIAFgEHAE8A+wFDAAWA7sAAXSk9U4AdiD7UnEi8DYwLC37Olk3zebAvcBlvbTPRcARhenHUgyD0/7O7GQfk4Dz0vQPgQeBYcAg4C7g8rTsIuAhYEdgG+A3wOJCPSuBA6vqngL8B7AU+E432/sq4Iou2mYeMDy9z1+n+DciSygXAJsAuwHPA1+s2nZ22nazGnXvmdp5hzTfAezek/2mf7sngSuBLVL7frbGv1/dejqLy6/mvHymYLVUvrl/AXgaeLm4MCJuiYilEfGniLgZeBY4KC27Ns0/SnbgPb839lnDVSmGV4F/IfvGXs9+wBxJ2wP/FfjriFgW2bf1W4HR6XrBfwNOj4iXI+J14D7SGUayDfBmVd0jgH2AORFxRfWOJd1SvKAsaVbhW/ubqc7O/CQiFqf3eTFZN9uBwNCIuCgiPoiI54FrgZOqtr0qbfsu61oDbArsLWnjiFgUEc/1cL8HkSXq70XE2xHxXkQ8XGPfndXTVVxWMicFq+UG4K/J+tCnVC+UdGq6EPy6pNeBfYHiheJrU9k/RurHl3SKpLfS657u7rOG3xem3wG2rLWSJKVY5gKfA+ZGxNLCKkOAZcDhwMKIWFhYNpi1k8JrwFaFujclO7v5K2Av1b4r6ZPAwrT+AGBNRKxOy7YCXu/8bbK4MP0i2UF3F2CHSvunf4PzgO072XYt6X2eDVwIrJD0C0k79HC/w4EXC++vnrr1NBCXlcxJwdYRES+SXfwdA9xeXKbsLptrgbOAj0bENmRdDUrLtwR+BFwHXChpcKrzxojYMr2O7s4+e2hX4E+p/qFkff5FY8muBQwh6wKqvM8BZNcQ5hTWnQN8ojC/L9k32/8Ajgf+h6TDCnVsAgyMDy9w70pKEMleZN0tnRlemN45xbgYeCEitim8toqIMVXbdvqwlIi4KSI+y4dddZf2cL+LgZ0buLDdaT1dxGUlc1Kwes4gu2D6dlX5FmT/UVcCSDqN7OBY8WPgiYj4GnA3cE0v7LMn8ovMwOPAp9MF8S0lXUT2Lfd6susSn5H0cUkfIevv3521zxSmA39RmB8JzIuI1RExC5gATJW0a1q+FzBU0oOSHgTuIEuglbOMTwEzuoh/gqSdUnI9D7iZ7HrKqnRBdjNJAyTtK+nARhtF0p6SDktxvAe8S5bgerLfx8jOui6RtIWkQZIOqbH7uvU0EJeVzEnBaoqI5yJiZo3yp4DLyS7CLifrr/81gKSxwFHAmWn17wCjJJ3Sk3320H6kb/up7ovJzgyWkB20D4uIdyLiAeDnwG+BmWQXad8hu75RMQUYI6lyu+yItF4l/inATcCdym5P3Re4OiIOjYhDganA/LT6ccCDVV1ZtdxEdm3j+fT6YUSsAY5N+38B+APwz8DWjTcLmwKXpG1/D2xHdvBf7/0Wln8ceImsjU+s3nEX9XQVl5VMEX4cp1k1SWcCx0TEsVXl/xNYERE/aqCOfwBmRcQtaf4W4NyIWCjpUeCMiJjXyfaLgK9FxP09eCvd1qr9Wnvwj1rMAEkHk3V9LCa76HwR2bf5tUREd7617gP8rDD/cbJv3UTEn61vrGZlclIwy4wkuwayMfA74G8j4pGeVBgRx1XNj+xJfWbN4O4jMzPL+UKzmZnl+nT30ZAhQ6Kjo6PVYZiZ9SlPPPHEHyJiaK1lfTopdHR0MHNmb9/BaGbWv0l6sd4ydx+ZmVnOScHMzHJOCmZmlnNSMDOznJOCmZnlnBTMzCznpGBmZjknBTMzyzkpmJlZrk//orknOibe3fC6iy45psRIzMzah88UzMws56RgZmY5JwUzM8s5KZiZWc5JwczMck4KZmaWc1IwM7Ock4KZmeWcFMzMLOekYGZmudKSgqThkn4paYGk+ZK+ncovlPSypNnpNaawzbmSFkp6RtIXy4rNzMxqK3Pso9XAdyNilqStgCckzUjLroyIy4orS9obOAnYB9gBuF/SJyJiTYkxmplZQWlnChGxLCJmpek3gQXAjp1sMhb4RUS8HxEvAAuBg8qKz8zM1tWUawqSOoCRwKOp6CxJcyRdL2nbVLYjsLiw2RJqJBFJ4yXNlDRz5cqVJUZtZrbhKT0pSNoSuA04OyJWAVcDuwMjgGXA5ZVVa2we6xRETIqI0RExeujQoSVFbWa2YSo1KUjamCwh3BgRtwNExPKIWBMRfwKu5cMuoiXA8MLmOwFLy4zPzMzWVubdRwKuAxZExBWF8mGF1Y4H5qXpacBJkjaVtCuwB/BYWfGZmdm6yrz76BDgq8BcSbNT2XnAyZJGkHUNLQK+ARAR8yVNBZ4iu3Npgu88MjNrrtKSQkQ8TO3rBNM72eZi4OKyYjIzs875F81mZpZzUjAzs5yTgpmZ5ZwUzMws56RgZmY5JwUzM8s5KZiZWc5JwczMck4KZmaWc1IwM7Ock4KZmeWcFMzMLOekYGZmOScFMzPLOSmYmVnOScHMzHJOCmZmlnNSMDOznJOCmZnlnBTMzCznpGBmZjknBTMzyzkpmJlZzknBzMxyTgpmZpZzUjAzs5yTgpmZ5ZwUzMws56RgZmY5JwUzM8uVlhQkDZf0S0kLJM2X9O1UPljSDEnPpr/bpnJJukrSQklzJI0qKzYzM6utzDOF1cB3I2Iv4GBggqS9gYnAAxGxB/BAmgc4GtgjvcYDV5cYm5mZ1VBaUoiIZRExK02/CSwAdgTGApPTapOBL6XpscCUyDwCbCNpWFnxmZnZuppyTUFSBzASeBTYPiKWQZY4gO3SajsCiwubLUll1XWNlzRT0syVK1eWGbaZ2Qan9KQgaUvgNuDsiFjV2ao1ymKdgohJETE6IkYPHTq0t8I0MzNKTgqSNiZLCDdGxO2peHmlWyj9XZHKlwDDC5vvBCwtMz4zM1tbmXcfCbgOWBARVxQWTQPGpelxwF2F8lPTXUgHA29UupnMzKw5BpZY9yHAV4G5kmansvOAS4Cpks4AXgJOSMumA2OAhcA7wGklxmZmZjWUlhQi4mFqXycAOLzG+gFMKCseMzPrmn/RbGZmOScFMzPLOSmYmVnOScHMzHJOCmZmlnNSMDOznJOCmZnlnBTMzCznpGBmZjknBTMzyzkpmJlZzknBzMxyTgpmZpZzUjAzs5yTgpmZ5RpKCpL2LTsQMzNrvUbPFK6R9Jikb0raptSIzMysZRpKChHxWeAUYDgwU9JNkr5QamRmZtZ0DV9TiIhnge8D5wB/AVwl6WlJXy4rODMza65GrynsL+lKYAFwGHBsROyVpq8sMT4zM2uigQ2u9xPgWuC8iHi3UhgRSyV9v5TIzMys6RpNCmOAdyNiDYCkjYBBEfFORNxQWnRmZtZUjV5TuB/YrDC/eSozM7N+pNGkMCgi3qrMpOnNywnJzMxapdGk8LakUZUZSZ8C3u1kfTMz64MavaZwNnCLpKVpfhhwYjkhmZlZqzSUFCLicUmfBPYEBDwdEX8sNTIzM2u6Rs8UAA4EOtI2IyUREVNKicrMzFqioaQg6QZgd2A2sCYVB+CkYGbWjzR6pjAa2DsiotGKJV0P/CWwIiL2TWUXAl8HVqbVzouI6WnZucAZZEnnWxFxb6P7MjOz3tHo3UfzgI91s+6fAUfVKL8yIkakVyUh7A2cBOyTtvk/kgZ0c39mZtZDjZ4pDAGekvQY8H6lMCKOq7dBRDwkqaPB+scCv4iI94EXJC0EDgJ+0+D2ZmbWCxpNChf24j7PknQqMBP4bkS8BuwIPFJYZ0kqW4ek8cB4gJ133rkXwzIzs0afp/ArYBGwcZp+HJi1Hvu7muyC9QhgGXB5Klet3daJZVJEjI6I0UOHDl2PEMzMrJ5Gh87+OnAr8NNUtCNwZ3d3FhHLI2JNRPyJbNTVg9KiJWQP8KnYCVhavb2ZmZWr0QvNE4BDgFWQP3Bnu+7uTNKwwuzxZBewAaYBJ0naVNKuwB7AY92t38zMeqbRawrvR8QHUtbLI2kgdbp3KiT9HDgUGCJpCfAD4FBJI9K2i4BvAETEfElTgaeA1cCEyjDdZmbWPI0mhV9JOg/YLD2b+ZvAv3S2QUScXKP4uk7Wvxi4uMF4zMysBI12H00k+8HZXLJv99PJntdsZmb9SKMD4lUuDF9bbjhmZtZKjY599AI1riFExG69HpGZmbVMd8Y+qhgEnAAM7v1wzMyslRr98dorhdfLEfEj4LCSYzMzsyZrtPtoVGF2I7Izh61KicjMzFqm0e6jywvTq8l+Y/CVXo+mTXVMvLuh9RZdckzJkZiZlavRu48+X3YgZmbWeo12H32ns+URcUXvhGNmZq3UnbuPDiQbowjgWOAhYHEZQZmZWWt05yE7oyLiTcgfq3lLRHytrMDMzKz5Gh3mYmfgg8L8B0BHr0djZmYt1eiZwg3AY5LuIPtl8/HAlNKiMjOzlmj07qOLJd0DfC4VnRYRvy0vLDMza4VGu48ANgdWRcSPgSXpYThmZtaPNPo4zh8A5wDnpqKNgf9XVlBmZtYajZ4pHA8cB7wNEBFL8TAXZmb9TqNJ4YOICNLw2ZK2KC8kMzNrlUaTwlRJPwW2kfR14H78wB0zs36n0buPLkvPZl4F7AlcEBEzSo3MzMyarsukIGkAcG9EHAE4EZiZ9WNddh9FxBrgHUlbNyEeMzNroUZ/0fweMFfSDNIdSAAR8a1SojIzs5ZoNCncnV5mZtaPdZoUJO0cES9FxORmBWRmZq3T1TWFOysTkm4rORYzM2uxrpKCCtO7lRmImZm1XldJIepMm5lZP9TVheYDJK0iO2PYLE2T5iMiPlJqdGZm1lSdJoWIGNCsQMzMrPW68zyFbpF0vaQVkuYVygZLmiHp2fR321QuSVdJWihpjqRRZcVlZmb1lZYUgJ8BR1WVTQQeiIg9gAfSPMDRwB7pNR64usS4zMysjtKSQkQ8BLxaVTwWqPzmYTLwpUL5lMg8QjYa67CyYjMzs9rKPFOoZfuIWAaQ/m6XyncEFhfWW5LK1iFpvKSZkmauXLmy1GDNzDY0zU4K9ahGWc1bYCNiUkSMjojRQ4cOLTksM7MNS7OTwvJKt1D6uyKVLwGGF9bbCVja5NjMzDZ4zU4K04BxaXoccFeh/NR0F9LBwBuVbiYzM2ueRkdJ7TZJPwcOBYZIWgL8ALiE7NGeZwAvASek1acDY4CFwDvAaWXFZWZm9ZWWFCLi5DqLDq+xbgATyorFzMwa0y4Xms3MrA04KZiZWc5JwczMck4KZmaWc1IwM7Ock4KZmeWcFMzMLOekYGZmOScFMzPLOSmYmVnOScHMzHJOCmZmlnNSMDOznJOCmZnlnBTMzCznpGBmZjknBTMzyzkpmJlZzknBzMxyTgpmZpZzUjAzs9zAVgfQn3RMvLuh9RZdckzJkZiZrR+fKZiZWc5JwczMck4KZmaWc1IwM7Ock4KZmeWcFMzMLOdbUlvAt66aWbtqSVKQtAh4E1gDrI6I0ZIGAzcDHcAi4CsR8Vor4jMz21C1svvo8xExIiJGp/mJwAMRsQfwQJo3M7MmaqdrCmOByWl6MvClFsZiZrZBalVSCOA+SU9IGp/Kto+IZQDp73a1NpQ0XtJMSTNXrlzZpHDNzDYMrbrQfEhELJW0HTBD0tONbhgRk4BJAKNHj46yAmwHviBtZs3WkjOFiFia/q4A7gAOApZLGgaQ/q5oRWxmZhuypicFSVtI2qoyDRwJzAOmAePSauOAu5odm5nZhq4V3UfbA3dIquz/poj4N0mPA1MlnQG8BJzQgtjMzDZoTU8KEfE8cECN8leAw5sdj5mZfaidbkk1M7MWc1IwM7Ock4KZmeWcFMzMLOekYGZmOScFMzPLOSmYmVnOScHMzHJOCmZmlvPjOPsBj6ZqZr3FZwpmZpZzUjAzs5y7jzYg7mYys674TMHMzHJOCmZmlnNSMDOznJOCmZnlnBTMzCznpGBmZjknBTMzyzkpmJlZzknBzMxyTgpmZpZzUjAzs5yTgpmZ5ZwUzMws56RgZmY5D51t62h0iG3o/WG2231471a2jVkzOClYj7T7QbxR3TnYm/VnTgrWr/lgb9Y9bZcUJB0F/BgYAPxzRFzS4pCsF/T2wdkHe7NytFVSkDQA+CfgC8AS4HFJ0yLiqdZGZlaeVnXB9ZeuP+tdbZUUgIOAhRHxPICkXwBjAScF63N8dtQ8/SnBtfq9tFtS2BFYXJhfAvxZcQVJ44HxafYtSc+s576GAH9Yz21boS/F25dihb4Vb9Nj1aU92ryt2raL99JWsXZFl/Yo3l3qLWi3pKAaZbHWTMQkYFKPdyTNjIjRPa2nWfpSvH0pVuhb8falWKFvxduXYoXy4m23H68tAYYX5ncClrYoFjOzDU67JYXHgT0k7SppE+AkYFqLYzIz22C0VfdRRKyWdBZwL9ktqddHxPySdtfjLqgm60vx9qVYoW/F25dihb4Vb1+KFUqKVxHR9VpmZrZBaLfuIzMzayEnBTMzy/X7pCDpKEnPSFooaWKN5ZtKujktf1RSR/OjzGMZLumXkhZImi/p2zXWOVTSG5Jmp9cFrYg1xbJI0twUx8wayyXpqtS2cySNakWcKZY9C202W9IqSWdXrdOytpV0vaQVkuYVygZLmiHp2fR32zrbjkvrPCtpXAvj/d+Snk7/1ndI2qbOtp1+bpoU64WSXi78W4+ps22nx48mxntzIdZFkmbX2bbnbRsR/fZFdrH6OWA3YBPgSWDvqnW+CVyTpk8Cbm5hvMOAUWl6K+B3NeI9FPjXVrdtimURMKST5WOAe8h+f3Iw8GirYy58Ln4P7NIubQv8OTAKmFco+1/AxDQ9Ebi0xnaDgefT323T9LYtivdIYGCavrRWvI18bpoU64XA3zXwOen0+NGseKuWXw5cUFbb9vczhXzYjIj4AKgMm1E0Fpicpm8FDpdU60d0pYuIZRExK02/CSwg+5V3XzUWmBKZR4BtJA1rdVDA4cBzEfFiqwOpiIiHgFerioufzcnAl2ps+kVgRkS8GhGvATOAo0oLNKkVb0TcFxGr0+wjZL8zark6bduIRo4fva6zeNOx6SvAz8vaf39PCrWGzag+yObrpA/0G8BHmxJdJ1I31kjg0RqLPy3pSUn3SNqnqYGtLYD7JD2Rhh+p1kj7t8JJ1P9P1S5tC7B9RCyD7AsDsF2Nddq1jU8nO0uspavPTbOclbq6rq/TNdeObfs5YHlEPFtneY/btr8nhS6HzWhwnaaStCVwG3B2RKyqWjyLrNvjAOAfgTubHV/BIRExCjgamCDpz6uWt2PbbgIcB9xSY3E7tW2j2rGNzwdWAzfWWaWrz00zXA3sDowAlpF1yVRru7YFTqbzs4Qet21/TwqNDJuRryNpILA163eq2SskbUyWEG6MiNurl0fEqoh4K01PBzaWNKTJYVZiWZr+rgDuIDvdLmrHYUuOBmZFxPLqBe3UtsnySndb+ruixjpt1cbpQvdfAqdE6uSu1sDnpnQRsTwi1kTEn4Br68TQbm07EPgycHO9dXqjbft7Umhk2IxpQOWOjf8C/Hu9D3PZUn/hdcCCiLiizjofq1zzkHQQ2b/hK82LMo9jC0lbVabJLjLOq1ptGnBqugvpYOCNSndIC9X9ptUubVtQ/GyOA+6qsc69wJGStk1dIEemsqZT9oCsc4DjIuKdOus08rkpXdW1rePrxNBuw+4cATwdEUtqLey1ti37SnqrX2R3wPyO7C6C81PZRWQfXIBBZF0JC4HHgN1aGOtnyU5P5wCz02sMcCZwZlrnLGA+2Z0QjwCfaVGsu6UYnkzxVNq2GKvIHpr0HDAXGN3iz8LmZAf5rQtlbdG2ZIlqGfBHsm+oZ5Bd23oAeDb9HZzWHU32VMLKtqenz+9C4LQWxruQrA++8tmt3NW3AzC9s89NC2K9IX0m55Ad6IdVx5rm1zl+tCLeVP6zyme1sG6vt62HuTAzs1x/7z4yM7NucFIwM7Ock4KZmeWcFMzMLOekYGZmOScFMzPLOSmYmVnOScFKJykk3VCYHyhppaR/7YW6jyjWvR7bbybpV5IGSPqGpGtS+caSbpA0OQ09UmvbTSQ9lIYf6M4+u2wPSf/ZzTp7rR3SfMNtsb7tYO3JScGa4W1gX0mbpfkvAC/3Ut0HkP2Cc32dDtweEWuA/YE5kj5CNsLnSxExLiL+WGvDyIZTfgA4sZv77LI9IuIz3ayzN9sButEWPWgHa0NOCtYs9wDHpOm1xh+SdGca6nd+ZbhfSQemYY0HpTFd5kvat0a9B5ANqYCkT6ZvrPMl3V8ZzE7SXql8jqTvSVpY2P4UPhxTaD+ywRAfBG6JiPMbeF93pjq6q257pJjfSn87lD2J79r0vu4rJJOi3mwH6H5brG87WLtpxlgefm3YL+Atsm+et5KNNTWbwlPO+HBMn83IBvD6aJr/IXAZ2fhJ59ap+0lgKLAp2XgvI1P5OcDFwECyIbEr5VcDd6bpTYDfF+p6DVgJHNuN9zYAWNmb7VFZJ/3tIBuGekSanwr8TZntsD5tsT7t4Fd7vnymYE0REXPIDnAnA9OrFn9LUmUQuuHAHqn8IrKuldFkj6ZcS+rf/khErCR7KtnDEfHbtPgpsofSfBl4sqq80s0yBHg91TWc7GA9l+yxqNX7uqswfUul7z2y7pYPKqNTNqqL9qj2QkRUnsn7RNquGFuvtUOqr25b9HY7WPtxUrBmmkb2zb/YdXQo2ZDAn47s4Ta/Jfv2DNlzh7cke171INa1N9kjSyvTcwvL9iM78O1P6lZJ9i3Mv1uod3+yg+QpwPmSRhZiHE42amXFRvFh3ztk387fKwYmaYI+fND6DjVihxrtUcf7hek1ZN/6i3qzHaBOW6xPO1jf46RgzXQ9cFFEFA9aWwOvRcQ7kj4JHFxYNgn4e7IneF1ao768H53sQu3eAJJ2A74KTCEbKvsTqXwE8Dekb8iRPdN4gKRBZAfPuZE97+FrwM2Stk51fwrYW9I1kiZTeNCKpI+SdZusdQE2Iv4pIkakV70Hs9Rqj/XRm+0A9dui2+1gfY9vIbOmiezhID+uKv434ExJc4BnyLqQkHQqsDoibkpdFP8p6bCI+PfCtgeQPQMDsvHxx0iaS/bN9/SIeCXdpnm3pMeB3wCLIuL5Qh33kT3HYj/g7hTnDElTyQ7af0V2MPxuRDwu6RiyvvuKz9N190932mN99GY73E/9tniKEtrB2oufp2D9mqQtIz1iU9L3yB6w8/3C8pHAdyLiq53UcQdwYkR8IOm/A7dWvt1Lup3sIvgzpb6RHnI7WKOcFKxfk/T3ZI9R/CPwa7ID3/tV65wOTK7qH2+k7k2AkyJiSm/FWxa3gzXKScHMzHK+0GxmZjknBTMzyzkpmJlZzknBzMxyTgpmZpZzUjAzs5yTgpmZ5f4/xvprNS5p4e0AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.hist(const_list_df[const_list_df_noVPAS['max-min']!=0]['max-min'], bins=30)\n",
    "plt.title('Max-Min $log(K_H)$ per species')\n",
    "plt.xlabel('Max $log(K_H)$ - Min $log(K_H)$')\n",
    "plt.ylabel('Frequency')\n",
    "hist = plt.gcf()\n",
    "\n",
    "hist.savefig('range_per_species.png')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T17:41:31.487616Z",
     "start_time": "2020-01-26T17:41:30.955200Z"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEcCAYAAAAoSqjDAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3dfbxVZZ338c9XEPEpAUFDRNFiTPMB6OhQVmNijeIo1uRk4ySjFPkamnKcukVrGsepe/S+K8tpXhamE3hnCeYDk1QiZU1NPqAiqGgcFeUIwUkFVHxI+91/rGsvF5t9ztmHc9be5xy+79drv/Za17rWtX77Ovvs317XWnstRQRmZmYAOzU7ADMz6zucFMzMLOekYGZmOScFMzPLOSmYmVnOScHMzHJOCmZmlnNSMDOznJOC9QuSHpJ0XA/WXyNpQi/E8W+SzuuFdu6W9PaettNf9PTvZ43jpGA5SaslvSppZFX5MkkhaVyzthkRb4+IO7ZzG8OBNwMrexjrKOAs4NvFtlOcB6Z5Sbokva6jOmnuK8AlPYmnP+nJ388ay0nBqj0BfLQyI+kIYNd+vs0jgN9GxCs9bOdvgUUR8VKhbALwXEQ8KWl34AbgeOCYiHigk7YWAu+TNLqHMXVJ0uCyt2EDh5OCVbuW7NtwxXRgXrGCpNmSHpP0vKSHJX0wlb9F0rOSJqX5/ST9vo5hg3q2uVrSCYXpz0paLmmTpOslDe2k/SOBFWldSbpA0pOSNkqaL2mvtGwnSV9MQ01rJZ2S9mKGp3ZOAn5R1fYEYJmkA4BfAZuA4yNiQ2cvOCJeBu4FPtBRnfQ6L0x9/Jyk/6y8ztS3P5TULukJSZ+use4FkpYDL9ZKDGn50+nv+KikKb2w3bGSbkzLn5H0zUKbJxTqddhOR3FZg0SEH34QEQCrgROAR4FDgUHAGuBAIIBxqd7pwH5kXyo+ArwIjE7LPkE2TLMb8FPgK720zdXACYXpu1MMI9L2zu1kG3OAi9L0l4A7gNHAUOAW4Ktp2SXAL4ExwDDgN8CaQjvtwNFVbc8D/htYC5zfzf6+AvhaF33zIDA2vc5fp/h3IksoXwSGAAcDjwN/XrXusrTurjXaPiT1835pfhzwlp5sN/3tHgAuB3ZP/fvuGn+/DtvpLC4/GvPwnoLVUvnm/n7gEeDp4sKIWBARayPijxFxPbAKOCYtuyrN30X2wfv53thmDVekGJ4F/ovsG3tHjgCWS9oX+HvgryNiXWTf1m8AWtLxgn8AzomIpyNiI3AbaQ8jGQY8X9X2BODtwPKI+Fr1hiUtKB5QlnRf4Vv786nNznwzItak1/llsmG2o4FREXFJRLwaEY8DVwFnVK17RVr3Jbb1OrALcJiknSNidUQ81sPtHkOWqD8XES9GxMsR8asa2+6sna7ispI5KVgt1wJ/TTaGPq96oaSz0oHgjZI2AocDxQPFV6Wyf480ji/pTEkvpMePu7vNGn5XmN4C7FGrkiSlWFYA7wFWRMTaQpWRwDpgCtAaEa2FZSPYOik8B+xZaHsXsr2bvwQOVe2zkt4GtKb6g4DXI+K1tGxPYGPnL5M1heknyT50DwT2q/R/+htcBOzbybpbSa/zPOBiYIOkH0jar4fbHQs8WXh9HemwnTrispI5Kdg2IuJJsoO/U4Ebi8uUnWVzFfApYO+IGEY21KC0fA/g68DVwMWSRqQ2vxcRe6THSd3ZZg8dBPwxtT+KbMy/aBrZsYCRZENAldc5iOwYwvJC3eXAnxTmDyf7ZvvfwAeBf5V0fKGNIcDgeOMA90GkBJEcSjbc0pmxhekDUoxrgCciYljhsWdETK1at9ObpUTEdRHxbt4Yqrush9tdAxxQx4HtTtvpIi4rmZOCdWQG2QHTF6vKdyf7R20HkHQ22YdjxTeAeyPi48CtwLd6YZs9kR9kBu4B3pkOiO8h6RKyb7nXkB2XeJekt0p6E9l4/1vYek9hEfBnhfmJwIMR8VpE3AfMAuZLOigtPxQYJekOSXcAN5El0MpexjuAxV3EP0vS/im5XgRcT3Y8ZXM6ILurpEGSDpd0dL2dIukQScenOF4GXiJLcD3Z7t1ke12XStpd0lBJx9bYfIft1BGXlcxJwWqKiMciYmmN8oeBr5IdhF1PNl7/awBJ04ATgXNT9fOBSZLO7Mk2e+gI0rf91PaXyfYM2sg+tI+PiC0RsQT4PnA/sJTsIO0WsuMbFfOAqZIqp8tOSPUq8c8DrgNuVnZ66uHAlRFxXEQcB8wHHkrVTwXuqBrKquU6smMbj6fHlyLideCUtP0ngN8D3wH2qr9b2AW4NK37O2Afsg//7d5uYflbgafI+vgj1Rvuop2u4rKSKcK34zSrJulc4OSIOKWq/H8DGyLi63W08W/AfRGxIM0vAC6MiFZJdwEzIuLBTtZfDXw8Im7vwUvptmZt1/oG/6jFDJA0mWzoYw3ZQedLyL7NbyUiuvOt9e3AdwvzbyX71k1E/On2xmpWJicFs8xEsmMgOwO/Bf42Iu7sSYMRcWrV/MSetGfWCB4+MjOznA80m5lZrl8PH40cOTLGjRvX7DDMzPqVe++99/cRMarWsn6dFMaNG8fSpb19BqOZ2cAm6cmOlnn4yMzMck4KZmaWc1IwM7Ock4KZmeWcFMzMLOekYGZmOScFMzPLOSmYmVnOScHMzHL9+hfNPTFu9q111Vt96cklR2Jm1nd4T8HMzHJOCmZmlnNSMDOznJOCmZnlnBTMzCznpGBmZjknBTMzyzkpmJlZzknBzMxyTgpmZpZzUjAzs1xpSUHSIZKWFR6bJZ0naYSkxZJWpefhqb4kXSGpVdJySZPKis3MzGorLSlExKMRMSEiJgDvALYANwGzgSURMR5YkuYBTgLGp8dM4MqyYjMzs9oaNXw0BXgsIp4EpgFzU/lc4LQ0PQ2YF5k7gWGSRjcoPjMzo3FJ4Qzg+2l634hYB5Ce90nlY4A1hXXaUpmZmTVI6UlB0hDgVGBBV1VrlEWN9mZKWippaXt7e2+EaGZmSSP2FE4C7ouI9Wl+fWVYKD1vSOVtwNjCevsDa6sbi4g5EdESES2jRo0qMWwzsx1PI5LCR3lj6AhgITA9TU8HbimUn5XOQpoMbKoMM5mZWWOUejtOSbsB7wc+WSi+FJgvaQbwFHB6Kl8ETAVayc5UOrvM2MzMbFulJoWI2ALsXVX2DNnZSNV1A5hVZjxmZtY5/6LZzMxyTgpmZpZzUjAzs5yTgpmZ5ZwUzMws56RgZmY5JwUzM8s5KZiZWc5JwczMck4KZmaWc1IwM7Ock4KZmeWcFMzMLOekYGZmOScFMzPLOSmYmVnOScHMzHJOCmZmlnNSMDOzXKlJQdIwSTdIekTSSknvlDRC0mJJq9Lz8FRXkq6Q1CppuaRJZcZmZmbbKntP4RvATyLibcBRwEpgNrAkIsYDS9I8wEnA+PSYCVxZcmxmZlaltKQg6U3Ae4GrASLi1YjYCEwD5qZqc4HT0vQ0YF5k7gSGSRpdVnxmZratMvcUDgbagf+UdL+k70jaHdg3ItYBpOd9Uv0xwJrC+m2pbCuSZkpaKmlpe3t7ieGbme14ykwKg4FJwJURMRF4kTeGimpRjbLYpiBiTkS0RETLqFGjeidSMzMDyk0KbUBbRNyV5m8gSxLrK8NC6XlDof7Ywvr7A2tLjM/MzKqUlhQi4nfAGkmHpKIpwMPAQmB6KpsO3JKmFwJnpbOQJgObKsNMZmbWGINLbv/vge9JGgI8DpxNlojmS5oBPAWcnuouAqYCrcCWVNfMzBqo1KQQEcuAlhqLptSoG8CsMuMxM7PO+RfNZmaWc1IwM7Ock4KZmeWcFMzMLOekYGZmOScFMzPLOSmYmVnOScHMzHJOCmZmlnNSMDOznJOCmZnlnBTMzCznpGBmZjknBTMzyzkpmJlZzknBzMxyTgpmZpZzUjAzs1ypSUHSakkrJC2TtDSVjZC0WNKq9Dw8lUvSFZJaJS2XNKnM2MzMbFuN2FN4X0RMiIjKvZpnA0siYjywJM0DnASMT4+ZwJUNiM3MzAqaMXw0DZibpucCpxXK50XmTmCYpNFNiM/MbIdVdlII4DZJ90qamcr2jYh1AOl5n1Q+BlhTWLctlW1F0kxJSyUtbW9vLzF0M7Mdz+CS2z82ItZK2gdYLOmRTuqqRllsUxAxB5gD0NLSss1yMzPbfqXuKUTE2vS8AbgJOAZYXxkWSs8bUvU2YGxh9f2BtWXGZ2ZmWystKUjaXdKelWngA8CDwEJgeqo2HbglTS8EzkpnIU0GNlWGmczMrDHqGj6SdHhEPNjNtvcFbpJU2c51EfETSfcA8yXNAJ4CTk/1FwFTgVZgC3B2N7dnZmY9VO8xhW9JGgJ8l+zDfWNXK0TE48BRNcqfAabUKA9gVp3xmJlZCeoaPoqIdwNnko35L5V0naT3lxqZmZk1XN3HFCJiFfAF4ALgz4ArJD0i6UNlBWdmZo1VV1KQdKSky4GVwPHAKRFxaJq+vMT4zMysgeo9pvBN4Crgooh4qVKYfoPwhVIiMzOzhqs3KUwFXoqI1wEk7QQMjYgtEXFtadGZmVlD1XtM4XZg18L8bqnMzMwGkHqTwtCIeKEyk6Z3KyckMzNrlnqTwovF+xtIegfwUif1zcysH6r3mMJ5wAJJlWsRjQY+Uk5IZmbWLHUlhYi4R9LbgEPIrmb6SET8odTIzMys4bpz6eyjgXFpnYmSiIh5pURlZmZNUe8F8a4F3gIsA15PxQE4KZiZDSD17im0AIeli9aZmdkAVe/ZRw8Cby4zEDMza7569xRGAg9Luht4pVIYEaeWEpWZmTVFvUnh4jKDMDOzvqHeU1J/IelAYHxE3C5pN2BQuaGZmVmj1Xvp7E8ANwDfTkVjgJvLCsrMzJqj3gPNs4Bjgc2Q33Bnn3pWlDRI0v2SfpTmD5J0l6RVkq5Pt/lE0i5pvjUtH9fdF2NmZj1Tb1J4JSJercxIGkz2O4V6fIbs5jwVlwGXR8R44DlgRiqfATwXEW8lu3HPZXW2b2ZmvaTepPALSRcBu6Z7My8A/qurlSTtD5wMfCfNi+xubTekKnOB09L0tDRPWj4l1TczswapNynMBtqBFcAngUVk92vuyteB/wX8Mc3vDWyMiNfSfBvZ8QnS8xqAtHxTqm9mZg1S79lHfyS7HedV9TYs6S+ADRFxr6TjKsW1mq9jWbHdmcBMgAMOOKDecMzMrA71XvvoCWp8QEfEwZ2sdixwqqSpwFDgTWR7DsMkDU57A/sDlctxtwFjgbZ0zGIv4Nka25wDzAFoaWnxZTfMzHpRd659VDEUOB0Y0dkKEXEhcCFA2lP4bEScKWkB8GHgB8B04Ja0ysI0/5u0/Ge+1pKZWWPVdUwhIp4pPJ6OiK+THTDeHhcA50tqJTtmcHUqvxrYO5WfT3Ycw8zMGqje4aNJhdmdyPYc9qx3IxFxB3BHmn4cOKZGnZfJ9kDMzKxJ6h0++mph+jVgNfBXvR6NmZk1Vb1nH72v7EDMzKz56h0+Or+z5RHxtd4Jx8zMmqk7Zx8dTXaGEMApwC9JPzYzM7OBoTs32ZkUEc8DSLoYWBARHy8rMDMza7x6L3NxAPBqYf5VYFyvR2NmZk1V757CtcDdkm4i+2XzB4F5pUVlZmZNUe/ZR1+W9GPgPano7Ii4v7ywzMysGeodPgLYDdgcEd8guz7RQSXFZGZmTVLv7Tj/mezyFBemop2B/1dWUGZm1hz17il8EDgVeBEgItbSjctcmJlZ/1BvUng1XbE0ACTtXl5IZmbWLPUmhfmSvk12L4RPALfTjRvumJlZ/1Dv2UdfSfdm3gwcAnwxIhaXGpmZmTVcl0lB0iDgpxFxAuBEYGY2gHU5fBQRrwNbJO3VgHjMzKyJ6v1F88vACkmLSWcgAUTEp0uJyszMmqLepHBrepiZ2QDWaVKQdEBEPBURcxsVkJmZNU9XxxRurkxI+mF3GpY0VNLdkh6Q9JCkf0nlB0m6S9IqSddLGpLKd0nzrWn5uG6+FjMz66GukoIK0wd3s+1XgOMj4ihgAnCipMnAZcDlETEeeA6YkerPAJ6LiLcCl6d6ZmbWQF0lhehgukuReSHN7pweARwP3JDK5wKnpelpaZ60fIqkYlIyM7OSdZUUjpK0WdLzwJFperOk5yVt7qpxSYMkLQM2kP3G4TFgY0S8lqq0AWPS9BjS7T3T8k3A3jXanClpqaSl7e3t9bxGMzOrU6cHmiNiUE8aT79xmCBpGHATcGitaum51l7BNnsnETEHmAPQ0tLSrb0XMzPrXHfup7DdImIjcAcwmez6SZVktD+wNk23AWMB0vK9gGcbEZ+ZmWVKSwqSRqU9BCTtCpwArAR+Dnw4VZsO3JKmF6Z50vKfpSuzmplZg9T747XtMRqYm66dtBMwPyJ+JOlh4AeSvgTcD1yd6l8NXCuplWwP4YwSYzMzsxpKSwoRsRyYWKP8ceCYGuUvA6eXFY+ZmXWtIccUzMysf3BSMDOznJOCmZnlnBTMzCznpGBmZjknBTMzyzkpmJlZzknBzMxyTgpmZpZzUjAzs5yTgpmZ5ZwUzMws56RgZmY5JwUzM8s5KZiZWc5JwczMck4KZmaWc1IwM7NcaUlB0lhJP5e0UtJDkj6TykdIWixpVXoensol6QpJrZKWS5pUVmxmZlZbmXsKrwH/GBGHApOBWZIOA2YDSyJiPLAkzQOcBIxPj5nAlSXGZmZmNZSWFCJiXUTcl6afB1YCY4BpwNxUbS5wWpqeBsyLzJ3AMEmjy4rPzMy21ZBjCpLGAROBu4B9I2IdZIkD2CdVGwOsKazWlsqq25opaamkpe3t7WWGbWa2wyk9KUjaA/ghcF5EbO6sao2y2KYgYk5EtEREy6hRo3orTDMzo+SkIGlnsoTwvYi4MRWvrwwLpecNqbwNGFtYfX9gbZnxmZnZ1so8+0jA1cDKiPhaYdFCYHqang7cUig/K52FNBnYVBlmMjOzxhhcYtvHAh8DVkhalsouAi4F5kuaATwFnJ6WLQKmAq3AFuDsEmMzM7MaSksKEfErah8nAJhSo34As8qKx8zMuuZfNJuZWc5JwczMck4KZmaWc1IwM7Ock4KZmeWcFMzMLOekYGZmOScFMzPLOSmYmVnOScHMzHJOCmZmlnNSMDOzXJlXSR0Qxs2+te66qy89ucRIzMzK5z0FMzPLOSmYmVnOScHMzHJOCmZmlnNSMDOzXGlJQdI1kjZIerBQNkLSYkmr0vPwVC5JV0hqlbRc0qSy4jIzs46VuafwXeDEqrLZwJKIGA8sSfMAJwHj02MmcGWJcZmZWQdKSwoR8Uvg2ariacDcND0XOK1QPi8ydwLDJI0uKzYzM6ut0ccU9o2IdQDpeZ9UPgZYU6jXlsq2IWmmpKWSlra3t5carJnZjqavHGhWjbKoVTEi5kRES0S0jBo1quSwzMx2LI1OCusrw0LpeUMqbwPGFurtD6xtcGxmZju8RieFhcD0ND0duKVQflY6C2kysKkyzGRmZo1T2gXxJH0fOA4YKakN+GfgUmC+pBnAU8DpqfoiYCrQCmwBzi4rLjMz61hpSSEiPtrBoik16gYwq6xYzMysPn3lQLOZmfUBTgpmZpZzUjAzs5yTgpmZ5Xw7zl5U7607fdtOM+urvKdgZmY57yk0gfcozKyv8p6CmZnlvKfQh3mPwswazXsKZmaWc1IwM7Ock4KZmeWcFMzMLOcDzQOAD0ibWW/xnoKZmeW8p7AD8R6FmXXFewpmZpZzUjAzs5yTgpmZ5fpUUpB0oqRHJbVKmt3seMzMdjR95kCzpEHAfwDvB9qAeyQtjIiHmxvZjqfeA9LdUe/B694+GF7Ga6mHD9Zbf9VnkgJwDNAaEY8DSPoBMA1wUhgAevvDuVkf9mVo1llhPhvNaulLSWEMsKYw3wb8aXUlSTOBmWn2BUmPbuf2RgK/3851m6E/xdufYoUS4tVlvdnaVm02vG97+Fr603uhP8UKPYv3wI4W9KWkoBplsU1BxBxgTo83Ji2NiJaettMo/Sne/hQr9K94+1Os0L/i7U+xQnnx9qUDzW3A2ML8/sDaJsViZrZD6ktJ4R5gvKSDJA0BzgAWNjkmM7MdSp8ZPoqI1yR9CvgpMAi4JiIeKnGTPR6CarD+FG9/ihX6V7z9KVboX/H2p1ihpHgVsc2wvZmZ7aD60vCRmZk1mZOCmZnlBnxS6OrSGZJ2kXR9Wn6XpHGNjzKPZaykn0taKekhSZ+pUec4SZskLUuPLzYj1hTLakkrUhxLayyXpCtS3y6XNKkZcaZYDin02TJJmyWdV1WnaX0r6RpJGyQ9WCgbIWmxpFXpeXgH605PdVZJmt7EeP+vpEfS3/omScM6WLfT902DYr1Y0tOFv/XUDtZt+KV3Ooj3+kKsqyUt62DdnvdtRAzYB9kB68eAg4EhwAPAYVV1/g74Vpo+A7i+ifGOBial6T2B39aI9zjgR83u2xTLamBkJ8unAj8m+w3KZOCuZsdceF/8Djiwr/Qt8F5gEvBgoez/ALPT9GzgshrrjQAeT8/D0/TwJsX7AWBwmr6sVrz1vG8aFOvFwGfreJ90+vnRqHirln8V+GJZfTvQ9xTyS2dExKtA5dIZRdOAuWn6BmCKpFo/pCtdRKyLiPvS9PPASrJfevdX04B5kbkTGCZpdLODAqYAj0XEk80OpCIifgk8W1VcfG/OBU6rseqfA4sj4tmIeA5YDJxYWqBJrXgj4raIeC3N3kn2W6Om66Bv61HP50ev6yze9Nn0V8D3y9r+QE8KtS6dUf0hm9dJb+hNwN4Nia4TaRhrInBXjcXvlPSApB9LentDA9taALdJujddfqRaPf3fDGfQ8T9VX+lbgH0jYh1kXxiAfWrU6at9fA7ZXmItXb1vGuVTaajrmg6G5vpi374HWB8RqzpY3uO+HehJoZ5LZ9R1eY1GkrQH8EPgvIjYXLX4PrJhj6OAfwdubnR8BcdGxCTgJGCWpPdWLe+LfTsEOBVYUGNxX+rbevXFPv488BrwvQ6qdPW+aYQrgbcAE4B1ZEMy1fpc3wIfpfO9hB737UBPCvVcOiOvI2kwsBfbt6vZKyTtTJYQvhcRN1Yvj4jNEfFCml4E7CxpZIPDrMSyNj1vAG4i290u6ouXLjkJuC8i1lcv6Et9m6yvDLel5w016vSpPk4Huv8CODPSIHe1Ot43pYuI9RHxekT8Ebiqgxj6Wt8OBj4EXN9Rnd7o24GeFOq5dMZCoHLGxoeBn3X0Zi5bGi+8GlgZEV/roM6bK8c8JB1D9jd8pnFR5nHsLmnPyjTZQcYHq6otBM5KZyFNBjZVhkOaqMNvWn2lbwuK783pwC016vwU+ICk4WkI5AOprOEknQhcAJwaEVs6qFPP+6Z0Vce2PthBDH3t0jsnAI9ERFuthb3Wt2UfSW/2g+wMmN+SnUXw+VR2CdkbF2Ao2VBCK3A3cHATY3032e7pcmBZekwFzgXOTXU+BTxEdibEncC7mhTrwSmGB1I8lb4txiqyGyc9BqwAWpr8XtiN7EN+r0JZn+hbskS1DvgD2TfUGWTHtpYAq9LziFS3BfhOYd1z0vu3FTi7ifG2ko3BV967lbP69gMWdfa+aUKs16b35HKyD/rR1bGm+W0+P5oRbyr/buW9Wqjb633ry1yYmVluoA8fmZlZNzgpmJlZzknBzMxyTgpmZpZzUjAzs5yTgpmZ5ZwUzMws56RgpZMUkq4tzA+W1C7pR73Q9gnFtrdj/V0l/ULSIEmflPStVL6zpGslzU2XHqm17hBJv0yXH+jONrvsD0n/0802e60f0nzdfbG9/WB9k5OCNcKLwOGSdk3z7wee7qW2jyL7Bef2Oge4MSJeB44Elkt6E9kVPp+KiOkR8YdaK0Z2OeUlwEe6uc0u+yMi3tXNNnuzH6AbfdGDfrA+yEnBGuXHwMlpeqvrD0m6OV3q96HK5X4lHZ0uazw0XdPlIUmH12j3KLJLKiDpbekb60OSbq9czE7Soal8uaTPSWotrH8mb1xT6AiyiyHeASyIiM/X8bpuTm10V4f9kWJ+IT2PU3YnvqvS67qtkEyKerMfoPt9sb39YH1NI67l4ceO/QBeIPvmeQPZtaaWUbjLGW9c02dXsgt47Z3mvwR8hez6SRd20PYDwChgF7LrvUxM5RcAXwYGk10Su1J+JXBzmh4C/K7Q1nNAO3BKN17bIKC9N/ujUic9jyO7DPWEND8f+Jsy+2F7+mJ7+sGPvvnwnoI1REQsJ/uA+yiwqGrxpyVVLkI3Fhifyi8hG1ppIbs15VbS+PabIqKd7K5kv4qI+9Pih8luSvMh4IGq8sowy0hgY2prLNmH9Qqy26JWb+uWwvSCyth7ZMMtr1auTlmvLvqj2hMRUbkn771pvWJsvdYPqb0O+6K3+8H6HicFa6SFZN/8i0NHx5FdEvidkd3c5n6yb8+Q3Xd4D7L7VQ9lW4eR3bK0Mr2isOwIsg++I0nDKsnhhfmXCu0eSfYheSbweUkTCzGOJbtqZcVO8cbYO2Tfzl8uBiZplt640fp+NWKHGv3RgVcK06+Tfesv6s1+gA76Ynv6wfofJwVrpGuASyKi+KG1F/BcRGyR9DZgcmHZHOCfyO7gdVmN9vJxdLIDtYcBSDoY+Bgwj+xS2X+SyicAf0P6hhzZPY0HSRpK9uG5IrL7PXwcuF7SXqntdwCHSfqWpLkUbrQiaW+yYZOtDsBGxH9ExIT06OjGLLX6Y3v0Zj9Ax33R7X6w/senkFnDRHZzkG9UFf8EOFfScuBRsiEkJJ0FvBYR16Uhiv+RdHxE/Kyw7lFk98CA7Pr4UyWtIPvme05EPJNO07xV0j3Ab4DVEfF4oY3byO5jcQRwa4pzsaT5ZB/af0n2YfiPEXGPpJPJxu4r3kfXwz/d6Y/t0Zv9cDsd98XDlNAP1rf4fgo2oEnaI9ItNiV9juwGO18oLJ8InB8RH+ukjZuAj0TEq5L+Bbih8u1e0o1kB8EfLfWF9JD7werlpGADmqR/IruN4h+AX5N98L1SVeccYG7V+Hg9bQ8BzoiIeb0Vb1ncD1YvJwUzM8v5QLOZmeWcFMzMLOekYGZmOScFMzPLOSmYmVnOScHMzHJOCmZmlrPZzSQAAAAGSURBVPv/0j0/NStE3lgAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.hist(const_list_df[const_list_df_VPAS['max-min']!=0]['max-min'], bins=30)\n",
    "plt.title('Max-Min $log(K_H)$ per species')\n",
    "plt.xlabel('Max $log(K_H)$ - Min $log(K_H)$')\n",
    "plt.ylabel('Frequency')\n",
    "hist = plt.gcf()\n",
    "\n",
    "hist.savefig('range_per_species.png')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T17:47:25.807667Z",
     "start_time": "2020-01-26T17:47:24.388897Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "VPAS\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Num</th>\n",
       "      <th>min</th>\n",
       "      <th>max</th>\n",
       "      <th>mean</th>\n",
       "      <th>median</th>\n",
       "      <th>max-min</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Species</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>hexanedioic acid</th>\n",
       "      <td>4</td>\n",
       "      <td>-18.005165</td>\n",
       "      <td>-2.397895</td>\n",
       "      <td>-11.765305</td>\n",
       "      <td>-13.329080</td>\n",
       "      <td>15.607270</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>(2,4-dichlorophenoxy)-ethanoic acid</th>\n",
       "      <td>11</td>\n",
       "      <td>-15.520259</td>\n",
       "      <td>1.966113</td>\n",
       "      <td>-6.163701</td>\n",
       "      <td>-6.877296</td>\n",
       "      <td>17.486372</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>pebulate</th>\n",
       "      <td>4</td>\n",
       "      <td>-23.858760</td>\n",
       "      <td>3.194183</td>\n",
       "      <td>-4.310896</td>\n",
       "      <td>1.710496</td>\n",
       "      <td>27.052943</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>iodofenphos</th>\n",
       "      <td>2</td>\n",
       "      <td>-23.858760</td>\n",
       "      <td>-0.788457</td>\n",
       "      <td>-12.323609</td>\n",
       "      <td>-12.323609</td>\n",
       "      <td>23.070303</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                                     Num        min       max       mean  \\\n",
       "Species                                                                    \n",
       "hexanedioic acid                       4 -18.005165 -2.397895 -11.765305   \n",
       "(2,4-dichlorophenoxy)-ethanoic acid   11 -15.520259  1.966113  -6.163701   \n",
       "pebulate                               4 -23.858760  3.194183  -4.310896   \n",
       "iodofenphos                            2 -23.858760 -0.788457 -12.323609   \n",
       "\n",
       "                                        median    max-min  \n",
       "Species                                                    \n",
       "hexanedioic acid                    -13.329080  15.607270  \n",
       "(2,4-dichlorophenoxy)-ethanoic acid  -6.877296  17.486372  \n",
       "pebulate                              1.710496  27.052943  \n",
       "iodofenphos                         -12.323609  23.070303  "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3018 species with 0 constant(s)\n",
      "852 species with 1 constant(s)\n",
      "267 species with 2 constant(s)\n",
      "113 species with 3 constant(s)\n",
      "85 species with 4 constant(s)\n",
      "65 species with 5 constant(s)\n",
      "38 species with 6 constant(s)\n",
      "26 species with 7 constant(s)\n",
      "21 species with 8 constant(s)\n",
      "24 species with 9 constant(s)\n",
      "22 species with 10 constant(s)\n",
      "12 species with 11 constant(s)\n",
      "10 species with 12 constant(s)\n",
      "6 species with 13 constant(s)\n",
      "9 species with 14 constant(s)\n",
      "4 species with 15 constant(s)\n",
      "6 species with 16 constant(s)\n",
      "7 species with 17 constant(s)\n",
      "3 species with 18 constant(s)\n",
      "4 species with 19 constant(s)\n",
      "5 species with 20 constant(s)\n",
      "3 species with 21 constant(s)\n",
      "1 species with 22 constant(s)\n",
      "5 species with 23 constant(s)\n",
      "5 species with 24 constant(s)\n",
      "2 species with 25 constant(s)\n",
      "2 species with 26 constant(s)\n",
      "1 species with 28 constant(s)\n",
      "3 species with 29 constant(s)\n",
      "1 species with 37 constant(s)\n",
      "1 species with 39 constant(s)\n",
      "2 species with 40 constant(s)\n",
      "2 species with 42 constant(s)\n",
      "1 species with 45 constant(s)\n",
      "2 species with 49 constant(s)\n",
      "1 species with 60 constant(s)\n",
      "1 species with 61 constant(s)\n",
      "1 species with 74 constant(s)\n",
      "1 species with 83 constant(s)\n",
      "This accounts for all species\n",
      "noVPAS\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Num</th>\n",
       "      <th>min</th>\n",
       "      <th>max</th>\n",
       "      <th>mean</th>\n",
       "      <th>median</th>\n",
       "      <th>max-min</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Species</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>(2,4-dichlorophenoxy)-ethanoic acid</th>\n",
       "      <td>5</td>\n",
       "      <td>-15.520259</td>\n",
       "      <td>1.966113</td>\n",
       "      <td>-6.359637</td>\n",
       "      <td>-6.877296</td>\n",
       "      <td>17.486372</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                                     Num        min       max      mean  \\\n",
       "Species                                                                   \n",
       "(2,4-dichlorophenoxy)-ethanoic acid    5 -15.520259  1.966113 -6.359637   \n",
       "\n",
       "                                       median    max-min  \n",
       "Species                                                   \n",
       "(2,4-dichlorophenoxy)-ethanoic acid -6.877296  17.486372  "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3018 species with 0 constant(s)\n",
      "852 species with 1 constant(s)\n",
      "267 species with 2 constant(s)\n",
      "113 species with 3 constant(s)\n",
      "85 species with 4 constant(s)\n",
      "65 species with 5 constant(s)\n",
      "38 species with 6 constant(s)\n",
      "26 species with 7 constant(s)\n",
      "21 species with 8 constant(s)\n",
      "24 species with 9 constant(s)\n",
      "22 species with 10 constant(s)\n",
      "12 species with 11 constant(s)\n",
      "10 species with 12 constant(s)\n",
      "6 species with 13 constant(s)\n",
      "9 species with 14 constant(s)\n",
      "4 species with 15 constant(s)\n",
      "6 species with 16 constant(s)\n",
      "7 species with 17 constant(s)\n",
      "3 species with 18 constant(s)\n",
      "4 species with 19 constant(s)\n",
      "5 species with 20 constant(s)\n",
      "3 species with 21 constant(s)\n",
      "1 species with 22 constant(s)\n",
      "5 species with 23 constant(s)\n",
      "5 species with 24 constant(s)\n",
      "2 species with 25 constant(s)\n",
      "2 species with 26 constant(s)\n",
      "1 species with 28 constant(s)\n",
      "3 species with 29 constant(s)\n",
      "1 species with 37 constant(s)\n",
      "1 species with 39 constant(s)\n",
      "2 species with 40 constant(s)\n",
      "2 species with 42 constant(s)\n",
      "1 species with 45 constant(s)\n",
      "2 species with 49 constant(s)\n",
      "1 species with 60 constant(s)\n",
      "1 species with 61 constant(s)\n",
      "1 species with 74 constant(s)\n",
      "1 species with 83 constant(s)\n",
      "This accounts for all species\n"
     ]
    }
   ],
   "source": [
    "for dataset in ['VPAS', 'noVPAS']:\n",
    "    print(dataset)\n",
    "    display(eval('const_list_df_%s[const_list_df_%s[\"max-min\"] > 15]'%(dataset, dataset))) #these vary by more than 15 log units\n",
    "    num_species = 0\n",
    "    for num in range(0, 100):\n",
    "        if len(const_list_df[const_list_df['Num'] == num] != 0):\n",
    "            print(len(const_list_df[const_list_df['Num'] == num]), \"species with %i constant(s)\" % num)\n",
    "        num_species += len(const_list_df[const_list_df['Num'] == num])\n",
    "    if num_species == len(const_list_df):\n",
    "        print(\"This accounts for all species\")\n",
    "    else:\n",
    "        print(\"Not all species accounted for - there are %i species missing from this count\" % (len(const_list_df) - num_species))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Marking species with varying HLCs\n",
    "One-hot encoding to mark these species - not to be used as a feature to train models with."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T18:09:14.951360Z",
     "start_time": "2020-01-26T18:08:47.296834Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "vary_list_0_VPAS 4632 constants total. 1179 varying constants\n",
      "vary_list_1_VPAS 4632 constants total. 665 varying constants\n",
      "vary_list_2_VPAS 4632 constants total. 465 varying constants\n",
      "vary_list_3_VPAS 4632 constants total. 296 varying constants\n",
      "vary_list_4_VPAS 4632 constants total. 182 varying constants\n",
      "vary_list_5_VPAS 4632 constants total. 110 varying constants\n",
      "vary_list_6_VPAS 4632 constants total. 82 varying constants\n",
      "vary_list_7_VPAS 4632 constants total. 60 varying constants\n",
      "vary_list_8_VPAS 4632 constants total. 40 varying constants\n",
      "vary_list_9_VPAS 4632 constants total. 29 varying constants\n",
      "vary_list_0_noVPAS 4632 constants total. 687 varying constants\n",
      "vary_list_1_noVPAS 4632 constants total. 318 varying constants\n",
      "vary_list_2_noVPAS 4632 constants total. 200 varying constants\n",
      "vary_list_3_noVPAS 4632 constants total. 110 varying constants\n",
      "vary_list_4_noVPAS 4632 constants total. 60 varying constants\n",
      "vary_list_5_noVPAS 4632 constants total. 44 varying constants\n",
      "vary_list_6_noVPAS 4632 constants total. 34 varying constants\n",
      "vary_list_7_noVPAS 4632 constants total. 27 varying constants\n",
      "vary_list_8_noVPAS 4632 constants total. 17 varying constants\n",
      "vary_list_9_noVPAS 4632 constants total. 9 varying constants\n"
     ]
    }
   ],
   "source": [
    "for dataset in ['VPAS', 'noVPAS']:\n",
    "    const_list_df = eval('const_list_df_%s'%dataset)\n",
    "    for vary_by in range(0, 10):\n",
    "        exec('vary_list_%s_%s = []'% (str(vary_by), dataset)) # create many empty lists at once\n",
    "        for i in range(0, len(eval('const_list_df_%s'%dataset))):\n",
    "            exec('max_min = const_list_df_%s.iloc[i][\"max-min\"]'%dataset)\n",
    "            # one hot encoding - for vary_list_x, append 1 if range is larger than x and 0 if not\n",
    "            if max_min > vary_by:\n",
    "                exec('vary_list_%i_%s.append(1)'%(vary_by, dataset))\n",
    "            else:\n",
    "                exec('vary_list_%i_%s.append(0)'%(vary_by, dataset))\n",
    "        exec('print(\"vary_list_%i_%s\", len(vary_list_%i_%s), \"constants total.\", vary_list_%i_%s.count(1), \"varying constants\")'\n",
    "             %(vary_by, dataset, vary_by, dataset, vary_by, dataset)) #ensure all constants are processed properly"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Picking constants to use for each species\n",
    "Henry's Law constants vary by a fair amount - values between species can differ by orders of magnitude. Perhaps its an idea to take the log of the constants as we get them and then work out how to split from there?\n",
    "> Changed the code in `Getting data from the supplement` to use the log value from the start\n",
    "\n",
    "Now to figure out which values to use. I think a difference of one log unit is a reasonable threshold to drop, but will check this with Jeremy.\n",
    "\n",
    "> For now, I will not drop any constants, and will instead use the mean for each species.\n",
    "\n",
    "I'm going to create several columns in this dataframe - each holding a different choice of constant."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T18:23:44.694797Z",
     "start_time": "2020-01-26T18:23:44.689784Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "4632"
      ]
     },
     "execution_count": 71,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "len(constant_list_VPAS)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T18:25:11.805402Z",
     "start_time": "2020-01-26T18:25:10.957148Z"
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\jtd1g16\\AppData\\Local\\Continuum\\anaconda3\\lib\\site-packages\\numpy\\core\\fromnumeric.py:3257: RuntimeWarning: Mean of empty slice.\n",
      "  out=out, **kwargs)\n",
      "C:\\Users\\jtd1g16\\AppData\\Local\\Continuum\\anaconda3\\lib\\site-packages\\numpy\\core\\_methods.py:161: RuntimeWarning: invalid value encountered in double_scalars\n",
      "  ret = ret.dtype.type(ret / rcount)\n"
     ]
    }
   ],
   "source": [
    "for dataset in ['VPAS', 'noVPAS']:\n",
    "    exec('mean_constant_list_%s = []'%dataset) #take mean for each list of constants\n",
    "    for i in range(0, len(eval('constant_list_%s'%dataset))):\n",
    "        if type(constant_list[i]) == list:\n",
    "            eval('mean_constant_list_%s.append(np.mean(constant_list_%s[i]))'%(dataset, dataset))\n",
    "        else:\n",
    "            eval('mean_constant_list_%s.append(np.nan)'%dataset)\n",
    "\n",
    "    exec('first_constant_list_%s = []'%dataset) #use first constant in each case\n",
    "    for i in range(0, len(eval('constant_list_%s'%dataset))):\n",
    "        if type(eval('constant_list_%s'%dataset)[i]) == list:\n",
    "            if len(eval('constant_list_%s'%dataset)[i]) > 0:\n",
    "                eval('first_constant_list_%s.append(np.mean(constant_list_%s[i]))'%(dataset, dataset))\n",
    "            else:\n",
    "                eval('first_constant_list_%s.append(np.nan)'%dataset)\n",
    "        else:\n",
    "            eval('first_constant_list_%s.append(np.nan)'%dataset)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T18:25:17.916651Z",
     "start_time": "2020-01-26T18:25:17.909632Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "4632"
      ]
     },
     "execution_count": 73,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "len(first_constant_list_VPAS)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Construct dataframe of SMILES, names and constants for all species"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-30T09:44:09.123358Z",
     "start_time": "2020-01-30T09:44:09.062244Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0         O=O\n",
       "1     O=[O]=O\n",
       "2      [H][H]\n",
       "3         [H]\n",
       "4    [2H][2H]\n",
       "5       [OH-]\n",
       "6       O[O-]\n",
       "7          OO\n",
       "8         N#N\n",
       "9           N\n",
       "Name: 1, dtype: object"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "if not 'smiles_list' in dir():\n",
    "    smiles_list = pd.read_csv('../smiles_constants.csv').iloc[:, 2] #already constructed smiles list from CAS numbers\n",
    "smiles_list[0:10]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 88,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T18:39:45.937506Z",
     "start_time": "2020-01-26T18:39:45.645732Z"
    }
   },
   "outputs": [],
   "source": [
    "for dataset in ['VPAS', 'noVPAS']:\n",
    "    smiles_constant_species = pd.concat([pd.Series(species_list), #creating dataframe to hold all info found so far\n",
    "                                         pd.Series(smiles_list),\n",
    "                                         pd.Series(eval('mean_constant_list_%s'%dataset)), #kept first Kh and mean in each case\n",
    "                                         pd.Series(eval('first_constant_list_%s'%dataset)), #will train on both sets and see which is better\n",
    "                                        pd.Series(eval('vary_list_1_%s'%dataset)), #one-hot encoded variation in constants\n",
    "                                        pd.Series(eval('vary_list_2_%s'%dataset)), #to be used for highlighting after training\n",
    "                                        pd.Series(eval('vary_list_3_%s'%dataset)), # /not/ to be used as a feature\n",
    "                                        pd.Series(eval('vary_list_4_%s'%dataset)), #perhaps could be used to weight certain species?\n",
    "                                        pd.Series(eval('vary_list_5_%s'%dataset)),\n",
    "                                        pd.Series(eval('vary_list_6_%s'%dataset)),\n",
    "                                        pd.Series(eval('vary_list_7_%s'%dataset)),\n",
    "                                        pd.Series(eval('vary_list_8_%s'%dataset)),\n",
    "                                        pd.Series(eval('vary_list_9_%s'%dataset))],\n",
    "                                        axis=1)\n",
    "    smiles_constant_species.to_csv('smiles_constants_kh_noVPAS.csv') #export so all the work isn't lost\n",
    "    smiles_constant_species.set_index(0, inplace=True)\n",
    "    smiles_constant_species.rename({0+1:'SMILES', 1+1:'Kh_mean', 2+1:'Kh_first', 3+1:'Varying_1',\n",
    "                                   4+1:'Varying_2', 5+1:'Varying_3', 6+1:'Varying_4', 7+1:'Varying_5',\n",
    "                                   8+1:'Varying_6', 9+1:'Varying_7', 10+1:'Varying_8', 11+1:'Varying_9'}, axis=1, inplace=True)\n",
    "    exec('smiles_constant_species_%s = smiles_constant_species'%dataset)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 90,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T18:40:59.294546Z",
     "start_time": "2020-01-26T18:40:59.277500Z"
    }
   },
   "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>SMILES</th>\n",
       "      <th>Kh_mean</th>\n",
       "      <th>Kh_first</th>\n",
       "      <th>Varying_1</th>\n",
       "      <th>Varying_2</th>\n",
       "      <th>Varying_3</th>\n",
       "      <th>Varying_4</th>\n",
       "      <th>Varying_5</th>\n",
       "      <th>Varying_6</th>\n",
       "      <th>Varying_7</th>\n",
       "      <th>Varying_8</th>\n",
       "      <th>Varying_9</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>oxygen</th>\n",
       "      <td>O=O</td>\n",
       "      <td>11.267713</td>\n",
       "      <td>11.267713</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>ozone</th>\n",
       "      <td>O=[O]=O</td>\n",
       "      <td>9.508258</td>\n",
       "      <td>9.508258</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>hydrogen atom</th>\n",
       "      <td>[H][H]</td>\n",
       "      <td>12.859999</td>\n",
       "      <td>12.859999</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>hydrogen</th>\n",
       "      <td>[H]</td>\n",
       "      <td>11.768784</td>\n",
       "      <td>11.768784</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>deuterium</th>\n",
       "      <td>[2H][2H]</td>\n",
       "      <td>11.748648</td>\n",
       "      <td>11.748648</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                 SMILES    Kh_mean   Kh_first  Varying_1  Varying_2  \\\n",
       "0                                                                     \n",
       "oxygen              O=O  11.267713  11.267713          0          0   \n",
       "ozone           O=[O]=O   9.508258   9.508258          1          1   \n",
       "hydrogen atom    [H][H]  12.859999  12.859999          0          0   \n",
       "hydrogen            [H]  11.768784  11.768784          0          0   \n",
       "deuterium      [2H][2H]  11.748648  11.748648          0          0   \n",
       "\n",
       "               Varying_3  Varying_4  Varying_5  Varying_6  Varying_7  \\\n",
       "0                                                                      \n",
       "oxygen                 0          0          0          0          0   \n",
       "ozone                  1          1          0          0          0   \n",
       "hydrogen atom          0          0          0          0          0   \n",
       "hydrogen               0          0          0          0          0   \n",
       "deuterium              0          0          0          0          0   \n",
       "\n",
       "               Varying_8  Varying_9  \n",
       "0                                    \n",
       "oxygen                 0          0  \n",
       "ozone                  0          0  \n",
       "hydrogen atom          0          0  \n",
       "hydrogen               0          0  \n",
       "deuterium              0          0  "
      ]
     },
     "execution_count": 90,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "smiles_constant_species_VPAS.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Dropping species with no experimental constants"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 97,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T18:50:24.344877Z",
     "start_time": "2020-01-26T18:50:24.326830Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Species to drop in VPAS data:  2100\n",
      "Species to drop in noVPAS data:  3018\n"
     ]
    }
   ],
   "source": [
    "for dataset in ['VPAS', 'noVPAS']:\n",
    "    print('Species to drop in %s data: '%dataset, len(eval('const_list_df_%s[const_list_df_%s[\"Num\"] == 0]'%(dataset, dataset))))\n",
    "    no_constants = eval('const_list_df_%s[const_list_df_%s[\"Num\"] == 0]'%(dataset, dataset))\n",
    "    exec('to_drop_%s = list(no_constants.index.values)'%dataset)\n",
    "    constants = eval('const_list_df_%s[const_list_df_%s[\"Num\"] != 0]'%(dataset, dataset))\n",
    "    exec('to_keep_%s = list(constants.index.values)'%dataset)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 98,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T18:50:47.013146Z",
     "start_time": "2020-01-26T18:50:46.963013Z"
    }
   },
   "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>SMILES</th>\n",
       "      <th>Kh_mean</th>\n",
       "      <th>Kh_first</th>\n",
       "      <th>Varying_1</th>\n",
       "      <th>Varying_2</th>\n",
       "      <th>Varying_3</th>\n",
       "      <th>Varying_4</th>\n",
       "      <th>Varying_5</th>\n",
       "      <th>Varying_6</th>\n",
       "      <th>Varying_7</th>\n",
       "      <th>Varying_8</th>\n",
       "      <th>Varying_9</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>oxygen</th>\n",
       "      <td>O=O</td>\n",
       "      <td>11.267713</td>\n",
       "      <td>11.267713</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>ozone</th>\n",
       "      <td>O=[O]=O</td>\n",
       "      <td>9.508258</td>\n",
       "      <td>9.508258</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>hydrogen atom</th>\n",
       "      <td>[H][H]</td>\n",
       "      <td>12.859999</td>\n",
       "      <td>12.859999</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>hydrogen</th>\n",
       "      <td>[H]</td>\n",
       "      <td>11.768784</td>\n",
       "      <td>11.768784</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>deuterium</th>\n",
       "      <td>[2H][2H]</td>\n",
       "      <td>11.748648</td>\n",
       "      <td>11.748648</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                 SMILES    Kh_mean   Kh_first  Varying_1  Varying_2  \\\n",
       "0                                                                     \n",
       "oxygen              O=O  11.267713  11.267713          0          0   \n",
       "ozone           O=[O]=O   9.508258   9.508258          1          1   \n",
       "hydrogen atom    [H][H]  12.859999  12.859999          0          0   \n",
       "hydrogen            [H]  11.768784  11.768784          0          0   \n",
       "deuterium      [2H][2H]  11.748648  11.748648          0          0   \n",
       "\n",
       "               Varying_3  Varying_4  Varying_5  Varying_6  Varying_7  \\\n",
       "0                                                                      \n",
       "oxygen                 0          0          0          0          0   \n",
       "ozone                  1          1          0          0          0   \n",
       "hydrogen atom          0          0          0          0          0   \n",
       "hydrogen               0          0          0          0          0   \n",
       "deuterium              0          0          0          0          0   \n",
       "\n",
       "               Varying_8  Varying_9  \n",
       "0                                    \n",
       "oxygen                 0          0  \n",
       "ozone                  0          0  \n",
       "hydrogen atom          0          0  \n",
       "hydrogen               0          0  \n",
       "deuterium              0          0  "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "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>SMILES</th>\n",
       "      <th>Kh_mean</th>\n",
       "      <th>Kh_first</th>\n",
       "      <th>Varying_1</th>\n",
       "      <th>Varying_2</th>\n",
       "      <th>Varying_3</th>\n",
       "      <th>Varying_4</th>\n",
       "      <th>Varying_5</th>\n",
       "      <th>Varying_6</th>\n",
       "      <th>Varying_7</th>\n",
       "      <th>Varying_8</th>\n",
       "      <th>Varying_9</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>oxygen</th>\n",
       "      <td>O=O</td>\n",
       "      <td>11.267713</td>\n",
       "      <td>11.267713</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>ozone</th>\n",
       "      <td>O=[O]=O</td>\n",
       "      <td>9.508258</td>\n",
       "      <td>9.508258</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>hydrogen atom</th>\n",
       "      <td>[H][H]</td>\n",
       "      <td>12.859999</td>\n",
       "      <td>12.859999</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>hydrogen</th>\n",
       "      <td>[H]</td>\n",
       "      <td>11.768784</td>\n",
       "      <td>11.768784</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>deuterium</th>\n",
       "      <td>[2H][2H]</td>\n",
       "      <td>11.748648</td>\n",
       "      <td>11.748648</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                 SMILES    Kh_mean   Kh_first  Varying_1  Varying_2  \\\n",
       "0                                                                     \n",
       "oxygen              O=O  11.267713  11.267713          0          0   \n",
       "ozone           O=[O]=O   9.508258   9.508258          1          1   \n",
       "hydrogen atom    [H][H]  12.859999  12.859999          0          0   \n",
       "hydrogen            [H]  11.768784  11.768784          0          0   \n",
       "deuterium      [2H][2H]  11.748648  11.748648          0          0   \n",
       "\n",
       "               Varying_3  Varying_4  Varying_5  Varying_6  Varying_7  \\\n",
       "0                                                                      \n",
       "oxygen                 0          0          0          0          0   \n",
       "ozone                  1          1          0          0          0   \n",
       "hydrogen atom          0          0          0          0          0   \n",
       "hydrogen               0          0          0          0          0   \n",
       "deuterium              0          0          0          0          0   \n",
       "\n",
       "               Varying_8  Varying_9  \n",
       "0                                    \n",
       "oxygen                 0          0  \n",
       "ozone                  0          0  \n",
       "hydrogen atom          0          0  \n",
       "hydrogen               0          0  \n",
       "deuterium              0          0  "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "for dataset in ['VPAS', 'noVPAS']:\n",
    "    \n",
    "    exec('filtered_smiles_constant_species_%s = smiles_constant_species_%s.loc[to_keep_%s]'%(dataset, dataset, dataset))\n",
    "    eval('filtered_smiles_constant_species_%s.rename({1:\"SMILES\", 2:\"Constant\", 3:\"Varying\"}, axis=1, inplace=True)'%dataset)\n",
    "    display(eval('filtered_smiles_constant_species_%s.head()'%dataset))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Organic/inorganic split"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 99,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T18:55:32.203396Z",
     "start_time": "2020-01-26T18:55:30.724465Z"
    }
   },
   "outputs": [],
   "source": [
    "for dataset in ['VPAS', 'noVPAS']:\n",
    "    exec('carbons_list_%s = [] ; no_carbons_list_%s = []'%(dataset, dataset))\n",
    "    for i in range(0, len(eval('filtered_smiles_constant_species_%s'%dataset))):\n",
    "        smiles = eval('filtered_smiles_constant_species_%s.SMILES.iloc[i]'%dataset)\n",
    "        try:\n",
    "            mol = Chem.MolFromSmiles(smiles) # create rdkit mol object\n",
    "            smarts = Chem.MolToSmarts(mol) # use mol object to create SMARTS\n",
    "        except:\n",
    "            smarts = 'x' #mark failed species in the list\n",
    "        name = eval('filtered_smiles_constant_species_%s.index[i]'%dataset)\n",
    "        if type(smarts) == str:\n",
    "            if '#6' in smarts:\n",
    "                eval('carbons_list_%s.append(i)'%dataset) # RDKit SMARTS strings work by atomic number (hence carbon = #6)\n",
    "            else:\n",
    "                eval('no_carbons_list_%s.append(i)'%dataset) # append non-carbon-containing species to a seperate list"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 105,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T18:57:43.606767Z",
     "start_time": "2020-01-26T18:57:43.599749Z"
    }
   },
   "outputs": [],
   "source": [
    "organics_VPAS = filtered_smiles_constant_species_VPAS.iloc[carbons_list_VPAS] # new df containing carbon compounds only\n",
    "organics_noVPAS = filtered_smiles_constant_species_noVPAS.iloc[carbons_list_noVPAS] # new df containing carbon compounds only"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Removing metal-containing species"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 108,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T19:03:34.693041Z",
     "start_time": "2020-01-26T19:03:33.210099Z"
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\jtd1g16\\AppData\\Local\\Continuum\\anaconda3\\lib\\site-packages\\pandas\\core\\frame.py:4117: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  errors=errors,\n"
     ]
    }
   ],
   "source": [
    "for dataset in ['VPAS', 'noVPAS']:\n",
    "    drop_list = []\n",
    "    # removing metals using list of common metals.\n",
    "    # aniline compound at beginnning contains (corrupted?) SMILES string full of tungsten.\n",
    "    metals_removal_list = [\"4,4',4''-methanetriyltris(N,N-diethylaniline)\", 'Hg', 'Se', 'As', 'Zn', 'Si', 'Na', 'Al', 'Pb', 'Sn']\n",
    "    for i in range(0, len(eval('organics_%s'%dataset))):\n",
    "        for metal in metals_removal_list:\n",
    "            if metal in eval('organics_%s.SMILES[i]'%dataset):\n",
    "                drop_list.append(eval('organics_%s.index[i]'%dataset))\n",
    "                #create drop list in loop, then remove all metal-containing species at once\n",
    "\n",
    "    eval('organics_%s.drop(labels=drop_list, inplace=True, axis=0)'%dataset)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Export organics as csv"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 111,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T19:05:58.516433Z",
     "start_time": "2020-01-26T19:05:58.416168Z"
    }
   },
   "outputs": [],
   "source": [
    "organics_VPAS.to_csv('organics_constants_VPAS.csv')\n",
    "organics_noVPAS.to_csv('organics_constants_noVPAS.csv')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Create RDKit descriptors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 114,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T19:15:49.331266Z",
     "start_time": "2020-01-26T19:15:49.326254Z"
    }
   },
   "outputs": [],
   "source": [
    "def try_func(list_name, function, mol):\n",
    "    \"\"\"Tries function on an rdkit mol object; if function throws error then np.nan\n",
    "    is appended to the list - if the function succeeds then the output is appended\"\"\"\n",
    "    try:    \n",
    "        var = (function(mol))\n",
    "        list_name.append(var)\n",
    "    except:\n",
    "        list_name.append(np.nan)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 115,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T19:16:23.241425Z",
     "start_time": "2020-01-26T19:16:17.941334Z"
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\jtd1g16\\AppData\\Local\\Continuum\\anaconda3\\lib\\site-packages\\ipykernel_launcher.py:1: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
      "Try using .loc[row_indexer,col_indexer] = value instead\n",
      "\n",
      "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  \"\"\"Entry point for launching an IPython kernel.\n"
     ]
    }
   ],
   "source": [
    "for dataset in ['VPAS', 'noVPAS']:\n",
    "    wt_list = [] ; logp_list = [] ; fraction_sp3_list = []\n",
    "    ali_homo_rings_list = [] ; ali_hetero_rings_list = [] \n",
    "    aro_homo_rings_list = [] ; aro_hetero_rings_list = []\n",
    "    amide_bonds_list = [] ; hba_list = [] ; hbd_list = []\n",
    "    tpsa_list = []\n",
    "\n",
    "    for i in range(0, len(eval('organics_%s'%dataset))):\n",
    "        exec('smiles = (organics_%s.SMILES[i])'%dataset)\n",
    "        try:\n",
    "            mol = Chem.AddHs(Chem.MolFromSmiles(smiles))\n",
    "        except Exception as e:\n",
    "            pass\n",
    "\n",
    "\n",
    "        if type(mol) != str and type(mol) != None:\n",
    "            try_func(wt_list, descriptor.CalcExactMolWt, mol)\n",
    "            try_func(logp_list, Crippen.MolLogP, mol)\n",
    "            try_func(fraction_sp3_list, descriptor.CalcFractionCSP3, mol)\n",
    "            try_func(ali_homo_rings_list, descriptor.CalcNumAliphaticCarbocycles, mol)\n",
    "            try_func(ali_hetero_rings_list, descriptor.CalcNumAliphaticHeterocycles, mol)\n",
    "            try_func(aro_homo_rings_list, descriptor.CalcNumAromaticCarbocycles, mol)\n",
    "            try_func(aro_hetero_rings_list, descriptor.CalcNumAromaticHeterocycles, mol)\n",
    "            try_func(amide_bonds_list, descriptor.CalcNumAmideBonds, mol)\n",
    "            try_func(hba_list, descriptor.CalcNumHBA, mol)\n",
    "            try_func(hbd_list, descriptor.CalcNumHBD, mol)\n",
    "            try_func(tpsa_list, descriptor.CalcTPSA, mol)\n",
    "\n",
    "    exec(\"organics_%s['wt'] = wt_list\"%dataset)\n",
    "    exec(\"organics_%s['logp'] = logp_list\"%dataset)\n",
    "    exec(\"organics_%s['fraction_sp3'] = fraction_sp3_list\"%dataset)\n",
    "    exec(\"organics_%s['ali_homo_rings'] = ali_homo_rings_list\"%dataset)\n",
    "    exec(\"organics_%s['ali_hetero_rings'] = ali_hetero_rings_list\"%dataset)\n",
    "    exec(\"organics_%s['amide_bonds'] = amide_bonds_list\"%dataset)\n",
    "    exec(\"organics_%s['aro_homo_rings'] = aro_homo_rings_list\"%dataset)\n",
    "    exec(\"organics_%s['aro_hetero_rings'] = aro_hetero_rings_list\"%dataset)\n",
    "    exec(\"organics_%s['hba'] = hba_list\"%dataset)\n",
    "    exec(\"organics_%s['hbd'] = hbd_list\"%dataset)\n",
    "    exec(\"organics_%s['tpsa'] = tpsa_list\"%dataset)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 123,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T19:30:12.521278Z",
     "start_time": "2020-01-26T19:26:26.325879Z"
    }
   },
   "outputs": [
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "ce06b64c21754cc99e85be5640107ab8",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "IntProgress(value=0, description='Progress', max=2196)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\jtd1g16\\AppData\\Local\\Continuum\\anaconda3\\lib\\site-packages\\ipykernel_launcher.py:1: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
      "Try using .loc[row_indexer,col_indexer] = value instead\n",
      "\n",
      "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  \"\"\"Entry point for launching an IPython kernel.\n"
     ]
    },
    {
     "data": {
      "application/vnd.jupyter.widget-view+json": {
       "model_id": "b302895951814ac2bf60eb619e1b6f6f",
       "version_major": 2,
       "version_minor": 0
      },
      "text/plain": [
       "IntProgress(value=0, description='Progress', max=1404)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "for dataset in ['VPAS', 'noVPAS']:\n",
    "    asph_list = [] ; ecce_list = [] ; isf_list = [] ; npr1_list = [] ; npr2_list = [] \n",
    "    pmi1_list = [] ; pmi2_list = [] ; pmi3_list = [] ; rog_list = [] ; si_list = [] \n",
    "\n",
    "\n",
    "    max_count = len(eval('organics_%s'%dataset))\n",
    "    f = IntProgress(min=0, max=max_count, description='Progress') # instantiate the bar\n",
    "    display(f) # display the bar\n",
    "    count = 0\n",
    "\n",
    "\n",
    "    #Place this part inside the loop\n",
    "\n",
    "\n",
    "    from rdkit.Chem import Descriptors3D\n",
    "\n",
    "    for i in range(0, len(eval('organics_%s'%dataset))):\n",
    "        exec('smiles = (organics_%s.SMILES[i])'%dataset)\n",
    "        try:\n",
    "            mol = Chem.AddHs(Chem.MolFromSmiles(smiles))\n",
    "            AllChem.EmbedMolecule(mol) #embed to calc 3D descriptors\n",
    "            Chem.rdForceFieldHelpers.MMFFOptimizeMolecule(mol) #optimize\n",
    "        except Exception as e:\n",
    "            mol=None\n",
    "        try_func(asph_list, Descriptors3D.Asphericity, mol)\n",
    "        try_func(ecce_list, Descriptors3D.Eccentricity, mol)\n",
    "        try_func(isf_list, Descriptors3D.InertialShapeFactor, mol)\n",
    "        try_func(npr1_list, Descriptors3D.NPR1, mol)\n",
    "        try_func(npr2_list, Descriptors3D.NPR2, mol)\n",
    "        try_func(pmi1_list, Descriptors3D.PMI1, mol)\n",
    "        try_func(pmi2_list, Descriptors3D.PMI2, mol)\n",
    "        try_func(pmi3_list, Descriptors3D.PMI3, mol)\n",
    "        try_func(rog_list, Descriptors3D.RadiusOfGyration, mol)\n",
    "        try_func(si_list, Descriptors3D.SpherocityIndex, mol)\n",
    "        f.value += 1 # signal to increment the progress bar\n",
    "        count += 1\n",
    "\n",
    "    exec(\"organics_%s['asph'] = asph_list\"%dataset)\n",
    "    exec(\"organics_%s['ecce'] = ecce_list\"%dataset)\n",
    "    exec(\"organics_%s['isf'] = isf_list\"%dataset)\n",
    "    exec(\"organics_%s['npr1'] = npr1_list\"%dataset)\n",
    "    exec(\"organics_%s['npr2'] = npr2_list\"%dataset)\n",
    "    exec(\"organics_%s['pmi1'] = pmi1_list\"%dataset)\n",
    "    exec(\"organics_%s['pmi2'] = pmi2_list\"%dataset)\n",
    "    exec(\"organics_%s['pmi3'] = pmi3_list\"%dataset)\n",
    "    exec(\"organics_%s['rog'] = rog_list\"%dataset)\n",
    "    exec(\"organics_%s['si'] = si_list\"%dataset)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 126,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2020-01-26T21:42:10.480833Z",
     "start_time": "2020-01-26T21:42:10.163991Z"
    }
   },
   "outputs": [],
   "source": [
    "organics_VPAS.to_csv('rdkit_desc_VPAS.csv')\n",
    "organics_noVPAS.to_csv('rdkit_desc_noVPAS.csv')"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.3"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {
    "height": "calc(100% - 180px)",
    "left": "10px",
    "top": "150px",
    "width": "248.188px"
   },
   "toc_section_display": true,
   "toc_window_display": false
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}