ModeleLineaire1.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Modèle linéaire\n",
    "\n",
    "Considérons la cas classique d'une fonction affine :\n",
    "\n",
    "$$y=ax+b$$\n",
    "\n",
    "Ici, $a$ et $b$ sont des réels. Ces deux nombres définissent entièrement la courbe et permet donc d'obtenir une relation **affine** entre $x$ et $y$. En statistique, cette relation est à la base des modèles dit **linéaires**, où une variable réponse se définit comme une somme de variables explicatives où chacune de ces dernières sont multipliés par un coefficient.\n",
    "\n",
    "\n",
    "## Modèle linéaire simple\n",
    "\n",
    "![](https://upload.wikimedia.org/wikipedia/commons/thumb/3/3a/Linear_regression.svg/438px-Linear_regression.svg.png)\n",
    "\n",
    "Dans le modèle linéaire simple (une seule variable explicative), on suppose que la variable réponse suit le modèle suivant :\n",
    "\n",
    "$$y_i=\\beta_0 + \\beta_1 x_i + \\varepsilon_i$$\n",
    "\n",
    "On remarque la ressemblance avec la fonction affine présentée ci-dessus. La différence réside dans l'existence du terme aléatoire (appelé bruit) $\\varepsilon_i$. Afin de considérer le modèle, il est nécessaire de se placer sous les hypothèses suivantes\n",
    "\n",
    "$$(\\mathcal{H}) : \\left\\{\\begin{matrix}\n",
    "\\mathbb{E}[\\varepsilon_i]=0\\\\ \n",
    "\\text{Cov}(\\varepsilon_i, \\varepsilon_j)=\\delta_{ij} \\sigma^2\n",
    "\\end{matrix}\\right.$$\n",
    "Les différents éléments qui interviennent sont :\n",
    "\n",
    "- $\\beta_0$ : l'ordonnée à l'origine (nommée *intercept*)\n",
    "- $\\beta_1$ : le coefficient directeur\n",
    "- $x_i$ : l'observation $i$\n",
    "- $y_i$ : le $i$-ème prix\n",
    "- $\\varepsilon_i$ : le bruit aléatoire liée à la $i$-ème observation\n",
    "\n",
    "La solution peut se calculer facilement via les formules fermées suivantes :\n",
    "\n",
    "$$\\hat{\\beta}_1=\\frac{\\sum_{i=1}^n (x_i - \\bar{x})(y_i - \\bar{y})}{\\sum_{i=1}^n (x_i - \\bar{x})^2} \\qquad \\hat{\\beta}_0 = \\hat{y} - \\hat{\\beta}_1 \\bar{x}$$\n",
    "\n",
    "## Modèle linéaire multiple\n",
    "\n",
    "Dans le cas multiple (pour $p$ variables explicatives), pour la $i$-ème observation, le modèle s'écrit :\n",
    "\n",
    "$$y_i= \\beta_0 + \\sum_{j=1}^p \\beta_j x_{ij} + \\varepsilon_i$$\n",
    "\n",
    "Ainsi, une observation $x_i$ n'est plus une valeur, mais un **vecteur** $(x_{i1}, \\dots, x_{ip})$. Il est plus commode de regrouper ces prix $y_i$ et ces vecteurs d'observations $x_i$ dans des matrices :\n",
    "\n",
    "$$Y=X \\beta + \\varepsilon$$\n",
    "\n",
    "Sous les hypothèses équivalentes du modèle simple en plus grand dimension\n",
    "\n",
    "$$(\\mathcal{H}) : \\left\\{\\begin{matrix}\n",
    "\\text{rank}(X)=p\\\\ \n",
    "\\mathbb{E}[\\varepsilon]=0 \\text{ et }\\text{Var}(\\varepsilon)=\\sigma^2 I_p\n",
    "\\end{matrix}\\right.$$\n",
    "\n",
    "Les différents éléments qui interviennent sont :\n",
    "\n",
    "- $\\beta$ : le vecteur directeur\n",
    "- $X$ : la matrice des observations\n",
    "- $Y$ : le vecteur de prix\n",
    "- $\\varepsilon$ : le vecteur de bruit\n",
    "\n",
    "Avec $X=( \\mathbf{1}, X_1, \\dots, X_n)$, $Y=(y_1, \\dots, y_n)^\\top$ et $\\varepsilon=(\\varepsilon_1, \\dots, \\varepsilon_n)^\\top$. La solution des MCO (Moindres Carrés Ordinaires) est alors :\n",
    "\n",
    "$$\\hat{\\beta}= (X^\\top X)^{-1} X^\\top Y$$\n",
    "\n",
    "Vous pouvez d'ailleurs faire la démonstration de votre coté ! Pour plus d'information mathématiques, je vous conseil le portail de wikipédia qui est très bien fait : [lien ici](https://fr.wikipedia.org/wiki/Portail:Probabilit%C3%A9s_et_statistiques)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Implémenter une régression linéaire \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "#importer vos librairies \n",
    "import pandas as pd \n",
    "import numpy as np \n",
    "from sklearn import linear_model \n",
    "from sklearn.metrics import mean_squared_error, r2_score"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "data loaded\n"
     ]
    }
   ],
   "source": [
    "#charger les données \n",
    "#price_availability.csv\n",
    "#listings_final.csv\n",
    "#attention aux valeurs manquantes !!\n",
    "price = pd.read_csv('data/price_availability.csv', sep=';')\n",
    "listing = pd.read_csv('data/listings_final.csv', sep=';')\n",
    "listing = listing.drop(589)\n",
    "print('data loaded')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Données d'entrée"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "L'objectif ici est de charger les données pour créer les matrices $X$ et $Y$ du modèle linéaire. **Attention**, il n'est pas nécessaire de rajouter le vecteur colonne $\\mathbf{1}$ en première colonne, car *scikit-learn* le fait automatiquement !"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "#définir vos variables de travail X et Y\n",
    "list(listing.columns)\n",
    "X = listing.loc[:, ['listing_id', 'person_capacity', 'bedrooms', 'bathrooms']]\n",
    "Y = []"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "#construire l'array des prix \n",
    " \n",
    "for i, row in X.iterrows():\n",
    "    y = 0 \n",
    "    ID = int(row['listing_id'])\n",
    "    subset = price[price['listing_id'] == ID]\n",
    "    y = subset['local_price'].mean()\n",
    "    Y.append(y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "Y = np.asarray(Y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "En *Machine Learning*, on a l'habitude de couper l'ensemble de données en deux sous-ensembles :\n",
    "\n",
    "- Un ensemble d'entraînement (*train set*), sur lequel le modèle va être calibré.\n",
    "- Un ensemble de test (*test set*), qui ne sera pas utilisé pendant le calibrage mais permettra de vérifier l'aptitude du modèle à généraliser sur de nouvelles observations inconnues.\n",
    "\n",
    "En général, on découpe l'ensemble de données (*split*) en prenant $\\alpha \\%$ de l'ensemble pour entraînement et $1-\\alpha \\%$ comme test. Dans la plus part des cas, on considère que $\\alpha=10,20\\ ou \\ 30\\%$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 63,
   "metadata": {},
   "outputs": [],
   "source": [
    "#utiliser la méthode split de sklearn en splitant avec un alpha=30 et un random state=42 \n",
    "#afficher la shape de vos données \n",
    "from sklearn.model_selection import train_test_split \n",
    "X_train, X_test, y_train, y_test = train_test_split(X,Y, test_size=0.3, random_state=42)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(699,)"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Entraînement"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Pour information, *scikit-learn* utilise le solveur OLS (Ordinary Least Squares) de *numpy*."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "LinearRegression(copy_X=True, fit_intercept=True, n_jobs=None,\n",
       "         normalize=False)"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#créer l'objet de régression et entrainer le sur notre ensemble d'entraînement\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "LinearRegression(copy_X=True, fit_intercept=True, n_jobs=None, normalize=False)"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from sklearn.linear_model import LinearRegression\n",
    "reg = LinearRegression().fit(X_train, y_train)\n",
    "reg"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "On affiche le vecteur des coefficients pour interpréter rapidement le modèle."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Coefficients beta_j : \n",
      " [2.47834609e-06 3.23451085e+01 1.43887639e+01 7.75318453e+01]\n",
      "Coefficients INTERCEPT beta_0 : \n",
      " -79.79953133672862\n"
     ]
    }
   ],
   "source": [
    "#afficher les coefficients\n",
    "#que remarquez vous ? \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([2.47834609e-06, 3.23451085e+01, 1.43887639e+01, 7.75318453e+01])"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "reg.coef_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-79.79953133672873"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "reg.intercept_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Validation du modèle\n",
    "\n",
    "### Le coefficient de détermination $R^2$\n",
    "\n",
    "Par la suite, nous ferons l'hypothèse de gaussianité sur les bruits. Dans l'idée, nous aimerions obtenir une valeur numérique qui nous indique à quel point la régression linéaire a un sens sur nos données. Pour cela, introduisons les notations suivantes :\n",
    "\n",
    "- $SCT=\\|Y-\\hat{y} \\mathbf{1}\\|^2$ est la somme des carrés totaux\n",
    "- $SCE=\\|\\hat{Y}-\\hat{y} \\mathbf{1}\\|^2$ est la somme des carrés expliqués\n",
    "- $SCR=\\|\\hat{\\varepsilon}\\|^2$ est la somme des carrés résiduels\n",
    "\n",
    "L'idée est de décomposer la somme des carrés totaux comme la somme des carrés que le modèle explique, en plus de la somme des carrés qui sont liés aux résidus (et donc que le modèle ne peut pas expliquer). On voit donc ici l'intérêt de calculer un coefficient à partir du $SCE$. Puisque l'on a la relation suivante :\n",
    "\n",
    "$$SCT=SCE+SCR \\text{ alors } 1=\\frac{SCE}{SCT}+\\frac{SCR}{SCT}$$\n",
    "\n",
    "Plus les résidus sont petits (et donc la régression est \"bonne\"), plus $SCR$ devient petit et donc $SCE$ devient grand. Le schéma inverse s'opère de la même façon. Dans le meilleur des cas, on obtient $SCR=0$ et donc $SCE=SCT$ d'où le premier membre vaut $1$. Dans le cas contraite, $SCE=0$ et automatiquement, le premier membre est nul. C'est ainsi que l'on définit le coefficient de détermination $R^2$ comme \n",
    "$$R^2=\\frac{SCE}{SCT}=1-\\frac{SCR}{SCT}$$\n",
    "Ainsi, $R^2 \\in [0,1]$. Plus $R^2$ est proche de $1$, plus la régression linéaire a du sens. Au contraire, si $R^2$ est proche de $0$, le modèle linéaire possède un faible pouvoir explicatif."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "taille du vecteur Y_pred : 999 \n"
     ]
    }
   ],
   "source": [
    "#faire une prediction sur X_test\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([113.96717254, 137.68977495,  70.97927863,  64.89744358,\n",
       "       245.81666527, 170.4788043 , 561.0179989 , 106.9418946 ,\n",
       "       383.27653253, 233.75294378, 162.22364018, 612.08191179,\n",
       "       159.32890963, 251.58804631, 118.1210298 ,  78.77744367,\n",
       "       124.68010657, 194.48881109, 434.03826221, 790.83639952,\n",
       "       586.83384388, 171.78524314, 221.95377062, 128.08854869,\n",
       "       128.62279089, 175.3227481 , 121.6775822 , 154.48438242,\n",
       "       156.15601944,  69.49655503,  92.87029058, 216.30267547,\n",
       "       112.59683353, 754.37308581, 455.87067501, 287.66991756,\n",
       "        67.79340934, 359.1764754 , 129.24823819, 214.40674172,\n",
       "       169.73401224,  60.79931294, 209.06281997, 236.21743727,\n",
       "       320.30632589, 123.62663053,  77.81376409, 198.81582295,\n",
       "        78.94863791, 145.15976783, 127.84082014, 152.26503044,\n",
       "       224.58908772,  69.26447379, 298.66775754, 164.39766518,\n",
       "       214.32584501, 185.47409534, 218.11339591, 119.70295068,\n",
       "        77.9772229 , 169.06272431, 141.63669477, 263.0876178 ,\n",
       "       188.41291802, 458.44228337, 118.66693955, 131.70905247,\n",
       "       209.49147903, 345.54429022, 197.34660175, 394.61109274,\n",
       "       198.20161495, 162.16239973, 184.34533509, 487.48883024,\n",
       "       109.84428074, 117.72615492, 119.68566371, 131.93004214,\n",
       "       228.88831719, 108.60269627, 228.89023047, 142.17461235,\n",
       "       310.3515368 , 112.37928876, 193.37097094,  72.55334513,\n",
       "       178.34262967, 195.67568608, 326.78427698,  81.19456421,\n",
       "       282.83619968,  48.30382614, 123.43324817, 133.67848061,\n",
       "        78.61651227, 189.43117688, 158.84300955, 258.02266785,\n",
       "       190.22068078, 177.0186278 , 243.49686317,  96.25899615,\n",
       "       184.43946763, 145.623853  , 104.86118574, 200.66149635,\n",
       "       130.88459397, 195.60082318, 245.68001918, 184.35430409,\n",
       "       300.57480297,  61.76247208,  96.58031817, 129.12396698,\n",
       "        73.5802751 , 605.29950136, 522.13449942,  87.89552183,\n",
       "       246.08249439, 443.31574508, 110.39853062, 128.23899421,\n",
       "       217.80494243, 194.38683559, 159.71566246,  62.33204806,\n",
       "       371.82149099, 140.72636603, 156.3426255 , 132.64695579,\n",
       "       355.70442418, 360.4112272 ,  80.0623546 , 449.36984161,\n",
       "        77.23409826, 111.19595335,  77.6648125 , 248.65426981,\n",
       "       124.60112613,  78.88402793,  80.54872309, 230.48283117,\n",
       "        83.76296151, 516.74216441, 102.85410261, 107.52542266,\n",
       "       316.53231831, 135.74406757, 166.48472937, 194.96877108,\n",
       "       272.31559052, 136.56450917, 283.99387404, 142.9966163 ,\n",
       "       239.95300002, 125.02999896, 580.47720356, 108.50639963,\n",
       "       290.95116352, 189.5267414 , 199.48103037, 116.71567126,\n",
       "        52.56102782, 143.29622505, 110.55652768, 242.82208018,\n",
       "       261.52554153, 126.60530512, 240.90803018, 364.26128633,\n",
       "       130.44362134, 129.45285591, 116.50295729, 137.25854273,\n",
       "       175.478702  ,  96.99450484, 167.73686678, 210.91415491,\n",
       "       211.48043616,  48.52797011, 170.24032595, 169.76904663,\n",
       "        71.6627817 , 221.00065572, 204.75764283, 435.69333157,\n",
       "       261.64813874, 204.62187904, 143.16409946, 130.38279828,\n",
       "       139.26692736,  73.99393882, 370.493038  ,  25.79682375,\n",
       "       156.30663297, 172.08365533, 288.60255699, 111.2906802 ,\n",
       "       104.55680519, 284.23603646, 147.36196272, 342.68064868,\n",
       "       103.01736664, 121.54807439,  78.68921951, 145.74257173,\n",
       "        74.5235936 , 218.24553537, 138.19297339, 129.64358941,\n",
       "       139.66352483, 275.01623397, 225.12515944,  90.95076246,\n",
       "       179.05092066, 172.22820537, 166.65315395, 110.26473959,\n",
       "       130.22909869,  90.54430874, 106.88804757, 134.95567594,\n",
       "       222.18263102, 106.04156599, 143.29648528, 147.03195508,\n",
       "       116.37520548, 429.14785525, 144.10921427, 140.67050906,\n",
       "       191.2243786 , 552.93793883, 148.99563146, 165.30650939,\n",
       "       126.81369438, 164.12509865, 294.81832755, 123.76619116,\n",
       "       240.24658787, 144.54008564, 102.13022779, 125.68324885,\n",
       "       117.88869165, 112.02526198, 143.75120329, 420.34991199,\n",
       "       200.95992144, 248.05276839, 161.49042941, 144.21318832,\n",
       "       119.63303603, 147.00699367,  47.59312805, 282.94835845,\n",
       "       134.18015783, 125.31443924, 146.8484965 , 173.7342611 ,\n",
       "       196.49038834, 282.88726226, 156.74506034, 160.42507417,\n",
       "       270.76315033, 100.82498451, 488.88494145, 144.5937686 ,\n",
       "       117.66548202, 235.14905066, 348.04575284, 206.22541358,\n",
       "       248.06182872, 224.21884017, 141.87210889,  84.56346977,\n",
       "       111.36733099, 222.79156907, 178.79353303, 114.36098968,\n",
       "       302.57041013, 153.20687432,  79.59677301,  38.89872862,\n",
       "       196.97786414, 135.91952257, 178.15343521, 152.17990966,\n",
       "       141.14245555, 215.36975971, 123.38849915, 117.98105974,\n",
       "       187.06464035, 255.98983285, 147.70609248, 437.08630447,\n",
       "       377.00558382,  90.84526918,  67.77306459, 562.83614767])"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pred = reg.predict(X_test)\n",
    "pred"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Mean squared error: 21320.77\n",
      "R2 : 0.37\n"
     ]
    }
   ],
   "source": [
    "#afficher l'erreur des moindres carrées sur l'ensemble d'entrainement ainsi que le R2\n",
    "#que remarquez vous ? "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.3516173774799306"
      ]
     },
     "execution_count": 34,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "reg.score(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 69,
   "metadata": {},
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "fit() missing 1 required positional argument: 'y'",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-69-4f6a15b71412>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;32mfrom\u001b[0m \u001b[0msklearn\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtree\u001b[0m \u001b[0;32mimport\u001b[0m \u001b[0mDecisionTreeRegressor\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      2\u001b[0m \u001b[0mregr\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mDecisionTreeRegressor\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m \u001b[0mregr\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mfit\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mX_train\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0my_train\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m: fit() missing 1 required positional argument: 'y'"
     ]
    }
   ],
   "source": [
    "from sklearn.tree import DecisionTreeRegressor\n",
    "regr = DecisionTreeRegressor\n",
    "regr.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}