等分散の検定

Ftest.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "from sklearn import datasets\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy.stats import f"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "iris = datasets.load_iris()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['sepal length (cm)',\n",
       " 'sepal width (cm)',\n",
       " 'petal length (cm)',\n",
       " 'petal width (cm)']"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "iris[\"feature_names\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(['setosa', 'versicolor', 'virginica'], \n",
       "      dtype='<U10')"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "iris[\"target_names\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "X = iris[\"data\"]\n",
    "ys = iris[\"target\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "# A: setosa, B: versicolor, sepal length\n",
    "A = X[ys == 0][:, 0]\n",
    "B = X[ys == 1][:, 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.text.Text at 0x7f1b9b5b32b0>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAD8CAYAAACYebj1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAD7dJREFUeJzt3XuMZ2V9x/H3RxAvIJeRrVK5rLYEoqkoTrhUSgVSWxTB\nVlrBtAq13WIRtE1stI0V/aPGXpKqJGy3Uiv1WqkotkC1tiq21XQWV7xbChp2izLQdUEwIPDtH3P2\ncXYyu3Nm2TOHmXm/kl9+v/Oc53fOd5LJfOac85zzpKqQJAngUWMXIEl65DAUJEmNoSBJagwFSVJj\nKEiSGkNBktQYCpKkxlCQJDWGgiSp2XvsAhbr4IMPrrVr145dhiQtKxs3bryjqtYs1G/ZhcLatWuZ\nmpoauwxJWlaSfKdPP08fSZIaQ0GS1BgKkqTGUJAkNYaCJKkZLBSSHJVk06zXXUleO6dPkrwjyU1J\nbkxy7FD1SJIWNtiQ1Kr6JvAsgCR7AVuAq+Z0Ox04snsdD1zWvUuSRrBUp49OA/6nquaOkz0LuKJm\nfB44MMkhS1STJGmOpQqFc4APzNP+FODWWcubu7YdJFmXZCrJ1PT09EAlShpDkt16aRiDh0KSfYAz\ngQ/v7jaqakNVTVbV5Jo1C96lLWkZqap5X7tat3299rylOFI4Hbihqr43z7otwGGzlg/t2iRJI1iK\nUDiX+U8dAVwNvLwbhXQCsK2qbluCmiRJ8xj0gXhJ9gV+AfidWW0XAFTVeuAa4AXATcC9wPlD1iNJ\n2rVBQ6Gq7gGeOKdt/azPBVw4ZA2SpP68o1mS1BgKkqTGUJAkNYaCJKkxFCRJjaEgSWoMBUlSYyhI\nkhpDQZLUGAqSpMZQkCQ1hoIkqTEUJEmNoSBJagwFSVJjKEiSGkNBktQYCpKkxlCQJDWGgiSpMRQk\nSY2hIElqDAVJUmMoSJIaQ0GS1BgKkqTGUJAkNYOGQpIDk1yZ5BtJvp7kxDnrn5dkW5JN3euPh6xH\nkrRrew+8/bcD11XV2Un2AR4/T5/rq+qMgeuQJPUwWCgkOQA4GTgPoKruB+4fan+SpIdvyNNHTwWm\ngXcn+WKSdyXZd55+Jyb5UpJrkzxjwHokSQsYMhT2Bo4FLquqZwP3AK+f0+cG4IiqOgZ4J/DR+TaU\nZF2SqSRT09PTA5YsSavbkKGwGdhcVV/olq9kJiSaqrqrqn7Qfb4GeHSSg+duqKo2VNVkVU2uWbNm\nwJIlaXUbLBSq6rvArUmO6ppOA742u0+SJydJ9/m4rp47h6pJkrRrQ48+ugh4Xzfy6Gbg/CQXAFTV\neuBs4FVJHgB+CJxTVTVwTZKknRg0FKpqEzA5p3n9rPWXApcOWYMkqT/vaJYkNYaCJKkxFCRJjaEg\nSWoMBUlSYyhIkhpDQZLUGAqSpMZQkCQ1hoKkJTExMUGS3i9gUf2TMDExMfJPufwN/ewjSQJg69at\nDP1os+1hot3nkYIkqTEUJEmNoSBJagwFSVJjKEiSGkNBktQYCpKkxlCQJDWGgiSpMRQkSY2hIElq\nDAVJUmMoSJIaQ0GS1BgKkqTGUJAkNYaCJKlZMBSSPDfJJ5N8K8nNSW5JcnOfjSc5MMmVSb6R5OtJ\nTpyzPknekeSmJDcmOXZ3fxBJ0sPXZzrOy4HfAzYCDy5y+28Hrquqs5PsAzx+zvrTgSO71/HAZd27\nJGkEfUJhW1Vdu9gNJzkAOBk4D6Cq7gfun9PtLOCKmpm49fPdkcUhVXXbYvcnSXr4dhoKs07l/FuS\nPwM+Aty3fX1V3bDAtp8KTAPvTnIMM0car6mqe2b1eQpw66zlzV3bDqGQZB2wDuDwww9fYLeSpN21\nqyOFv5izPDnrcwGn9tj2scBFVfWFJG8HXg+8cbFFVtUGYAPA5ORkLfb7kqR+dhoKVXUKQJKnVdUO\nF5aTPK3HtjcDm6vqC93ylcyEwmxbgMNmLR/atUmSRtBnSOqV87R9eKEvVdV3gVuTHNU1nQZ8bU63\nq4GXd6OQTmDm+oXXEyRpJLu6pnA08AzggCS/MmvV/sBje27/IuB93cijm4Hzk1wAUFXrgWuAFwA3\nAfcC5y/6J5Ak7TG7uqZwFHAGcCDwolntdwO/3WfjVbWJHa9FAKyftb6AC3tVKkka3K6uKXwM+FiS\nE6vqP5ewJknSSPrcp/CyJOfOadsGTHXBIUlaIfqEwmOAo/nxxeWXALcAxyQ5papeO1RxklaOetP+\ncMkBw+9DD0ufUHgm8NyqehAgyWXA9cBJwJcHrE3SCpI338XMZcQB95FQlwy6ixWvz5DUg4D9Zi3v\nC0x0IXHf/F+RJC1HfY4U/hTYlOTTQJh5ntGfJNkX+JcBa5MkLbEFQ6GqLk9yDXBc1/SHVfW/3efX\nDVaZJGnJ9Z1k51HMPNxuK/DTSU4eriRJ0lgWPFJI8jbgpcBXgYe65gI+O2BdkqQR9Lmm8GLgqKry\novIylmTR3xl6pIikR54+oXAz8GgcabSs7ewPfBL/+Etq+oTCvcyMPvoUO06yc/FgVUmSRtEnFK7u\nXpKkFa7PkNT3JHkccHhVfXMJapIkjWTBIalJXgRsAq7rlp+VxCMHSVqB+tyncAkzN659H9ocCX2m\n45QkLTN9QuFHVbVtTttD8/aUJC1rfS40fzXJy4C9khwJXAz8x7BlSZLG0OdI4SJm5mq+D/gAcBfg\nHAqStAL1GX10L/BH3UuStILtNBSSfJyZZxzNq6rOHKQiSdJodnWk8OdLVoUk6RFhp6FQVZ9ZykIk\nSePrO5+CJGkV6DMkVZL2iN15hPtiHHTQQYNufzUwFCQticU+ot3Huo/D0UeSpGbQ0UdJvg3cDTwI\nPFBVk3PWPw/4GHBL1/SRqnrLw92vJGn3LMXoo1Oq6o5drL++qs7YQ/uSJD0MC15T6J539Fbg6cBj\nt7dXlU9KlaQVps+Q1HcDlwEPAKcAVwDv7bn9Aj6RZGOSdTvpc2KSLyW5Nskzem5XkjSAPqHwuKr6\nFJCq+k5VXQK8sOf2T6qqY4HTgQuTnDxn/Q3AEVV1DPBO4KPzbSTJuiRTSaamp6d77lqStFh9QuG+\nJI8C/jvJq5P8MrBfn41X1Zbu/XbgKmYm65m9/q6q+kH3+Rrg0UkOnmc7G6pqsqom16xZ02fXkqTd\n0CcUXgM8npl5FJ4D/AbwioW+lGTfJE/Y/hl4PvCVOX2enO5uliTHdfXcuZgfQJK05/R5dPZ/AXRH\nCxdX1d09t/0k4Krub/7ewPur6rokF3TbXQ+cDbwqyQPAD4FzyrtVJGk0fUYfTTJzsXn7f/3bgN+s\nqo27+l5V3QwcM0/7+lmfLwUuXWTNkqSB9HnMxd8Av1tV1wMkOYmZkHjmkIVJkpZen2sKD24PBICq\n+hwzw1MlSStMnyOFzyT5K2bmZy7gpcCnkxwLUFU3DFifJGkJ9QmF7dcF3jSn/dnMhMSpe7QiSdJo\n+ow+OmUpCpEkjW/BawpJnpTk8iTXdstPT/LK4UuTJC21Phea/xb4Z+Anu+VvAa8dqiBJ0nj6hMLB\nVfX3wEMAVfUAM/MjSJJWmD6hcE+SJ9LNwpbkBGDboFVJkkbRZ/TR7wNXAz+V5N+BNcw8nkKStML0\nGX10Q5KfB44CAnyzqn40eGWSpCXXZ/TRrzIzp8JXgRcDH9p+45okaWXpc03hjVV1d/fMo9OAy5mZ\niU2StML0evZR9/5C4K+r6p+AfYYrSZI0lj6hsKV79tFLgWuSPKbn9yRJy0yfP+6/xszNa79YVd8H\nJoDXDVqVJGkUfUYf3Qt8ZNbybcBtQxYlSRqHp4FWkImJCZIs6gUsqv/ExMTIP6WkIfW5eU3LxNat\nWxl6iuvtQSJpZfJIQZLUGAqSpMZQkCQ1hoIkqTEUJEmNoSBJagwFSVJjKEiSGkNBktQMGgpJvp3k\ny0k2JZmaZ32SvCPJTUludPIeSRrXUjzm4pSqumMn604HjuxexzMzec/xS1CTJGkeY58+Ogu4omZ8\nHjgwySEj1yRJq9bQoVDAJ5JsTLJunvVPAW6dtby5a5MkjWDo00cnVdWWJD8BfDLJN6rqs4vdSBco\n6wAOP/zwPV3jilFv2h8uOWD4fUhasQYNhara0r3fnuQq4DhgdihsAQ6btXxo1zZ3OxuADQCTk5PD\nPht6Gcub71qSR2fXJYPuQtKIBjt9lGTfJE/Y/hl4PvCVOd2uBl7ejUI6AdjWzewmSRrBkEcKTwKu\n6iZl2Rt4f1Vdl+QCgKpaD1wDvAC4CbgXOH/AeiRJCxgsFKrqZuCYedrXz/pcwIVD1SBJWpyxh6RK\nkh5BDAVJUrMUdzRL0k511x0XvW7okXarlaEgaVT+cX9k8fSRJKkxFCRJjaEgSWoMBUlSYyhIkhpH\nH60wuxrCtyccdNBBg25f0rgMhRVkd4b2JXFIoKTG00eSpMZQkCQ1hoIkqTEUJEmNoSBJagwFSVJj\nKEiSGkNBktQYCpKkxlCQJDWGgiSpMRQkSY2hIElqDAVJUmMoSJIaQ0GS1BgKkqRm8FBIsleSLyb5\nx3nWnZdkOsmm7vVbQ9cjSdq5pZiO8zXA14H9d7L+Q1X16iWoQ5K0gEGPFJIcCrwQeNeQ+5Ek7RlD\nnz76S+APgId20eclSW5McmWSw+brkGRdkqkkU9PT04MUKkkaMBSSnAHcXlUbd9Ht48Daqnom8Eng\nPfN1qqoNVTVZVZNr1qwZoFpJEgx7pPBc4Mwk3wY+CJya5L2zO1TVnVV1X7f4LuA5A9YjSVrAYKFQ\nVW+oqkOrai1wDvCvVfXrs/skOWTW4pnMXJCWJI1kKUYf7SDJW4CpqroauDjJmcADwP8B5y11PZKk\nH0tVjV3DokxOTtbU1NTYZawYSVhuvwOSFi/JxqqaXKifdzRLkhpDQZLUGAqSpMZQkCQ1hoIkqTEU\nJEmNoSBJagwFSVJjKEiSmiV/zIXGkWTR67zTWVp9DIVVwj/wkvrw9JEkqTEUJEmNoSBJagwFSVJj\nKEiSGkNBktQYCpKkxlCQJDXLbo7mJNPAd8auYwU5GLhj7CKkefi7uWcdUVVrFuq07EJBe1aSqT6T\neUtLzd/NcXj6SJLUGAqSpMZQ0IaxC5B2wt/NEXhNQZLUeKQgSWoMhVUsyYuTVJKjx65F2i7Jg0k2\nJflSkhuS/OzYNa0mhsLqdi7wue5deqT4YVU9q6qOAd4AvHXsglYTQ2GVSrIfcBLwSuCckcuRdmZ/\nYOvYRawmTse5ep0FXFdV30pyZ5LnVNXGsYuSgMcl2QQ8FjgEOHXkelYVjxRWr3OBD3afP4inkPTI\nsf300dHALwFXJMnYRa0WDkldhZJMAJuBaaCAvbr3I8pfCI0syQ+qar9Zy98Dfqaqbh+xrFXDI4XV\n6Wzg76rqiKpaW1WHAbcAPzdyXdIOupFxewF3jl3LauE1hdXpXOBtc9r+oWv/7NKXI+1g+zUFgACv\nqKoHxyxoNfH0kSSp8fSRJKkxFCRJjaEgSWoMBUlSYyhIkhpDQZLUGAqSpMZQkCQ1/w/1JR+aCTgO\ndQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f1b9da2bcf8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.boxplot(np.c_[A, B])\n",
    "plt.xticks([1, 2], [\"A\", \"B\"])\n",
    "plt.ylabel(\"sepal length\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "s_A: 0.35248968721345136, s_B: 0.5161711470638634\n"
     ]
    }
   ],
   "source": [
    "s_A = np.sqrt(np.sum((A - np.mean(A)) ** 2) / (A.shape[0] - 1))\n",
    "s_B = np.sqrt(np.sum((B - np.mean(B)) ** 2) / (B.shape[0] - 1))\n",
    "print(\"s_A: {}, s_B: {}\".format(s_A, s_B))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "F: 0.46634291316869936\n"
     ]
    }
   ],
   "source": [
    "F = s_A ** 2 / s_B ** 2\n",
    "print(\"F: {}\".format(F))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "p1 = (1 / 49 - 1 / 50) / (1 / 40 - 1 / 50)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "F_40_40 = 2.2958\n",
    "F_50_40 = 2.2295\n",
    "F_40_50 = 2.1644\n",
    "F_50_50 = 2.0967"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "F_49_40: 2.234912244897959, F_49_50: 2.1022265306122447\n"
     ]
    }
   ],
   "source": [
    "F_49_40 = p1 * F_40_40 + (1 - p1) * F_50_40\n",
    "F_49_50 = p1 * F_40_50 + (1 - p1) * F_50_50\n",
    "print(\"F_49_40: {}, F_49_50: {}\".format(F_49_40, F_49_50))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "F_49_49: 2.113058017492711\n"
     ]
    }
   ],
   "source": [
    "F_49_49 = p1 * F_49_40 + (1 - p1) * F_49_50\n",
    "print(\"F_49_49: {}\".format(F_49_49))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "inv F_49_49: 0.47324777252759437\n"
     ]
    }
   ],
   "source": [
    "inv_F_49_49 = 1 / F_49_49\n",
    "print(\"inv F_49_49: {}\".format(inv_F_49_49))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.47325020730190659"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f.ppf(0.005, A.shape[0] - 1, B.shape[0] - 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2.113047146246799"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f.ppf(0.995, A.shape[0] - 1, B.shape[0] - 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.0043285941813499922"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f.cdf(F, A.shape[0] - 1, B.shape[0] - 1)"
   ]
  }
 ],
 "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.6.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}