chi2.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 52,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "from IPython.display import display\n",
    "from scipy.stats import chi2\n",
    "from scipy.stats import chisquare\n",
    "\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let $Z_1, Z_2,\\ldots Z_k$ be independent $\\mathcal{N} (0, 1)$\n",
    "distributed random variables. The sum of their squares\n",
    "\n",
    "$$\n",
    "Q = \\sum_{i=1}^{k}Z_i^2\n",
    "$$\n",
    "\n",
    "is chi-squared distributed $Q \\sim \\chi^2(k)$. Where $k$ is the\n",
    "number of degrees of freedom (the only parameter of the distribution)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "def chi2_(k, size):\n",
    "    \n",
    "    return np.array([(np.random.randn(k) ** 2).sum() for _ in range(size)])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Chi-2 Test (Pearson)\n",
    "### Fairnes of Dice\n",
    "https://en.wikipedia.org/wiki/Pearson%27s_chi-squared_test#Examples"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A dice is thrown 60 times. The null hypothesis is that the dice is fair\n",
    "and each result (side of the dice) $i$ is expected to occur 10 times\n",
    "(10 = 60 / 6). Observed values are the actual results. Is the dice fair?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "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>i</th>\n",
       "      <th>observed</th>\n",
       "      <th>expected</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1</td>\n",
       "      <td>5</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2</td>\n",
       "      <td>8</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>3</td>\n",
       "      <td>9</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>4</td>\n",
       "      <td>8</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>5</td>\n",
       "      <td>10</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>6</td>\n",
       "      <td>20</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   i  observed  expected\n",
       "0  1         5        10\n",
       "1  2         8        10\n",
       "2  3         9        10\n",
       "3  4         8        10\n",
       "4  5        10        10\n",
       "5  6        20        10"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data = pd.DataFrame({\n",
    "    'i': [1, 2, 3, 4, 5, 6],\n",
    "    'observed': [5, 8, 9, 8, 10, 20],\n",
    "    'expected': [10] * 6})\n",
    "\n",
    "data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The test statistics $T$ is computed as\n",
    "\n",
    "$$\n",
    "T = \\sum_{i=1}^{6}\\frac{(O_i-E_i)^2}{E_i}\n",
    "$$\n",
    "\n",
    "where $E_i$ is the expected result and $O_i$ is the observed value."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "13.4"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "T = np.sum((data['observed'] - data['expected']) ** 2 / data['expected'])\n",
    "T"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "11.070497693516351"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# PPF is the inverse CDF (we want to know the value at 95 % percentile)\n",
    "# k = 6 - 1 = 5 (number of degrees of freedom)\n",
    "\n",
    "chi2.ppf(0.95, 5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 70,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# We reject our null hypothesis at 95 % significance level\n",
    "T < 11.07"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.95"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "chi2.cdf(11.070497693516351, 5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.019905220334774376"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# alternatively, we can check the percentile for our statistics (p-value)\n",
    "# it is lower than 0.05, we reject the null hypothesis at 95 % significance level\n",
    "1 - chi2.cdf(13.4, 5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Power_divergenceResult(statistic=13.4, pvalue=0.01990522033477438)"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
     
    }
   ],
   "source": [
    "chisquare(data['observed'], data['expected'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Workers in Neighborhood\n",
    "https://en.wikipedia.org/wiki/Chi-squared_test#Example_chi-squared_test_for_categorical_data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "City of 1M people, sample of 650 workers. The null hypothesis\n",
    "is that each persons' place of residence (neighborhood A, B, C, D) is\n",
    "independent of their working class."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "metadata": {},
   "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>class</th>\n",
       "      <th>neighborhood</th>\n",
       "      <th>observed</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>white_collar</td>\n",
       "      <td>A</td>\n",
       "      <td>90</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>white_collar</td>\n",
       "      <td>B</td>\n",
       "      <td>60</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>white_collar</td>\n",
       "      <td>C</td>\n",
       "      <td>104</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>white_collar</td>\n",
       "      <td>D</td>\n",
       "      <td>95</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>blue_collar</td>\n",
       "      <td>A</td>\n",
       "      <td>30</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>blue_collar</td>\n",
       "      <td>B</td>\n",
       "      <td>50</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>blue_collar</td>\n",
       "      <td>C</td>\n",
       "      <td>51</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>blue_collar</td>\n",
       "      <td>D</td>\n",
       "      <td>20</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>no_colar</td>\n",
       "      <td>A</td>\n",
       "      <td>30</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>no_colar</td>\n",
       "      <td>B</td>\n",
       "      <td>40</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>no_colar</td>\n",
       "      <td>C</td>\n",
       "      <td>45</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>no_colar</td>\n",
       "      <td>D</td>\n",
       "      <td>35</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "           class neighborhood  observed\n",
       "0   white_collar            A        90\n",
       "1   white_collar            B        60\n",
       "2   white_collar            C       104\n",
       "3   white_collar            D        95\n",
       "4    blue_collar            A        30\n",
       "5    blue_collar            B        50\n",
       "6    blue_collar            C        51\n",
       "7    blue_collar            D        20\n",
       "8       no_colar            A        30\n",
       "9       no_colar            B        40\n",
       "10      no_colar            C        45\n",
       "11      no_colar            D        35"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "data = pd.DataFrame(\n",
    "    {\n",
    "        'A': [90, 30, 30],\n",
    "        'B': [60, 50, 40],\n",
    "        'C': [104, 51, 45],\n",
    "        'D': [95, 20, 35]\n",
    "    },\n",
    "    index=['white_collar', 'blue_collar', 'no_colar']\n",
    ")\n",
    "\n",
    "data = pd.DataFrame(data.stack()).reset_index()\n",
    "data.columns = ['class', 'neighborhood', 'observed']\n",
    "\n",
    "n = data['observed'].sum()\n",
    "\n",
    "display(data)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is our initial data, we have the observed number of people from each\n",
    "class living in each neighborhood. The first step is to calculate the totals\n",
    "for each class, neighboord."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "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>class</th>\n",
       "      <th>class_total</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>blue_collar</td>\n",
       "      <td>151</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>no_colar</td>\n",
       "      <td>150</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>white_collar</td>\n",
       "      <td>349</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          class  class_total\n",
       "0   blue_collar          151\n",
       "1      no_colar          150\n",
       "2  white_collar          349"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>neighborhood</th>\n",
       "      <th>neighborhood_total</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>A</td>\n",
       "      <td>150</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>B</td>\n",
       "      <td>150</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>C</td>\n",
       "      <td>200</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>D</td>\n",
       "      <td>150</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "  neighborhood  neighborhood_total\n",
       "0            A                 150\n",
       "1            B                 150\n",
       "2            C                 200\n",
       "3            D                 150"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "class_totals = (\n",
    "    data.groupby(['class'])\n",
    "        .agg({'observed': sum})\n",
    "        .reset_index()\n",
    "        .rename(columns={'observed': 'class_total'})\n",
    ")\n",
    "\n",
    "neighborhood_totals = (\n",
    "    data.groupby(['neighborhood'])\n",
    "        .agg({'observed': sum})\n",
    "        .reset_index()\n",
    "        .rename(columns={'observed': 'neighborhood_total'})\n",
    ")\n",
    "\n",
    "display(class_totals)\n",
    "display(neighborhood_totals)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We have $n = 650$ observations with the following class share: $blue = \\frac{151}{650} \\approx 0.232,\n",
    "no \\approx 0.231, white \\approx 0.537$.\n",
    "\n",
    "Under the null hypothesis, share of each class in each neighborhood is\n",
    "the same as the global class share.\n",
    "For example, the expected number of $white$ collar workers in neighborhood $C$\n",
    "is $200 \\times \\frac{349}{650} \\approx 0.537 \\times 200 \\approx 107$.\n",
    "\n",
    "We get all the expected values by joining the `class_totals` and\n",
    "`neighborhood totals` with the original data and computing\n",
    "`(class_total / 650) * neighborhood_total` for each row."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 63,
   "metadata": {},
   "outputs": [],
   "source": [
    "final = (\n",
    "    data.merge(class_totals, on='class')\n",
    "        .merge(neighborhood_totals, on='neighborhood')    \n",
    ")\n",
    "\n",
    "final['expected'] = final.apply(\n",
    "    lambda r: (r['class_total'] / n) * r['neighborhood_total'], axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "metadata": {},
   "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>class</th>\n",
       "      <th>neighborhood</th>\n",
       "      <th>observed</th>\n",
       "      <th>class_total</th>\n",
       "      <th>neighborhood_total</th>\n",
       "      <th>expected</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>white_collar</td>\n",
       "      <td>A</td>\n",
       "      <td>90</td>\n",
       "      <td>349</td>\n",
       "      <td>150</td>\n",
       "      <td>80.538462</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>blue_collar</td>\n",
       "      <td>A</td>\n",
       "      <td>30</td>\n",
       "      <td>151</td>\n",
       "      <td>150</td>\n",
       "      <td>34.846154</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>no_colar</td>\n",
       "      <td>A</td>\n",
       "      <td>30</td>\n",
       "      <td>150</td>\n",
       "      <td>150</td>\n",
       "      <td>34.615385</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>white_collar</td>\n",
       "      <td>B</td>\n",
       "      <td>60</td>\n",
       "      <td>349</td>\n",
       "      <td>150</td>\n",
       "      <td>80.538462</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>blue_collar</td>\n",
       "      <td>B</td>\n",
       "      <td>50</td>\n",
       "      <td>151</td>\n",
       "      <td>150</td>\n",
       "      <td>34.846154</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>no_colar</td>\n",
       "      <td>B</td>\n",
       "      <td>40</td>\n",
       "      <td>150</td>\n",
       "      <td>150</td>\n",
       "      <td>34.615385</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>white_collar</td>\n",
       "      <td>C</td>\n",
       "      <td>104</td>\n",
       "      <td>349</td>\n",
       "      <td>200</td>\n",
       "      <td>107.384615</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>blue_collar</td>\n",
       "      <td>C</td>\n",
       "      <td>51</td>\n",
       "      <td>151</td>\n",
       "      <td>200</td>\n",
       "      <td>46.461538</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>no_colar</td>\n",
       "      <td>C</td>\n",
       "      <td>45</td>\n",
       "      <td>150</td>\n",
       "      <td>200</td>\n",
       "      <td>46.153846</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>white_collar</td>\n",
       "      <td>D</td>\n",
       "      <td>95</td>\n",
       "      <td>349</td>\n",
       "      <td>150</td>\n",
       "      <td>80.538462</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>blue_collar</td>\n",
       "      <td>D</td>\n",
       "      <td>20</td>\n",
       "      <td>151</td>\n",
       "      <td>150</td>\n",
       "      <td>34.846154</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>no_colar</td>\n",
       "      <td>D</td>\n",
       "      <td>35</td>\n",
       "      <td>150</td>\n",
       "      <td>150</td>\n",
       "      <td>34.615385</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "           class neighborhood  observed  class_total  neighborhood_total  \\\n",
       "0   white_collar            A        90          349                 150   \n",
       "1    blue_collar            A        30          151                 150   \n",
       "2       no_colar            A        30          150                 150   \n",
       "3   white_collar            B        60          349                 150   \n",
       "4    blue_collar            B        50          151                 150   \n",
       "5       no_colar            B        40          150                 150   \n",
       "6   white_collar            C       104          349                 200   \n",
       "7    blue_collar            C        51          151                 200   \n",
       "8       no_colar            C        45          150                 200   \n",
       "9   white_collar            D        95          349                 150   \n",
       "10   blue_collar            D        20          151                 150   \n",
       "11      no_colar            D        35          150                 150   \n",
       "\n",
       "      expected  \n",
       "0    80.538462  \n",
       "1    34.846154  \n",
       "2    34.615385  \n",
       "3    80.538462  \n",
       "4    34.846154  \n",
       "5    34.615385  \n",
       "6   107.384615  \n",
       "7    46.461538  \n",
       "8    46.153846  \n",
       "9    80.538462  \n",
       "10   34.846154  \n",
       "11   34.615385  "
      ]
     },
     "execution_count": 64,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "final"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We compute the same test statistics $T$ as in the previous example."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 65,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "24.571202858582595"
      ]
     },
     "execution_count": 65,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "T = np.sum((final['observed'] - final['expected']) ** 2 / final['expected'])\n",
    "T"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The number of degrees of freedom $k$ of our $\\chi^2(k)$ distribution is\n",
    "\n",
    "$$\n",
    "k = uv - 1 = 12 - 1 = 11\n",
    "$$\n",
    "\n",
    "where $u$ is the number of working classes and $v$ is the number of neighborhoods."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "19.67513757268249"
      ]
     },
     "execution_count": 72,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# PPF is the inverse CDF (we want to know the value at 95 % percentile)\n",
    "# k = 6\n",
    "\n",
    "chi2.ppf(0.95, 11)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 73,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# we reject our null hypothesis at 95 % significance level\n",
    "T < 19.68"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 74,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.010532857023036746"
      ]
     },
     "execution_count": 74,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# alternatively, we can check the percentile for our statistics (p-value)\n",
    "# it is lower than 0.05, we reject the null hypothesis at 95 % significance level\n",
    "1 - chi2.cdf(24.57, 11)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 75,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Power_divergenceResult(statistic=24.571202858582595, pvalue=0.010528619368354175)"
      ]
     },
     "execution_count": 75,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "chisquare(final['observed'], final['expected'])"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}