.. _td2acorrectionsession3Arst:

================================================================
2A.ml - Statistiques descriptives avec scikit-learn - correction
================================================================


.. only:: html

    **Links:** :download:`notebook <td2a_correction_session_3A.ipynb>`, :downloadlink:`html <td2a_correction_session_3A2html.html>`, :download:`python <td2a_correction_session_3A.py>`, :downloadlink:`slides <td2a_correction_session_3A.slides.html>`, :githublink:`GitHub|_doc/notebooks/td2a_ml/td2a_correction_session_3A.ipynb|*`


ACP, CAH, régression linéaire, correction.

.. code:: ipython3

    %matplotlib inline

.. code:: ipython3

    import matplotlib.pyplot as plt

.. code:: ipython3

    from jyquickhelper import add_notebook_menu
    add_notebook_menu()






.. contents::
    :local:





.. code:: ipython3

    # Répare une incompatibilité entre scipy 1.0 et statsmodels 0.8.
    from pymyinstall.fix import fix_scipy10_for_statsmodels08
    fix_scipy10_for_statsmodels08()

Prérequis de l’énoncé
~~~~~~~~~~~~~~~~~~~~~

.. code:: ipython3

    import pyensae.datasource
    f = pyensae.datasource.download_data("dads2011_gf_salaries11_dbase.zip",
                              website="https://www.insee.fr/fr/statistiques/fichier/2011542/")
    
    import pandas
    try:
        from dbfread import DBF
        use_dbfread = True
    except ImportError as e :
        use_dbfread = False
        
    if use_dbfread:
        import os
        from pyensae.sql.database_exception import ExceptionSQL
        from pyensae.datasource import dBase2sqllite
        print("convert dbase into sqllite")
        try:
            dBase2sqllite("salaries2011.db3", "varlist_salaries11.dbf", overwrite_table="varlist")
            dBase2sqllite("salaries2011.db3", "varmod_salaries11.dbf", overwrite_table="varmod")
            dBase2sqllite("salaries2011.db3", 'salaries11.dbf', overwrite_table="salaries", fLOG = print)
        except ExceptionSQL:
            print("La base de données est déjà renseignée.")
    else :
        print("use of zipped version")
        import pyensae.datasource
        db3 = pyensae.datasource.download_data("salaries2011.zip")
        # pour aller plus vite, données à télécharger au
        # http://www.xavierdupre.fr/enseignement/complements/salaries2011.zip


.. parsed-literal::
    convert dbase into sqllite
    unable to execute a SQL request (1)(file salaries2011.db3)
    create table varlist (VARIABLE TEXT, LIBELLE TEXT, TYPE TEXT, LONGUEUR REAL)
    La base de données est déjà renseignée.


Exercice 1 : CAH (classification ascendante hiérarchique)
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

Le point commun de ces méthodes est qu’elles ne sont pas supervisées.
L’objectif est de réduire la complexité des données. Réduire le nombre
de dimensions pour l’ACP ou segmenter les observations pour les k-means
et la CAH.

.. code:: ipython3

    import pandas, numpy, matplotlib.pyplot as plt
    from sklearn.decomposition import PCA
    import pyensae.datasource
    
    pyensae.datasource.download_data("eleve_region.txt")
    df = pandas.read_csv("eleve_region.txt", sep="\t", encoding="utf8", index_col=0)
    print(df.shape)
    df.head(n=5)


.. parsed-literal::
    (27, 21)






.. raw:: html

    <div>
    <style scoped>
        .dataframe tbody tr th:only-of-type {
            vertical-align: middle;
        }

        .dataframe tbody tr th {
            vertical-align: top;
        }

        .dataframe thead th {
            text-align: right;
        }
    </style>
    <table border="1" class="dataframe">
      <thead>
        <tr style="text-align: right;">
          <th></th>
          <th>1993</th>
          <th>1994</th>
          <th>1995</th>
          <th>1996</th>
          <th>1997</th>
          <th>1998</th>
          <th>1999</th>
          <th>2000</th>
          <th>2001</th>
          <th>2002</th>
          <th>...</th>
          <th>2004</th>
          <th>2005</th>
          <th>2006</th>
          <th>2007</th>
          <th>2008</th>
          <th>2009</th>
          <th>2010</th>
          <th>2011</th>
          <th>2012</th>
          <th>2013</th>
        </tr>
        <tr>
          <th>académie</th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
        </tr>
      </thead>
      <tbody>
        <tr>
          <th>Aix-Marseille</th>
          <td>241357</td>
          <td>242298</td>
          <td>242096</td>
          <td>242295</td>
          <td>243660</td>
          <td>244608</td>
          <td>245536</td>
          <td>247288</td>
          <td>249331</td>
          <td>250871</td>
          <td>...</td>
          <td>250622</td>
          <td>248208</td>
          <td>245755</td>
          <td>243832</td>
          <td>242309</td>
          <td>240664</td>
          <td>240432</td>
          <td>241336</td>
          <td>239051</td>
          <td>240115</td>
        </tr>
        <tr>
          <th>Amiens</th>
          <td>198281</td>
          <td>196871</td>
          <td>195709</td>
          <td>194055</td>
          <td>192893</td>
          <td>191862</td>
          <td>189636</td>
          <td>185977</td>
          <td>183357</td>
          <td>180973</td>
          <td>...</td>
          <td>175610</td>
          <td>172110</td>
          <td>168718</td>
          <td>165295</td>
          <td>163116</td>
          <td>162548</td>
          <td>163270</td>
          <td>164422</td>
          <td>165275</td>
          <td>166345</td>
        </tr>
        <tr>
          <th>Besançon</th>
          <td>116373</td>
          <td>115600</td>
          <td>114282</td>
          <td>113312</td>
          <td>112076</td>
          <td>110261</td>
          <td>108106</td>
          <td>105463</td>
          <td>103336</td>
          <td>102264</td>
          <td>...</td>
          <td>100117</td>
          <td>98611</td>
          <td>97038</td>
          <td>95779</td>
          <td>95074</td>
          <td>94501</td>
          <td>94599</td>
          <td>94745</td>
          <td>94351</td>
          <td>94613</td>
        </tr>
        <tr>
          <th>Bordeaux</th>
          <td>253551</td>
          <td>252644</td>
          <td>249658</td>
          <td>247708</td>
          <td>247499</td>
          <td>245757</td>
          <td>244992</td>
          <td>243047</td>
          <td>243592</td>
          <td>245198</td>
          <td>...</td>
          <td>244805</td>
          <td>244343</td>
          <td>242602</td>
          <td>242933</td>
          <td>243146</td>
          <td>244336</td>
          <td>246806</td>
          <td>250626</td>
          <td>252085</td>
          <td>255761</td>
        </tr>
        <tr>
          <th>Caen</th>
          <td>145435</td>
          <td>144369</td>
          <td>141883</td>
          <td>140658</td>
          <td>139585</td>
          <td>137704</td>
          <td>135613</td>
          <td>133255</td>
          <td>131206</td>
          <td>129271</td>
          <td>...</td>
          <td>125552</td>
          <td>123889</td>
          <td>122550</td>
          <td>121002</td>
          <td>119857</td>
          <td>119426</td>
          <td>119184</td>
          <td>119764</td>
          <td>119010</td>
          <td>119238</td>
        </tr>
      </tbody>
    </table>
    <p>5 rows × 21 columns</p>
    </div>



.. code:: ipython3

    for c in df.columns:
        if c != "1993":
            df[c] /= df ["1993"]
    df["1993"] /= df["1993"]
    df.head()






.. raw:: html

    <div>
    <style scoped>
        .dataframe tbody tr th:only-of-type {
            vertical-align: middle;
        }

        .dataframe tbody tr th {
            vertical-align: top;
        }

        .dataframe thead th {
            text-align: right;
        }
    </style>
    <table border="1" class="dataframe">
      <thead>
        <tr style="text-align: right;">
          <th></th>
          <th>1993</th>
          <th>1994</th>
          <th>1995</th>
          <th>1996</th>
          <th>1997</th>
          <th>1998</th>
          <th>1999</th>
          <th>2000</th>
          <th>2001</th>
          <th>2002</th>
          <th>...</th>
          <th>2004</th>
          <th>2005</th>
          <th>2006</th>
          <th>2007</th>
          <th>2008</th>
          <th>2009</th>
          <th>2010</th>
          <th>2011</th>
          <th>2012</th>
          <th>2013</th>
        </tr>
        <tr>
          <th>académie</th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
          <th></th>
        </tr>
      </thead>
      <tbody>
        <tr>
          <th>Aix-Marseille</th>
          <td>1.0</td>
          <td>1.003899</td>
          <td>1.003062</td>
          <td>1.003886</td>
          <td>1.009542</td>
          <td>1.013470</td>
          <td>1.017315</td>
          <td>1.024574</td>
          <td>1.033038</td>
          <td>1.039419</td>
          <td>...</td>
          <td>1.038387</td>
          <td>1.028385</td>
          <td>1.018222</td>
          <td>1.010255</td>
          <td>1.003944</td>
          <td>0.997129</td>
          <td>0.996168</td>
          <td>0.999913</td>
          <td>0.990446</td>
          <td>0.994854</td>
        </tr>
        <tr>
          <th>Amiens</th>
          <td>1.0</td>
          <td>0.992889</td>
          <td>0.987029</td>
          <td>0.978687</td>
          <td>0.972826</td>
          <td>0.967627</td>
          <td>0.956400</td>
          <td>0.937947</td>
          <td>0.924733</td>
          <td>0.912710</td>
          <td>...</td>
          <td>0.885662</td>
          <td>0.868011</td>
          <td>0.850904</td>
          <td>0.833640</td>
          <td>0.822651</td>
          <td>0.819786</td>
          <td>0.823427</td>
          <td>0.829237</td>
          <td>0.833539</td>
          <td>0.838936</td>
        </tr>
        <tr>
          <th>Besançon</th>
          <td>1.0</td>
          <td>0.993358</td>
          <td>0.982032</td>
          <td>0.973697</td>
          <td>0.963076</td>
          <td>0.947479</td>
          <td>0.928961</td>
          <td>0.906250</td>
          <td>0.887972</td>
          <td>0.878761</td>
          <td>...</td>
          <td>0.860311</td>
          <td>0.847370</td>
          <td>0.833853</td>
          <td>0.823035</td>
          <td>0.816976</td>
          <td>0.812053</td>
          <td>0.812895</td>
          <td>0.814149</td>
          <td>0.810764</td>
          <td>0.813015</td>
        </tr>
        <tr>
          <th>Bordeaux</th>
          <td>1.0</td>
          <td>0.996423</td>
          <td>0.984646</td>
          <td>0.976955</td>
          <td>0.976131</td>
          <td>0.969261</td>
          <td>0.966243</td>
          <td>0.958572</td>
          <td>0.960722</td>
          <td>0.967056</td>
          <td>...</td>
          <td>0.965506</td>
          <td>0.963684</td>
          <td>0.956817</td>
          <td>0.958123</td>
          <td>0.958963</td>
          <td>0.963656</td>
          <td>0.973398</td>
          <td>0.988464</td>
          <td>0.994218</td>
          <td>1.008716</td>
        </tr>
        <tr>
          <th>Caen</th>
          <td>1.0</td>
          <td>0.992670</td>
          <td>0.975577</td>
          <td>0.967154</td>
          <td>0.959776</td>
          <td>0.946842</td>
          <td>0.932465</td>
          <td>0.916251</td>
          <td>0.902162</td>
          <td>0.888858</td>
          <td>...</td>
          <td>0.863286</td>
          <td>0.851851</td>
          <td>0.842644</td>
          <td>0.832001</td>
          <td>0.824128</td>
          <td>0.821164</td>
          <td>0.819500</td>
          <td>0.823488</td>
          <td>0.818304</td>
          <td>0.819871</td>
        </tr>
      </tbody>
    </table>
    <p>5 rows × 21 columns</p>
    </div>



.. code:: ipython3

    from sklearn.cluster import AgglomerativeClustering
    
    ward = AgglomerativeClustering(linkage='ward', compute_full_tree=True).fit(df)
    ward




.. parsed-literal::
    AgglomerativeClustering(compute_full_tree=True)



.. code:: ipython3

    from scipy.cluster.hierarchy import dendrogram
    import matplotlib.pyplot as plt
    
    dendro = [ ]
    for a,b in ward.children_:
        dendro.append([a, b, float(len(dendro)+1), len(dendro)+1])
        # le dernier coefficient devrait contenir le nombre de feuilles dépendant de ce noeud
        # et non le dernier indice
        # de même, le niveau (3ème colonne) ne devrait pas être le nombre de noeud
        # mais la distance de Ward
        
    fig = plt.figure( figsize=(8,8) )
    ax = fig.add_subplot(1,1,1) 
    r = dendrogram(dendro, color_threshold=1, labels=list(df.index),
                   show_leaf_counts=True, ax=ax, orientation="left")



.. image:: td2a_correction_session_3A_11_0.png


Je reprends également le graphique montrant la matrice de corrélations
qu’on peut également obtenir avec `seaborn :
clustermap <https://seaborn.pydata.org/examples/structured_heatmap.html>`__.

.. code:: ipython3

    from scipy.spatial.distance import pdist, squareform
    
    data_dist = pdist(df)
    
    fig = plt.figure(figsize=(8,8))
    
    # x ywidth height
    ax1 = fig.add_axes([0.05,0.1,0.2,0.6])
    Z1 = dendrogram(dendro, orientation='right',labels=list(df.index))
    ax1.set_xticks([])
    
    # Compute and plot second dendrogram.
    ax2 = fig.add_axes([0.3,0.71,0.6,0.2])
    Z2 = dendrogram(dendro)
    ax2.set_xticks([])
    ax2.set_yticks([])
    
    # Compute and plot the heatmap
    axmatrix = fig.add_axes([0.3,0.1,0.6,0.6])
    idx1 = Z1['leaves']
    idx2 = Z2['leaves']
    
    D = squareform(data_dist)
    D = D[idx1,:]
    D = D[:,idx2]
    im = axmatrix.matshow(D, aspect='auto', origin='lower', cmap=plt.cm.YlGnBu)
    axmatrix.set_xticks([])
    axmatrix.set_yticks([])
    
    # Plot colorbar.
    axcolor = fig.add_axes([0.91,0.1,0.02,0.6])
    plt.colorbar(im, cax=axcolor)
    plt.title("Matrice de corrélation et dendogramme.");



.. image:: td2a_correction_session_3A_13_0.png


Exercice 2 : régression linéaire
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

Ce sont trois méthodes supervisées : on s’en sert pour expliquer prédire
le lien entre deux variables :math:`X` et :math:`Y` (ou ensemble de
variables) ou prédire :math:`Y` en fonction de :math:`X`. On suppose que
les données ont déjà été téléchargées (voir l’énoncé de cet exercice).
On continue avec l’extraction de l’âge, du sexe et du salaire comme
indiqué dans l’énoncé.

.. code:: ipython3

    import sqlite3, pandas
    con = sqlite3.connect("salaries2011.db3")
    df = pandas.io.sql.read_sql("select * from varmod", con)
    con.close()
    
    values = df[ df.VARIABLE == "TRNNETO"].copy()
    
    def process_intervalle(s):
        # [14 000 ; 16 000[ euros
        acc = "0123456789;+"
        s0 = "".join(c for c in s if c in acc)
        spl = s0.split(';')
        if len(spl) != 2:
            raise ValueError("Unable to process '{0}'".format(s0))
        try:
            a = float(spl[0])
        except Exception as e:
            raise ValueError("Cannot interpret '{0}' - {1}".format(s, spl))
        b = float(spl[1]) if "+" not in spl[1] else None
        if b is None:
            return a
        return (a+b) / 2.0
    
    values["montant"] = values.apply(lambda r : process_intervalle(r["MODLIBELLE"]), axis=1)
    values.head()






.. raw:: html

    <div>
    <style scoped>
        .dataframe tbody tr th:only-of-type {
            vertical-align: middle;
        }

        .dataframe tbody tr th {
            vertical-align: top;
        }

        .dataframe thead th {
            text-align: right;
        }
    </style>
    <table border="1" class="dataframe">
      <thead>
        <tr style="text-align: right;">
          <th></th>
          <th>VARIABLE</th>
          <th>MODALITE</th>
          <th>MODLIBELLE</th>
          <th>montant</th>
        </tr>
      </thead>
      <tbody>
        <tr>
          <th>8957</th>
          <td>TRNNETO</td>
          <td>00</td>
          <td>[0 ; 200[ euros</td>
          <td>100.0</td>
        </tr>
        <tr>
          <th>8958</th>
          <td>TRNNETO</td>
          <td>01</td>
          <td>[200 ; 500[ euros</td>
          <td>350.0</td>
        </tr>
        <tr>
          <th>8959</th>
          <td>TRNNETO</td>
          <td>02</td>
          <td>[500 ; 1 000[ euros</td>
          <td>750.0</td>
        </tr>
        <tr>
          <th>8960</th>
          <td>TRNNETO</td>
          <td>03</td>
          <td>[1 000 ; 1 500[ euros</td>
          <td>1250.0</td>
        </tr>
        <tr>
          <th>8961</th>
          <td>TRNNETO</td>
          <td>04</td>
          <td>[1 500 ; 2 000[ euros</td>
          <td>1750.0</td>
        </tr>
      </tbody>
    </table>
    </div>



.. code:: ipython3

    import sqlite3, pandas
    con = sqlite3.connect("salaries2011.db3")
    data = pandas.io.sql.read_sql("select TRNNETO,AGE,SEXE from salaries", con)
    con.close()
    
    salaires = data.merge (values, left_on="TRNNETO", right_on="MODALITE" )
    salaires["M"] = salaires.apply(lambda r: 1 if r["SEXE"] == "1" else 0, axis=1)
    salaires["F"] = salaires.apply(lambda r: 1 if r["SEXE"] == "2" else 0, axis=1)
    data = salaires[["AGE","M","F","montant"]]
    data = data [data.M + data.F > 0]
    data.head()






.. raw:: html

    <div>
    <style scoped>
        .dataframe tbody tr th:only-of-type {
            vertical-align: middle;
        }

        .dataframe tbody tr th {
            vertical-align: top;
        }

        .dataframe thead th {
            text-align: right;
        }
    </style>
    <table border="1" class="dataframe">
      <thead>
        <tr style="text-align: right;">
          <th></th>
          <th>AGE</th>
          <th>M</th>
          <th>F</th>
          <th>montant</th>
        </tr>
      </thead>
      <tbody>
        <tr>
          <th>0</th>
          <td>49.0</td>
          <td>1</td>
          <td>0</td>
          <td>750.0</td>
        </tr>
        <tr>
          <th>1</th>
          <td>27.0</td>
          <td>1</td>
          <td>0</td>
          <td>750.0</td>
        </tr>
        <tr>
          <th>2</th>
          <td>22.0</td>
          <td>1</td>
          <td>0</td>
          <td>750.0</td>
        </tr>
        <tr>
          <th>3</th>
          <td>26.0</td>
          <td>1</td>
          <td>0</td>
          <td>750.0</td>
        </tr>
        <tr>
          <th>4</th>
          <td>29.0</td>
          <td>0</td>
          <td>1</td>
          <td>750.0</td>
        </tr>
      </tbody>
    </table>
    </div>



On supprime les valeurs manquantes :

.. code:: ipython3

    nonull = data.dropna().copy()
    nonull.shape




.. parsed-literal::
    (2240691, 4)



version scikit-learn
^^^^^^^^^^^^^^^^^^^^

La `régression
linéraire <http://scikit-learn.org/stable/modules/linear_model.html>`__
suit :

.. code:: ipython3

    nonull[["AGE","M"]].dropna().shape




.. parsed-literal::
    (2240691, 2)



.. code:: ipython3

    from sklearn import linear_model
    clf = linear_model.LinearRegression()
    clf.fit (nonull[["AGE","M"]].values, nonull.montant.values)
    clf.coef_, clf.intercept_, "R^2=", clf.score(
        nonull[["AGE","M"]], nonull.montant)


.. parsed-literal::
    C:\Python395_x64\lib\site-packages\sklearn\base.py:443: UserWarning: X has feature names, but LinearRegression was fitted without feature names
      warnings.warn(




.. parsed-literal::

    (array([ 310.98873096, 4710.02901965]),
     4927.219006697269,
     'R^2=',
     0.13957345814666222)



On prend un échantillon aléatoire :

.. code:: ipython3

    import random
    val = nonull.copy()
    val["rnd"] = val.apply(lambda r: random.randint(0, 1000), axis=1)
    ech = val[val["rnd"] == 1]
    ech.shape




.. parsed-literal::
    (2278, 5)



On sépare homme et femmes :

.. code:: ipython3

    homme = ech[ech.M == 1]
    femme = ech[ech.M == 0]
    
    predh = clf.predict(homme[["AGE","M"]])
    predf = clf.predict(femme[["AGE","M"]])
    
    import matplotlib.pyplot as plt
    plt.figure(figsize=(16,6))
    plt.plot(homme.AGE, homme.montant, "r.")
    plt.plot(femme.AGE + 0.2, femme.montant, "b.")
    plt.plot(homme.AGE, predh, "ro-", label="hommes")
    plt.plot(femme.AGE, predf, "bo-", label="femmes")
    plt.legend()
    plt.title("Montant moyen par âge et genre");


.. parsed-literal::
    C:\Python395_x64\lib\site-packages\sklearn\base.py:443: UserWarning: X has feature names, but LinearRegression was fitted without feature names
      warnings.warn(
    C:\Python395_x64\lib\site-packages\sklearn\base.py:443: UserWarning: X has feature names, but LinearRegression was fitted without feature names
      warnings.warn(



.. image:: td2a_correction_session_3A_25_1.png


version statsmodels
^^^^^^^^^^^^^^^^^^^

.. code:: ipython3

    import statsmodels.api as sm
    nonull["one"] = 1.0  # on ajoute la constante
    model = sm.OLS(nonull.montant, nonull [["AGE","M", "one"]])
    results = model.fit()
    print("coefficients",results.params)
    results.summary()


.. parsed-literal::
    C:\Python395_x64\lib\site-packages\statsmodels\tsa\base\tsa_model.py:7: FutureWarning: pandas.Int64Index is deprecated and will be removed from pandas in a future version. Use pandas.Index with the appropriate dtype instead.
      from pandas import (to_datetime, Int64Index, DatetimeIndex, Period,
    C:\Python395_x64\lib\site-packages\statsmodels\tsa\base\tsa_model.py:7: FutureWarning: pandas.Float64Index is deprecated and will be removed from pandas in a future version. Use pandas.Index with the appropriate dtype instead.
      from pandas import (to_datetime, Int64Index, DatetimeIndex, Period,


.. parsed-literal::

    coefficients AGE     310.988731
    M      4710.029020
    one    4927.219007
    dtype: float64






.. raw:: html

    <table class="simpletable">
    <caption>OLS Regression Results</caption>
    <tr>
      <th>Dep. Variable:</th>         <td>montant</td>     <th>  R-squared:         </th>  <td>   0.140</td>  
    </tr>
    <tr>
      <th>Model:</th>                   <td>OLS</td>       <th>  Adj. R-squared:    </th>  <td>   0.140</td>  
    </tr>
    <tr>
      <th>Method:</th>             <td>Least Squares</td>  <th>  F-statistic:       </th>  <td>1.817e+05</td> 
    </tr>
    <tr>
      <th>Date:</th>             <td>Sat, 12 Feb 2022</td> <th>  Prob (F-statistic):</th>   <td>  0.00</td>   
    </tr>
    <tr>
      <th>Time:</th>                 <td>18:57:08</td>     <th>  Log-Likelihood:    </th> <td>-2.4023e+07</td>
    </tr>
    <tr>
      <th>No. Observations:</th>      <td>2240691</td>     <th>  AIC:               </th>  <td>4.805e+07</td> 
    </tr>
    <tr>
      <th>Df Residuals:</th>          <td>2240688</td>     <th>  BIC:               </th>  <td>4.805e+07</td> 
    </tr>
    <tr>
      <th>Df Model:</th>              <td>     2</td>      <th>                     </th>      <td> </td>     
    </tr>
    <tr>
      <th>Covariance Type:</th>      <td>nonrobust</td>    <th>                     </th>      <td> </td>     
    </tr>
    </table>
    <table class="simpletable">
    <tr>
       <td></td>      <th>coef</th>     <th>std err</th>      <th>t</th>      <th>P>|t|</th>  <th>[0.025</th>    <th>0.975]</th>  
    </tr>
    <tr>
      <th>AGE</th> <td>  310.9887</td> <td>    0.600</td> <td>  518.682</td> <td> 0.000</td> <td>  309.814</td> <td>  312.164</td>
    </tr>
    <tr>
      <th>M</th>   <td> 4710.0290</td> <td>   14.658</td> <td>  321.323</td> <td> 0.000</td> <td> 4681.299</td> <td> 4738.759</td>
    </tr>
    <tr>
      <th>one</th> <td> 4927.2190</td> <td>   26.032</td> <td>  189.272</td> <td> 0.000</td> <td> 4876.196</td> <td> 4978.242</td>
    </tr>
    </table>
    <table class="simpletable">
    <tr>
      <th>Omnibus:</th>       <td>132298.944</td> <th>  Durbin-Watson:     </th>  <td>   0.183</td> 
    </tr>
    <tr>
      <th>Prob(Omnibus):</th>   <td> 0.000</td>   <th>  Jarque-Bera (JB):  </th> <td>158649.997</td>
    </tr>
    <tr>
      <th>Skew:</th>            <td> 0.614</td>   <th>  Prob(JB):          </th>  <td>    0.00</td> 
    </tr>
    <tr>
      <th>Kurtosis:</th>        <td> 3.440</td>   <th>  Cond. No.          </th>  <td>    150.</td> 
    </tr>
    </table><br/><br/>Notes:<br/>[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.



On reproduit le même dessin :

.. code:: ipython3

    import random
    val = nonull.copy()
    val["rnd"] = val.apply(lambda r: random.randint(0,1000), axis=1)
    ech = val[val["rnd"] == 1]

.. code:: ipython3

    homme = ech [ ech.M == 1]
    femme = ech [ ech.M == 0]
    
    predh = results.predict(homme[["AGE","M","one"]])
    predf = results.predict(femme[["AGE","M","one"]])
    
    import matplotlib.pyplot as plt
    def graph(homme, femme, predh, predf):
        fig, ax = plt.subplots(1, 1, figsize=(16, 6))
        ax.plot(homme.AGE, homme.montant, "r.")
        ax.plot(femme.AGE + 0.2, femme.montant, "b.")
        ax.plot(homme.AGE, predh, "ro-", label="hommes")
        ax.plot(femme.AGE, predf, "bo-", label="femmes")
        ax.legend()
        ax.set_title("Montant moyen par âge et genre");
        return ax
    
    graph(homme, femme, predh, predf);



.. image:: td2a_correction_session_3A_30_0.png


On ajoute l’intervalle de confiance sur un échantillon :

.. code:: ipython3

    from statsmodels.sandbox.regression.predstd import wls_prediction_std
    
    prstd, iv_l, iv_u = wls_prediction_std(results)
    
    val = nonull.copy()
    val["rnd"] = val.apply(lambda r: random.randint(0, 1000), axis=1)
    val["pred"] = prstd
    val["up"] = iv_u
    val["down"] = iv_l
    
    ech = val[val["rnd"] == 1]
    ech.head()






.. raw:: html

    <div>
    <style scoped>
        .dataframe tbody tr th:only-of-type {
            vertical-align: middle;
        }

        .dataframe tbody tr th {
            vertical-align: top;
        }

        .dataframe thead th {
            text-align: right;
        }
    </style>
    <table border="1" class="dataframe">
      <thead>
        <tr style="text-align: right;">
          <th></th>
          <th>AGE</th>
          <th>M</th>
          <th>F</th>
          <th>montant</th>
          <th>one</th>
          <th>rnd</th>
          <th>pred</th>
          <th>up</th>
          <th>down</th>
        </tr>
      </thead>
      <tbody>
        <tr>
          <th>86</th>
          <td>34.0</td>
          <td>0</td>
          <td>1</td>
          <td>750.0</td>
          <td>1.0</td>
          <td>1</td>
          <td>10964.546920</td>
          <td>36990.964539</td>
          <td>-5989.292820</td>
        </tr>
        <tr>
          <th>355</th>
          <td>48.0</td>
          <td>0</td>
          <td>1</td>
          <td>750.0</td>
          <td>1.0</td>
          <td>1</td>
          <td>10964.547492</td>
          <td>41344.807893</td>
          <td>-1635.451707</td>
        </tr>
        <tr>
          <th>1046</th>
          <td>68.0</td>
          <td>0</td>
          <td>1</td>
          <td>750.0</td>
          <td>1.0</td>
          <td>1</td>
          <td>10964.559457</td>
          <td>47564.605962</td>
          <td>4584.299462</td>
        </tr>
        <tr>
          <th>2008</th>
          <td>21.0</td>
          <td>0</td>
          <td>1</td>
          <td>750.0</td>
          <td>1.0</td>
          <td>1</td>
          <td>10964.552144</td>
          <td>32948.121273</td>
          <td>-10032.156559</td>
        </tr>
        <tr>
          <th>3761</th>
          <td>18.0</td>
          <td>1</td>
          <td>0</td>
          <td>750.0</td>
          <td>1.0</td>
          <td>1</td>
          <td>10964.553468</td>
          <td>36725.186695</td>
          <td>-6255.096328</td>
        </tr>
      </tbody>
    </table>
    </div>



Puis on l’ajoute au graphe précédent :

.. code:: ipython3

    homme = ech[ech.M == 1]
    femme = ech[ech.M == 0]
    
    predh = results.predict(homme[["AGE","M","one"]])
    predf = results.predict(femme[["AGE","M","one"]])
    
    ax = graph(homme, femme, predh, predf)
    ax.plot(homme.AGE, homme.up, 'r-')
    ax.plot(homme.AGE, homme.down, 'r-')
    ax.plot(femme.AGE, femme.up, 'b-')
    ax.plot(femme.AGE, femme.down, 'b-')
    ax.set_title("Montant moyen par âge et genre avec écart-type");



.. image:: td2a_correction_session_3A_34_0.png