Les vélos sur le pont de Fremont

URL = "https://data.seattle.gov/api/views/65db-xm6k/rows.csv?accessType=DOWNLOAD"

import pandas as pd

# on a déja le fichier en local
local_file = "data/fremont.csv"

# ceci nécessite alors
# %pip install requests

from pathlib import Path

if Path(local_file).exists():
    print(f"le fichier {local_file} est déjà là")
else:
    print(f"allons chercher le fichier {local_file}")

    import requests
    req = requests.get(URL)

    # doit afficher 200
    print(req.status_code)

    # on sauve tel quel dans le fichier local
    with open(local_file, 'w') as writer:
        writer.write(req.text)

!head $local_file

%pip install seaborn

import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns

# version naïve
data = pd.read_csv(local_file); data.shape

data.head()

!grep '01/01/2014 12:00:00 AM' $local_file

data.drop_duplicates(inplace=True)

data.shape

# intéressant aussi, pour voir notamment les points manquants
data.info();

# ou tout simplement
data.dtypes

# la version lente

# avec cette forme, on demanderait à read_csv:
# de mettre la date comme index,
# et de parser les dates
# mais on ne va pas le faire car
# 1. c'est peu sûr
# 2. c'est trop lent

# data = pd.read_csv(local_file, index_col='Date', parse_dates=True); data.head()

# une meilleure idée, pour améliorer ces deux points d'un coup
# est de fournir nous-même le format des dates

data.index = pd.to_datetime(data.Date, format="%m/%d/%Y %I:%M:%S %p")
del data['Date']
data.head()

# les noms de colonne ne sont pas pratiques du tout
data.columns = ['Total', 'West', 'East']

data[data['Total'].isna()]

# ou encore si on préfère

data[data.isna().any(axis=1)].shape

# ceci va nous permettre d'avoir des colonnes d'entiers - avec des nan

data = data.convert_dtypes(convert_integer=True)
data.head()

# ça ressemble à ceci

data.dtypes

# on a toujours les n/a, mais ce n'est pas grave
data[data.isna().any(axis=1)].shape

%matplotlib inline

# no longer working ?
# plt.style.use('seaborn')

sns.set(rc={'figure.figsize': (12, 4)})
#plt.rcParams["figure.figsize"] = (12, 4)

# un premier jet, pas terrible du tout

data[['East', 'West']].plot();

# c'est plus lisible avec seulement un point par semaine
# on pourrait faire la moyenne aussi bien sûr,
# ça donnerait le même dessin mais avec les Y divisés par 7

# le point c'est qu'on a quelques années de plus que sur la vidéo

data.resample('W').sum().plot();

# c'est facile de couper, la date correpond à l'index de la df
# et grâce au type 'datetime' on peut simplement faire une comparaison

data = data[data.index.year <= 2017]

data.tail(3)

data.resample('W').sum().plot();

data.shape

# la forme du resample() est de:
data.resample('1W').sum().shape

# on vérifie que la version resamplée a bien 
# 7 * 24 = 168 fois moins d'entrées que la version brute
# puisqu'on a une mesure par heure et qu'on ré-échatillonne sur une semaine

(full, _), (resampled, _) = data.shape, data.resample('1W').sum().shape

full / resampled , 7 * 24

data.resample('1W').sum().plot();

# la somme sur une fenêtre tournante d'un an

# mais : méfiez-vous de l'échelle des Y

data.resample('1D').sum().rolling(365).sum().plot();

# on fait en sorte que le bas de l'échelle des Y soit bien 0
# pour eviter l'effet de loupe

ax = data.resample('1D').sum().rolling(365).sum().plot()
ax.set_ylim(0, None);

# regardons la tendance des profils journaliers en moyenne

data.groupby(data.index.time).mean().plot();

# mais pour y voir un peu mieux on veut afficher les jours 
# individuellement les uns des autres

# on veut donc dessiner autant de courbes que de jours
# et chaque courbe a en X l'heure de la journée et en Y le nombre (total) de passages

# pour ça on calcule une pivot table
# une courbe par jour: les colonnes sont les jours
# en X les heures: les lignes sont les heures

# ça pourrait se faire à coups de groupby/unstack
# mais c'est quand même plus simple comme ceci

pivoted = data.pivot_table('Total', index=data.index.time, columns=data.index.date)

# regardons le coin en haut à gauche

pivoted.iloc[:5, :5]

# on confirme les dimensions 

pivoted.shape

# du coup on n'a plus qu'à dessiner
# l'astuce qui tue c'est alpha=0.01 pour éviter les gros patés

pivoted.plot(legend=False, alpha=0.01);

pivoted.T.shape

# ! pip install sklearn
from sklearn.decomposition import PCA

# on transpose, et on remplit les trous avec des 0
pca_input = pivoted.T.fillna(0)

# on utilise le PCA de scikit-learn comme une boite noire
pca_output = PCA(2, svd_solver='full').fit_transform(pca_input)

# on obtient un tableau numpy avec deux colonnes seulement
# car on a demandé les deux premières composantes principales

type(pca_output), pca_output.shape

# on voit effectivement que cet ACP semble bien séparer deux clusters

plt.scatter(pca_output[:, 0], pca_output[:, 1], marker='.', color='green');

# pour les trouver ces deux clusters, Jake utilise une GaussianMixture

from sklearn.mixture import GaussianMixture
gmm = GaussianMixture(2)

# c'est ici que tout se passe
labels = gmm.fit(pca_input).predict(pca_input)


# la sortie est une association jour -> type
labels.shape, labels

# cette prédiction est bien en phase avec les deux clusters de tout à l'heure

plt.scatter(pca_output[:, 0], pca_output[:, 1], c=labels, cmap='rainbow')
plt.colorbar();

# pour vérifier notre classification on peut redessiner
# les jours classés label==0 

# ça correspond donc aux jours de la semaine (à moins que ce soit l'inverse..)

pivoted.loc[:, labels==0].plot(legend=False, alpha=0.01);

# et les jours classés label==1

pivoted.loc[:, labels==1].plot(legend=False, alpha=0.01);

# notre index horizontal n'est pas de type DatetimeIndex
pivoted.columns

# un index qui contient toutes nos dates et de type DatetimeIndex
dates = pd.DatetimeIndex(pivoted.columns)

# ceci nous calcule un index sur les jours
# mais avec comme valeur 0 pour le lundi, ... et 6 pour le dimanche

dayofweek = pd.DatetimeIndex(pivoted.columns).dayofweek

dayofweek.shape, dayofweek

# qu'on va utiliser pour mettre les jours en couleur
# les jours de weekend sont en orange et rouge

plt.scatter(pca_output[:, 0], pca_output[:, 1], c=dayofweek, cmap='rainbow')
# pour la légende
plt.colorbar();

# on remarque dans le cluster rouge-orange
# des jours d'une couleur qui jure

# pour comprendre à quoi ils correspondent 

odd_index = (labels == 1) & (dayofweek < 5)
odd_index.shape, odd_index

# afficher les 48 jours qui sont dans cette catégorie

# comme on peut s'y attendre
# on y retrouve les jours fériés (4 juillet, Noel, ...)

odd_dates = dates[odd_index]
odd_dates, len(odd_dates)

# pour rafficher seulement ces jours-là

pivoted[odd_dates].plot(legend=False, alpha=0.3);