.. _velibtrajectoriesrst:

===============================================
2A.ml - Déterminer la vitesse moyenne des vélib
===============================================


.. only:: html

    **Links:** :download:`notebook <velib_trajectories.ipynb>`, :downloadlink:`html <velib_trajectories2html.html>`, :download:`PDF <velib_trajectories.pdf>`, :download:`python <velib_trajectories.py>`, :downloadlink:`slides <velib_trajectories.slides.html>`, :githublink:`GitHub|_doc/notebooks/challenges/velib/velib_trajectories.ipynb|*`


Ce notebook explicite une solution pour calculer la vitesse moyenne des
velib sachant qu’on ne connaît que l’état des stations à intervalle
réguliers.

.. code:: ipython3

    %matplotlib inline

Même si je propose quelques jeux de données, il est possible de créer le
sien en s’inspirant du code suivant : `Récupérer les données Velib et
les
visualiser <http://www.xavierdupre.fr/app/manydataapi/helpsphinx/notebooks/api_velib_jcdecaux.html#apivelibjcdecauxrst>`__.
Le premier lien décrit les données plus en détail : elles sont
constituées des décomptes des vélos et places disponibles pour chaque
station d’une même ville et pour chaque minute. La méthode proposée ici
est celle des **appariements** décrites par le premier lien.
L’algorithme peut être décrit en deux étapes :

1. Construction d’une base d’événements : vélos reposés et retirés.
2. Appariement des vélos retirés avec les vélos reposés.

La première partie ne pose pas de difficulté particulière. On peut juste
penser à retirer les premiers vélos réposés qui ne pourront pas être
appariés de toutes façons. On fait de même pour les derniers vélos
retirés.

La seconde partie est constitutée de deux éléments :


-  un coût d’appariement
-  la minimisation de l’appariement

Pour le coût, on peut y mettre à peu près toutes les contraintes
imaginables (vitesse trop grande, durée trop grande ou trop petite).
Pour la seconde, le code optimise de façon très naïve : on part d’un
appariement aléatoire, on tire deux appariements
:math:`a_1 \rightarrow b_1` et :math:`a_2 \rightarrow b_2`. On inverse :
:math:`a_1 \rightarrow b_2` et :math:`a_2 \rightarrow b_1`. Si
l’appariement est moins coûteux, on garde. Le code complet est
disponible dans le module *velib_trajectories*.

.. code:: ipython3

    from pyquickhelper.loghelper import str2datetime
    from ensae_projects.datainc import besancon_df
    import pandas
    
    jeu = besancon_df()
    df = pandas.read_csv(jeu, sep="\t", encoding="utf8")
    df ["collect_date"] = df.apply( lambda r: str2datetime(r["collect_date"]),axis=1)
    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>address</th>
          <th>available_bike_stands</th>
          <th>available_bikes</th>
          <th>banking</th>
          <th>bike_stands</th>
          <th>bonus</th>
          <th>collect_date</th>
          <th>contract_name</th>
          <th>file</th>
          <th>last_update</th>
          <th>lat</th>
          <th>lng</th>
          <th>name</th>
          <th>number</th>
          <th>status</th>
        </tr>
      </thead>
      <tbody>
        <tr>
          <th>0</th>
          <td>Place Jacobins</td>
          <td>2</td>
          <td>8</td>
          <td>0</td>
          <td>10</td>
          <td>0</td>
          <td>2014-05-22 02:40:03.954829</td>
          <td>Besancon</td>
          <td>velib_data.2014-05-22_02-40-02.681782.txt</td>
          <td>2014-05-22 00:09:26</td>
          <td>47.234142</td>
          <td>6.033465</td>
          <td>14 - JACOBINS</td>
          <td>14</td>
          <td>OPEN</td>
        </tr>
        <tr>
          <th>1</th>
          <td>Rue Charles NODIER Sortie du tunnel, direction...</td>
          <td>1</td>
          <td>11</td>
          <td>0</td>
          <td>12</td>
          <td>0</td>
          <td>2014-05-22 02:40:03.954829</td>
          <td>Besancon</td>
          <td>velib_data.2014-05-22_02-40-02.681782.txt</td>
          <td>2014-05-22 02:38:48</td>
          <td>47.232098</td>
          <td>6.035136</td>
          <td>23 - RIVOTTE</td>
          <td>23</td>
          <td>OPEN</td>
        </tr>
        <tr>
          <th>2</th>
          <td>Avenue de la Gare d'Eau face Conseil Général</td>
          <td>4</td>
          <td>6</td>
          <td>0</td>
          <td>10</td>
          <td>0</td>
          <td>2014-05-22 02:40:03.954829</td>
          <td>Besancon</td>
          <td>velib_data.2014-05-22_02-40-02.681782.txt</td>
          <td>2014-05-21 13:24:02</td>
          <td>47.232140</td>
          <td>6.022843</td>
          <td>16 - GARE D'EAU</td>
          <td>16</td>
          <td>OPEN</td>
        </tr>
        <tr>
          <th>3</th>
          <td>18 rue Xavier Marmier</td>
          <td>6</td>
          <td>6</td>
          <td>0</td>
          <td>12</td>
          <td>0</td>
          <td>2014-05-22 02:40:03.954829</td>
          <td>Besancon</td>
          <td>velib_data.2014-05-22_02-40-02.681782.txt</td>
          <td>2014-05-22 00:00:44</td>
          <td>47.241256</td>
          <td>6.006690</td>
          <td>30 - XAVIER MARMIER</td>
          <td>30</td>
          <td>OPEN</td>
        </tr>
        <tr>
          <th>4</th>
          <td>Rue Battant (au niveau du Square Bouchot)</td>
          <td>1</td>
          <td>9</td>
          <td>0</td>
          <td>10</td>
          <td>0</td>
          <td>2014-05-22 02:40:03.954829</td>
          <td>Besancon</td>
          <td>velib_data.2014-05-22_02-40-02.681782.txt</td>
          <td>2014-05-22 02:31:01</td>
          <td>47.242990</td>
          <td>6.022070</td>
          <td>02 - SQUARE BOUCHOT</td>
          <td>2</td>
          <td>OPEN</td>
        </tr>
      </tbody>
    </table>
    </div>



**On vérifie les données pour la première date :**

.. code:: ipython3

    from manydataapi.velib import DataCollectJCDecaux as DataVelibCollect
    dt = df["file"][0]
    subset = df [ df["file"] == dt ]
    fig,ax,plt = DataVelibCollect.draw(subset, figsize=(16,6))
    ax.set_title("Besançon - {0} - {1} stations".format(dt.replace("besancon.","") \
                                                          .replace(".txt","").replace("_", " "),
                                                        len(subset)));



.. image:: velib_trajectories_5_0.png


**On calcule les événements (1 vélo apparu, 1 vélo disparu) :**

.. code:: ipython3

    import ensae_projects.challenge.velib_trajectories as velib
    events = list(sorted(velib.enumerate_events(df)))
    print("nb events",len(events))
    events [:2]


.. parsed-literal::
    nb events 589




.. parsed-literal::

    [('velib_data.2014-05-22_04-17-04.635259.txt',
      Timestamp('2014-05-22 04:17:04.857260'),
      '24 - MEDIATHEQUE',
      47.24010992877723,
      6.027358757731063,
      1),
     ('velib_data.2014-05-22_05-29-04.949889.txt',
      Timestamp('2014-05-22 05:29:05.262905'),
      '27 - MAIRIE',
      47.23687885054159,
      6.0222303583977155,
      1)]



**On calcule le meilleur appariement (ça prend un peu de temps).** On ne
prend que 200 itérations mais il en faudrait plus.

.. code:: ipython3

    params = velib.ParemetreCoutTrajet()
    print(params)
    mindist, moyenne, appariement, positif, negatif = velib.appariement(events, iter=50, params=params)
    print("vitesse moyenne", moyenne)


.. parsed-literal::
    {'max_speed': 50, 'high_speed': 25, 'low_speed': 5, 'high_time': 0.75, 'low_time': 0.1}
    iteration  0 : app- 144 / 301 min 21743314023.08393 vitesse  5.185466030632195  nbchange 0
    iteration  10 : app- 32 / 301 min 4217360521.9227724 vitesse  4.644576073119561  nbchange 272
    iteration  20 : app- 26 / 301 min 3217656643.672325 vitesse  5.27657037459076  nbchange 104
    iteration  30 : app- 19 / 301 min 2250635665.002002 vitesse  5.872166710865633  nbchange 88
    iteration  40 : app- 16 / 301 min 1848259504.5709043 vitesse  6.277284851644164  nbchange 56
    vitesse moyenne 6.407462005952014


La vitesse évolue encore et il faudrait faire tourner l’algorithme plus
longtemps et le recoder de façon plus efficace. Voyons déjà ce que cela
donne en terme de déplacement.

.. code:: ipython3

    import matplotlib.pyplot as plt
    app = [(positif[a], negatif[b] ) for a,b in appariement]
    plt.figure(figsize=(16,6))
    for deb,fin in app :
        x = [deb[3], fin[3]]
        y = [deb[4], fin[4]]
        if x[0] > 0 and y[0] and x[1] > 0 and y[1] > 0 :
            # on enlève les trajets aberrants
            plt.plot(x,y,color="b")



.. image:: velib_trajectories_11_0.png


Difficile de dire sur ce simple dessin si cela est sensé. Il faudrait
faire tourner l’algorithme plus longtemps et le recoder de façon plus
efficace.