2_relative_model_fit.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import os\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "from bambi import Model\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "plt.rc(\"axes.spines\", top=False, right=False)\n",
    "\n",
    "from plotting_functions import show_img, cm2inch, plot_correlation\n",
    "from pymc3 import plot_posterior"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load subject summary of behavioral measures"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "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 tr th {\n",
       "        text-align: left;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr>\n",
       "      <th></th>\n",
       "      <th>subject</th>\n",
       "      <th colspan=\"8\" halign=\"left\">rt</th>\n",
       "      <th>best_chosen</th>\n",
       "      <th>gaze_influence</th>\n",
       "      <th>dataset</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th></th>\n",
       "      <th>Unnamed: 1_level_1</th>\n",
       "      <th>mean</th>\n",
       "      <th>std</th>\n",
       "      <th>min</th>\n",
       "      <th>max</th>\n",
       "      <th>se</th>\n",
       "      <th>q1</th>\n",
       "      <th>q3</th>\n",
       "      <th>iqr</th>\n",
       "      <th>mean</th>\n",
       "      <th>Unnamed: 11_level_1</th>\n",
       "      <th>Unnamed: 12_level_1</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0</td>\n",
       "      <td>1762.538462</td>\n",
       "      <td>1221.112889</td>\n",
       "      <td>682.0</td>\n",
       "      <td>7866.0</td>\n",
       "      <td>128.716600</td>\n",
       "      <td>1107.00</td>\n",
       "      <td>1970.50</td>\n",
       "      <td>863.50</td>\n",
       "      <td>0.857143</td>\n",
       "      <td>0.040579</td>\n",
       "      <td>krajbich2010</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1</td>\n",
       "      <td>1672.050505</td>\n",
       "      <td>868.963929</td>\n",
       "      <td>667.0</td>\n",
       "      <td>4826.0</td>\n",
       "      <td>87.778612</td>\n",
       "      <td>1108.50</td>\n",
       "      <td>1862.00</td>\n",
       "      <td>753.50</td>\n",
       "      <td>0.737374</td>\n",
       "      <td>0.202390</td>\n",
       "      <td>krajbich2010</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2</td>\n",
       "      <td>4477.180000</td>\n",
       "      <td>3088.070399</td>\n",
       "      <td>943.0</td>\n",
       "      <td>15825.0</td>\n",
       "      <td>310.362753</td>\n",
       "      <td>2146.25</td>\n",
       "      <td>5728.50</td>\n",
       "      <td>3582.25</td>\n",
       "      <td>0.850000</td>\n",
       "      <td>0.032181</td>\n",
       "      <td>krajbich2010</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3</td>\n",
       "      <td>1882.613636</td>\n",
       "      <td>541.973214</td>\n",
       "      <td>1003.0</td>\n",
       "      <td>4742.0</td>\n",
       "      <td>58.105628</td>\n",
       "      <td>1519.50</td>\n",
       "      <td>2180.25</td>\n",
       "      <td>660.75</td>\n",
       "      <td>0.693182</td>\n",
       "      <td>0.071577</td>\n",
       "      <td>krajbich2010</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>4</td>\n",
       "      <td>1727.370000</td>\n",
       "      <td>1007.248844</td>\n",
       "      <td>673.0</td>\n",
       "      <td>5891.0</td>\n",
       "      <td>101.232318</td>\n",
       "      <td>1160.50</td>\n",
       "      <td>1868.00</td>\n",
       "      <td>707.50</td>\n",
       "      <td>0.840000</td>\n",
       "      <td>0.141787</td>\n",
       "      <td>krajbich2010</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "             subject           rt                                            \\\n",
       "  Unnamed: 1_level_1         mean          std     min      max          se   \n",
       "0                  0  1762.538462  1221.112889   682.0   7866.0  128.716600   \n",
       "1                  1  1672.050505   868.963929   667.0   4826.0   87.778612   \n",
       "2                  2  4477.180000  3088.070399   943.0  15825.0  310.362753   \n",
       "3                  3  1882.613636   541.973214  1003.0   4742.0   58.105628   \n",
       "4                  4  1727.370000  1007.248844   673.0   5891.0  101.232318   \n",
       "\n",
       "                             best_chosen      gaze_influence  \\\n",
       "        q1       q3      iqr        mean Unnamed: 11_level_1   \n",
       "0  1107.00  1970.50   863.50    0.857143            0.040579   \n",
       "1  1108.50  1862.00   753.50    0.737374            0.202390   \n",
       "2  2146.25  5728.50  3582.25    0.850000            0.032181   \n",
       "3  1519.50  2180.25   660.75    0.693182            0.071577   \n",
       "4  1160.50  1868.00   707.50    0.840000            0.141787   \n",
       "\n",
       "              dataset  \n",
       "  Unnamed: 12_level_1  \n",
       "0        krajbich2010  \n",
       "1        krajbich2010  \n",
       "2        krajbich2010  \n",
       "3        krajbich2010  \n",
       "4        krajbich2010  "
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "subject_summary = pd.read_csv('results/descriptives/subject_summary.csv', index_col=[0], header=[0, 1])\n",
    "subject_summary.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 2.2 Relative model fit (individual level model comparison)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load WAICs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "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>subject</th>\n",
       "      <th>dWAIC</th>\n",
       "      <th>dSE</th>\n",
       "      <th>dataset</th>\n",
       "      <th>warning</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0</td>\n",
       "      <td>-0.74</td>\n",
       "      <td>2.48</td>\n",
       "      <td>krajbich2010</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1</td>\n",
       "      <td>-21.53</td>\n",
       "      <td>9.18</td>\n",
       "      <td>krajbich2010</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2</td>\n",
       "      <td>0.98</td>\n",
       "      <td>1.48</td>\n",
       "      <td>krajbich2010</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3</td>\n",
       "      <td>0.81</td>\n",
       "      <td>3.01</td>\n",
       "      <td>krajbich2010</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>4</td>\n",
       "      <td>-9.55</td>\n",
       "      <td>5.75</td>\n",
       "      <td>krajbich2010</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   subject  dWAIC   dSE       dataset  warning\n",
       "0        0  -0.74  2.48  krajbich2010      0.0\n",
       "1        1 -21.53  9.18  krajbich2010      0.0\n",
       "2        2   0.98  1.48  krajbich2010      0.0\n",
       "3        3   0.81  3.01  krajbich2010      0.0\n",
       "4        4  -9.55  5.75  krajbich2010      0.0"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "waic_summaries = [file for file in os.listdir('results/model_comparison/multiplicative_vs_nobias/')\n",
    "                  if file.endswith('waic.csv')]\n",
    "\n",
    "n_krajbich2010 = 39\n",
    "n_krajbich2011 = 30\n",
    "n_folke2016 = 24\n",
    "n_tavares2017 = 25\n",
    "n = n_krajbich2010 + n_krajbich2011 + n_folke2016 + n_tavares2017\n",
    "\n",
    "waic = pd.DataFrame(dict(subject=np.arange(n)))\n",
    "waic['dWAIC'] = np.nan\n",
    "waic['dSE'] = np.nan\n",
    "waic['dataset'] = np.nan\n",
    "waic['warning'] = np.nan\n",
    "\n",
    "for file in waic_summaries:\n",
    "    _, _, _, subject, _ = file.split('_')\n",
    "    subject = int(subject)\n",
    "    dataset = subject_summary.loc[subject, 'dataset'].values[0]\n",
    "\n",
    "    waic_summary = pd.read_csv(os.path.join('results', 'model_comparison', 'multiplicative_vs_nobias', file), index_col='Unnamed: 0')\n",
    "    waic.loc[waic['subject'] == subject, 'dWAIC'] = waic_summary.loc['GLAM', 'WAIC'] - waic_summary.loc['nobias', 'WAIC']\n",
    "    waic.loc[waic['subject'] == subject, 'dSE'] = waic_summary['dSE'].values[1]\n",
    "    waic.loc[waic['subject'] == subject, 'dataset'] = dataset\n",
    "    waic.loc[waic['subject'] == subject, 'warning'] = waic_summary['var_warn'].sum()\n",
    "    \n",
    "waic.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "count    118.000000\n",
       "mean     -62.081949\n",
       "std       80.795577\n",
       "min     -400.640000\n",
       "25%      -82.307500\n",
       "50%      -33.605000\n",
       "75%       -7.430000\n",
       "max        1.600000\n",
       "Name: dWAIC, dtype: float64"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "waic['dWAIC'].describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "109"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.sum(waic['dWAIC'] < 0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.923728813559322"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.mean(waic['dWAIC'] < 0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "dataset\n",
       "folke2016       1.000000\n",
       "krajbich2010    0.794872\n",
       "krajbich2011    0.966667\n",
       "tavares2017     1.000000\n",
       "dtype: float64"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "waic.groupby('dataset').apply(lambda x: np.mean(x['dWAIC'] < 0))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Run mixed effects regressions, accounting for dataset level"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "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 tr th {\n",
       "        text-align: left;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr>\n",
       "      <th></th>\n",
       "      <th>subject</th>\n",
       "      <th colspan=\"8\" halign=\"left\">rt</th>\n",
       "      <th>best_chosen</th>\n",
       "      <th>gaze_influence</th>\n",
       "      <th>dataset</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th></th>\n",
       "      <th>Unnamed: 1_level_1</th>\n",
       "      <th>mean</th>\n",
       "      <th>std</th>\n",
       "      <th>min</th>\n",
       "      <th>max</th>\n",
       "      <th>se</th>\n",
       "      <th>q1</th>\n",
       "      <th>q3</th>\n",
       "      <th>iqr</th>\n",
       "      <th>mean</th>\n",
       "      <th>Unnamed: 11_level_1</th>\n",
       "      <th>Unnamed: 12_level_1</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0</td>\n",
       "      <td>1762.538462</td>\n",
       "      <td>1221.112889</td>\n",
       "      <td>682.0</td>\n",
       "      <td>7866.0</td>\n",
       "      <td>128.716600</td>\n",
       "      <td>1107.00</td>\n",
       "      <td>1970.50</td>\n",
       "      <td>863.50</td>\n",
       "      <td>0.857143</td>\n",
       "      <td>0.040579</td>\n",
       "      <td>krajbich2010</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1</td>\n",
       "      <td>1672.050505</td>\n",
       "      <td>868.963929</td>\n",
       "      <td>667.0</td>\n",
       "      <td>4826.0</td>\n",
       "      <td>87.778612</td>\n",
       "      <td>1108.50</td>\n",
       "      <td>1862.00</td>\n",
       "      <td>753.50</td>\n",
       "      <td>0.737374</td>\n",
       "      <td>0.202390</td>\n",
       "      <td>krajbich2010</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2</td>\n",
       "      <td>4477.180000</td>\n",
       "      <td>3088.070399</td>\n",
       "      <td>943.0</td>\n",
       "      <td>15825.0</td>\n",
       "      <td>310.362753</td>\n",
       "      <td>2146.25</td>\n",
       "      <td>5728.50</td>\n",
       "      <td>3582.25</td>\n",
       "      <td>0.850000</td>\n",
       "      <td>0.032181</td>\n",
       "      <td>krajbich2010</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>3</td>\n",
       "      <td>1882.613636</td>\n",
       "      <td>541.973214</td>\n",
       "      <td>1003.0</td>\n",
       "      <td>4742.0</td>\n",
       "      <td>58.105628</td>\n",
       "      <td>1519.50</td>\n",
       "      <td>2180.25</td>\n",
       "      <td>660.75</td>\n",
       "      <td>0.693182</td>\n",
       "      <td>0.071577</td>\n",
       "      <td>krajbich2010</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>4</td>\n",
       "      <td>1727.370000</td>\n",
       "      <td>1007.248844</td>\n",
       "      <td>673.0</td>\n",
       "      <td>5891.0</td>\n",
       "      <td>101.232318</td>\n",
       "      <td>1160.50</td>\n",
       "      <td>1868.00</td>\n",
       "      <td>707.50</td>\n",
       "      <td>0.840000</td>\n",
       "      <td>0.141787</td>\n",
       "      <td>krajbich2010</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "             subject           rt                                            \\\n",
       "  Unnamed: 1_level_1         mean          std     min      max          se   \n",
       "0                  0  1762.538462  1221.112889   682.0   7866.0  128.716600   \n",
       "1                  1  1672.050505   868.963929   667.0   4826.0   87.778612   \n",
       "2                  2  4477.180000  3088.070399   943.0  15825.0  310.362753   \n",
       "3                  3  1882.613636   541.973214  1003.0   4742.0   58.105628   \n",
       "4                  4  1727.370000  1007.248844   673.0   5891.0  101.232318   \n",
       "\n",
       "                             best_chosen      gaze_influence  \\\n",
       "        q1       q3      iqr        mean Unnamed: 11_level_1   \n",
       "0  1107.00  1970.50   863.50    0.857143            0.040579   \n",
       "1  1108.50  1862.00   753.50    0.737374            0.202390   \n",
       "2  2146.25  5728.50  3582.25    0.850000            0.032181   \n",
       "3  1519.50  2180.25   660.75    0.693182            0.071577   \n",
       "4  1160.50  1868.00   707.50    0.840000            0.141787   \n",
       "\n",
       "              dataset  \n",
       "  Unnamed: 12_level_1  \n",
       "0        krajbich2010  \n",
       "1        krajbich2010  \n",
       "2        krajbich2010  \n",
       "3        krajbich2010  \n",
       "4        krajbich2010  "
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "subject_summary.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def prepare_df_for_bambi(waic, subject_summary):\n",
    "    df = pd.DataFrame(dict(subject=waic['subject'],\n",
    "                           dWAIC=waic['dWAIC'],\n",
    "                           dataset=waic['dataset'],\n",
    "                           mean_rt=subject_summary['rt']['mean'],\n",
    "                           best_chosen=100 * subject_summary['best_chosen']['mean'],\n",
    "                           gaze_influence=100 * subject_summary['gaze_influence']['Unnamed: 11_level_1']))\n",
    "    return df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "mfx_df = prepare_df_for_bambi(waic, subject_summary)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "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>best_chosen</th>\n",
       "      <th>dWAIC</th>\n",
       "      <th>dataset</th>\n",
       "      <th>gaze_influence</th>\n",
       "      <th>mean_rt</th>\n",
       "      <th>subject</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>85.714286</td>\n",
       "      <td>-0.74</td>\n",
       "      <td>krajbich2010</td>\n",
       "      <td>4.057934</td>\n",
       "      <td>1762.538462</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>73.737374</td>\n",
       "      <td>-21.53</td>\n",
       "      <td>krajbich2010</td>\n",
       "      <td>20.238963</td>\n",
       "      <td>1672.050505</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>85.000000</td>\n",
       "      <td>0.98</td>\n",
       "      <td>krajbich2010</td>\n",
       "      <td>3.218073</td>\n",
       "      <td>4477.180000</td>\n",
       "      <td>2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>69.318182</td>\n",
       "      <td>0.81</td>\n",
       "      <td>krajbich2010</td>\n",
       "      <td>7.157652</td>\n",
       "      <td>1882.613636</td>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>84.000000</td>\n",
       "      <td>-9.55</td>\n",
       "      <td>krajbich2010</td>\n",
       "      <td>14.178662</td>\n",
       "      <td>1727.370000</td>\n",
       "      <td>4</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   best_chosen  dWAIC       dataset  gaze_influence      mean_rt  subject\n",
       "0    85.714286  -0.74  krajbich2010        4.057934  1762.538462        0\n",
       "1    73.737374 -21.53  krajbich2010       20.238963  1672.050505        1\n",
       "2    85.000000   0.98  krajbich2010        3.218073  4477.180000        2\n",
       "3    69.318182   0.81  krajbich2010        7.157652  1882.613636        3\n",
       "4    84.000000  -9.55  krajbich2010       14.178662  1727.370000        4"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mfx_df.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "n_samples = 20000"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Previous results loaded from \"results/mixed_effects_models/2_relative_model_comparison_b_mean_rt.csv\"...\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>mean</th>\n",
       "      <th>sd</th>\n",
       "      <th>hpd0.95_lower</th>\n",
       "      <th>hpd0.95_upper</th>\n",
       "      <th>effective_n</th>\n",
       "      <th>gelman_rubin</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>Intercept</th>\n",
       "      <td>3018.413606</td>\n",
       "      <td>740.492000</td>\n",
       "      <td>1697.350531</td>\n",
       "      <td>4664.071945</td>\n",
       "      <td>236</td>\n",
       "      <td>1.012394</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset[folke2016]</th>\n",
       "      <td>-2.264231</td>\n",
       "      <td>5.428171</td>\n",
       "      <td>-15.265071</td>\n",
       "      <td>7.562233</td>\n",
       "      <td>138</td>\n",
       "      <td>1.007388</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset[krajbich2010]</th>\n",
       "      <td>3.929004</td>\n",
       "      <td>7.591286</td>\n",
       "      <td>-7.126642</td>\n",
       "      <td>21.757592</td>\n",
       "      <td>84</td>\n",
       "      <td>1.012804</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset[krajbich2011]</th>\n",
       "      <td>0.766852</td>\n",
       "      <td>5.568722</td>\n",
       "      <td>-9.497250</td>\n",
       "      <td>15.329559</td>\n",
       "      <td>65</td>\n",
       "      <td>1.024069</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset[tavares2017]</th>\n",
       "      <td>-0.316185</td>\n",
       "      <td>5.330837</td>\n",
       "      <td>-11.273393</td>\n",
       "      <td>15.452821</td>\n",
       "      <td>47</td>\n",
       "      <td>1.029550</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset_sd</th>\n",
       "      <td>6.098573</td>\n",
       "      <td>5.438404</td>\n",
       "      <td>0.001251</td>\n",
       "      <td>17.039704</td>\n",
       "      <td>123</td>\n",
       "      <td>1.003957</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset[folke2016]</th>\n",
       "      <td>2253.237608</td>\n",
       "      <td>779.963533</td>\n",
       "      <td>693.567832</td>\n",
       "      <td>3799.449096</td>\n",
       "      <td>337</td>\n",
       "      <td>1.002170</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset[krajbich2010]</th>\n",
       "      <td>-758.520704</td>\n",
       "      <td>747.559973</td>\n",
       "      <td>-2340.662162</td>\n",
       "      <td>671.110617</td>\n",
       "      <td>268</td>\n",
       "      <td>1.011362</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset[krajbich2011]</th>\n",
       "      <td>-480.632014</td>\n",
       "      <td>750.650099</td>\n",
       "      <td>-2112.649017</td>\n",
       "      <td>975.235279</td>\n",
       "      <td>293</td>\n",
       "      <td>1.008871</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset[tavares2017]</th>\n",
       "      <td>-1004.691123</td>\n",
       "      <td>762.858145</td>\n",
       "      <td>-2508.961853</td>\n",
       "      <td>507.980766</td>\n",
       "      <td>430</td>\n",
       "      <td>1.006865</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset_sd</th>\n",
       "      <td>1521.538302</td>\n",
       "      <td>480.055906</td>\n",
       "      <td>704.049022</td>\n",
       "      <td>2452.825535</td>\n",
       "      <td>851</td>\n",
       "      <td>1.003293</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>mean_rt_sd</th>\n",
       "      <td>1073.580659</td>\n",
       "      <td>72.016383</td>\n",
       "      <td>940.315575</td>\n",
       "      <td>1219.393132</td>\n",
       "      <td>2757</td>\n",
       "      <td>1.000033</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC</th>\n",
       "      <td>1.139426</td>\n",
       "      <td>5.610384</td>\n",
       "      <td>-17.189608</td>\n",
       "      <td>9.696600</td>\n",
       "      <td>46</td>\n",
       "      <td>1.032335</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                                    mean          sd  hpd0.95_lower  \\\n",
       "Intercept                    3018.413606  740.492000    1697.350531   \n",
       "dWAIC|dataset[folke2016]       -2.264231    5.428171     -15.265071   \n",
       "dWAIC|dataset[krajbich2010]     3.929004    7.591286      -7.126642   \n",
       "dWAIC|dataset[krajbich2011]     0.766852    5.568722      -9.497250   \n",
       "dWAIC|dataset[tavares2017]     -0.316185    5.330837     -11.273393   \n",
       "dWAIC|dataset_sd                6.098573    5.438404       0.001251   \n",
       "1|dataset[folke2016]         2253.237608  779.963533     693.567832   \n",
       "1|dataset[krajbich2010]      -758.520704  747.559973   -2340.662162   \n",
       "1|dataset[krajbich2011]      -480.632014  750.650099   -2112.649017   \n",
       "1|dataset[tavares2017]      -1004.691123  762.858145   -2508.961853   \n",
       "1|dataset_sd                 1521.538302  480.055906     704.049022   \n",
       "mean_rt_sd                   1073.580659   72.016383     940.315575   \n",
       "dWAIC                           1.139426    5.610384     -17.189608   \n",
       "\n",
       "                             hpd0.95_upper  effective_n  gelman_rubin  \n",
       "Intercept                      4664.071945          236      1.012394  \n",
       "dWAIC|dataset[folke2016]          7.562233          138      1.007388  \n",
       "dWAIC|dataset[krajbich2010]      21.757592           84      1.012804  \n",
       "dWAIC|dataset[krajbich2011]      15.329559           65      1.024069  \n",
       "dWAIC|dataset[tavares2017]       15.452821           47      1.029550  \n",
       "dWAIC|dataset_sd                 17.039704          123      1.003957  \n",
       "1|dataset[folke2016]           3799.449096          337      1.002170  \n",
       "1|dataset[krajbich2010]         671.110617          268      1.011362  \n",
       "1|dataset[krajbich2011]         975.235279          293      1.008871  \n",
       "1|dataset[tavares2017]          507.980766          430      1.006865  \n",
       "1|dataset_sd                   2452.825535          851      1.003293  \n",
       "mean_rt_sd                     1219.393132         2757      1.000033  \n",
       "dWAIC                             9.696600           46      1.032335  "
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "result_path = os.path.join('results', 'mixed_effects_models', '2_relative_model_comparison_b_mean_rt.csv')\n",
    "if not os.path.isfile(result_path):\n",
    "    m = Model(mfx_df)\n",
    "    b_mean_rt_fit = m.fit('mean_rt ~ dWAIC',\n",
    "                          random=['dWAIC|dataset'],\n",
    "                          categorical=['dataset'],\n",
    "                          samples=n_samples)\n",
    "    b_summary = b_mean_rt_fit.summary(ranefs=True)\n",
    "    b_summary.to_csv(result_path)\n",
    "else:\n",
    "    print('Previous results loaded from \"{}\"...'.format(result_path))\n",
    "    b_summary = pd.read_csv(result_path, index_col=0)\n",
    "b_summary"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Previous results loaded from \"results/mixed_effects_models/2_relative_model_comparison_c_best_chosen.csv\"...\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>mean</th>\n",
       "      <th>sd</th>\n",
       "      <th>hpd0.95_lower</th>\n",
       "      <th>hpd0.95_upper</th>\n",
       "      <th>effective_n</th>\n",
       "      <th>gelman_rubin</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>1|dataset[folke2016]</th>\n",
       "      <td>-3.037726</td>\n",
       "      <td>3.692019</td>\n",
       "      <td>-10.930851</td>\n",
       "      <td>3.819971</td>\n",
       "      <td>15732</td>\n",
       "      <td>0.999990</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset[krajbich2010]</th>\n",
       "      <td>1.492997</td>\n",
       "      <td>3.258928</td>\n",
       "      <td>-5.091004</td>\n",
       "      <td>8.290051</td>\n",
       "      <td>11386</td>\n",
       "      <td>1.000025</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset[krajbich2011]</th>\n",
       "      <td>-2.438621</td>\n",
       "      <td>3.349371</td>\n",
       "      <td>-9.703574</td>\n",
       "      <td>4.007553</td>\n",
       "      <td>12798</td>\n",
       "      <td>1.000005</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset[tavares2017]</th>\n",
       "      <td>4.928456</td>\n",
       "      <td>3.863594</td>\n",
       "      <td>-1.695214</td>\n",
       "      <td>13.199760</td>\n",
       "      <td>12516</td>\n",
       "      <td>0.999999</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset_sd</th>\n",
       "      <td>5.416193</td>\n",
       "      <td>2.985705</td>\n",
       "      <td>0.239002</td>\n",
       "      <td>11.240824</td>\n",
       "      <td>11559</td>\n",
       "      <td>0.999977</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Intercept</th>\n",
       "      <td>81.193962</td>\n",
       "      <td>3.156818</td>\n",
       "      <td>74.836674</td>\n",
       "      <td>87.899151</td>\n",
       "      <td>10785</td>\n",
       "      <td>1.000003</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>best_chosen_sd</th>\n",
       "      <td>7.379101</td>\n",
       "      <td>0.509492</td>\n",
       "      <td>6.441272</td>\n",
       "      <td>8.433300</td>\n",
       "      <td>33344</td>\n",
       "      <td>0.999987</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC</th>\n",
       "      <td>0.103594</td>\n",
       "      <td>0.056698</td>\n",
       "      <td>-0.013960</td>\n",
       "      <td>0.216177</td>\n",
       "      <td>5169</td>\n",
       "      <td>1.000144</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset[folke2016]</th>\n",
       "      <td>0.083300</td>\n",
       "      <td>0.063005</td>\n",
       "      <td>-0.037683</td>\n",
       "      <td>0.213422</td>\n",
       "      <td>7895</td>\n",
       "      <td>1.000177</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset[krajbich2010]</th>\n",
       "      <td>0.023669</td>\n",
       "      <td>0.074038</td>\n",
       "      <td>-0.122391</td>\n",
       "      <td>0.174967</td>\n",
       "      <td>12197</td>\n",
       "      <td>1.000027</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset[krajbich2011]</th>\n",
       "      <td>0.077069</td>\n",
       "      <td>0.062611</td>\n",
       "      <td>-0.043336</td>\n",
       "      <td>0.207148</td>\n",
       "      <td>8186</td>\n",
       "      <td>1.000194</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset[tavares2017]</th>\n",
       "      <td>-0.101397</td>\n",
       "      <td>0.057438</td>\n",
       "      <td>-0.213677</td>\n",
       "      <td>0.017578</td>\n",
       "      <td>4940</td>\n",
       "      <td>1.000161</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset_sd</th>\n",
       "      <td>0.113445</td>\n",
       "      <td>0.046964</td>\n",
       "      <td>0.039819</td>\n",
       "      <td>0.206234</td>\n",
       "      <td>11103</td>\n",
       "      <td>0.999998</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                                  mean        sd  hpd0.95_lower  \\\n",
       "1|dataset[folke2016]         -3.037726  3.692019     -10.930851   \n",
       "1|dataset[krajbich2010]       1.492997  3.258928      -5.091004   \n",
       "1|dataset[krajbich2011]      -2.438621  3.349371      -9.703574   \n",
       "1|dataset[tavares2017]        4.928456  3.863594      -1.695214   \n",
       "1|dataset_sd                  5.416193  2.985705       0.239002   \n",
       "Intercept                    81.193962  3.156818      74.836674   \n",
       "best_chosen_sd                7.379101  0.509492       6.441272   \n",
       "dWAIC                         0.103594  0.056698      -0.013960   \n",
       "dWAIC|dataset[folke2016]      0.083300  0.063005      -0.037683   \n",
       "dWAIC|dataset[krajbich2010]   0.023669  0.074038      -0.122391   \n",
       "dWAIC|dataset[krajbich2011]   0.077069  0.062611      -0.043336   \n",
       "dWAIC|dataset[tavares2017]   -0.101397  0.057438      -0.213677   \n",
       "dWAIC|dataset_sd              0.113445  0.046964       0.039819   \n",
       "\n",
       "                             hpd0.95_upper  effective_n  gelman_rubin  \n",
       "1|dataset[folke2016]              3.819971        15732      0.999990  \n",
       "1|dataset[krajbich2010]           8.290051        11386      1.000025  \n",
       "1|dataset[krajbich2011]           4.007553        12798      1.000005  \n",
       "1|dataset[tavares2017]           13.199760        12516      0.999999  \n",
       "1|dataset_sd                     11.240824        11559      0.999977  \n",
       "Intercept                        87.899151        10785      1.000003  \n",
       "best_chosen_sd                    8.433300        33344      0.999987  \n",
       "dWAIC                             0.216177         5169      1.000144  \n",
       "dWAIC|dataset[folke2016]          0.213422         7895      1.000177  \n",
       "dWAIC|dataset[krajbich2010]       0.174967        12197      1.000027  \n",
       "dWAIC|dataset[krajbich2011]       0.207148         8186      1.000194  \n",
       "dWAIC|dataset[tavares2017]        0.017578         4940      1.000161  \n",
       "dWAIC|dataset_sd                  0.206234        11103      0.999998  "
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "result_path = os.path.join('results', 'mixed_effects_models', '2_relative_model_comparison_c_best_chosen.csv')\n",
    "if not os.path.isfile(result_path):\n",
    "    m = Model(mfx_df)\n",
    "    c_best_chosen_fit = m.fit('best_chosen ~ dWAIC',\n",
    "                              random=['dWAIC|dataset'],\n",
    "                              categorical=['dataset'],\n",
    "                              samples=n_samples)\n",
    "    c_summary = c_best_chosen_fit.summary(ranefs=True)\n",
    "    c_summary.to_csv(result_path)\n",
    "else:\n",
    "    print('Previous results loaded from \"{}\"...'.format(result_path))\n",
    "    c_summary = pd.read_csv(result_path, index_col=0)\n",
    "c_summary"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Results loaded from: results/mixed_effects_models/2_relative_model_comparison_c_best_chosen_posteriorplot.png\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZ0AAAEaCAYAAADKaorCAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAE2hJREFUeJzt3X9w1PWdx/HXJwSmll8XL1FDsKDM\ncNWBk6SMcuOAnOXHiaNcD4saroNVODxPThjHS1YsoFUSWh209LiCd6WcI1zB3GG1AsEftP4Yq0fA\nisjkDrHUEHE5OQiRlvz43h/cfpvN7ma/m2Tfu8k+HzMZdr/fz/eb93747veVz/f73e86z/MEAICF\nvEwXAADIHYQOAMAMoQMAMEPoAADMEDoAADOEDgDADKEDADBD6AAAzBA6AAAz+Sm25/YFAIBEXLIG\njHQAAGYIHQCAGUIHAGCG0AEAmCF0AABmCB0AgBlCBwBghtABAJghdJCz9uzZI0lyLvrzbD19Pnv2\nbDnnNHDgwKQ1OOc0dOjQoCUDfR6hg5x39OhRPfnkk/7zUCjkPw6Hw1q5cqX/vHPAdOac0/PPPy/P\n89TS0qL8/MQ3/XDOyfM8NTU1JV0v0F8QOsg5zrmonfyll16qpUuXSpLWrl2rVatW6ZprrpEkXXTR\nRVqxYkXU8p7nqbKyMma9O3bs0JIlS6Kmtba2JqyjY5gBucJ5Xkq3U+Pea+jTIqMLSRowYIDa2tqi\npif6V5JmzJih2tramPXEa5vM8uXL9cgjj8StC+jDkg7ZU73hJ9Bv7Nq1K6X2u3fv9kc4Q4YMiZk/\nefJkNTU1BTpHM3nyZL399tuaNGlSSjUAfR0jHeSURCOd9vZ23XfffRo3bpwWLVqk5uZmHTx4UI2N\njbr55pslSYMHD1Zzc7O/rsrKSlVXV8cd9aRaCyMd9BNJRzqEDnLKgQMHNH78eEnSK6+8ouuvv96f\n13nHnywU4h1Wa2ho0MiRI/02eXl5frB1lpeX5y9H4KCfIHQAAGY4pwNkWrwr3aqrqzNQCZB5jHQA\nAL2Fbw4FAGQPQgcAYIbQAQCYIXQAAGYIHQCAGUIHAGCG0AEAmCF0AABmCB0AgBlCBwBghtABAJgh\ndAAAZggdAIAZvtoAyEKNzS3aVH8qUNvK0sI0VwP0HkY6AAAzhA4AwAyhAwAwQ+gAAMwQOgAAM4QO\nAMAMoQMAMEPoAADMEDoAADOEDgDADKEDADDDvdeAPo77tKEvYaQDADBD6AAAzBA6AAAzhA4AwAwX\nEgBGqvedCNx2/tjhaawEyBxGOgAAM4QOAMAMoQMAMEPoAADMEDoAADOEDgDADKEDADBD6ABZbs+P\nn+zWcr96bpPa29p6uRqgZwgdIIs9dE2Jpt65RKGyorjzQ2VFevYnP9aGhbNjpl9zy3y9//LP9MO/\nnh41Lz8/X865tNUMdIXQAeJwzqmmpiZq5+ycU0VFhf/80KFDcs7pqquu8qc1NjbKORezU3fO6em/\n+UbC8IgnVFakR3/VIEmaMOubcdtU1YU17447Y6aX3jhXknTVzG+o4eD+qHnLly8PXAPQ27gNDpDA\nnDlz5HmenHM6ePCg/3j16tWSzoeO53mSpLy8PLW3t6u4uNif5pzzH0vSwg3/Iel8mFTVhf3HnZXd\ndKvmb3kmatqtj67T7n+q1vS/rUxad+26Ks397j/GnRepacWKFUnXA6QDoQMEcMUVV0iSNm7cKEm6\n4YYbtHPnzph2O3bs0KxZswKvNxI+FjYsnB0VgkAmEDpANxQWFsbdgc+aNStqpJNM0JHOTx+6R7c+\nui5QbTPuCWnr8r/T3EeiRztH9r4Vc7iQEII1zukA3fDMM89E7cDHjBkTNT/ozryqLhzz882Hf+jP\n+86kkZKk/S9t85cJcl5o34tbJUm/rt2ukisn+OvzPM+vjcBBJrgUNzy2UqCbUv1qg031p3q9bWVp\nYeAagG5IOrxnpAMAMMM5HSDHBB1xMSpCOjDSAQCYIXQAAGYIHQCAGUIHAGCG0AEAmCF0AABmuGQa\nvaa8pjzq+eY5mzNUCYBsxUgHAGCG0AEAmCF0AABmOKcD9AC3lAFSw0gHAGCG0AEAmCF0AABmCB0A\ngBlCBwBghtABAJghdAAAZggdAIAZQgcAYIbQAQCYIXQAAGYIHQCAGUIHAGCG0AEAmCF0AABmCB0A\ngBlCBwBghtABAJghdAAAZggdAIAZQgcAYIbQAQCYIXSQkz788EM55+LOW7t2rZxzMfOHDRsWM+3l\n9d/vdg2hsiKFyoq6vfyDEy+OWT6yzq7W7XmeQmVFemz6lTHLNtZ/4D/ftm1bt2sDEiF0kJOuvPLK\nhPMWL14sz/PkeZ4fMi0tLTp9+rQ8z9OMGTMkSc45TVv0QLd+f6isSFV1YVXVhbsVPKGyIq36z+Oq\nqgtrxJBBUfMi662qC8dd9sGvXaSqurCW7T7o/+5IPT+4barfbu7cuSnXBSRD6CDnOOfkeV5Ky7z5\n5pv+41dffVWSUlpHcXFxyuFy7mxzoNHQTX91i1p+dzbQOrdu3aqlNW8lnD/84hGSzodQqn0EBEHo\nIOfU1tYGatfU1OQ/njp1qn/IrbW1NeGhuY7e3rZRobIifff6P1FjY2PCkceIr/5p1PNQWZFGDBmk\n1nPnuhyxRLzw78+pdl2VJGlA/kCdO9usqpnj44bV3Llztbb8+pjpkRHXqePHJEnfWvNM0tcHdEd+\npgsALKUyyhk2bFhU28jjkydP+ofeIjvrzsGw+sYy/W/jb5MGhiSd+fwPbUJlRRp77de1Z9cObao/\nlXCZjofl6v7rYz2y8TlJ0qPvnA+N0K73E46Q8gcOijsvUmvk9eTl5am9vV0jRozQsWPHkr4OIAhC\nBzmn4yglUQh1FU4XXnhh0uCq+HmdJOnE0cN64i8nKSQlDKDTnzX6jyNtRhUMVktLi8pX/7PGT58d\nd7lI26+NvUir9n7WZT0drfjlYf9x5/D5zXvv+uttaWmRJDU2NgroLRxeQ06JXCAQCY3Iv52DqKKi\nQpWVlaqsrIxavqSkJKVzHYVfGaOqurA8z9PmigX+9AuG/ZE23nurXlqzIu5yvznZrKq6sMZPnx13\nVBIqK9Kpzxq1YcHNGjToDxcSrJxyub44dVKPfv2KmPYR/7r0Wzod/jTuCO1H357lP87PP/836XXX\nXRf49QLJuBRPFnJmEQmV15RHPd88Z3OGKrFTve9EoHaVpYWB20rS/LHDuzy81t22qdRRWVoYqB3Q\nQdKTnYx0AABmCB0AgBkuJAA64fATkD6MdAAAZhjpAEiIUR96GyMdAIAZQgcAYIbQAQCYIXQAAGYI\nHQCAGUIHAGCG0AEAmCF0AABmCB0AgBlCBwBghtABAJghdAAAZggdAIAZQgcAYIbQAQCYIXQAAGYI\nHQCAGUIHAGCG0AEAmCF0AABmCB0AgBlCBwBghtABAJghdAAAZggdAIAZQgcAYIbQAQCYIXQAAGYI\nHQCAGUIHAGCG0AEAmMnPdAFAujU2t2hT/alAbStLC9NcDZDbCB0AvaKxuSVw2+LBA9NYCbIZoQOg\nVwQdTUqMKHMZ53QAAGYIHQCAGUIHAGCG0AEAmCF0AABmCB0AgBlCBwBghtABAJghdAAAZggdAIAZ\nQgcAYIbQAQCYIXQAAGYIHQCAGUIHAGCG0AEAmCF0AABmCB0AgBlCBwBghtABAJghdAAAZggdAIAZ\nQgcAYCY/0wUA3dHY3KJN9acCtZ0/dniaqwEQFCMdAIAZQgcAYIbQAQCYIXQAAGYIHQCAGUIHAGCG\n0AEAmOFzOgDMVe87EbhtZWlhGiuBNUY6AAAzhA4AwAyhAwAwQ+gAAMwQOgAAM4QOAMAMoQMAMEPo\nACk6+v5ehcqKdOjQIX/a6hvLFCor8n8ifl27Peq5pJjnXc2LPH/5R9/z11334k+j5o8YMkirZ5Xq\nzOfhHr0uwAKhA6Tg3NlmbalYoKq6sKZNmxY1r6ou7P9EbF/1gKrqwnLOSZLa29uj5qcisu4LS0ZF\nhdOxM+dU8dI+XTCsoMtAA7IBoQOkYMW1o1Xx0j5J0ieffJJ0J//VKTOini+beHGPaxhdOinu9AH5\n+Zqz4im98P1lPf4dQLoQOkAKZt67TB/+Ypckafv27VHz/ue3R/Tm5vVRQbTvxa0KlRXp3Xff1fq7\nbgo0yokcMksUaP/24CJ9aWj8r+CeOLtcb23ZEPTlAOa49xqyRmNzS6ZLSGrqnUv8MHj88cf96RU/\nr5MkXVu+SC2/O6sTRw+r8Ctj/JCZWFqo3zc3qbH+gC4ec4V2r6vSzMUPxf0dVXVhzR87XJvqT0UF\nz84fPCI5p9tWrU/XywPSjtBB1thUfypw2/lj4/+lbyESJPeXFuqX4daY+V57uzwveppzTlV1YYXK\nilRVF9bMxQ9p/44aTbhhTuDf+xd/vzxpm1BZkVbt/SzwOgFrHF4DUvTav6zR8cOH5JzTn916l6Tz\nO/sTRw/rFz9Zq9p1VSoaNSZqmXPnzkmSvvNavcIf/7fWzb8hpcDpymPLH9SzD3xbobIiFY4a41+0\nAGQjRjpAiv78rqWSJM/z/Fv0R0Y/192xWNfdsThmmYEDB0qSvjy8QF8eXqB7Nu2Iu+7O53wiz6fd\n/Q8J288fO1wjUxglApnESAcAYIbQAQCYIXQAAGYIHQCAGUIHAGCGq9cAZLWWdk8nzsZ+Hiqe4sED\n01wNeorQAZDVTpxtDfzB4crSwjRXg57i8BoAwAyhAwAwQ+gAAMwQOgAAM1xIgLSK3JssGU4AA7mB\nkQ4AwAyhAwAwQ+gAAMwQOgAAM4QOAMAMoQMAMEPoAADMEDoAADOETj/hnJNzTo899ljSNh0NGTIk\nZlrEBx98kHAekK2q952I+ols95dPvDZqekRkfudtfcyYMXLOaeXKlTG/I9EySI47EvQDzjl5nuc/\nXrZsWdI2kcdnzpxJ+MYZN25c3OmJ7jJw8OTvU64dSKdQWZGq6sL+43gi7wUp+r1x+PBhSYoJnY5t\nkDpGOv3Mxx9/rJqamh6vhzcW+puqurA23TevR+tobm7Wtm3btHPnTr3zzju9VFluYaTTz4waNUoL\nFy7UnDlzoqaXlJSktJ7a2treLAvICodeT7xdX3311Vq1alWXy99///1av369PM+T53n8cdYNhE4f\nFTkk1nmDb21t1ZQpU2LaNzQ0pLRu3kjojy4cOTru9KefflpHjhxRKBTqcvkpU6bowIEDksT5nG4i\ndPqoRKEwaNAgtbe3p2XdQF+2YcHNeuBn78ZM37t3r+6++261tbUlXUd5ebnmzevZIbpcR+j0A5Fh\nfmcdRyzV1dV+my1btkS16fgvgYP+pKouHPcCglBZkSr/f1ufOHGipNj3QMf31MMPP+xPHzBggD/v\nhRdeSF/x/RSh00/EC4uO0yoqKlRRURFouVTmA9kucvVaommJtvFE01tbW3unsBxF6ECNzS3aVH8q\nUFu+bA1ATxA6AHJWY3NLoHbFgwemuZLcQegAyFmM8O3x4VAAgBlCBwBghsNr/VQqFwfMHzs8zdUA\nwHmEDgAEkOhGt53NHzucCw+6wOE1AIAZl+KH//ikYAYF/UtLOv/XViqH17pqm+iW8DHmdHre85td\nAz0W78Oh0vkr0oK+p1Jpm+MjnaQ3pGOkAwAwwzmdDEvlLy0A2e+PvzQg8IdOI3Lp80KEDgD0sqAh\nIqV+9Whf/0OV0EmDvr5RdNbxmHhKx8GfS88xc6l3z1l1t450tZX6/+tLRbpeHzKjX4VOOoe0/S1I\nACAT0nb12pYtW3Ty5MnUK8pSJ0+eVEFBQabL6BPC4bCKigJe8Zbj6Kvg6KvgequvCgoKdPvtt6ey\nSPKvU41813fAn5y1YcOGTJfQZ6xcuTLTJfQZ9FVw9FVwGeyrpDnCJdMAADN8OBQA0Fv4cCgAIHsQ\nOgAAM4ROEs45Oec0ePDguPNLS0vlnNPUqVNtC8sCt9xyi98/iUTmd9UmF5w5cyZwX+W6ZH317LPP\nsl11EOmH48ePdzk/W/qK0OnChAkT1NbWJs/zNH36dLW0xH4OaN++fUrxvFi/UVNT41+REm+DXrJk\niY4dOybP87RmzZoMVJg9hg4d2mVfScrZ7aizZH112223+fM///xzzZw5MwNVZgfnnN8Xl1xySdw2\nHa8cy4bgIXS68N577ykv73wXbd++XYMGDcpwRdnj3nvvTbqTfOqpp1RcXCzpfAA98cQTFqVlpcsu\nuyzTJfQZyfpqwIAB/uOCggLV19enu6Q+oa8cbSF00C0bNmxIeZn169enoZK+oby8PNMl9Bmp9JVz\nTkeOHEljNX3HvHnzdOzYsbjzIofXsmE03a9ug9MdEyZMiDt9//79xpX0LZMnT055mSlTpqShkr7h\njTfeyHQJfUbQvsqWnWi2eP3117VgwYK48yL9lA19xud0utDxP+i1115TW1ubpk2bFrft1KlTtWfP\nHsPqMu/yyy/XRx99JCn+xlxSUqKGhgZJUlNTk4YOHWpeY7bo2D9dvfGzYaeQaUH6in46L+h2Fa99\nukpK1iDnRzpd6Xji7YILLtAXX3whKfY/OiLX3ggTJkzwX3/H1x0J4IaGhqj+yaW+6ez06dNx+yre\nthSvXS5J1lejR4/2n0fkal913Ed9+umn/vRE+6hs6CdGOgCA3sIdCQAA2YPQAQCYIXQAAGYIHQCA\nGUIHAGCG0AEAmCF0AABmCB0AgBlCBwBghtABAJghdAAAZlK94Wfmv3YOANBnMdIBAJghdAAAZggd\nAIAZQgcAYIbQAQCYIXQAAGYIHQCAGUIHAGCG0AEAmCF0AABm/g9WfIb+gv2qyAAAAABJRU5ErkJg\ngg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1343aef60>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot_path = os.path.join('results', 'mixed_effects_models', '2_relative_model_comparison_c_best_chosen_posteriorplot.png')\n",
    "if not os.path.isfile(plot_path):\n",
    "    plot_posterior(m.backend.trace, varnames=['dWAIC'], ref_val=0)\n",
    "    plt.savefig(plot_path, dpi = 330)\n",
    "else:\n",
    "    print('Results loaded from: {}'.format(plot_path))\n",
    "    show_img(plot_path)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": false,
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Previous results loaded from \"results/mixed_effects_models/2_relative_model_comparison_d_gaze_influence.csv\"...\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>mean</th>\n",
       "      <th>sd</th>\n",
       "      <th>hpd0.95_lower</th>\n",
       "      <th>hpd0.95_upper</th>\n",
       "      <th>effective_n</th>\n",
       "      <th>gelman_rubin</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>1|dataset[folke2016]</th>\n",
       "      <td>2.514997</td>\n",
       "      <td>3.523715</td>\n",
       "      <td>-3.850334</td>\n",
       "      <td>10.189233</td>\n",
       "      <td>15389</td>\n",
       "      <td>0.999987</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset[krajbich2010]</th>\n",
       "      <td>-1.871083</td>\n",
       "      <td>3.129804</td>\n",
       "      <td>-8.439211</td>\n",
       "      <td>4.084390</td>\n",
       "      <td>11903</td>\n",
       "      <td>1.000017</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset[krajbich2011]</th>\n",
       "      <td>1.149952</td>\n",
       "      <td>3.104042</td>\n",
       "      <td>-5.255793</td>\n",
       "      <td>7.622787</td>\n",
       "      <td>14329</td>\n",
       "      <td>1.000006</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset[tavares2017]</th>\n",
       "      <td>-3.049778</td>\n",
       "      <td>3.596354</td>\n",
       "      <td>-10.863824</td>\n",
       "      <td>2.893496</td>\n",
       "      <td>12576</td>\n",
       "      <td>1.000105</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1|dataset_sd</th>\n",
       "      <td>4.569710</td>\n",
       "      <td>3.274987</td>\n",
       "      <td>0.001060</td>\n",
       "      <td>10.837882</td>\n",
       "      <td>8581</td>\n",
       "      <td>0.999981</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Intercept</th>\n",
       "      <td>11.267927</td>\n",
       "      <td>2.976148</td>\n",
       "      <td>5.648708</td>\n",
       "      <td>17.665784</td>\n",
       "      <td>12336</td>\n",
       "      <td>1.000009</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC</th>\n",
       "      <td>-0.208423</td>\n",
       "      <td>0.109235</td>\n",
       "      <td>-0.415993</td>\n",
       "      <td>0.010731</td>\n",
       "      <td>5731</td>\n",
       "      <td>0.999987</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset[folke2016]</th>\n",
       "      <td>-0.085523</td>\n",
       "      <td>0.113878</td>\n",
       "      <td>-0.314477</td>\n",
       "      <td>0.132726</td>\n",
       "      <td>6219</td>\n",
       "      <td>0.999989</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset[krajbich2010]</th>\n",
       "      <td>-0.622183</td>\n",
       "      <td>0.133612</td>\n",
       "      <td>-0.884679</td>\n",
       "      <td>-0.361419</td>\n",
       "      <td>7433</td>\n",
       "      <td>0.999984</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset[krajbich2011]</th>\n",
       "      <td>-0.124041</td>\n",
       "      <td>0.113100</td>\n",
       "      <td>-0.351546</td>\n",
       "      <td>0.088626</td>\n",
       "      <td>6239</td>\n",
       "      <td>0.999977</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset[tavares2017]</th>\n",
       "      <td>0.136744</td>\n",
       "      <td>0.109802</td>\n",
       "      <td>-0.085081</td>\n",
       "      <td>0.345225</td>\n",
       "      <td>5767</td>\n",
       "      <td>0.999980</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dWAIC|dataset_sd</th>\n",
       "      <td>0.281091</td>\n",
       "      <td>0.071636</td>\n",
       "      <td>0.154776</td>\n",
       "      <td>0.425652</td>\n",
       "      <td>14851</td>\n",
       "      <td>0.999980</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>gaze_influence_sd</th>\n",
       "      <td>7.414348</td>\n",
       "      <td>0.511515</td>\n",
       "      <td>6.445630</td>\n",
       "      <td>8.429494</td>\n",
       "      <td>35581</td>\n",
       "      <td>0.999975</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                                  mean        sd  hpd0.95_lower  \\\n",
       "1|dataset[folke2016]          2.514997  3.523715      -3.850334   \n",
       "1|dataset[krajbich2010]      -1.871083  3.129804      -8.439211   \n",
       "1|dataset[krajbich2011]       1.149952  3.104042      -5.255793   \n",
       "1|dataset[tavares2017]       -3.049778  3.596354     -10.863824   \n",
       "1|dataset_sd                  4.569710  3.274987       0.001060   \n",
       "Intercept                    11.267927  2.976148       5.648708   \n",
       "dWAIC                        -0.208423  0.109235      -0.415993   \n",
       "dWAIC|dataset[folke2016]     -0.085523  0.113878      -0.314477   \n",
       "dWAIC|dataset[krajbich2010]  -0.622183  0.133612      -0.884679   \n",
       "dWAIC|dataset[krajbich2011]  -0.124041  0.113100      -0.351546   \n",
       "dWAIC|dataset[tavares2017]    0.136744  0.109802      -0.085081   \n",
       "dWAIC|dataset_sd              0.281091  0.071636       0.154776   \n",
       "gaze_influence_sd             7.414348  0.511515       6.445630   \n",
       "\n",
       "                             hpd0.95_upper  effective_n  gelman_rubin  \n",
       "1|dataset[folke2016]             10.189233        15389      0.999987  \n",
       "1|dataset[krajbich2010]           4.084390        11903      1.000017  \n",
       "1|dataset[krajbich2011]           7.622787        14329      1.000006  \n",
       "1|dataset[tavares2017]            2.893496        12576      1.000105  \n",
       "1|dataset_sd                     10.837882         8581      0.999981  \n",
       "Intercept                        17.665784        12336      1.000009  \n",
       "dWAIC                             0.010731         5731      0.999987  \n",
       "dWAIC|dataset[folke2016]          0.132726         6219      0.999989  \n",
       "dWAIC|dataset[krajbich2010]      -0.361419         7433      0.999984  \n",
       "dWAIC|dataset[krajbich2011]       0.088626         6239      0.999977  \n",
       "dWAIC|dataset[tavares2017]        0.345225         5767      0.999980  \n",
       "dWAIC|dataset_sd                  0.425652        14851      0.999980  \n",
       "gaze_influence_sd                 8.429494        35581      0.999975  "
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "result_path = os.path.join('results', 'mixed_effects_models', '2_relative_model_comparison_d_gaze_influence.csv')\n",
    "if not os.path.isfile(result_path):\n",
    "    m = Model(mfx_df)\n",
    "    d_gaze_influence_fit = m.fit('gaze_influence ~ dWAIC', \n",
    "                                 random=['dWAIC|dataset'],\n",
    "                                 categorical=['dataset'],\n",
    "                                 samples=n_samples)\n",
    "    d_summary = d_gaze_influence_fit.summary(ranefs=True)\n",
    "    d_summary.to_csv(result_path)\n",
    "else:\n",
    "    print('Previous results loaded from \"{}\"...'.format(result_path))\n",
    "    d_summary = pd.read_csv(result_path, index_col=0)\n",
    "d_summary"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Results loaded from: results/mixed_effects_models/2_relative_model_comparison_d_gaze_influence_posteriorplot.png\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZ0AAAEaCAYAAADKaorCAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAFFtJREFUeJzt3X10VPWdx/HPTUiUkwBSGSpY5EEb\nbAGF8Ykqwcej1K0KPmyRuoctVeyxKuDWNVMsRFQSW7QetXtKwId0hW5FKlFbRMvKCipWCWDR2mwL\n1AqxTgrlIaImYfYPdq7zmJlJMt95er/OycnM7/7une9cZubD73dv7jiBQEAAAFgoynQBAIDCQegA\nAMwQOgAAM4QOAMAMoQMAMEPoAADMEDoAADOEDgDADKEDADDTK8X+XL4AABCPk6gDIx0AgBlCBwBg\nhtABAJghdAAAZggdAIAZQgcAYIbQAQCYIXQAAGYIHRSsdevWSZIcJ/zv2bp7/4orrpDjOCopKUlY\ng+M46tOnT7IlAzmP0EHBe//99/Xggw+6930+n3vb7/erurravR8ZMJEcx1FDQ4MCgYDa2trUq1f8\ni344jqNAIKADBw4k3C6QLwgdFBzHccI+5IcMGaI5c+ZIkh5++GEtXLhQZ511liRp4MCBmj9/ftj6\ngUBAVVVVUdtdvXq1Zs+eHdbW3t4et47QMAMKhRMIpHQ5Na69hpwWHF1IUnFxsTo6OsLa4/2WpIsv\nvlgvvvhi1HZi9U1k3rx5WrBgQcy6gByWcMie6gU/gbyxZs2alPq/9NJL7ginvLw8anllZaUOHDiQ\n1DGayspKbdy4UePHj0+pBiDXMdJBQYk30jl8+LBmzZql0aNH68Ybb1Rra6veffddNTc36/LLL5ck\nlZWVqbW11d1WVVWVamtrY456Uq2FkQ7yRMKRDqGDgrJt2zaNGTNGkrR27VpdcMEF7rLID/5EoRBr\nWm3Xrl360pe+5PYpKipygy1SUVGRux6BgzxB6AAAzHBMB8i0WGe61dbWZqASIPMY6QAAegrfHAoA\nyB6EDgDADKEDADBD6AAAzBA6AAAzhA4AwAyhAwAwQ+gAAMwQOgAAM4QOAMAMoQMAMEPoAADMEDoA\nADN8tQHyXnNrm+qb9iXVt2rcgDRXAxQ2RjoAADOEDgDADKEDADDDMR3kpFSO00yv6JfmagAki5EO\nAMAMoQMAMEPoAADMEDoAADOEDgDADGevIWvUbm5Jui9npAG5iZEOAMAMoQMAMEPoAADMEDoAADOE\nDgDADKEDADBD6AAAzBA6AAAzhA4AwAxXJEBaJXuVgapxA9JcCYBswEgHAGCG0AEAmCF0AABmCB0A\ngBlCBwXnvfUv6qPtTXGX7z74mX54b62aW9uifqrmL3BvBzU1NWnJkiUWpQM5j9BBQfF5PTq58mJ5\nhn9ZPq8nanlVVZV+/r/7VXbp9RpcXqqfbvyT6pv2qb5pnwaXl+qYyTdpcHmp6pv2SZIcx1FFRYVu\nuOEGOY5j/XSAnEPoIG2GDRsm6cgH/X3/5HVvzz9naFg/n9cT9YFdPXGEfF6PXql/JKxf8PfN35ne\nrdriBURtba17u6bRr3sv+qr7mDWNfrc9VmABSIy/00Fa/eTqCapp9Gtt3SLNO3uoahr9+vObG9Tc\ntE2DKka7H+ZV4wbIcRz3g736le2SpH0fNat1799V1v9YSdInB/arptGvByaN1hl3HHmMeAEQ3Fao\nkRMuSrr2O88crHve2N1pH7/f7waY1+tNettAoSJ0kFZznt4gSbpw5vd14czvS5JOPGOC6m64QjOX\nNEg6Ehq+iPVCg+TYE0bo+6vekCQd3aevJGnr9r/qwXW/V//BJ8QMl6AXHlqgD75wlP7c3lvn/uut\nSdf95qondXR5XxWXlHTaz+PxKBAISIo/egLwOUIHGTX6wm/oWz9+XFXjBrhXLwidypLkTs3F09lI\nZ9Kt8zS9op97DOaPG36bsKbmpm165u7btHDTR8k+DUnSmjVrtHHjRo0fPz6l9YBCQuggo7atfd69\nfXCPX+Vf8GjgiJH68E9/0HEnfUU+r0fHDBrS6TY6G+nEExydBAWD7sCBA3po6vlR2wwexwn9HemS\nSy6J2i6AcE6KbxLeUUhJKtdeS7avpLDRS0/2TaUOrhcHREk4x8zZawAAM4QOAMAMoQMAMEPoAADM\nEDoAADOcMg10A2e6AalhpAMAMEPoAADMML0GoEumrZwWdn/5VcszVAlyCSMdAIAZQgcAYIbQAQCY\nIXQAAGYIHQCAGUIHAGCG0AEAmCF0AABmCB0AgBlCBwBghtABAJghdAAAZggdAIAZQgcAYIavNkDK\n+LZMAF3FSAcAYIbQAQCYIXQAAGYIHQCAGUIHAGCG0AEAmCF0AABmCB0AgBlCBwBghtABAJghdAAA\nZrj2Grrl57Ov0x9eWaM5K1/TwOFfDlvmOI4kqabRH3d9n9cTttzn9ajPgIH6wYvvuG1Pz79FV9/1\ncJfq83k9urO4WPe8+WHK6772X0v13I98kqR73tit4pKSmP0WTT5Lf39/u06ddKWmLlzsPu6I08/R\nDXWr3H4zlzZoxGnnJPXYXLcO+YrQQZfdfcFI3d7wpo7u01c1k06R74W33WW1Xz9VgUCg04uDblv7\nnLyXfdO9/+PLz4gKqMhQSkVw3ekV/TS4vDTl7Zw99XqdPfX6TuuI1R5s++u2RrfNcZwuPw8gnzC9\nhi77+B97dHSfvpIk3wtvy+f1HLnt9ahq9daE6y+7fUbY/T0f7HRv7939V7V9ciilD+ofXXa6W4Mk\nDRlzesJ1tq19Xj6vRxuW/SzpxwlV1Cv+/9uGjPZKkv746m8VCAS6tH0g3xA66BHvvPybsPs+r0ef\nffZZWAhELo8MlJpGv3xej3xej/oPHqJ5Z5+Q8HEbn39KPq9H8ycM078/95a7zeambTr/O3Pirhd8\nnNEXfkM1jX5N+NZ3O+139xu7opYdOnRIvfsco40rHtcLDy1wn2vo85CkJ265NuHzAAoFoYMuC/1w\nDR2lBJeVlpa6fUK1fXJIs375P3G3WdPo1wNXfi1s3VjhVXfDFVox73uqafTrrg3hj19+7EAd3PNR\nzMfweT2aMvf+pEZRwXrmnzM0alnv3r3VurdF46/5tibdOi/mesFwDR7fihfCQKHgmA66JfSDu1fp\nUUmt8/gtU7Vj02vu/cbnfqmqiOmnYwYNSbidmUsaJEn7PmpW7aRTwurpc+xA/WrBHJ0x+bq4NQcD\noHr9Dh1VVt7pY8145CntaHxdw71fS1hXqKk1dSn1B/IdIx10S9Nr/60tq1fK5/Xoa9/8jqTPR0D7\n9+8Pm0YLfsjPXNLgjgS8l30z5oH4GT99Kuka+g0c5G5v5V2zwpb9/qUGjR46WJX/8r2o9YLrtH/2\nacxjOj6vR3t2/UVvPvOfWvrdK93AiRytbFj2Mz2/6M6Y6596yZSknwdQCJwUD3ByNBQpfV11OvpK\n0vSKfqpv2tfjfdNVc6rPLxdOmZ62clrY/eVXLc9QJcgiTqIOjHQAAGYIHQCAGUIHAGCG0AEAmCF0\nAABmCB0AgBlCBwBghtABAJghdAAAZggdAIAZQgcAYIbQAQCY4asNoObWtpQuiAkAXcVIBwBghtAB\nAJghdAAAZjimA2ShtsMBtRxqT6rvoLKSNFcD9BxCB8hCLYfaObkDeYnpNQCAGUIHAGCG0AEAmCF0\nAABmCB0AgBlCBwBghtABAJghdAAAZggdAIAZQgcAYIbQAQCYIXQAAGYIHQCAGUIHAGCG0AEAmCF0\nAABmCB0AgBlCBwBghtABAJghdAAAZggdAIAZQgcAYKZXpgtAerQdDqjlUHumy4CBVP6tB5WVpLka\noHOETp5qOdSu+qZ9SfWdXtEvzdUgnVL5t64aNyDN1QCdY3oNAGCG0AEAmCF0AABmCB0AgBlCBwBg\nhtABAJghdAAAZggdAIAZQgcAYIbQAQCYIXQAAGYIHQCAGUIHAGCG0AEAmCF0AABmCB0AgBlCBwBg\nhtABAJghdAAAZggdAIAZQgcAYKZXpgsAYKt2c0tS/arGDUhzJShEjHQAAGYIHQCAGUIHAGCG0AEA\nmCF0AABmOHsthyR71pEkTa/ol8ZKAKBrGOkAAMwQOgAAM4QOAMAMoQMAMEPoAADMEDoAADOEDgDA\nDKEDADBD6CDr+Lwe+byemG0+r0d3XzAyrP2uiScmvd1494PbHjJkSMzHXXb7t1N9GgBi4IoEyCo+\nr0c1jf6o25LCbkvS0huvdNs62ttV3KtX1DqpqGn0q2rcADmO425j+Glna+aSBklSdeXwI7/X7+jS\n9gEw0kEO277pVff2X7a+oUAgoIWbPur2dk/9+lXa7/8wqr16/Q592nqw29sHChkjHeSE4748Sh3t\n7brzzEE67qSvavqWLVr41t/cKbKaRn9So5zIKbZYtq5eqan3/qxH6gYQjtBBVrntV6+7wXDsCSPc\n9lm/XCfp83AJCoZMcKotGDzxAii0LTKAXnhogcpPq+jy9ByAxAgdZBXPsJPCjukk6+Ce7gfFpFvn\n6eZxA1K6mne+62xfvLv3U/f2V/sfZVEO8gDHdJBVfnD6F/XJgf2688zBmnTrPLf9se/9sz79+KB8\nXo8GVYwKW8fn9Wj2ivVpqWfPrr/ohYcWqLpyeLdOUgBwBCMdZJWFb/1NknTP73aHtc/46VOSos9g\ni2wL3k7UL956kap+s0WSwgIQQNcx0gEAmCF0AABmCB0AgBlCBwBghhMJMqy5tS2pfoPKStJcCQCk\nH6GTYfVN+5LqVzVuQJorAYD0Y3oNAGCG0AEAmCF0AABmCB0AgBlOJADQI5I9E1PibMxCRugA6BHJ\nnokpcTZmIWN6DQBghtABAJghdAAAZggdAIAZQgcAYIbQAQCYIXQAAGYIHQCAGUInhOM4chxH9957\nb8J+sdaL9Oyzz8ZdBuQ7n9fj/kSq3dyi2s0tn7/nftfsth1z3PFyHEfX3PWwaje3uOvwXsoPXJHg\n/zmOo0Ag4N6eO3du3H6RAoFAWHvwjfKIb75qGv1hbaH4q2zkq3de/o2uu79eo86/VLvfe1uvLl+s\nc6bdGNbH5/W474/Q21Wrt2rv7ve1/a1Xw/pHvs+QmxjpxLBz506tXLky5rL33nsvqW34vB7d/ORL\nennpAz1ZGpATnvy36Rp1/qWSpMEnn6LnF93Zaf+i4mKLspAFCJ0Yhg4dqsWLF0e1O46jkSNHJr0d\nn9ej86+/Tf8x/etatfD2niwRyCunXDwl0yXASEGHTrw54vb2dk2cODGs7Z133lFra2tK2w9OF9xU\nv1pvPP1El+sE8t3OzRszXQKMFPQxneAxnEilpaU6fPhwWNuoUaNS2vYZU67rcl1Arjt10pVhx2CG\njj2z0/7/+PADi7KQBQo6dELFO0gZeoJBPMH1gr9rGv268oc/CTtrJzjqAQrB1IWLY77+Q08YmL1i\nvdvnK+dOcvuGrvd09a2qCjnBJ/R3ovclshOhEyLWiziZtsj7wTPVCBoUsliv/9C2L554csI+oQiZ\n/FDQx3QAALYY6QDIas2tbUl/Kyl/+5b9CJ0kxfrjzlh40QNAfEyvAQDMMNIBkFeYlchujHQAAGYI\nHQCAGUIHAGDGSfEPrgr2r7NSmSeO7Pu3P7+nB6+pTEdZQOZcFXE/4sLss1es1xdPPLnbDzO9ol9K\np0xzTCejEn73RMGeSMC5/wAIKHtMrwEAzBA6AAAzBTu9ZinehQ1TnX9Otq+U2jx4OufMeX6p95Xy\n//nlomSf3/SKfhpUVpLmanJXXoVOc2tbpksAUOCOPbo4pc+iQguotJ299otf/EJ79+5NvaIesnfv\nXvXv3z9jj58qv98vj8eTuGOWyKX9m2v7NpfqzaXXgZRb+1bK/P7t37+/rr322lRWSXj2mgKBQCo/\nOaOuri7TJaSkuro60yWkJJf2b67t21yqN5deB4FAbu3bQCD39m8giRzhRAIAgBn+OBQA0FMSTq8x\n0gEAmCF0AABm8jJ0HMeR4zh67LHHOl0e/Mm0ZOvIplrLyso6Xe44jhYtWmRcXfx6YtmzZ09WvQ5+\n/etfd1rLM888k9Fak9lXweW33HKLYWWd19LV5ZauvvrqTusZM2aMu/z11183rq5n5d0xHcdxFHxO\nobfj9cm0ZOoNLpOU0brHjh2rxsZGFRUVafLkyVqxYoVKSuL/jUGm93OifdvR0aHi4uKY/TMhldeC\ndZ2rVq2SJE2ePFmbN2/WK6+8olmzZsWti32ZmlT2XbbUHEdhH9NZt26dNm3aFHPZfffdp46ODuOK\nuuakk07KihfZ1q1bVVR05CWzatUqlZaWZrii7gkNnExraWlRfX19psuIa8qUKZo8ebIkady4cZo9\ne3an/dvb23X//fdblJbzbr755qx4f1vJ69A599xztXjx4qj2QCCgO+64Q+vXr8+a4bUkHX/88Zku\noUcEpwEeffTRTJfiSrRvM/2/x7q6Ok2bNi1jj9/TiouLY773EK2uri7TJZjKycvgjB07Nqrtmmuu\n0dy5c8PaduzYoYkTJ8bdznnnndfTpUXZvXu3Lr300qj2LVu2RLXt2rUrqs36wzDWvpVi1xtP6DTB\njBkzeqSueGLVm+y+DbLax53t28rKSm3YsMHkNWmls/cePldZmfx3bWX6P0c9oSCP6cTrnwm5NJcb\n+vgvv/yyOjo6dNFFFyXVPxOSPb63b98+9e3b17q8mLVk63GIadOmadmyZe7MwIQJE7Rhw4a4dZWX\nl+vgwYOmNcarJdv2ZSwjRozQ9u3bJeXGsehOFN6XuAUCAfeNsXTpUrc98kUY2j+TQusNrSUbX2Ch\ntfbu3Vsff/yxpPj7tqGhwb7IEIn27RNPPCFJ6tevX9g6mfLss892Wu/OnTs1fPhwt23Hjh0aNmyY\nSW3Lly+P+b4J/bcP3d833XSTSV3xJPO+Ci7PlpN0YtVx3nnnad26dVlVa3fl3UgHAJAxhX32GgAg\nuxA6AAAzhA4AwAyhAwAwQ+gAAMwQOgAAM4QOAMAMoQMAMEPoAADMEDoAADOEDgDATKoX/MyeL58B\nAOQcRjoAADOEDgDADKEDADBD6AAAzBA6AAAzhA4AwAyhAwAwQ+gAAMwQOgAAM4QOAMDM/wENKB9Q\nu0koHwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x134416fd0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot_path = os.path.join('results', 'mixed_effects_models', '2_relative_model_comparison_d_gaze_influence_posteriorplot.png')\n",
    "if not os.path.isfile(plot_path):\n",
    "    plot_posterior(m.backend.trace, varnames=['dWAIC'], ref_val=0)\n",
    "    plt.savefig(plot_path, dpi=330)\n",
    "else:\n",
    "    print('Results loaded from: {}'.format(plot_path))\n",
    "    show_img(plot_path)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Plot DIC differences between full and no-gaze bias GLAM variant"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def figure_4(waic, colors=None, fontsize=7, figsize=cm2inch(18,6), plot_type='bar', log_scale=False, ylim=(-100, 5)):\n",
    "\n",
    "    # Log transform: add minimum dWAIC + 1, then log-transform, new indifference then at log(minimum dWAIC + 1)\n",
    "    waic = waic.copy()\n",
    "    if log_scale:\n",
    "        min_dWAIC = np.min(waic['dWAIC'])\n",
    "        indifference_level = np.log(np.abs(min_dWAIC) + 1)\n",
    "        waic['plotted_dWAIC'] = np.log(waic['dWAIC'] + np.abs(min_dWAIC) + 1)\n",
    "    else:\n",
    "        indifference_level = 0\n",
    "        waic['plotted_dWAIC'] = waic['dWAIC']\n",
    "    \n",
    "    if colors is None:\n",
    "        colors = np.array(['C{}'.format(i) for i in range(10)])\n",
    "\n",
    "    datasets = ['krajbich2010', 'krajbich2011', 'folke2016', 'tavares2017']\n",
    "    dataset_labels = ['Krajbich 2010', 'Krajbich 2011', 'Folke 2016', 'Tavares 2017']\n",
    "    \n",
    "    fig = plt.figure(figsize=figsize)\n",
    "\n",
    "    ax0 = plt.subplot2grid((1, 3), (0, 0))\n",
    "    ax1 = plt.subplot2grid((1, 3), (0, 1), colspan=2)\n",
    "\n",
    "    sorting_idx = np.argsort(waic['plotted_dWAIC'])\n",
    "    color_idx = pd.Categorical(waic['dataset'], categories=datasets).codes\n",
    "    patches = []\n",
    "    \n",
    "    for i, dataset in enumerate(datasets):\n",
    "        \n",
    "        if plot_type == 'bar':\n",
    "            ax1.bar(np.arange(len(waic))[waic['dataset'][sorting_idx] == dataset],\n",
    "                        waic['plotted_dWAIC'][sorting_idx][waic['dataset'][sorting_idx] == dataset],\n",
    "                        color=colors[i])\n",
    "        \n",
    "        else:\n",
    "            ax1.scatter(np.arange(len(waic))[waic['dataset'][sorting_idx] == dataset],\n",
    "                        waic['plotted_dWAIC'][sorting_idx][waic['dataset'][sorting_idx] == dataset],\n",
    "                        color=colors[i], marker='+')\n",
    "        \n",
    "        bars = ax0.bar(i, 100 * np.mean(waic.loc[waic['dataset'] == dataset, 'dWAIC'] < 0), color=colors[i])\n",
    "        patches.append(bars[0])        \n",
    "\n",
    "    # Legend\n",
    "    ax1.legend(patches,\n",
    "               ['Krajbich 2010', 'Krajbich 2011', 'Folke 2016', 'Tavares 2017'],\n",
    "               ncol=2,\n",
    "               loc='lower right',\n",
    "               handlelength=1, handleheight=1,\n",
    "               fontsize=7, columnspacing=1,\n",
    "               frameon=False)\n",
    "    \n",
    "    # Labels    \n",
    "    ax0.tick_params(axis='both', which='major', labelsize=fontsize)\n",
    "    ax0.set_xticks(np.arange(len(datasets)))\n",
    "    ax0.set_xticklabels(dataset_labels, fontsize=fontsize, rotation=45)\n",
    "    ax0.set_yticks(np.arange(0, 110, 20))\n",
    "    ax0.set_ylim(0, 100)\n",
    "    ax0.set_ylabel('% Participants\\nwith gaze bias', fontsize=fontsize)\n",
    "    \n",
    "    ax1.set_xlabel('Participant', fontsize=fontsize)\n",
    "    if log_scale:\n",
    "        ax1.set_ylabel('Log(dWAIC))', fontsize=fontsize)\n",
    "    else:\n",
    "        ax1.set_ylabel('WAIC(with bias)-\\nWAIC(without bias)',fontsize=fontsize)\n",
    "    ax1.set_xticks([])\n",
    "    ax1.tick_params(axis='both', which='major', labelsize=fontsize)\n",
    "    \n",
    "    \n",
    "    ax1.axhline(indifference_level, color='black', lw=1, zorder=9)\n",
    "    \n",
    "    ax1.set_ylim(ylim)\n",
    "    \n",
    "    # Labels\n",
    "    ax0.text(-0.4, 1.05, 'a', transform=ax0.transAxes,\n",
    "            fontsize=fontsize, fontweight='bold', va='top')\n",
    "    \n",
    "    ax1.text(-0.175, 1.05, 'b', transform=ax1.transAxes,\n",
    "            fontsize=fontsize, fontweight='bold', va='top')\n",
    "    \n",
    "    \n",
    "        \n",
    "    fig.tight_layout()\n",
    "    \n",
    "    return fig"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": false,
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfsAAACkCAYAAACQJIdEAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAIABJREFUeJzt3Xuc1nP+//HHs3SakCJfRa3TOifa\nvljn8yHKqeS0Cqlstp8dspViiK+EJI0OjsU6JG2Uld1lHRcrKhORopSECBEaev3+eH+u6Zprrjl0\nHWaumXndb7e5Xdf1Ob3f13Q1r+v9eb/fr7fMDOecc87VXQ1qugLOOeecyy4P9s4551wd58HeOeec\nq+M82DvnnHN1nAd755xzro7zYO+cc87VcR7sXa0hqbekVTVdD+ecq20qDPaSTpP0nqQVkiZWV6Wc\nq0AzSc9K+kzSgTVdGeecqw0qa9kb0BAQ0FfSNtmvknMVagx0B54DBtVwXZxzrlaoLNj3BV4Ezo9e\nb57JwiXtLOleSdOi1+dKulvSFEnNo5/J0bbzMlm2q9Us+nHOOVcFm1Wy/wlgNNA6G4Wb2UfAxbFg\nD5xuZj0knQKcEW2bZmYzJT0G/DUb9XC1ynrgb0AH4PQarotzztUKFQZ7M7sPuK+a6gIbW2vLCH/M\nAYqix1+TnSCpL+EOBHvttdfv3n333axWsKZ1mNyh8oM2UVGvojLbFu6xZ8bL2fP9hck2q6rnm9kD\nwAMZqk61O/HEE2327Nk1XQ3nXN1Spb+huToavz2wIvrZIdqWtK5mNsnMOptZ52bNmlVX/ZzbZKtX\nr67pKjjn6qnKbuNnlaStgRuB/SUNAWZIGg80AwZEh42TdDIws4aq6Vy5JDUH7iJ0L7xgZt7V5JzL\nOTUa7M3sK6B/wuaHE15fWE3VSU9Biyxc89vMX9Nl2hn4uBLnXJp2HPw0AEtHnpyV69dosHeuDtiB\nSsaVOOfql1jgLs/SkSez4+CnSz1mmwd759ITG1cyjyTjSuIHkEavq69mzrmcpJuTPyY+rwqzqs1C\n9mDvXHqmU8G4EjObBEwC6Ny5s82ZM2eTLl7erIg931/Iwj32LPNYFcnOjc3ymHrTL+GYs1cC0GGn\n9qW3v78wdFkVfEuHyR0o6lVUUu5ZQ0r/OYmdE9te1KuozLnlzS5JrEfiOSXXAhY+2rbC8uPfZ6z8\n2Lk7/hR6DZ+ZcWWp8mLXjP1uinoVbWyBJZRbUsdI7HcW/34T/41idV27cCQAW+w5uNTrpSNPLvN+\nE1uLsXOKPv6k1HspqWPc7yzx/S5tem6p1yXXinuf8Y/JzilTj7j3m/i7Siy3jNi/ySacG/9+y6tz\nTGLdS65Vj7pKPdg7lwYz+4HaMq7EuXqi7G3xENSXJryuT3J16p1zzjm3Saqj77u28mDvnHPO1XEe\n7J1zztVq3qKvnPfZO+ecyz0FFfezL8VtCg/2zjnnqleyUfCx4O6t9KzwYO+ccy6rNgbw+jcKPld4\nsHfOOZcZia3zejSPPdd5sHfOOZeyZEtku9zjwd4559wm87712iWngr2k9sBY4GtgEfAJcBTQBLg0\nylbmXJVI2gzYH9gSWGRmy2u4Ss7Veh7ka6ecCvZAB8JyoQ9Fy4X+zsx6SDqFsJTogzVbPVdbSLoe\naAt8DKwBukjaGrjDzObWaOWcq0V8cF3dkGvB/nVgmqSLCIH9pGj7MsIXAeeqqtDMPk/cKKl1TVTG\nudoiFty9BV+35FqwvxC41sxekjQN2BBtb09YSrSM+CVE27dvX2Z/ZesKp8L/E+Q+M/tcUjdgFlAI\nLDCzQjP7soar5lzO8UF2dV+uBfvZQIGkcwkJkt6WNB5oBgxIdkLiEqLVVE9XOxwCfAc8Bxyc7sUk\n9QEOAFoBN5jZPEkjgTxgnZkNTrcM55zLhpwK9ma2AOiesPnhZMc6F0/SqWb2ZMLmdsDlwAXA4emW\nYWb3APdI2h/oKulroJGZDZR0i6R2PgjQOZeLcirYO5eGI4HEYN8LaGxmP0i6cVMuJqkDcFPC5osI\nM0UGAtcAOwCx4P5Jwmvncp53SdYfHuxdXTE+ybYLgLMkNSKs8HhkVS9mZkXAKfHbouvcBYwxs+WS\nRAjwEO4izEi8TmVjSpyrEZ7Zrt7xYO9qrWgq3UWE/niTtAZ4xMz+FR2yH/CamRVIuioDRd4M7AZc\nKuk5M3tcUrGk0cDPyW7h+5gSlxPKrCDn6puUgr2kXsCjwAPAB2ZWUJV9zmVYPuGzdpuZbZC0BXCq\npBPM7FngG6CJpPOBTukWZmb5SbYNSfe6zjmXbam27PcBjgamRI9V3edcxpjZ1ZIaRYF+F+ArM3so\n7pBYP30XwEfKu/ojasn7lDoXk2qw3xLoB5xNQr9mJfucy7T/kzQGuAH4FTgfQNLlhAx6AAIOAjJx\nK9+5nOSD7VxFUg32U8zsVQBJj8bvMLN+seeSBqZRN+eqYkvgVMLI+fPjtj8GNAZi/eSq5no5Vz18\nsJ2rgk0O9pL2Bs6JBkOJMEDq5bj9g4HTgR+AbYB9M1NV55J6AdjZzN6R9GFso5l9JukPwGmETIxP\nEtIuO+dcvdMghXO6A7tHj6cTWlDxtgNmm9nRwOT0qudcpZ4BFkm6APg5Yd9BZnammfUg3MZ3zrl6\naZNb9mZ2naSmhGlNTdh4mzRmNdBA0nC8Ve+y73bC5/ifwLGEgaGxO1Am6bAarJtzzuWEVPvs7wFe\nA34hBPuXYjvM7AaAKKXomHQr6FwlvgBkZlMk7RC3vTvwJRtnhPgCOK5O2DgQzwfkuapLNdi/ZWaF\niRsljQJ2jL0kfBE4K8UynKuKF4FfJc0A3o1tNLPraq5KzmVObPqcj7Z36Ug12F8o6VhgHWBmdhbh\nyVUAURpRzMwzhrmskdSEMEAP4BXKdik555wjxWBvZuX2xUsaBPwvob90jpndkmrlnKvEUKA9pYP8\nRbEnkg42s/9Ezw8zs5dxrpbwhDguk1JNl9sNOA/YAlhvZqfF7d421tKXdOsmXrcBMIIwd3oOUAwc\nRRgIeKmZ/ZBKfV3dZGbXSsoDukWbSla9iwbonSvpG5JMEXUuF/mtepctqUy9AzgJeI8wh/n12EZJ\nJwNtJP0hmuO83SZe91TCKmLFwArgdDO7BJgKnJFiXV3ddjfh81IcPY+pbIqoczmjqFeRt+RdVqXa\nZ/8V0Bw4gPAHNWYbwhSo2JeIf27idXcH/mNmEyVNIyRDgZAMpUOKdXV12xIzewJA0n6xjdEU0SPi\njvP+fJdzvCXvqkuqwb6QMO3ufGBCbKOZpZtEZwWwPnr+KxtTnLaP9pXh64XXX5IeB3aXdAjhs7I1\nMDzukD2ix+bAocRNEU2jzDbAq0A3M1sgaSSQB6wzM19sxzmXk1IN9nsBfYCGwDuZqw7TgTujRCgv\nAWskjQeaAQOSneDrhddfUWa8ivZPjD2PzRDJgKuAx6NrtgcamdlASbdIapdsTXvnYvZ8fyEAfsPe\nVbdUg/2JwLmEYH8P8Fxsh6S2QCczmyVpPzObV9WLmtk64OKEzQ+nWEdXz0Utf4t+kt4ZquDcDoTF\ndeK9DEwDjotebw/EgvsnhPEmHuydczkn1WD/NdAmOn+FpLwoUENYanQNMIswAtpXvnM1orKWfyXn\nFpGwRHN0l6kNIc/+1sAthAAP0A6YkXgd72ZyzuWCVIP9bsCNca/HsXF+8zeAr7noqoWkfQhBuTGA\nmV0ft+8SQivcgH+b2YSkF6kiM7s0um4BMM3MPpFULGk08HOyW/jezeScywWpLHHbhNB/HusDtbhW\nPcAi4ILoj/B/06+icxUaARQAPyXZt29czoeMrdNgZgVxz4dk6rquDorWmo/11TtXU1Jp2Q8l3LKE\njfnvL4rbPxt4MNqXl1btnKvcfKDIzDYk2dcibvpdq2qsk3PO5ZRUlri9VtKeZrYQQNLuCYfMBu4y\ns7GSrgKuyUA9nSsjGoDXGHhK0o/ErdMQuRw4J+65c87VS6lm0IsfMd8rYd/jhEF7hYQ0t85lRTQA\nb5SZnRI9H5dwSFtCSudtgMuqu37OOZcrUh2g11pS0+h5Ykrc5WY2XdJiSqcvdS6jJO0FnCNpTbTp\nQkonzhkBXAv8XN11c/WXp711uSjVYD+GjYH89oR990nqTOivH5RqxZyrgh5szH9vlM1/Px9YUE5/\nvnOZVeCTkFzuSmU0/k6EltLIaFPidKJCwu3TNwlL3aadotS5cjxISHITU/JZrEJ/vnPO1RuptOyP\nB/6HuKl3wPVx+9cSRkdfH61t71y2HEf4LMLGmSHXQ+jPl3SImb0KIOnwmqmiq/O8Re9qgVRG40+U\ndK6ZPQwgKXGA3txo+5PAyvSr6Fxy0Wfx/xES5pRao6EK/fnOpcXnzrvaJJXb+AOAkyRtRWhNHQPE\nr3b3t+hxBvBj2jV0rmKvAMdKGgw0MbMzo+2V9ec751y9kcpt/HeBltHjBkr3mUJoPX0WPW8j6W9m\ndiPOZcfxhDz1S4E34raX25/vXFr8tr2rhVK5jf+CpAPM7MVyDnnGzIYDSLoe2DKdCjpXiYaEL50/\nU/pOUrn9+c45V9+kOvWum6RdgO8Io5yvitu3k6Q/xJ4Db2/KhSU1B14k5DvfPbpGI6C/mXnrzCX6\nG3AkcDDQCfgHVNyf75xz9U2qwf68uOeJAfhiNq733cfMNjWhyV+AqYTsfp3M7DxJlwGHEtYTdy7e\nicALwPgk8+nL6893zrl6JdVgvx3wBzYudFOyEE4U3GelclFJxwHvAU2BFsCX0a5lbFw33DkAJJ0N\n3GFmv8Rt2wVobWavU35/fqrlbQ5cR5i//w8zmynpCvzuU72wdOTJNV0F51KWarC/BFgFjKJsbvx0\nHAk0B/YCfgVi06baA0lvxUrqC/QFaN++fQar4mqBz4CJkrYm9Nv/ArxGtH485ffnp+oSwv+ZDcBy\nSY3xu0/OuVog1WD/OaH1vQHYNn6HpO2AI4gWwTGzKVW9qJldHV2jN7Aa2E3SHdG17irnnElEf9w7\nd+7srar65dX4gaKSGplZcdz+pP35VSGpA3BTwuYvCFP4nieM9v8zldx98i+jtVP8HHrPde/qglSD\n/UNAMXAV8FzCvknAvcA3qVbKzB5I9VxXrxRKutTMNkhqCIwFLo3bX1F/foXMrAg4JX6bpGHAGjMr\nliTgK8KKelDO3Sf/MlrL+LQ6V0elklTnfOBkYKyZDUxyyDwzezLtmjlXueJYEDezXyXF991X1p+f\ninuAmyVdDEw1s/WS3q7s7pNzztW0VFr2XQiD88YT+kdLRIuPbCPpWTZOy/PFR1y2bB67dS+pEWHt\n+pjK+vM3mZmtImGMipmNTvV6zjlXXVIJ9gbsBjSJ8o9jZu9Fjz0k7W1m70JJv6dz2TIaeETSr4Sp\nmvGZGivrz3fOuXqjQQrnLCLkHV8SPXaP7ZCUR1h8pFn0vF9GaulcEmb2jpl1N7OeZtbDzObF7S6U\n1AAgrj/fuaT2fH+hL2zj6rRU0uVeV8HuQmA/wnr2RsJtfucyKeo2ig16E6W7jcrtz3cuNhBvac3W\nwrlqk+po/PJcSpjbHOOjj13WmFmP2HNJW1I6s2NF/fmuvvLR9q6eSjnYR1OPOgEfmtl30eYhQLvY\nIYRgf1GS053LCEmHEgbNHUT4/MVU1J/v6hm/Re/qu1Sm3sUGOuUDHxGC+QAAM7tW0qnAi2aW8jx7\n56pC0huEnA9/AW40s5I0zdHiN93LO9c55+qTVFr2oyRNB34CdiS04BOvWSCpFfClmV2RXhWdK9fN\nQDdgZ8pmcqyoP9/VA96ad26jVAbo/VlST6ANMAVYnHBIU2A9IbXo8rRr6Fw5zGw6MF1SW+BCSbPM\n7JRoX0X9+a4O8/S2zpWVym38PsAuhJH2FwALgEfjDjkC+JYwRe+F9KvoXHKSRhPSNb9oZmX65Cvo\nz3d1SKwF7yHeufKlMs++tZkNAX5nZsOAFQn7ryQE+ROAx9OrnnMVGku4fT9W0hOSSgJ+1J+/P6E/\n/5X4/nxXRxR866PrnauiVPrsF0V/VB8FMLNXEvZfR1jms5+ZrU6zfs5V5BNgIfA/0U9e3L5y+/Od\nc66+SaXP/gngiQr2/zmtGjlXdSuBqcDNZvZp/I6K+vOdq212HPx0ufuWjjy5GmuSOzpMLj8bu4/b\nKCuV2/hZI+k0SXdLekzS8ZKukDRO0sRoXr9z8XYD/gkMkHR7tAQtEPrzJZ0MfGdmN3qgd/XNAw88\nwKxZs1i+fDk9evRg5cqVJfs2bCi74vOgQYMA6N699IzVp556ipdffrnM8QUFBSxYsCBp2YMGDaJf\nv3707NmTtWvX8sMPP9CrVy8uueQS/vrXvwIwevRoOnXqVHKNZMdk0gMPPECXLl3o378/M2bMKLN/\n6dKlXHnllSWPVXHPPffQt29funfvzrx5IVv34MGDGThwIIMHDwbg73//O126dGHcuHEl59133338\n6U9/YtiwYUmvmw2ZzqCXFjObAcyQ1JKQFKWxmZ0n6TLgUEL3gHMxDaKfZoRpoPG38ccCRxH681sA\n75vZ1ekUJul3wDDge+BfZjZZ0hXATkAjoL+ZedZIlzPef/99HnroISZNmkTLli055JBD6Nq1K0cf\nfTTPPvssq1evZrvttmPIkCF8/PHHAPz8889ce+21fPDBBwwdOpSvv/6aBg0asGzZMoYNG8a2227L\n6aefDsCkSZNYu3Ythx9+OBdeeGFJubfccgsAt99+O/PmzWPp0qV0796drl270rNnT8477zzy8/P5\n7rvvSs6ZPn16mWMy7Y9//COnnBK+9y9evJiCggKaNm1K165d6dixY6lj58yZw5QpU7jtttsYNWoU\nq1evZu3atdx+++1ssUVIyNmnTx/69OnD3LlzmTlzJq1ataK4uJixY8cyaNAgli9fTpcuXcjLyyv5\nUvPFF1/w2GOPsf/++9OmTZuMv8fy5FTLPs4wwtrhX0avlwE7JDtQUl9JcyTN+fLLL5Md4uquMUBL\nYIyZnW5m8QsvxfrzF0XP85Kcv6kOAm4DLgSOl9QY6GRmlxEGgx+agTKcy5hbb72VgQMH0rJlSwDy\n8vIYPHgwnTt3ZsOGDWy55ZZMnz691DnFxcUMHTqUO++8k8LCwpLthYWFXHPNNdx2220cemj4qPfs\n2ZP777+f2bNnlyl71apVzJkzh4MPPpgVK1bQrl1IrtqwYcMyxwJVOiZdd911F/3792f+/PlMmDCB\nESNGcM899/DQQw+VOu7111/nkUceYcyYMSxevJiXXnqJFi1a0KhRIxYuLJ2/4ZdffmHs2LH07t2b\nTz/9tOQ9tG/fnhUrEsevw0cffUSrVq0YOXIky5YtY8mSJVl5r4lyqmUf3aofCTwDvEnItQ/QHngn\n2TlmNolojfLOnTt7q6oeMbPeFewutz+/KqLlmW9K2DyC8AXjV8JA1K2p5AuppL5AXwj/+Z2rTmPH\njmX06NE0bdqUTp060aJFCwDmzp2LJEaMGMErrySOsQZJJPacmhkNGpRuH8aul3jsp59+yuDBgyks\nLKRhw4bssMMOrFixgv322y9pFwJQpWPSFd+ynzJlSkm9E+vfpk0bPv/8c9atW8eGDRvYe++9KSgo\nKHO94uJiBgwYwOWXX067du0ws5IAv3z5ck477bQy52y//fa0atUKgK222orvv/8+k2+xXDkV7IE/\nAccCLYBdgbcl3QE0Ae6qyYq5Wmc3Qs6HAZKaAV+Z2Q1VPdnMioBS/fySxgE9zeyTKEPfecA20e6k\nX0j9y2ia4qbWxQZd+eCrqsvLy+PBBx/kggsu4KqrrirZvuuuuzJ//nxuvfVWYndEYwGvUaNG3HDD\nDXz44YcMGTKEt99+GwiBsqCggDZt2tCtW7cKy+3WrRu77bYbgwcP5tJLL+WMM87gsssu4+mnn6Zr\n164ATJ48mVmzZrFw4UKGDRuW9Jhs6tevH8OHDycvL49zzjmn1L7f/OY3DBgwgH79+jFx4kQaNGhA\nfn4+P/74I0OHDi1pvf/lL39h0aJFjB8/nmOOOYYePXrQqFEj8vPzadKkCe3ateO1115j9OjRrFmz\nhjZt2nDmmWfSqlUr8vPzKS4uLtN9kC05FezNbCy+7rjLjIr681M1nZAuei3wppmtl+RfSF1O6t27\nd8nzxx8PKU+mTZsGhBb5E0+ESVVXXnklM2fOZOuttwZg5syZpa6z7777ljyfPHlyyfODDz645Pmj\njz5a6py33nqrTH3uv//+Uq979epFr169Kjwmk+J/HwC77bYbDz74YKltt956a6nH2EDB2OtEo0eP\nLrPtpptK3xD8/e9/z1NPPVVq24gRI6pe8QzJqWDvXAaNISR3GmNmyzJxQTN7Hng+YVvZ/+0ubZ7X\nvrRsT6/r2rVrtbSmM8nv8GwaD/auTqqkP9/lKP8D7lx2eLB3ztU4D/LOZZcHe+ecy3UFLSrYVz/X\nB1i4x57l7vNuoLJydZ69c865NNS2DHofffQRF198cZnyM+XHH3+kf//+dOvWjcMOO4z+/funPMd9\nU6YGTpw4sSTLXlFRuIN11VVXMXDgQIYOHQqEQZEnnXQSEyZMAGDlypX079+f/v370759+4xMz/OW\nvXPO1VG1KYPezjvvzL333pu1YN+sWTMmTJjACy+8wIIFC+jduzcjR47km2++Yf/99+eoo45i0qRJ\njBw5kjvvvJOOHTuydOlS5s+fz/fff09hYSEFBQV8//337LfffrRp04ZnnnmGdevW0bNnT1auXMkr\nr7xC8+bNGTVqFJttFsJrv34h19ebb77JrFmz2HzzzZHE2LFjyc/PZ+XKlXTt2pVmzZqxePFiANq2\nbcuECRNYtWoV69evZ/PNN0/7/XvL3jnn6qjalEGvukmiuLiYVq1aMXXqVHbeeWeWLl2KmfHqq69y\n2GGHYWY0btyYTz75hHfeCWk0zjnnHHr37s3YsWPZaqutaNu2LW+88QZLliyhY8eO5OfnlwT6mOLi\nYgoLC7ngggtK/S7atWvHp5+Wn/PrvvvuKzNlMFXesnfOuTqqNmXQq26zZs2iY8eOnHvuuRx99NEA\nHH744dx7773ssssuSGL69Ok8+eSTDB8+nHXr1gEb37OZMXz48FJfXubNm0d+fj4333wzO+20EwDr\n169nwIABXHHFFWy//fasX7+ep58OqxiuWLGC7bffPmn9zIyXX3655FZ/ujzYO+eqnWfDqx61KYNe\n27Ztufrqq5k7dy433XQTQ4YMycavpESnTp0YPnw4y5cvL/kC0rNnT37729+WfAFq3bo1o0aNYs6c\nOZxwwgmlzr/sssvo06cPLVu25MADD2TNmjUsWbKERo0aldxJgZC0aPHixRQWFnLcccdx5plnYmbk\n5+fTrFkz2rZty6uvvsodd9zBt99+y3bbbcdpp53Gv/71L4466qiMvV/VpUW6OnfubHPmzCm1raJ1\noFOVNMFFRaNlU5VklG1FazinKtkf3IpGuqaqnBGy9Wbp4mSfz8qU9++w5/sLWbjHnmUeqyLZubHP\n1dSbfgnHnB0Gc3XYqX3p7e8vDJ/1gm/pMLkDRb2KSso9a0jptkPsnNj2ol5FJee6TZTl0fgzZ87k\n6aefLhkgVhv4aPwSVfob6i1755zLdVn+glQbM+jVs4CeNh+g55xzztVxHuydc9XHb+E7VyP8Nr5z\nLqt8EJ5zNS/ng72k5oSlQ9cDL5jZX2u4Sq6ekHQAcAWw3MyujLZdAewENAL6A22AW4BfgfvN7N81\nVF3nnCtXbbiNfwYwzcwuASqe7+FcBpnZf4G/xF5Lagx0MrPLgCLgUOBiYCTQG7ikBqrpnHOVyvmW\nPbAD4Q8rhNZTKZL6An2jl99L+iDFcrYBVlflQN2cYgmbWA7XpTUrrervp3fas9+qVpaSljPbzE5M\ntwLpktQBuClh80UJr7cGvoyeLyN8NncgtPw3JCYWia6bqc9n4oWTP6Z57l6xJ9fFniwovT12bPTZ\nLPXZ6V26mDLnOOeyoUp/Q2tDsF9B+IM6jyR3IsxsEjAp3UIkzTGzzulep76VU91lZYuZFQGnJG6X\ntGPcy68IX2wA2gPvEH0+JX1HEpn6fDrnXDpqQ7CfDoyTdDIws6Yr4+oPSbsB1wJ7S+prZpMkvS3p\nDqAJYSzJEsJt/F+Ae2quts45V746lUEvHXWtxe0te+ecczG1YYBedamuW611rZzqLsu5nCWpt6RZ\nkiZLuqqC43aSdGv0vJukw5Ic01/SLmnUo0y3lKu/vGXvnHMZIqk3sNrMZkn6N/AfwsDO2WY2Q9K7\nwINAY+BE4GHge8IA1yLgBuAL4G/AscA0oDuwJbCW0G00DRgKbAXMN7O7Jb1C6PLsDFwJXA/kAQ+b\n2axqeOsux3nL3jnnMusSSXcD9xPGdnwOnBftW2lmI4EpwH/MbFzceQOA683sCjNLXHd2upldSxhE\naoTxVl8DZ0X715rZaOAR4AjgFTzQuzge7Os5SQ0rPyprZfucLFcX3R3lBWkNPAX8H7BFtC+WLzjZ\nou4qZztsHEzdCOgCvGdm18Rt/yF6LCZ8wciNReNdzqgNo/FzlqRdCV0hH2a5nEZmViypgZll7D+x\npBZAf0kPm9nyTF23CuVuCzQ1s08kybwvydVN/yFkWTyEkAE03ipgV0n5hBY6hNkdBZI+I3xJiHeW\npJ6E2/tzgZGS2gDlfVmfD1wtaTMzm5H+W3G1nffZpyj6j3cA0ILQHzctS+WcDRwHDDKzrzMV8CU1\nItw27AS8CXwMvGZmX6V77UrK3ZrQT7k5cLyZ/VDJKc7Va5IKCFlEF9R0XVzt5bfxU9eR0Cd3OXCK\npDMyXUAUGHckJG8ZLSkvgy37BsD2wHfAo8A+hMCfNZKaEuat3wqMA34bbffPoXPlMLMCD/QuXf5H\nNnVvAzuZ2feExVL2lbRFJedskqiVfYuZ3QE8Q0gudKWkdqleU1KepJZm9jNh9G5roCXwLnC8pM2y\n1ZduZj8BhWb2T8Jo5JOi7d6/6JxzWeTBfhNIOl3S76OX84F9JJ0PFBBG2a7NZDlRf9uvAGb2GLAz\n0CbV/vXoTsF4oLektmb2GjCB0K94BDDWzH7JdB+6pFaSjozezwcAZjYZ2EPSIZksyznnXFneZ19F\nkvYGxgJvAVPNbI6ktoS+5w43u0ONAAAIrUlEQVRm9kSGy5kTlfNWtH0r4IQo6JPKwLboi8rehDm9\nbYApZvaDpCZAMzP7JhPvIaHM5oQvGKuB1wlTiH6J9p0A/NfM1mS6XOeccxt5sK8CSc3M7Mfo+ZHA\nQcBsYFvg+bjgldbI8iTlHAA8SwjML5nZumjfJg3Si74o7Gxmb0evdwD2B04GFhG16FOtdwXlijBa\nuAMhIcgxhNXiPgGWAj/6SHznnMs+n3pXCUndgP0lrTOzW8zsBUl5wH2EVuo/YsemGejLK+f+qJzZ\nceVsSqBvCQwBmkk6CPinmX0o6X+BPQhjArIR6FsR1np/yczeiLY1IKwPnwf8MfblxTnnXHZ5n305\nFOxJWKX7LmBLScOj3Z8SboHfUAvKORT4xcz+RLiVfmy0fRFwrpktSePaSUW37q8GPgO6SOoeJe8p\nBn4FrvVb9845V3082JcjaqV/SBgFv87MhgONJZ1oZvPNbAykP20sW+VIah71xT8PrJZ0tJlNJSTy\nOMjM3jWzlenUvQKxOxz/AcYAewEHmdnHwBVmtihL5TrnnEvCb+MnIelMwhz01wnpKc+W9BKwjo2B\nDEhv2li2ypHUGriTMJhwJfAScEI08n0Dod8846Jy9yXcQXgE6EaY3vcOcJSk1zM1Y8E551zVecs+\nQdR3fgTwE3CAmd1FyGd9HNDAzJ7N5XKixDV9gCeAfxJG2c8htLDfAcaY2WcZeAuJ5bYiJMvZHTiS\nkJFvHiHp0FHAA7FphM4556qXt+zLagQMN7NvJV0fbXvdzFbEDshQPveslGNmP0maamZLogC8h6T9\ngIZm9mSadU4qGnV/PPComT0j6T7C9MH/Aq8Bzc3s64qu4ZxzLnt86l05ov7u8cBXwJtRf3dOlxP/\n5SDq498RmAy8D1yTjRZ9XNlNoqx8SPozsJjQRXGvmRVnq1znnHOVq/e38RMHvsVeR4Frb+DDTAT6\naiqnSexJ1MdvhNW1bshmoI/EB/QDgV7A0x7onXOu5tX7YG9mGyRtKekESdtHr2O/l4lmNgnSX3s9\n2+VIOgcYL+moaJob0ej3S80sKwPyonIPi8raEE0j3Iwwre/KVNP6Ouecy6x6fRtf0v7AlsA5hAB1\nIHC1mS1OOC6tZWWzXU6Ube8iwij4UYSBfj9ne4EZhWV+DwMGEdanXxNt3ywbiXqcc86lpl627ONa\nz7sAXYEiMxsNTALGSNo2/vhUg2Z1lCNpR+BLYDtC/3xLQkKbxqnUeRPK3Z6wLO4iwi37/5PUGMAD\nvXPO5ZZ6Nxo/NohNYTnaz4BXCavXdTSz5yStMrMvcrmc6JoHAu8RMu9dT5hqdxEhQ92dFpaTzQpJ\npxLWBSgGpgK/B5oT0uCuz1a5zjnnUlPvgn0UgJsTWtd3m9nfJP0C9JS0xszehfSn12WrHEnNgFbA\n0UBH4LOov/xZ4GlCatxVqda7CuUfTsh5Pw+YSPhysQMw2LKwap5zzrn01atgL6k9YUnaz4FXzOz5\naNdS4Dsz+yR2bJqL2mSlnGhAXx9CFrx7gBuB/SQ1An4HXGZmP6Ra7yqUvx8hac5nZnaNpH0Jt/L/\n6vPonXMud9WLPvtolHhDQirX30c/HaKR62MI69G/GDs2V8uJ+vRnAs2Awwl9888CC4FLshzozwU6\nExLl7BO9vpKwqp0Heuecy2F1fjR+bIR7lLzmfwhrqi8GfktI6do6U/Pos12OpCOALYAiwoA/A/4O\nrDWz1elcu5JyWwM3EzLk/UPS3sA2wGIz+zRb5TrnnMuM+nAbv5Ok94BxhH7mw4DnCAPJFpnZvyEj\nKXCzUk7cQL98oC1hsN8ZwAjgVOBrM/s2jXpXVv4phIGAM4DOktaa2WvZKs8551zm1Ydg3wUoAMZG\nrdLVwPfAqvhWaQZy3We8HEktgQMl/YvwpeFGM1sj6WrCmIBx2ZzmJukYYCDwIyGl70rgMElzszna\n3znnXGbV2T57SbtKGkZoAT8PnBDt2gf4ycz+m8vlaOMytasIt+s3B/4gqQewM/BplgP97kAbwvK4\nBUAPwoDD8R7onXOudqlzwT428C3KTvc5YfW124FvJD0HfGBm/8jlcqJr9wA+BAT0Be4FlgMNCS38\nrA2Kk7QPcKSZPQS0JozA/wp4z3w9euecq3Xq1AC9uP7trQjrws8CTgMOIqyrfqiZvRx/bK6WI2k7\nQp/8ScBDZjYtlbqmUG4HMyuSNI6QMGcVYbrf7dWwmI5zzrksqFPBHkDS5kA+cARhPvrZhMFyxwOX\nm9n6DAzGq5ZyJG0dXXcpMMfMPk+nzlUorzVwBeEOwlvAsWZ2g+e6d8652q1O3MaPks0gKY+QVS6P\nkLnuDuDfwALCbegtIPXBeNVVToyZfUUYBd8aWJfOtSoTTes7xswGA2sIYw/6RFkAf81m2c4557Kr\nzrTso0Vl/khY9W1nwm31lcCWZjZGUp6ZpR0wq6uchDKz0rKOWvLHEr70dSN8SfmYkInPJG2Tzfn7\nzjnnqketb9lL6qGwpvohhMBrZraA0MpuYWZjokN/rA3lJJPFQH830IKQF2AV0DPafVv06JnxnHOu\nDqgL8+wbE3K1vxytBjde0mAzmxE7IBN99NVYTtZJagWMAVYAawmJel4iLHCzDHgMUl/a1znnXG6p\ntS17SXtHc84PIaRuxcymAHMJqWRLpDlIrlrKqWZbAA8Cw6LXPQlz+T8FpprZspqqmHPOucyrtX32\nktoSgu0fgNeApoTW931m9kZtK6emSNqVsGDPZsAUM/PBeM45V8fU2mAPJelkBxJyty8lZJVbGe3L\n2C316iqnpkj6DbDGzL6r6bo455zLvFp7Gx/AzNYQpr5tBzTIVgCurnJqipkt80DvnHN1V61u2cdI\napHNld+quxznnHMuk+pEsI+prpZ2XWnRO+ecqx/qVLB3zjnnXFm1us/eOeecc5XzYO+cc87VcR7s\nnXPOuTrOg71zzjlXx3mwd8455+q4/w/rQn7wXidNIgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1343e53c8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig = figure_4(waic);\n",
    "fig.savefig('results/figures/figure_4_relative_model_comparison.png', dpi=330, bbox_inches='tight')\n",
    "fig.savefig('results/figures/figure_4_relative_model_comparison.pdf', dpi=330, bbox_inches='tight')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python [Root]",
   "language": "python",
   "name": "Python [Root]"
  },
  "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.5.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}