{ "cells": [ { "cell_type": "markdown", "source": [ "# Quickstart\n", "\n", "The partial and full confound tests can characterize whether a predictive model\n", "(e.g. a machine learning model) is biased by the confounder variable.\n", "\n", "In research using predictive modelling techniques, confounder-bias is often investigated\n", "(if investigated at all) by testing the association between the confounder variable and\n", "the predictive values.\n", "However, the significance of this association does not necessarily imply a significant\n", "confound-bias of the model, especially if the confounder is also associated to the\n", "true target variable. In this case, namely, the model still might not be directly\n", "driven by the confounder, i.e. the dependence of the predictions on the confounder\n", "can be explained solely by the confounder-target association.\n", "Put simply, this is what is tested by the proposed partial confounder test.\n", "\n", "Here we will apply the `partial confound test` on two simulated datasets:\n", "- H0: null-hypothesis dataset with no confounder bias,\n", " i.e. conditional independence between the predicted values\n", " and the confounder variable, given the observed target variable.\n", " Note that the (unconditional) association between the prediction and the target are significant\n", " but - according to the H0 of the partial confoudner test - can be fully explained by the association\n", " between the target and the confounder.\n", "- H1: alternative hypothesis with an explicit confounder bias. Here, the association between\n", " the predictions and the confounder is stronger than what could follow form the association\n", " between the target and the confounder.\n", "\n", "\n", "##### Import the necessary packages" ], "metadata": { "collapsed": false } }, { "cell_type": "code", "execution_count": 1, "outputs": [], "source": [ "from mlconfound.stats import partial_confound_test\n", "from mlconfound.simulate import simulate_y_c_yhat\n", "from mlconfound.plot import plot_null_dist, plot_graph\n", "\n", "import pandas as pd\n", "import seaborn as sns\n", "sns.set_style(\"whitegrid\")" ], "metadata": { "collapsed": false, "pycharm": { "name": "#%%\n" } } }, { "cell_type": "markdown", "source": [ "##### H0 simulations\n", "\n", "Next, we simulate some data from the null hypothesis.\n", "\n", "The simulation samples `y` and `c` (`n` datapoints each) from the multivariate normal distribution,\n", "so that their correlation is `cov_y_c`.\n", "Next, a 'prediction' of y is simulated like this:\n", "```\n", "yhat = y_ratio_yhat * y + c_ratio_yhat * c + e\n", "```\n", "where `e` is random standard-Gaussian noise.\n", "\n", "For the H0 simulation, `c_ratio_yhat` is set to zero." ], "metadata": { "collapsed": false } }, { "cell_type": "code", "execution_count": 2, "outputs": [], "source": [ "H0_y, H0_c, H0_yhat = simulate_y_c_yhat(cov_y_c=0.3,\n", " y_ratio_yhat=0.5, c_ratio_yhat=0,\n", " n=500, random_state=42)" ], "metadata": { "collapsed": false, "pycharm": { "name": "#%%\n" } } }, { "cell_type": "markdown", "source": [ "##### Partial confound test\n", "\n", "Now let's perform the partial confounder tests on this data.\n", "The number of permutations and the steps in the Markov-chain Monte Carlo process set to the default values\n", "(1000 and 50, respectively).\n", "Increase the number of permutations for more accurate p-value estimates.\n", "\n", "The random seed is set for reproducible results.\n", "The flag `return_null_dist` is set so that the full permutation-based null distribution is returned,\n", "e.g. for plotting purposes.\n", "\n", "The pandas dataframe is solely created for \"pretty-printing\" the results." ], "metadata": { "collapsed": false } }, { "cell_type": "code", "execution_count": 3, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "Permuting: 100%|██████████| 1000/1000 [00:02<00:00, 457.42it/s]\n" ] }, { "data": { "text/plain": " p ci lower ci upper R2(y,c) R2(y^,c) R2(y,y^)\n0 0.798 0.771754 0.822478 0.094132 0.031789 0.478693", "text/html": "
\n\n\n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n
pci lowerci upperR2(y,c)R2(y^,c)R2(y,y^)
00.7980.7717540.8224780.0941320.0317890.478693
\n
" }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ret=partial_confound_test(H0_y, H0_yhat, H0_c, return_null_dist=True,\n", " random_state=42)\n", "#pretty print results\n", "pd.DataFrame({\n", " 'p' : [ret.p],\n", " 'ci lower' : [ret.p_ci[0]],\n", " 'ci upper' : [ret.p_ci[1]],\n", " 'R2(y,c)' : [ret.r2_y_c],\n", " 'R2(y^,c)' : [ret.r2_yhat_c],\n", " 'R2(y,y^)' : [ret.r2_y_yhat],\n", "})\n" ], "metadata": { "collapsed": false, "pycharm": { "name": "#%%\n" } } }, { "cell_type": "markdown", "source": [ "Let's use the built-in plot functions of the package `mlconfound` for a graphical representation of the results." ], "metadata": { "collapsed": false } }, { "cell_type": "code", "execution_count": 4, "outputs": [ { "data": { "text/plain": "" }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" }, { "data": { "text/plain": "
", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEWCAYAAACJ0YulAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAAf0UlEQVR4nO3de1xUZf4H8M/ILcgL4cKgZq2Q2ai7iJd0FtMcEDVELsa6u2V42dBfiqHGppltr1JaW2tZrV6GWGu1bm0soiKmiRpYSCZeyiYva6WWM1PIRR1nYIbz+8P1eYlymfsMzOf9evWKmXnOc77PAH54zpzzHJkkSRKIiIgAdHF3AURE5DkYCkREJDAUiIhIYCgQEZHAUCAiIoGhQEREAkOBOq0lS5bgb3/7GwCgsrISY8aMsXjb6dOn48MPPwQAbN26FbNmzXJYXQkJCaisrAQArF27Fk899ZTD+l63bh2WLVvmsP7I+/i6uwAiTzdlyhRMmTKl3XZLliyBXC7HwoUL22y3fft2h9RVWVmJ7OxslJWViefmzp3rkL7Je3GmQOQiJpPJ3SUQtYuhQB5NpVJhw4YNSExMxLBhw5CVlQWj0QgAKCwsxO9///tm7QcMGIDvv//e6v18+umnmDhxIoYNG4YXXngBN17of+N+JElCTk4OlEolhg4disTERJw8eRIffPABtm3bhg0bNiA6Olr8xa5SqZCXl4fExEQMGTIEJpMJKpUKn332mei/oaEBWVlZiI6ORkpKCr755ptWx3P9kJher8fjjz8OnU6H6OhoREdHQ6vV3nI4qrS0FAkJCRg+fDimT5+O//73vxa9t+S9GArk8Xbs2IH8/HyUlpbixIkTKCwsdGj/Fy9exPz585GVlYUDBw7grrvuQlVVVYtt9+/fjy+++AI7d+7EoUOHkJubi+DgYEybNg2JiYmYPXs2Dh8+jHXr1olttm/fjry8PHzxxRfw9b31iG1paSkmTpyIzz//HJMnT8YTTzyBxsbGNmsOCgrC+vXrERYWhsOHD+Pw4cOQy+XN2nz77bdYvHgxnnnmGVRUVGDMmDGYO3cuGhoaRBtnv7fU8TAUyONNnz4dcrkcwcHBGDduHNRqtUP7LysrQ//+/TFx4kT4+fkhPT0dv/jFL1ps6+vriytXruDMmTOQJAmRkZEICwtrt/5evXrhtttua/H1QYMGiX3PnDkTDQ0NOHr0qN3jKikpwdixYxETEwM/Pz/Mnj0bBoMBhw8fblabM99b6ngYCuTxQkNDxdeBgYHQ6/UO7V+n0yE8PFw8lslk6NWrV4ttlUolHnnkEbzwwgtQKpVYvnw5Ll++3Gb/rfV13Y377tKlC+RyOXQ6nRUjaJlOp0Pv3r2b9d2rVy9otVrxnLPfW+p4GArUYQUGBsJgMIjHP/30k039hIaGQqPRiMeSJOHChQuttn/sscdQWFiIkpISfPfdd8jPzwdwLUxa0trz192476amJmi1WjH7CAwMxNWrV8XrN46xvX7DwsLw448/isfXx3XzYSaiGzEUqMO67777cOrUKajVahiNRqxdu9amfsaOHYtTp05h165dMJlMeOedd/Dzzz+32PbYsWM4evQoGhsbERgYCH9/f3Tpcu3XqGfPnjh//rzV+z9+/LjY98aNG+Hv74+oqCgxxuLiYpjNZpSVleHgwYNiu549e6K2thaXLl1qsd9Jkybhk08+QUVFBRobG/HWW2/B398f0dHRVtdI3oOhQB1Wv379MG/ePMyYMQPx8fEYNmyYTf2EhITg73//O1555RWMHDkS33//PYYOHdpi2ytXruDZZ5/F/fffj3HjxiE4OBizZ88GADz88MM4ffo0hg8fjieeeMLi/cfGxqKkpAQjRozAli1bsHbtWvj5+QEAli1bhr1792L48OHYtm0b4uLixHaRkZFISEhAXFwchg8f3uywEABERETgr3/9K1588UWMGjUKe/fuxbp16+Dv72/tW0ReRMab7BAR0XWcKRARkcBQICIigaFAREQCQ4GIiIQOvUrqkSNHEBAQ4O4ybGY0Gt1T/4kT1/4/YIDDu3bbmJyos42ps40H4Jhs6XvIkCEtvtahQyEgIAAKhcLdZdhMrVa7p/7/+79r/9+3z+Fdu21MTtTZxtTZxgNwTLb03RoePiIiIoGhQEREAkOBiIgEhgIREQkMBSIiEhgKREQkMBSIiEhgKBARkcBQICIigaFArWowNVm9zd39Ip1QCRG5Sode5oKcy9+3C36XV2HVNu9nKJ1UDRG5AmcKREQkMBSIiEhgKBARkcBQICIigaFAREQCQ4E6FGtPk7XltFoib8ZTUqlDsfY0WZ4iS2QdzhSIiEhgKBARkcBQICIigaFAREQCQ4GIiASnhcLSpUuhVCoxefJk8VxtbS1mzpyJ+Ph4zJw5E3V1dQAASZKwYsUKjB8/HomJiTh+/LizyiIiojY4LRRSU1ORn5/f7Lm8vDwolUrs2rULSqUSeXl5AICysjJ899132LVrF1588UU8//zzziqLnIzXERB1bE67TmHEiBE4f/58s+dKS0vx7rvvAgCSk5Mxffp0ZGdno7S0FMnJyZDJZBgyZAjq6+uh0+kQFhbmrPLISXgdAVHH5tKL16qrq8U/9KGhoaiurgYAaLVahIeHi3bh4eHQarXthoLRaIRarXZewU5mMBjcUv9dej0A4Gw7+1YoFND/r601rN3GmvfAlprsfY/d9X1yls42HoBjciS3XdEsk8kgk8ns6iMgIAAKhcJBFbmeWq12T/1BQQBg0b6D/tfWuu6t28ba98DZ/d/Mbd8nJ+ls4wE4Jlv6bo1Lzz7q2bMndDodAECn0yEkJAQAIJfLodFoRDuNRgO5XO7K0oiICC4OBZVKhaKiIgBAUVERYmNjmz0vSRKOHDmCbt268fMEIiI3cNrho0WLFuHzzz9HTU0NxowZg8zMTGRkZCArKwsFBQXo3bs3cnNzAQBjx47FJ598gvHjxyMwMBA5OTnOKouIiNrgtFB49dVXW3x+48aNtzwnk8nw5z//2VmlEBGRhXhFMxERCQwFIiISGApERCQwFIiISGAoEBGRwFDwIlx8joja47ZlLsj1ri9W99yFegDAC+0sXMfF6oi8D2cKREQkMBSIiEhgKBARkcBQICIigaFAREQCQ4GIiASGArkVr50g8iy8ToHc6vq1E5bitRNEzsWZQgfGv7KJyNE4U+jA+Fc2ETkaZwpERCQwFIiISGAoEBGRwFAgIiKBoUBERAJDgYiIBIYCEREJDAUiIhIYCkREJDAUiIhIYCgQEZHAUCAiIsEtC+L94x//wIcffgiZTIZ7770XL730EnQ6HRYtWoTa2loMGjQIL7/8Mvz9/d1RHhGR13L5TEGr1eKdd97Bf/7zHxQXF8NsNmP79u1YvXo1ZsyYgY8//hjdu3dHQUGBq0sjIvJ6bjl8ZDabYTAYYDKZYDAYEBoaigMHDmDChAkAgJSUFJSWlrqjNCIir+byw0dyuRyzZs3CuHHjEBAQgJiYGAwaNAjdu3eHr++1csLDw6HVatvty2g0Qq1WO7tkpzEYDHbVr1AooNfrrdpGr9fDbDaLry1pby1banJme3t/Ruz9PnmazjYegGNyJJeHQl1dHUpLS1FaWopu3brhySefRHl5uU19BQQEQKFQOLhC11Gr1XbXHxQUZHV7Hx8fi7e1tn9ba3Jme3vfY0d8nzxJZxsPwDHZ0ndrXH746LPPPsOdd96JkJAQ+Pn5IT4+HlVVVaivr4fJZAIAaDQayOVyV5dGnZC1tyzlLU7J27l8ptC7d28cPXoUV69exW233YaKigoMHjwYI0eOxM6dO5GQkIDNmzdDpVK5ujTqhHjLUiLruDwUoqKiMGHCBKSkpMDX1xcKhQLTpk3Dgw8+iIULFyI3NxcKhQJpaWmuLo2IyOu55TqFBQsWYMGCBc2e69u3L09DJSJyM17RTEREAkOBiIgEhgIREQkMBSIiEhgKREQkMBSIiEhgKBARkcBQICIigaFAREQCQ4GIiASGAhERCQwFIiISGApERCQwFIiISGAoEBGRwFAgIiKBoUBERAJDgYiIBIYCEREJDAUiIhIYCkREJDAUiIhIYCgQEZHAUCAiIoGhQHSDBlPTLc8pFAqr2hN1ZL7uLoDIk/j7dsHv8iqaPafX6xEUFNRi+/czlK4oi8hlLJoppKenW/QcERF1bG3OFIxGI65evYqamhrU1dVBkiQAwOXLl6HVal1SIBERuU6bofD+++9j48aN0Ol0SE1NFaHQtWtXPProoy4pkIiIXKfNUEhPT0d6ejreffddTJ8+3WE7ra+vx7PPPouTJ09CJpMhJycH/fr1w8KFC/HDDz+gT58+yM3NRY8ePRy2TyIiap9FHzRPnz4dVVVV+OGHH2A2m8XzycnJNu105cqVeOCBB7BmzRo0NDTAYDBg3bp1UCqVyMjIQF5eHvLy8pCdnW1T/0REZBuLPmjOzs7Gyy+/jEOHDuHLL7/El19+ia+++sqmHV66dAkHDx7Eww8/DADw9/dH9+7dUVpaKkImOTkZu3fvtql/IiKynUUzha+++golJSWQyWR27/D8+fMICQnB0qVL8c0332DQoEFYtmwZqqurERYWBgAIDQ1FdXV1u30ZjUao1Wq7a3IXg8FgV/0KhQJ6vd6qbfR6vZjtWbKttf3bso2nt29qamqzj472M2jvz50n4pgcx6JQ6N+/P3766Sfxj7Y9TCYTvv76ayxfvhxRUVFYsWIF8vLymrWRyWQWBVBAQECbFxZ5OrVabXf9rZ0/31Z7Hx8fi7e1tn9ba/Lk9m1dpwC0fXGbJ3LEz52n4Zis77s1FoVCTU0NEhIS8Otf/xp+fn7i+XXr1lldTHh4OMLDwxEVFQUAmDhxIvLy8tCzZ0/odDqEhYVBp9MhJCTE6r6JXK3B1AR/X+sWBrBlGyJXsSgUMjMzHbbD0NBQhIeH48yZM4iIiEBFRQUiIyMRGRmJoqIiZGRkoKioCLGxsQ7bJ5GztHQFdHt4FTR5MotC4f7773foTpcvX46nnnoKjY2N6Nu3L1566SU0NTUhKysLBQUF6N27N3Jzcx26TyIiap9FoRAdHS2O8Tc2NsJkMiEwMBBVVVU27VShUKCwsPCW5zdu3GhTf0RE5BgWhcLhw4fF15IkobS0FEeOHHFWTURE5CZWf9olk8kQFxeH/fv3O6MeIiJyI4tmCrt27RJfNzU14auvvkJAQIDTiiIiIvewKBT27t0rvvbx8UGfPn3wxhtvOK0oIiJyD4tC4aWXXnJ2HURE5AEs+kxBo9Fg3rx5UCqVUCqVyMzMhEajcXZtRETkYhaFwtKlS6FSqVBeXo7y8nKMGzcOS5cudXZtRETkYhaFwsWLFzF16lT4+vrC19cXqampuHjxorNrIyIiF7MoFIKDg7FlyxaYzWaYzWZs2bIFwcHBTi6NiIhczaJQyMnJwY4dOxATE4PRo0dj586d+Mtf/uLs2oiIyMUsOvtozZo1WLVqlbg9Zm1tLVatWsWzkoiIOhmLZgonTpxodr/k4ODgTndDCyIisjAUmpqaUFdXJx7X1tY2u1czERF1DhYdPpo1axamTZuGiRMnAgA++ugjzJ0716mFERGR61kUCsnJyRg8eDAOHDgAAHjttddwzz33OLUwIiJyPYtCAQDuueceBgERUSfHG8USEZHAUPAgDaYmd5dARF7O4sNH5HzW3gSeN4AnIkfjTIGIiASGAhERCQwFIiISGApERCQwFIiISGAoEBGRwFAgIiKBoUBERAJDgYiIBIYCEREJbgsFs9mM5ORkzJkzBwBw7tw5pKWlYfz48cjKykJDQ4O7SiMi8lpuC4V33nkHkZGR4vHq1asxY8YMfPzxx+jevTsKCgrcVRoRkddySyhoNBrs27cPDz/8MABAkiQcOHAAEyZMAACkpKSgtLTUHaUREXk1t6ySmpOTg+zsbFy5cgUAUFNTg+7du8PX91o54eHh0Gq17fZjNBqhVqudWqszGQyGZvUrFAro9Xqr+rCl/fX7a1uyrbX921qTJ7dvampqsw9b3iN3/tze/HPXGXBMjuPyUNi7dy9CQkIwePBgVFZW2tVXQEAAFAqFgypzPbVafUv9QUFBVvVhS3sfHx+Lt7W2f1tr8uT2er2+zT5seY/c+XPb0s9dR8cxWd93a1weClVVVdizZw/KyspgNBpx+fJlrFy5EvX19TCZTPD19YVGo4FcLnd1aUREXs/lnyksXrwYZWVl2LNnD1599VWMGjUKr7zyCkaOHImdO3cCADZv3gyVSuXq0oiIvJ7HXKeQnZ2Nt99+G+PHj0dtbS3S0tLcXRIRkddx6+04R44ciZEjRwIA+vbty9NQiYjczGNmCkRE5H4MBSIiEhgKREQkMBSIiEhgKBARkcBQIPJwDaYmp7YnupFbT0klovb5+3bB7/IqLG7/fobSidVQZ8eZAhERCQwFIiISGApERCQwFIiISGAoEBGRwFAgIiKBoUBERAJDgYiIBIYCEREJDAUiIhIYCkQuxrWJyJNx7SMiF+NaRuTJOFNwovb+IlQoFC6qhIjIMpwpOFF7fxHq9XoEBQWJx/yLkIjcjTMFIiISGApERCQwFIiISGAoEBGRwFAgIiKBoUBERAJDgYiIBIYCEREJLr947cKFC/jTn/6E6upqyGQy/Pa3v0V6ejpqa2uxcOFC/PDDD+jTpw9yc3PRo0cPV5dHROTVXD5T8PHxwZIlS1BSUoIPPvgAmzZtwunTp5GXlwelUoldu3ZBqVQiLy/P1aUREXk9l4dCWFgYBg0aBADo2rUrIiIioNVqUVpaiuTkZABAcnIydu/e7erSiIi8nlvXPjp//jzUajWioqJQXV2NsLAwAEBoaCiqq6vb3d5oNEKtVju7TJspFAro9fpWX29qarrl9bbat8SW9maz2eJtre3f1po8uX1L3yd7+ndETe1p6/fCYDB49O+NLTgmx3FbKFy5cgULFizAM888g65duzZ7TSaTQSaTtdtHQECAx680euOCdze7eUG89tpb239r7X18fCze1tr+ba3Jk9u39H2yp39H1NSetn4v1Gq1x//eWItjsr7v1rjl7KPGxkYsWLAAiYmJiI+PBwD07NkTOp0OAKDT6RASEuKO0oi8jrU3/eFNgjo3l88UJEnCsmXLEBERgZkzZ4rnVSoVioqKkJGRgaKiIsTGxrq6NCKvxJv+0I1cHgqHDh3Cli1bcO+99yIpKQkAsGjRImRkZCArKwsFBQXo3bs3cnNzXV0aEZHXc3koDB8+HCdOnGjxtY0bN7q4GiIiuhGvaCbqZHgbWLIHb8dJ1MnwNrBkD84UiIhIYCgQEZHAUCAiIoGhYCFesENE3oAfNFvI2gt8AH6AR0QdD2cKREQkMBSIiEhgKBARkcBQICIigaFAREQCQ4GIiASGAhERCQwFIiISGApERCQwFIiISGAoEBGRwFAgIiKBoUBERAJDgYiIBIYCEREJXhsKvGkOEdGtvPYmO9beNIc3zCEib+C1MwUiso21s2xbZuXWbnN3v0in9u9NvHamQES2sWWW7exb2Qbd5s+Zv4NwpkBERAJDgYiIBIYCEXkcHvN3H4/6TKGsrAwrV65EU1MT0tLSkJGR4e6SiMgNeHag+3jMTMFsNuOFF15Afn4+tm/fjuLiYpw+fdrdZRFRJ+TsM6gc0V6hUDh0H5bymJnCsWPHcPfdd6Nv374AgISEBJSWluKee+5xc2VE1Nk4+wwqR7TX6/UICgpqcxtnkEmSJDmlZyt99NFHKC8vx8qVKwEARUVFOHbsGJ577rlWtzly5AgCAgJcVSIRUadgNBoxZMiQFl/zmJmCLVobFBER2cZjPlOQy+XQaDTisVarhVwud2NFRETex2NC4Ve/+hW+++47nDt3Dg0NDdi+fTtUKpW7yyIi8ioec/jI19cXzz33HP74xz/CbDZj6tSp6N+/v7vLIiLyKh7zQTMREbmfxxw+IiIi92MoEBGRwFBwkrKyMkyYMAHjx49HXl7eLa83NDQgKysL48ePR1paGs6fPw8A+PTTT5GamorExESkpqaiosK6JYedydYxXffjjz8iOjoaGzZscFXJbbJnPN988w2mTZuGhIQEJCYmwmg0urL0Vtk6psbGRjz99NNITEzEpEmT8Oabb7q69Fa1N6aDBw8iJSUFAwcOxEcffdTstc2bNyM+Ph7x8fHYvHmzq0puk63jUavVzX7mSkpKnFOgRA5nMpmk2NhY6ezZs5LRaJQSExOlU6dONWvz3nvvScuXL5ckSZKKi4ulJ598UpIkSTp+/Lik0WgkSZKkEydOSKNHj3Zp7a2xZ0zXZWZmSpmZmVJ+fr6rym6VPeNpbGyUJk+eLKnVakmSJOnixYuSyWRyaf0tsWdMW7dulbKysiRJkiS9Xi+NGzdOOnfunEvrb4klYzp37pykVqul7OxsaceOHeL5mpoaSaVSSTU1NVJtba2kUqmk2tpaVw+hGXvGc+bMGenbb7+VJEmSNBqNFBMTI9XV1Tm8Rs4UnODGJTv8/f3Fkh032rNnD1JSUgAAEyZMQEVFBSRJwsCBA8X1Gf3794fRaERDQ4PLx3Aze8YEALt370afPn085owye8bz6aefYsCAAbjvvvsAAHfccQd8fHxcPoab2TMmmUyGq1evwmQywWAwwM/PD127dnXHMJqxZEx33nkn7rvvPnTp0vyfs/379yMmJgbBwcHo0aMHYmJiUF5e7sryb2HPePr164df/vKXAK5d1xUSEoKLFy86vEaGghNotVqEh4eLx3K5HFqt9pY2vXr1AnDtdNxu3bqhpqamWZudO3di4MCB8Pf3d37R7bBnTFeuXMH69esxf/58l9bcFnvG8+2330Imk2H27NlISUnB+vXrXVp7a+wZ04QJExAYGIjRo0dj3LhxmDVrFoKDg11ZfossGZMztnUWR9V07NgxNDY24q677nJkeQA86DoFau7UqVNYvXo13nrrLXeXYrfXXnsN6enpuP32291dikOYzWYcOnQIBQUFCAwMxIwZMzB48GAolR13+eZjx46hS5cuKC8vR319Pf7whz/gN7/5jVigkjyHTqdDdnY2Vq1adctswhE4U3ACS5bskMvluHDhAgDAZDLh0qVLuOOOOwAAGo0G8+fPx6pVq5zyl4At7BnT0aNHsXr1aqhUKmzcuBFvvvkm3nvvPZfWfzN7xhMeHo4RI0YgJCQEgYGBGDNmDI4fP+7S+ltiz5iKi4vxwAMPwM/PDz179sTQoUPx5ZdfurT+ltiz/I0nLp1jb02XL1/GnDlzsHDhQqet/cZQcAJLluxQqVTibIidO3di1KhRkMlkqK+vR0ZGBhYvXoxhw4a5o/wW2TOmTZs2Yc+ePdizZw/S09MxZ84cPProo+4YhmDPeEaPHo2TJ0+KY/AHDx70iCXe7RlTr169UFlZCeDaks1Hjx5FRESEy8dwM3uWvxk9ejT279+Puro61NXVYf/+/Rg9erSTK26bPeNpaGjAvHnzkJSUhIkTJzqvSId/dE2SJEnSvn37pPj4eCk2NlZ64403JEmSpNzcXGn37t2SJEmSwWCQMjMzpbi4OGnq1KnS2bNnJUmSpNdff12KioqSpkyZIv77+eef3TaOG9k6phutWbPGI84+kiT7xlNUVCQ99NBDUkJCgrRq1Sq31N8SW8d0+fJlKTMzU3rooYekSZMmSevXr3fbGG7W3piOHj0qPfDAA1JUVJR0//33Sw899JDY9sMPP5Ti4uKkuLg4qaCgwC3138zW8RQVFUkDBw5s9m/D119/7fD6uMwFEREJPHxEREQCQ4GIiASGAhERCQwFIiISGApERCQwFMhrKRQKJCUlYfLkyZg7dy7q6+sBtL8a5YIFC3Du3Ll2+zeZTMjIyMDIkSNx8uRJq2q7flLg2rVrxeOGhgY88sgjMJlMVvVFZA2GAnmt2267DVu2bEFxcTF69OiBf/7zn+L5VatWYfv27cjPz0dOTo4IjFOnTsFsNlu0/MPzzz+PiIgIvP7661i4cGGzK1nbs3XrVuTn58NoNGL9+vXYunUr/P39oVQqnbdkMhEYCkQAgCFDhoiFydpajXLbtm2IjY0FABQUFGDlypWij3//+9/IyckBcG29p27dumHJkiUYPnw4VqxYgUWLFuHSpUu37Pvnn3/GvHnzMGXKFEyZMgVVVVVISkpCeHg4NmzYgN69eyMpKQkAEBcXh23btjntfSBiKJDXM5vNqKioaHG5gZtXo6yqqsKgQYMAAJMmTcLevXvR2NgIACgsLMTUqVMBAPPnz8fTTz8t+omOjsamTZvQrVu3W/axYsUKjBgxAlu3bsXmzZvRv39/bNu2DRqNBrNnz8aPP/4ogqB///4esSYRdV5cJZW8lsFgQFJSErRaLSIjIxETE9Ps9ZZWo/zpp58QEhICALj99tsxatQo7Nu3DxEREWhsbMSAAQOsruPAgQN4+eWXAQA+Pj7o1q0bJk+eDJlMhrVr1+Lxxx8XnzH4+PjAz88Ply9f9oj7HVDnw5kCea3rnyns3bsXkiSJzxSA1lejDAgIaHbrzbS0NBQWFqKwsBCpqakOq00mkwEAMjMzmz0Gri2MFhAQ4LB9Ed2IoUBeLzAwEM8++yzefvttmEymNlejjIyMxNmzZ8XjqKgoaDQaFBcXY/LkyRbv87333hPLhyuVSmzatAnAtUNZLX3ucF1NTQ3uuOMO+Pn5WTNEIosxFIgADBw4EAMGDEBxcTF27NiBL774Aps3b0ZSUhKSkpKgVqsBAGPHjhVLTF83adIkDB06FD169LB4f2fOnBF3Nlu2bBkqKyuRmJiI1NRUnD59utXtKisr8eCDD1o9PiJLcZVUIisYDAY89thj+Ne//iXuyzxnzhzMmDHDqjuvzZkzB2vXrrX6Vqvz58/H4sWL0a9fP6u2I7IUQ4HISuXl5YiMjETXrl2RlpaGAQMGYM2aNU7fb0NDA0pKSpCcnOz0fZH3YigQEZHAzxSIiEhgKBARkcBQICIigaFAREQCQ4GIiIT/B9ao0h9oZy9NAAAAAElFTkSuQmCC\n" }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plot_null_dist(ret)" ], "metadata": { "collapsed": false, "pycharm": { "name": "#%%\n" } } }, { "cell_type": "markdown", "source": [ "The histogram shows the $R^2$ values between the predictions and the permuted confounder variable\n", "(conditional permutations). The red line indicates that the unpermuted $R^2$ is not \"extreme\",\n", "i.e. we have no evidence against the null ($p=0.8$)." ], "metadata": { "collapsed": false } }, { "cell_type": "code", "execution_count": 5, "outputs": [ { "data": { "text/plain": "", "image/svg+xml": "\n\n\n\n\n\n%3\n\n\n\nc\n\nc\n\n\n\ny\n\ny\n\n\n\nc--y\n\n0.094\n\n\n\nyhat\n\nŷ\n\n\n\nc--yhat\n\n0.032 (p=0.8)\n\n\n\ny--yhat\n\n0.479\n\n\n\n" }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "plot_graph(ret)" ], "metadata": { "collapsed": false, "pycharm": { "name": "#%%\n" } } }, { "cell_type": "markdown", "source": [ "The graph shows the unconditional $R^2$ values across the target $y$, confounder $c$ and predictions $\\hat{y}$." ], "metadata": { "collapsed": false } }, { "cell_type": "markdown", "source": [ "##### H1 simulations and test\n", "\n", "No let's apply the partial confounder test for H1, that is for a confounded model." ], "metadata": { "collapsed": false } }, { "cell_type": "code", "execution_count": 6, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "Permuting: 100%|██████████| 1000/1000 [00:01<00:00, 759.92it/s]\n" ] }, { "data": { "text/plain": " p ci lower ci upper R2(y,c) R2(y^,c) R2(y,y^)\n0 0.015 0.008419 0.02462 0.094132 0.079504 0.390694", "text/html": "
\n\n\n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n \n
pci lowerci upperR2(y,c)R2(y^,c)R2(y,y^)
00.0150.0084190.024620.0941320.0795040.390694
\n
" }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "H1_y, H1_c, H1_yhat = simulate_y_c_yhat(cov_y_c=0.3,\n", " y_ratio_yhat=0.4, c_ratio_yhat=0.1,\n", " n=500, random_state=42)\n", "ret=partial_confound_test(H1_y, H1_yhat, H1_c, num_perms=1000, return_null_dist=True,\n", " random_state=42, n_jobs=-1)\n", "\n", "#pretty print results\n", "pd.DataFrame({\n", " 'p' : [ret.p],\n", " 'ci lower' : [ret.p_ci[0]],\n", " 'ci upper' : [ret.p_ci[1]],\n", " 'R2(y,c)' : [ret.r2_y_c],\n", " 'R2(y^,c)' : [ret.r2_yhat_c],\n", " 'R2(y,y^)' : [ret.r2_y_yhat],\n", "})" ], "metadata": { "collapsed": false, "pycharm": { "name": "#%%\n" } } }, { "cell_type": "code", "execution_count": 7, "outputs": [ { "data": { "text/plain": "", "image/svg+xml": "\n\n\n\n\n\n%3\n\n\n\nc\n\nage\n\n\n\ny\n\nIQ\n\n\n\nc--y\n\n0.094\n\n\n\nyhat\n\nprediction\n\n\n\nc--yhat\n\n0.08 (p=0.02*)\n\n\n\ny--yhat\n\n0.391\n\n\n\n" }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" }, { "data": { "text/plain": "
", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEWCAYAAACJ0YulAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAAgNElEQVR4nO3de1xUZf4H8M8o9zQJl4uWVqKro22Il3RW0xyuhsjFWHcrwsuGbophxqaZbq9SW1trXa3WEGut1q2NRVTAdEMNNCQLL5s7mW6ZWDJTCCKOM8PA+f2hPj9RLjPAnLl93q+Xrxcz55nnfJ+ZkQ/n9hyFJEkSiIiIAHSzdwFEROQ4GApERCQwFIiISGAoEBGRwFAgIiKBoUBERAJDgVzW4sWL8ec//xkAUF5ejgkTJlj82tTUVHz44YcAgO3bt2PWrFldVldcXBzKy8sBAOvXr8fTTz/dZX1v2LABS5cu7bL+yP142LsAIkc3depUTJ06td12ixcvRnBwMBYuXNhmu8LCwi6pq7y8HFlZWSgpKRHPzZ07t0v6JvfFLQUimZjNZnuXQNQuhgI5NLVajU2bNiE+Ph4jR45EZmYmjEYjACAvLw+/+c1vmrUfPHgwvvvuO6vXc+DAAcTGxmLkyJF44YUXcP2F/tevR5IkrFq1CiqVCiNGjEB8fDy+/vprfPDBB9ixYwc2bdqE8PBw8Re7Wq1GdnY24uPjMXz4cJjNZqjVanz66aeif5PJhMzMTISHhyMpKQlfffVVq+O5tktMr9fj8ccfh06nQ3h4OMLDw6HVam/aHVVcXIy4uDiMGjUKqamp+N///mfRe0vui6FADm/nzp3IyclBcXExTpw4gby8vC7t//z585g/fz4yMzNx8OBB9O/fHxUVFS223b9/Pz7//HPs2rULX3zxBdauXQt/f39Mnz4d8fHxmD17Ng4fPowNGzaI1xQWFiI7Oxuff/45PDxu3mNbXFyM2NhYfPbZZ5gyZQqeeOIJNDQ0tFmzn58fNm7ciKCgIBw+fBiHDx9GcHBwszbffvstFi1ahGeffRZlZWWYMGEC5s6dC5PJJNrY+r0l58NQIIeXmpqK4OBg+Pv7Y9KkSdBoNF3af0lJCQYNGoTY2Fh4enoiLS0NP/vZz1ps6+HhgUuXLuGbb76BJEkIDQ1FUFBQu/X36dMHPj4+LS4fNmyYWPfMmTNhMplw9OjRTo+rqKgIEydOxLhx4+Dp6YnZs2fDYDDg8OHDzWqz5XtLzoehQA4vMDBQ/Ozr6wu9Xt+l/et0OoSEhIjHCoUCffr0abGtSqXCI488ghdeeAEqlQrLli1DfX19m/231tc116+7W7duCA4Ohk6ns2IELdPpdOjbt2+zvvv06QOtViues/V7S86HoUBOy9fXFwaDQTz+8ccfO9RPYGAgqqqqxGNJknDu3LlW2z/22GPIy8tDUVERTp8+jZycHABXwqQlrT1/zfXrbmpqglarFVsfvr6+uHz5slh+/Rjb6zcoKAg//PCDeHxtXDfuZiK6HkOBnNaQIUNw8uRJaDQaGI1GrF+/vkP9TJw4ESdPnsTu3bthNpvxzjvv4Keffmqx7bFjx3D06FE0NDTA19cXXl5e6Nbtyn+j3r174+zZs1av//jx42LdmzdvhpeXF8LCwsQYCwoK0NjYiJKSEhw6dEi8rnfv3qitrcXFixdb7Hfy5Mn45JNPUFZWhoaGBrz11lvw8vJCeHi41TWS+2AokNO6++67MW/ePMyYMQPR0dEYOXJkh/oJCAjAX/7yF7zyyisYM2YMvvvuO4wYMaLFtpcuXcJzzz2H++67D5MmTYK/vz9mz54NAHjooYdw6tQpjBo1Ck888YTF64+IiEBRURFGjx6Nbdu2Yf369fD09AQALF26FHv37sWoUaOwY8cOREZGiteFhoYiLi4OkZGRGDVqVLPdQgAwYMAA/OlPf8KLL76IsWPHYu/evdiwYQO8vLysfYvIjSh4kx0iIrqGWwpERCQwFIiISGAoEBGRwFAgIiLBqWdJPXLkCLy9vVtcZjQaW13mKjhG1+DSYzxxAgBgvOsu1x3jVc70ORqNRgwfPrzFZU4dCt7e3lAqlS0u02g0rS5zFRyja3DpMf7udwAAzV//6rpjvMqZPse2pjPh7iMiIhIYCkREJDAUiIhIYCgQEZHAUCAiIoGhQEREAkOBiIgEhgIREQkMBSIiEhgKBAAwmZts0paInItTT3NBXcfLoxt+nV1mUdv301U2roaI7IVbCkREJDAUiIhIYCgQEZHAUCAiIoGhQEREAkOBiIgEhgIREQkMBSIiEhgKREQkMBSIiEhgKBARkcBQICIigaFAREQCQ4GIiASGAhERCQwFIiISGApERCQwFIiISLBZKCxZsgQqlQpTpkwRz9XW1mLmzJmIjo7GzJkzceHCBQCAJElYsWIFoqKiEB8fj+PHj9uqLCIiaoPNQiE5ORk5OTnNnsvOzoZKpcLu3buhUqmQnZ0NACgpKcHp06exe/duvPjii3j++edtVRYREbXBZqEwevRo9OrVq9lzxcXFSExMBAAkJibi448/bva8QqHA8OHDUVdXB51OZ6vSiIioFbIeU6iurkZQUBAAIDAwENXV1QAArVaLkJAQ0S4kJARarVbO0oiICICHvVasUCigUCg61YfRaIRGo2lxmcFgaHWZq+jKMSqVSuj1eovby/Xe8nN0bv2vfqdceYzXuMoYZQ2F3r17Q6fTISgoCDqdDgEBAQCA4OBgVFVViXZVVVUIDg5utz9vb28olcoWl2k0mlaXuYquHqOfn5/FbeV6b/k5Ormr3ykfHx/XHeNVzvQ5thVesu4+UqvVyM/PBwDk5+cjIiKi2fOSJOHIkSPo2bOn2M1ERETysdmWwlNPPYXPPvsMNTU1mDBhAjIyMpCeno7MzEzk5uaib9++WLt2LQBg4sSJ+OSTTxAVFQVfX1+sWrXKVmVRFzCZm+DlYdnfE9a0JSL7s1kovPrqqy0+v3nz5pueUygU+MMf/mCrUqiLeXl0w6+zyyxq+366ysbVEFFX4p9wREQkMBSIiEhgKBARkcBQICIigaFAREQCQ4GIiASGAhERCQwFIiISGApERCQwFMihmMxNzR63NcHYjW2JqPPsNnU2UUtunEJDr9e3Onsrp9Ag6nrcUiAiIoGhQEREAkOBiIgEhgIREQkMBSIiEhgKREQkMBTIpngtAZFz4XUKZFPW3LoTsO7aA94rmqjrMRTIafFe0URdj386ERGRwFAgIiKBoUBERAJDgYiIBIYCEREJDAUiIhIYCi6MF44RkbV4nYIL43n8RGQtu4TC3/72N3z44YdQKBT4+c9/jpdeegk6nQ5PPfUUamtrMWzYMLz88svw8vKyR3lERG5L9t1HWq0W77zzDv71r3+hoKAAjY2NKCwsxJo1azBjxgz8+9//xq233orc3Fy5SyMicnt2OabQ2NgIg8EAs9kMg8GAwMBAHDx4EDExMQCApKQkFBcX26M0IiK3Jvvuo+DgYMyaNQuTJk2Ct7c3xo0bh2HDhuHWW2+Fh8eVckJCQqDVatvty2g0QqPRtLjMYDC0usxVtDdGpVIJvV5vcX+O0PbG9k1NTW2+3pq+HfX74Mrf1f5XPx9XHuM1rjJG2UPhwoULKC4uRnFxMXr27Iknn3wSpaWlHerL29sbSqWyxWUajabVZa7CkjH6+flZ3J8jtL2xvV6vb/P11vTtqN8Hl/6uXv18fHx8XHeMVznT59hWeMkeCp9++inuuOMOBAQEAACio6NRUVGBuro6mM1meHh4oKqqCsHBwXKXRkTk9mQ/ptC3b18cPXoUly9fhiRJKCsrw8CBAzFmzBjs2rULALB161ao1Wq5SyMicnuybymEhYUhJiYGSUlJ8PDwgFKpxPTp0/HAAw9g4cKFWLt2LZRKJVJSUuQujYjI7dnlOoUFCxZgwYIFzZ7r168fT0MlIrIzTnPhRG6ctsJZDmoRkfPgNBdO5MZpK9o7M4dTVxCRtbilQEREAkOBiIgEhgIREQkMBSIiEhgKREQkMBSIiEhgKBARkcBQICIigaFAREQCQ4GIiASGAhERCQwFIiISGApERCQwFIiISGAoEBGRwFAgIiKBoUBERAJDgYiIBIYCuYUb72/d1e2JXIVF92hOS0vD5s2b232OyFHdeH/r9vD+1uSu2gwFo9GIy5cvo6amBhcuXIAkSQCA+vp6aLVaWQokIiL5tBkK77//PjZv3gydTofk5GQRCj169MCjjz4qS4FERCSfNkMhLS0NaWlpePfdd5GamipXTUR2ZzI3wcvDskNu1rQlcnQWHVNITU1FRUUFvv/+ezQ2NornExMTbVUXkV1ZcwyCxx/IlVgUCllZWaisrMSQIUPQvXt3AIBCoWAoEBG5GItC4csvv0RRUREUCoWt6yEiIjuyaEfooEGD8OOPP3bZSuvq6rBgwQLExsZi8uTJOHz4MGprazFz5kxER0dj5syZuHDhQpetj4iILGPRlkJNTQ3i4uJw7733wtPTUzy/YcOGDq105cqVuP/++7Fu3TqYTCYYDAZs2LABKpUK6enpyM7ORnZ2NrKysjrUPxERdYxFoZCRkdFlK7x48SIOHTqEP/7xjwAALy8veHl5obi4GO+++y6AKwewU1NTGQpERDKzKBTuu+++Llvh2bNnERAQgCVLluCrr77CsGHDsHTpUlRXVyMoKAgAEBgYiOrq6nb7MhqN0Gg0LS4zGAytLnNWSqUSer1ePG5qamr2uCXtLXe0tje2b2+MctXRns5811zxu3pN/6vvoSuP8RpXGaNFoRAeHi4OMjc0NMBsNsPX1xcVFRVWr9BsNuO///0vli1bhrCwMKxYsQLZ2dnN2igUCosOant7e0OpVLa4TKPRtLrMmfn5+Ymf9Xp9s8fttbemb3u1vbF9e2OUq472dOa75qrfVQDA1ffQx8fHdcd4lTN9jm2Fl0WhcPjwYfGzJEkoLi7GkSNHOlRMSEgIQkJCEBYWBgCIjY1FdnY2evfuDZ1Oh6CgIOh0OgQEBHSofyIi6jirL8NUKBSIjIzE/v37O7TCwMBAhISE4JtvvgEAlJWVITQ0FGq1Gvn5+QCA/Px8REREdKh/IiLqOIu2FHbv3i1+bmpqwpdffglvb+8Or3TZsmV4+umn0dDQgH79+uGll15CU1MTMjMzkZubi759+2Lt2rUd7p+IiDrGolDYu3ev+Ll79+64/fbb8cYbb3R4pUqlEnl5eTc9z6m4iYjsy6JQeOmll2xdBxEROQCLjilUVVVh3rx5UKlUUKlUyMjIQFVVla1rIyIimVkUCkuWLIFarUZpaSlKS0sxadIkLFmyxNa1ERGRzCwKhfPnz2PatGnw8PCAh4cHkpOTcf78eVvXRkREMrMoFPz9/bFt2zY0NjaisbER27Ztg7+/v41LIyIiuVkUCqtWrcLOnTsxbtw4jB8/Hrt27RJzFxERkeuw6OyjdevWYfXq1ejVqxcAoLa2FqtXr+ZZSURELsaiLYUTJ06IQACu7E5yhYmfiIioOYtCoampqdlNb2pra5vdq5mIiFyDRbuPZs2ahenTpyM2NhYA8NFHH2Hu3Lk2LYyIiORnUSgkJibinnvuwcGDBwEAr732GgYOHGjTwoiISH4WhQIADBw4kEFAROTirJ46m4iIXBdDgYiIBIYCEREJDAUiIhIYCkREJDAUiIhIYCgQEZHAUCAiIoGhQEREAkOBiIgEhgIREQkMBSIiEhgKREQkMBSIiEhgKBARkWC3UGhsbERiYiLmzJkDAKisrERKSgqioqKQmZkJk8lkr9KIrGIyN9mkLZE9WHyTna72zjvvIDQ0FPX19QCANWvWYMaMGYiLi8Py5cuRm5uLhx9+2F7lEVnMy6Mbfp1dZlHb99NVNq6GqHPssqVQVVWFffv24aGHHgIASJKEgwcPIiYmBgCQlJSE4uJie5RGROTW7LKlsGrVKmRlZeHSpUsAgJqaGtx6663w8LhSTkhICLRabbv9GI1GaDSaFpcZDIZWlzkrpVIJvV4vHjc1NTV73JL2ljta2xvbtzdGueroyrY3fi9d8bt6Tf+r74srj/EaVxmj7KGwd+9eBAQE4J577kF5eXmn+vL29oZSqWxxmUajaXWZM/Pz8xM/6/X6Zo/ba29N3/Zqe2P79sYoVx1d2fbG76WrflcBAFffFx8fH9cd41XO9Dm2FV6yh0JFRQX27NmDkpISGI1G1NfXY+XKlairq4PZbIaHhweqqqoQHBwsd2lERG5P9mMKixYtQklJCfbs2YNXX30VY8eOxSuvvIIxY8Zg165dAICtW7dCrVbLXRoRkdtzmOsUsrKy8PbbbyMqKgq1tbVISUmxd0lERG7HbqekAsCYMWMwZswYAEC/fv2Qm5trz3KIiNyew2wpEBGR/TEUiIhIYCjYGac9ICJHYtdjCsQpEojIsXBLgYiIBIYCEREJDAUiIhIYCkQyaunEgtbmy+FJCGQPPNBMJKOWTixobdI/nlhA9sAtBSIiEhgKREQkMBSIiEhgKBARkcBQICIigaFAREQCQ4GIiASGAhERCQwFIiISGApERCQwFIiISGAo2AAnMiMiZ8UJ8WyAd1MjImfFLQUiIhIYCkREJDAUiIhIYCgQEZHAUCAiIoGhQEREguynpJ47dw6///3vUV1dDYVCgV/96ldIS0tDbW0tFi5ciO+//x6333471q5di169esldHhGRW5N9S6F79+5YvHgxioqK8MEHH2DLli04deoUsrOzoVKpsHv3bqhUKmRnZ8tdGpFDsfYiSF40SV1B9i2FoKAgBAUFAQB69OiBAQMGQKvVori4GO+++y4AIDExEampqcjKypK7PCKHYc1FkAAvhKSuYdcrms+ePQuNRoOwsDBUV1eLsAgMDER1dXW7rzcajdBoNC0uMxgMrS6zNaVSCb1eb3H7jrZtampq97Vy1NGVbW9s394Y5arDlm3bGqO1NdvrO9+a/lfrt+f/R7m4yhjtFgqXLl3CggUL8Oyzz6JHjx7NlikUCigUinb78Pb2hlKpbHGZRqNpdZkc/Pz8bN5Wr9e3+1o56ujKtje2b2+MctVhy7ZtjdHamu35nW/R1fp9fHwcr7YuZu/fOdZoK7zscvZRQ0MDFixYgPj4eERHRwMAevfuDZ1OBwDQ6XQICAiwR2lERG5N9lCQJAlLly7FgAEDMHPmTPG8Wq1Gfn4+ACA/Px8RERFyl0ZE5PZkD4UvvvgC27Ztw8GDB5GQkICEhAR88sknSE9Px4EDBxAdHY1PP/0U6enpcpdG5NSsOfuIZypRa2Q/pjBq1CicOHGixWWbN2+WuRoi18Ep26kr8IpmIiISGApERCQwFIiISGAoWIgH5ojIHfAezRbiQTwicgfcUiAiIoGhQEREAkOBiIgEhgIREQkMBSI3xCkxqDU8+4jIDfFsOmoNtxSIiEhgKBARkcBQICIigaFAREQCQ4GIiASGAhERCQwFIiISGApE1CZe6OZeePEaEbWJF7q5F24pEJHdcCvE8XBLgYjshlshjsdttxT4VwcR0c3cdkvBmr9QAP6VQkTuwW23FIiI6GYMBSLqMq3tllUqlTbru7NtqTm33X1ERF3vxt2yy8/VAQAWry2Gn5/fTe2t2S3Lg9Ly4JYCEbkca7cUuGXx/xxqS6GkpAQrV65EU1MTUlJSkJ6ebu+SiMgJ8USSjnOYLYXGxka88MILyMnJQWFhIQoKCnDq1Cl7l0VEbqArjle0dNzElsdBbLV14zBbCseOHcOdd96Jfv36AQDi4uJQXFyMgQMH2rkyInJ11h6vaKmtXq+/6biJrY6ZWNu3NRSSJEk26dlKH330EUpLS7Fy5UoAQH5+Po4dO4bly5e3+pojR47A29tbrhKJiFyC0WjE8OHDW1zmMFsKHdHaoIiIqGMc5phCcHAwqqqqxGOtVovg4GA7VkRE5H4cJhR+8Ytf4PTp06isrITJZEJhYSHUarW9yyIicisOs/vIw8MDy5cvx29/+1s0NjZi2rRpGDRokL3LIiJyKw5zoJmIiOzPYXYfERGR/TEUiIhIcMpQKCkpQUxMDKKiopCdnX3TcpPJhMzMTERFRSElJQVnz54Vy958801ERUUhJiYGpaWlcpZtlY6O8cCBA0hOTkZ8fDySk5NRVmb5xTBy68znCAA//PADwsPDsWnTJrlKtlpnxvjVV19h+vTpiIuLQ3x8PIxGo5ylW6yjY2xoaMAzzzyD+Ph4TJ48GW+++abcpVusvTEeOnQISUlJGDp0KD766KNmy7Zu3Yro6GhER0dj69atcpXccZKTMZvNUkREhHTmzBnJaDRK8fHx0smTJ5u1ee+996Rly5ZJkiRJBQUF0pNPPilJkiSdPHlSio+Pl4xGo3TmzBkpIiJCMpvNcg+hXZ0Z4/Hjx6WqqipJkiTpxIkT0vjx42Wt3VKdGeM1GRkZUkZGhpSTkyNX2VbpzBgbGhqkKVOmSBqNRpIkSTp//rzLfVe3b98uZWZmSpIkSXq9Xpo0aZJUWVkpa/2WsGSMlZWVkkajkbKysqSdO3eK52tqaiS1Wi3V1NRItbW1klqtlmpra+UeglWcbkvh+ukwvLy8xHQY19uzZw+SkpIAADExMSgrK4MkSSguLkZcXBy8vLzQr18/3HnnnTh27Jg9htGmzoxx6NCh4vqOQYMGwWg0wmQyyT6G9nRmjADw8ccf4/bbb3foM9Q6M8YDBw5g8ODBGDJkCADgtttuQ/fu3WUfQ3s6M0aFQoHLly/DbDbDYDDA09MTPXr0sMcw2mTJGO+44w4MGTIE3bo1/5W6f/9+jBs3Dv7+/ujVqxfGjRvn0HsoACfcfaTVahESEiIeBwcHQ6vV3tSmT58+AK6c6tqzZ0/U1NRY9FpH0JkxXm/Xrl0YOnQovLy8bF+0lTozxkuXLmHjxo2YP3++rDVbqzNj/Pbbb6FQKDB79mwkJSVh48aNstZuqc6MMSYmBr6+vhg/fjwmTZqEWbNmwd/fX87yLdKZ3xvO8jvneg5znQJ1rZMnT2LNmjV466237F1Kl3vttdeQlpaGW265xd6l2ExjYyO++OIL5ObmwtfXFzNmzMA999wDlcp1png+duwYunXrhtLSUtTV1eHhhx/GL3/5SzEpJtmH020pWDIdRnBwMM6dOwcAMJvNuHjxIm677TanmUqjM2MEgKqqKsyfPx+rV69G//795SvcCp0Z49GjR7FmzRqo1Wps3rwZb775Jt577z1Z67dEZ8YYEhKC0aNHIyAgAL6+vpgwYQKOHz8ua/2W6MwYCwoKcP/998PT0xO9e/fGiBEj8J///EfW+i3Rmd8bzvI753pOFwqWTIehVqvFUf5du3Zh7NixUCgUUKvVKCwshMlkQmVlJU6fPo17773XHsNoU2fGWFdXh/T0dCxatAgjR460R/kW6cwYt2zZgj179mDPnj1IS0vDnDlz8Oijj9pjGG3qzBjHjx+Pr7/+WuxzP3TokENOI9+ZMfbp0wfl5eUArkw7ffToUQwYMED2MbSnM1PwjB8/Hvv378eFCxdw4cIF7N+/H+PHj7dxxZ1k18PcHbRv3z4pOjpaioiIkN544w1JkiRp7dq10scffyxJkiQZDAYpIyNDioyMlKZNmyadOXNGvPaNN96QIiIipOjoaGnfvn12qd8SHR3j66+/LoWFhUlTp04V/3766Se7jaMtnfkcr1m3bp3Dnn0kSZ0bY35+vvTggw9KcXFx0urVq+1SvyU6Osb6+nopIyNDevDBB6XJkydLGzdutNsY2tPeGI8ePSrdf//9UlhYmHTfffdJDz74oHjthx9+KEVGRkqRkZFSbm6uXeq3Bqe5ICIiwel2HxERke0wFIiISGAoEBGRwFAgIiKBoUBERAJDgdyWUqlEQkICpkyZgrlz56Kurg4AoNFoms1OWlRU1Ox1CxYsQGVlZbv9m81mpKenY8yYMfj666+tqu3aSYHr168Xj00mEx555BGYzWar+iKyBkOB3JaPjw+2bduGgoIC9OrVC3//+9/F86tXr0ZhYSFycnKwatUqERgnT55EY2OjRVMxPP/88xgwYABef/11LFy4sNmVre3Zvn07cnJyYDQasXHjRmzfvh1eXl5QqVQ3hRRRV2IoEAEYPny4mKjs7rvvxl133QXgyjQFAQEBOH/+PABgx44diIiIAADk5uZi5cqVoo9//vOfWLVqFYAr8zP17NkTixcvxqhRo7BixQo89dRTuHjx4k3r/umnnzBv3jxMnToVU6dORUVFBRISEhASEoJNmzahb9++SEhIAABERkZix44dNnsfiBgK5PYaGxtRVlbW4tQFx44dQ0NDg5hDqqKiAsOGDQMATJ48GXv37kVDQwMAIC8vD9OmTQMAzJ8/H88884zoJzw8HFu2bEHPnj1vWseKFSswevRobN++HVu3bsWgQYOwY8cOVFVVYfbs2fjhhx9EEAwaNMgh5wci18FZUsltGQwGJCQkQKvVIjQ0FOPGjWu2XKfTISsrC6tXrxbz5P/4448ICAgAANxyyy0YO3Ys9u3bhwEDBqChoQGDBw+2uo6DBw/i5ZdfBgB0794dPXv2xJQpU6BQKLB+/Xo8/vjj4hhD9+7d4enpifr6eoe89wA5P24pkNu6dkxh7969kCRJHFMAgPr6esyZMwcLFy7E8OHDxfPe3t7NbouZkpKCvLw85OXlITk5uctqUygUAICMjIxmj4Ert7f09vbusnURXY+hQG7P19cXzz33HN5++22YzWaYTCbMmzcPCQkJiI2NbdY2NDQUZ86cEY/DwsJQVVWFgoICTJkyxeJ1vvfee2K6b5VKhS1btgC4siurpeMO19TU1OC2226Dp6enNUMkshhDgQjA0KFDMXjwYBQUFGDnzp34/PPPsXXrViQkJCAhIQEajQYAMHHiRDHd8zWTJ0/GiBEj0KtXL4vX980334i7jC1duhTl5eWIj49HcnIyTp061errysvL8cADD1g9PiJLcZZUIisYDAY89thj+Mc//iHumTxnzhzMmDHDqruizZkzB+vXr7f6Vqnz58/HokWLcPfdd1v1OiJLMRSIrFRaWorQ0FD06NEDKSkpGDx4MNatW2fz9ZpMJhQVFSExMdHm6yL3xVAgIiKBxxSIiEhgKBARkcBQICIigaFAREQCQ4GIiIT/A9fuHzC7G78rAAAAAElFTkSuQmCC\n" }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plot_null_dist(ret)\n", "# The labels on the graph plot can be customized:\n", "plot_graph(ret, y_name='IQ', yhat_name='prediction', c_name='age', outfile_base='example')" ], "metadata": { "collapsed": false, "pycharm": { "name": "#%%\n" } } }, { "cell_type": "markdown", "source": [ "The low p-value provides evidence against the null hypothesis of\n", "$y$ being conditionally independent on $c$ given $y$ and indicates that the model predictions are biased.\n", "\n", "-----------------------------------------------------------------------------\n", "Note\n", "\n", "For parametric corrections for multiple comparisons\n", "(e.g. false discover##### Partial confound testy rate in case of testing many confounders),\n", "permutation based p-values must be adjusted if they are zero.\n", "A decent option could be in this case to use the upper binomial confidence limit (`p_ci[1]`), instead.\n", "\n", "-----------------------------------------------------------------------------\n", "\n", "### References\n", "*Tamas Spisak, A conditional permutation-based approach to test confounder effect and center-bias in\n", "machine learning models, in prep, 2021.*" ], "metadata": { "collapsed": false } } ], "metadata": { "hide_input": false, "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.4" } }, "nbformat": 4, "nbformat_minor": 2 }