Introduction pratique à la statistique bayésienne

Applications en python

Thierry Phénix
thierry-phenix@imad-ge.ch

IMAD – Institution genevoise de maintien à domicile, Genève Suisse

October 23, 2023

Préambule

  • Ce cours repose sur l’idée que vous avez une connaissance minimale de ce que sont les statistiques (au moins fréquentistes)
  • Ce cours n’est pas un cours de mathématique mais de méthodologie. Cependant …
  • Ce cours se veut interactif. On ne se verra qu’une fois donc… n’hésitez pas à poser des questions !
  • Ce cours et ses sources seront disponibles en ligne à l’adresse suivante :

Objectifs du cours

  • Vous aider à faire le point sur vos connaissances en statistique
  • Vous présenter les principes de l’analyse Bayésienne
  • Vous donner les points d’entrée pour vous y mettre !

Plan du cours

  • Introduction
  • La démarche statistique
  • Cas pratiques
  • Les outils
  • Conclusion
  • Définition
  • Application : comptage de photon
    • Collecte des données
    • Analyse descriptive
    • Inférence frequentiste
    • Inférence bayesienne
  • Bilan

Introduction

Définition

La statistique est la discipline qui étudie des phénomènes à travers la collecte de données, leur traitement, leur analyse, l’interprétation des résultats et leur présentation afin de rendre ces données compréhensibles par tous. (wikipedia)

\(\Longrightarrow\) Trois professions

  • Data engineer
  • Data analyst
  • Data scientist

Introduction

Définition

La statistique est la discipline qui étudie des phénomènes à travers la collecte de données, leur traitement, leur analyse, l’interprétation des résultats et leur présentation afin de rendre ces données compréhensibles par tous. (wikipedia)

\(\Longrightarrow\) Trois professions

  • Data engineer
  • Data analyst
  • Data scientist

\(\Longrightarrow\) Deux approches

  • Statistiques descriptives (EDA)
    • Obtenir de l’information sur l’échantillon considéré
  • Statistiques inférentielles
    • Obtenir de l’information sur la population dont a été extrait l’échantillon

Introduction

Exemple 1 : Comptage simple des photons

Problème:

Des chercheurs pointent leur télescope vers le ciel et observent la lumière provenant d’une seule étoile.

  • Nous supposerons que le flux réel de l’étoile est constant dans le temps, c’est-à-dire qu’il a une valeur fixe \(F_{true}\)
  • Nous négligeons le bruit engendré par l’atmosphère terrestre et les autres sources lumineuses.
  • Pour chaque mesure nous enregistrons la valeur du flux de photons (nb photons /s) et l’erreur de mesure du téléscope \(\to (F_{i}, e_{i})\)

Question: Sur la base de ces mesures/erreurs,
quelle est notre meilleure estimation du flux réel ?

Introduction

Example 1 : Collecte de données

Génèration des données

Generate data
# Generating some simple photon count data
import numpy as np
from scipy import stats
np.random.seed(42)  # for repeatability

F_true = 1000  # true flux, say number of photons measured in 1 second
N = 50 # number of measurements
F = stats.poisson(F_true).rvs(N)  # N measurements of the flux
e = np.sqrt(F)  # errors on Poisson counts estimated via square root

Analyse descriptive

Plot data
%matplotlib inline
import matplotlib.pyplot as plt

fig, ax = plt.subplots()
ax.errorbar(F, np.arange(N), xerr=e, fmt='ok', ecolor='gray', alpha=0.5)
ax.vlines([F_true], 0, N, linewidth=5, alpha=0.2)
ax.set_xlabel("Flux");ax.set_ylabel("measurement number");

Introduction

Example 1 : Exploratory Data Analysis (EDA)

C’est une première forme d’analyse qui ne renseigne que sur l’échantillon

Introduction

Example 1 : Approche fréquentiste

Estimation par maximum de vraisemblance

Illustration du maximum de vraisemblance

Introduction

Example 1 : Approche frequentiste

Estimation par maximum de vraisemblance

  • On considère que les données sont aléatoires mais pas les paramètres.

  • On cherche à estimer la valeur de flux \(F\) qui maximise la probabilité \(P(D~|~F)\).

    • On calcule la distribution de probabilité des données \(D_i = (F_i, e_i)\)
    • En considérantt le flux réel \(F_{\rm true}\) connu
    • En faisant l’hypothèse de distribution gaussienne de l’erreur

    \[ P(D_i~|~F_{\rm true}) = \frac{1}{\sqrt{2\pi e_i^2}} \exp{\left[\frac{-(F_i - F_{\rm true})^2}{2 e_i^2}\right]} \]

    • En considérant les \(D_i\) indépendantes et identiquement distribuées : \[P(D~|~F) = \prod_{i=1}^N P(D_i~|~F_{\rm true})\]

Introduction

Example 1 : Approche frequentiste

Estimation par maximum de vraisemblance

  • La fonction de vraisemblance

    \[P(D~|~F_{\rm true}) = \prod_{i=1}^N P(D_i~|~F_{\rm true})\]

    est ici une fonction de \(F_{\rm true}\) donc définie sur l’espace du paramètre.

  • On construit sa version log (numérique)

    \[\log\mathcal{L} = -\frac{1}{2} \sum_{i=1}^N \left[ \log(2\pi e_i^2) + \frac{(F_i - F_{\rm true})^2}{e_i^2} \right]\]

  • Pour trouver la valeur de flux qui maximise cette fonction, on cherche la valeur \(F\) du flux telle que

    \[d\log\mathcal{L}/dF_{\rm true} = 0\]

Introduction

Example 1 : Approche frequentiste

Estimation par maximum de vraisemblance

  • Le résultat est le suivant \[ F_{\rm est} = \frac{\sum w_i F_i}{\sum w_i};~~w_i = 1/e_i^2 \]

  • Remarquez que dans le cas particulier où toutes les erreurs \(e_i\) sont égaux, cela se réduit à \[ F_{\rm est} = \frac{1}{N}\sum_{i=1}^N F_i \]

  • Conformément à l’intuition, \(F_{\rm est}\) est simplement la moyenne des données observées lorsque les erreurs sont égales.

Introduction

Example 1 : Approche frequentiste

Estimation par maximum de vraisemblance

  • Quelle est l’erreur de cette estimation ?
    En ajustant une approximation gaussienne à la courbe de vraisemblance à son maximum, on obtient comme estimation de son écart-type :

    \[ \sigma_{\rm est} = \left(\sum_{i=1}^N w_i \right)^{-1/2} \]

  • Résultat numérique

w = 1. / e ** 2
f_est = (w * F).sum() / w.sum()
e_est = w.sum() ** -0.5
print(f"""
    F_true = {F_true}
    F_est  = {f_est:.0f} +/- {e_est:.0f} (based on {N} measurements)
    """)

    F_true = 1000
    F_est  = 999 +/- 4 (based on 50 measurements)

Introduction

Example 1 : Approche bayésienne

Estimation de la distribution du flux

  • On considère que les données ET les paramètres sont probabilistes.
    On peut donc écrire la distribution conjointe (impossible en fréquentiste) \[ P(F_{\rm true},~D)\]

  • On utilise la règle de Bayes pour exprimer la distribution recherchée \[ P(F_{\rm true}~|~D) = \frac{P(D~|~F_{\rm true})~P(F_{\rm true})}{P(D)} \]

    • \(P(F_{\rm true})\) : Le prior ou la distribution du flux avant d’avoir des données.
    • \(P(F_{\rm true}~|~D)\) : La postérieure, ou la nouvelle probabilité du flux connaissant des données : c’est le résultat que nous voulons calculer.
    • \(P(D~|~F_{\rm true})\) : La vraisemblance, qui est proportionnelle à la \(\mathcal{L}(D~|~F_{\rm true})\) dans l’approche fréquentiste.
    • \(P(D)\) : Un simple terme de normalisation.

Introduction

Example 1 : Approche bayésienne

Estimation de la distribution du flux

Dans la grande majorité des cas on ne sait pas calculer cette probabilité !
Deux solutions :

Introduction

Example 1 : Approche bayésienne

Estimation de la distribution du flux

  1. Stochastic Variational Inference 1

Introduction

Example 1 : Approche bayésienne

Estimation de la distribution du flux

  1. Markov Chain Monte Carlo (MCMC)

Introduction

Example 1 : Approche bayésienne

Estimation de la distribution du flux

#| echo: true
# import packages
from jax import random
import numpyro
import numpyro.distributions as dist
from numpyro.infer import MCMC, NUTS, Predictive
# define model
def model(e_obs, F_obs=None):
    F = numpyro.sample("F", dist.Uniform(0, 10000))             # prior uniform
    numpyro.sample("obs", dist.Normal(F, e_obs), obs=F_obs)     # échantillon distribué ~ loi normal  
# run sampling
mcmc = MCMC(NUTS(model), num_warmup=500, num_samples=2000, progress_bar=False)
mcmc.run(random.PRNGKey(42), e_obs=e, F_obs=F)
# display result
mcmc.print_summary(0.89)

Introduction

Example 1 : Approche bayésienne

Estimation de la distribution du flux

# import packages
from jax import random
import numpyro
import numpyro.distributions as dist
from numpyro.infer import MCMC, NUTS, Predictive
# define model
def model(e_obs, F_obs=None):
    F = numpyro.sample("F", dist.Uniform(0, 10000))             # prior uniform
    numpyro.sample("obs", dist.Normal(F, e_obs), obs=F_obs)     # échantillon distribué ~ loi normal  
# run sampling
mcmc = MCMC(NUTS(model), num_warmup=500, num_samples=2000, progress_bar=False)
mcmc.run(random.PRNGKey(42), e_obs=e, F_obs=F)
# display result
mcmc.print_summary(0.89)

                mean       std    median      5.5%     94.5%     n_eff     r_hat
         F    998.98      4.56    998.85    992.00   1006.73    526.64      1.00

Number of divergences: 0

Introduction

Bilan

  • Les données sont tirées aléatoirement d’une même distribution (i.i.d)
  • L’inférence fréquentiste :
    • Cadre : les données sont aléatoires mais les paramètres sont déterministes
    • Pas de distribution sur le paramètre (résultat ponctuel)
    • Fait appel à un autre modèle pour avoir un niveau de confiance
    • Facile à calculer !
  • L’inférence bayésienne
    • Cadre : les données et les paramètres sont aléatoires
    • Répond directement à la question par une distribution sur le paramètre
    • L’incertitude ne dépend que des données
    • Compliqué à calculer (nécessite un échantillonneur)

Plan du cours

  • Introduction
  • La démarche statistique
  • Cas pratiques
  • Les outils
  • Conclusion
  • Juste pour être sûr !
  • La démarche
  • Point de vue fréquentiste
  • Point de vue bayésien

La démarche statistique

Juste pour être sûr !

Revoyons le concept Bayesian avec un jeu : Monty Hall

La démarche statistique

Juste pour être sûr !

Posons le problème de manière bayésienne :

  1. L’espace des probabilité est conditionné par le fait que le joueur a choisi la porte A
  2. L’hypothèse \(H\) que l’on considère est : La voiture se trouve derière la porte A, B ou C
  3. La données \(D\) est : l’animateur a ouvert la porte C

Question : Une fois que l’animateur a ouvert la porte C, le joueur a-t-il intérêt à changer son choix ?

\[P(H|D \land A) ?\]

La démarche statistique

Juste pour être sûr !

Utilisons le théorème de bayes :

\[P(H~|~D \land A) = \frac{P(D~|~H \land A)~P(H~|~A)}{P(D~|~A)}\]

  • \(P(H~|~A)\) : Probabilité que la voiture soit derrière les portes A, B ou C.
    Cette probabilité est indépendante du fait que le joueur est choisi A.
Monty Hall.
H \(P(H~|~A)\) \(P(D~|~H~\land~A)\) \(P(H~|~D~\land~A)\)
porte A 1/3
porte B 1/3
porte C 1/3

La démarche statistique

Juste pour être sûr !

Posons le problème de manière bayésienne :

\[P(H~|~D \land A) = \frac{P(D~|~H \land A)~P(H~|~A)}{P(D~|~A)}\]

  • \(P(H~|~A)\) : Probabilité que la voiture soit derrière les portes A, B ou C.
    Cette probabilité est indépendante du fait que le joueur est choisi A.
  • \(P(D~|~H~\land~A)\) : Probabilité que l’animateur ouvre la porte C sachant que le joueur a choisi la porte A et que la voiture se trouve derrière la porte H.
Monty Hall.
H \(P(H~|~A)\) \(P(D~|~H~\land~A)\) \(P(H~|~D~\land~A)\)
porte A 1/3 
   1/2
porte B 1/3 
     1
porte C 1/3 
     0

La démarche statistique

Juste pour être sûr !

Posons le problème de manière bayésienne :

\[P(H~|~D \land A) = \frac{P(D~|~H \land A)~P(H~|~A)}{P(D~|~A)}\]

  • \(P(H~|~A)\) : Probabilité que la voiture soit derrière les portes A, B ou C.
    Cette probabilité est indépendante du fait que le joueur est choisi A.
  • \(P(D~|~H~\land~A)\) : Probabilité que l’animateur ouvre la porte C sachant que le joueur a choisi la porte A et que la voiture se trouve derrière la porte H.
  • \(P(H~|~D~\land~A)\) : Probabilité postérieure de chaque porte sachant que l’animateur ouvre la porte C sachant que le joueur a choisi la porte A.
Monty Hall.
H \(P(H|A)\) \(P(D|H \land A)\) \(P(H|D \land A)\)
porte A 1/3 
 1/2
 
 1/3
porte B 1/3 
  1
 
 2/3
porte C 1/3 
  0
 
   0

La démarche statistique

En quoi ça consiste

Plan d’analyse

  1. Poser la question d’intérêt clairement
    • Identifier les variables prédictrices (variables indépendantes)
    • Identifier les variables prédites (variables dépendantes)
    • Identifier les variables non mesurées (mais influentes)

La démarche statistique

En quoi ça consiste

Plan d’analyse

  1. Poser la question d’intérêt clairement
  1. Concevoir le graph de dépendances causales
    • Définir les relations de causalité entre ces variables
    • Utiliser le lien suivant : https://www.dagitty.net

La démarche statistique

En quoi ça consiste

Plan d’analyse

  1. Poser la question d’intérêt clairement
  2. Concevoir le graph de dépendances causales
  1. Définir le modèle génératif
    • Comment sont générées les données ?
    • La VD évolue avec les VIs \(\to\) modèle linéaire
    • La VD est binaire \(\to\) modèle logistiques
    • La VD est un temps d’attente \(\to\) modèle de poisson

La démarche statistique

En quoi ça consiste

Exemple famille des distributions exponentielles

statistical rethinking 2nd edition page 283

La démarche statistique

En quoi ça consiste

Plan d’analyse

  1. Poser la question d’intérêt clairement
  2. Concevoir le graph de dépendances causales
  3. Définir le modèle génératif
  1. Implémenter le modèle et calculer les estimateurs
    • Construire le modèle
    • Définir les priors
    • Faire tourner l’inférence
    • Vérifier l’algorithme

La démarche statistique

En quoi ça consiste

Plan d’analyse

  1. Poser la question d’intérêt clairement
  2. Concevoir le graph de dépendances causales
  3. Définir le modèle génératif
  4. Implémenter le modèle et calculer les estimateurs
  1. Confronter le résultat à la réalité

    • Visualiser les résultats

    • Boucler entre 3 et 5 jusqu’à obtenir une réponse satisfaisante

Plan du cours

  • Introduction
  • La démarche statistique
  • Cas pratiques
  • les outils
  • Conclusion
  • EDA
  • Estimation de paramètres
  • Test d’hypothèse
  • Comparaison de modèles

Cas pratiques

Exploratory data analysis

Les données

Chargement du jeu de données 1

#| echo: false
import pandas as pd
df_howell1 = pd.read_csv("./Data/Howell1.csv", sep=";")
df_howell1.head(5)
height weight age male
0 151.765 47.825606 63.0 1
1 139.700 36.485807 63.0 0
2 136.525 31.864838 65.0 0
3 156.845 53.041914 41.0 1
4 145.415 41.276872 51.0 0

Cas pratiques

Exploratory data analysis

Les données

import seaborn as sns
sns.pairplot(data= df_howell1, hue="male", height=1.5, aspect=1.5)

Cas pratiques

Estimation de paramètres

Premier problème

Question : Peut-on prédire la taille à partir du poids après 18 ans ?

  • Préparation des données :

    # select observation based on age >= 18
data_sup18 = df_howell1[df_howell1.age>=18].copy()
# Center height -> intercept = weight at mean height
center = lambda x: x - x.mean() 
data_sup18["height_c"] =  center(data_sup18.height)

  • Pour Numpyro

    # create dict for inference with numpyro
dat1 = {
    "height" : data_sup18.height_c.values,
    "weight" : data_sup18.weight.values
}

Cas pratiques

Estimation de paramètres

Peut-on prédire le poids à partir de la taille après 18 ans ?

Définition du modèle

\[ \begin{align} W_i \sim Normal(\mu_i, \sigma) \\ mu_i = \alpha + \beta H_i \\ \alpha \sim Normal(0, 10) \\ \beta \sim HalfNormal(1) \\ \sigma \sim HalfNormal(1) \\ \end{align} \]

Cas pratiques

Estimation de paramètres

Peut-on prédire la taille à partir du poids après 18 ans ?

Implémentation du modèle

# Define the model
def mdl1(height, weight=None):
    # Define prior
    a = numpyro.sample('a', dist.Normal(0, 1))
    b = numpyro.sample('b', dist.Normal(0, 1))
    s = numpyro.sample('s', dist.HalfNormal(1))
    # Define likelihood
    numpyro.sample("obs", dist.Normal(a+b*height, s), obs=weight)
# plot model
numpyro.render_model(mdl1, model_args=dat1.values(), render_distributions=True)

Cas pratiques

Estimation de paramètres

Peut-on prédire la taille à partir du poids après 18 ans ?

  • Calcul des résultats
mcmc1 = MCMC(NUTS(mdl1), num_warmup=500, num_samples=1000, num_chains=2)
mcmc1.run(random.PRNGKey(42), **dat1)

Cas pratiques

Estimation de paramètres

Peut-on prédire la taille à partir du poids après 18 ans ?

  • Vérification de la convergence de l’algorithme

Cas pratiques

Estimation de paramètres

Peut-on prédire la taille à partir du poids après 18 ans ?

  • Visualisation des résultats
az.plot_posterior(idata1, figsize=(10,3));

Cas pratiques

Estimation de paramètres

Peut-on prédire la taille à partir du poids après 18 ans ?

  • Conclusion
mcmc1.print_summary(0.89)

                mean       std    median      5.5%     94.5%     n_eff     r_hat
         a     41.91      0.42     41.93     41.25     42.59   1014.03      1.00
         b      0.63      0.03      0.63      0.57      0.68   1297.06      1.00
         s      5.09      0.29      5.07      4.63      5.54    935.34      1.00

Number of divergences: 0
  • On peut dire :
    • Pour la taille moyenne, le poids est de 41.89kg et à 89% dans [41.25, 42.59]
    • Pour chaque centimètre supplémentaire, le poids augmente de 630g

Cas pratiques

Test d’hypothèse

Deuxième problème

Question : Quelle est la différence de poids entre les hommes et les femmes après 18 ans ?

Cas pratiques

Test d’hypothèse

Différence de poids hommes/femmes après 18 ans ?

Définition du modèle

\[ \begin{align} W_i \sim Normal(\mu_i, \sigma) \\ \mu_i = \alpha[male] \\ \alpha \sim HalfNormal(0, 10) \\ \sigma \sim HalfNormal(1) \\ \end{align} \]

Cas pratiques

Test d’hypothèse

Différence de poids hommes/femmes après 18 ans ?

Implémentation du modèle

# create dict for inference with numpyro
dat2 = {
    "male" : data_sup18.male.values,
    "weight" : data_sup18.weight.values
}
# Define the model
def mdl2(male, weight=None):
    # Define prior
    with numpyro.plate("male", 2):
        a = numpyro.sample('a', dist.HalfNormal(10))
    s = numpyro.sample('s', dist.HalfNormal(1))
    diff_a = numpyro.deterministic("diff_a", a[1]- a[0])
    # Define likelihood
    numpyro.sample("obs", dist.Normal(a[male], s), obs=weight)
# plot model
numpyro.render_model(mdl2, model_args=dat2.values(), render_distributions=True)

Cas pratiques

Test d’hypothèse

Différence de poids hommes/femmes après 18 ans ?

Graph du modèle

Cas pratiques

Test d’hypothèse

Différence de poids hommes/femmes après 18 ans ?

  • Calcul des résultats
mcmc2 = MCMC(NUTS(mdl2), num_warmup=500, num_samples=1000)
mcmc2.run(random.PRNGKey(42), **dat2)

Cas pratiques

Test d’hypothèse

Différence de poids hommes/femmes après 18 ans ?

  • Vérification de la convergence de l’algorithme

Cas pratiques

Test d’hypothèse

Différence de poids hommes/femmes après 18 ans ?

  • Visualisation des résultats
az.plot_posterior(idata2, figsize=(10,3));

Cas pratiques

Test d’hypothèse

Différence de poids hommes/femmes après 18 ans ?

  • Conclusion
az.summary(idata2, hdi_prob=0.89)
mean sd hdi_5.5% hdi_94.5% mcse_mean mcse_sd ess_bulk ess_tail r_hat
a[0] 41.755 0.387 41.159 42.353 0.013 0.009 856.0 552.0 NaN
a[1] 48.511 0.439 47.904 49.291 0.013 0.009 1080.0 555.0 NaN
diff_a 6.756 0.607 5.707 7.607 0.020 0.015 874.0 658.0 NaN
s 5.315 0.196 5.028 5.645 0.007 0.005 772.0 638.0 NaN
  • On peut dire que en moyenne la différence de poids entre les hommes et les femmes âgés de plus de 18 ans est de 6.7kg en faveur des homme et se touve à 89% dans [5.7, 7.6].

Cas pratiques

Comparaison de modèles

Troisième problème

Question : Quelle est le meilleur modèle pour prédire le poids en fonction de la taille ?

Cas pratiques

Comparaison de modèles

Meilleur modèle pour prédire le poids en fonction de la taille ?

  • Modèle linéaire
# This time with all the data
df_howell1["height_c"] =  center(df_howell1.height)
dat3 = {
    "height" : df_howell1.height_c.values,
    "weight" : df_howell1.weight.values
}
# Define the model
def mdl3a(height, weight=None):
    # Define prior
    a = numpyro.sample('a', dist.Normal(0, 10))
    b = numpyro.sample('b', dist.Normal(0, 1))
    s = numpyro.sample('s', dist.HalfNormal(1))
    # Define likelihood
    numpyro.sample("obs", dist.Normal(a+b*height, s), obs=weight)

Cas pratiques

Comparaison de modèles

Meilleur modèle pour prédire le poids en fonction de la taille ?

  • Modèle quadratic
# Define the model
def mdl3b(height, weight=None):
    # Define prior
    a = numpyro.sample('a', dist.Normal(0, 10))
    b1 = numpyro.sample('b1', dist.Normal(0, 1))
    b2 = numpyro.sample('b2', dist.Normal(0, 1))
    s = numpyro.sample('s', dist.HalfNormal(1))
    # Define likelihood
    numpyro.sample("obs", 
        dist.Normal(a+b1*height+b2*height**2, s), 
        obs=weight)
  • Modèle cubic
# Define the model
def mdl3c(height, weight=None):
    # Define prior
    a = numpyro.sample('a', dist.Normal(0, 10))
    b1 = numpyro.sample('b1', dist.Normal(0, 1))
    b2 = numpyro.sample('b2', dist.Normal(0, 1))
    b3 = numpyro.sample('b3', dist.Normal(0, 1))
    s = numpyro.sample('s', dist.HalfNormal(1))
    # Define likelihood
    numpyro.sample("obs", 
        dist.Normal(a+b1*height+b2*height**2+b3*height**3, s), 
        obs=weight)

Cas pratiques

Comparaison de modèles

Meilleur modèle pour prédire le poids en fonction de la taille ?

  • Résultats modèle 1

                mean       std    median      5.5%     94.5%     n_eff     r_hat
         a     35.59      0.20     35.59     35.28     35.91    577.64      1.00
         b      0.50      0.01      0.50      0.49      0.51   1082.96      1.00
         s      4.91      0.15      4.91      4.68      5.14   1011.70      1.00

Number of divergences: 0
  • Résultats modèle 2

                mean       std    median      5.5%     94.5%     n_eff     r_hat
         a     32.51      0.25     32.52     32.08     32.90    542.73      1.00
        b1      0.64      0.01      0.64      0.63      0.66    314.48      1.00
        b2      0.00      0.00      0.00      0.00      0.00    417.83      1.00
         s      3.91      0.12      3.91      3.73      4.10    389.55      1.00

Number of divergences: 0
  • Résultats modèle 3

                mean       std    median      5.5%     94.5%     n_eff     r_hat
         a      1.76      0.02      1.75      1.73      1.79      3.75      1.53
        b1      1.02      0.02      1.02      0.98      1.05      9.46      1.03
        b2      0.04      0.00      0.04      0.04      0.04    741.49      1.00
        b3      0.00      0.00      0.00      0.00      0.00     75.08      1.00
         s     16.30      0.38     16.30     15.74     16.83    117.14      1.00

Number of divergences: 0

Cas pratiques

Comparaison de modèles

Meilleur modèle pour prédire le poids en fonction de la taille ?

  • Comparaison des trois modèles
rank elpd_waic p_waic elpd_diff weight se dse warning scale
mdl quardatic 0 -1520.630761 3.796908 0.000000 0.945619 20.899299 0.000000 False log
mdl linear 1 -1648.960918 3.240282 128.330157 0.054381 14.881861 16.194413 False log
mdl cubic 2 -2422.694490 3.614403 902.063729 0.000000 19.844980 22.873266 True log

Plan du cours

  • Introduction
  • La démarche statistique
  • Cas pratiques
  • Les outils
  • Conclusion
  • Les livres
  • Les packages par langage

Les outils

Les livres

Plutôt pour les psychos

Livre complet qui décrit l’ensemble des cas auxquels vous serez confrontés. A lire !

Livre excellent accompagné de viséos très claires. Je recommande surtout pour l’analyse causale qui est fondamentale.

Les outils

Les livres

Plutôt pour python

Livre facile à comprendre qui met l’accent sur le package pymc. A lire si vous voulez utiliser ce package !

Livre excellent qui apporte une autre perspective et des graphisme chatoyants ;-).

Les outils

Les livres

Plutôt théoriques

Livre théorique support de Stan, le package de référence tous langages confondus. Ardu à lire mais vraiment intéressant

Livre théorique qui étend la cadre probabiliste du psychologue en l’intégrant dans une démarche plus généraliste de machine learning.

Les outils

Les packages

En python

C’est le package que je recommande pour commencer. Il est complet et permet de faire toutes les analyses que l’on veut. Il est un peu lent mais peut-être couplé avec jax qui est plus puissant.

Bambi est un package qui s’utilise conjointement avec Pymc et qui apporte la syntaxe que l’on retrouve dans R.

Les outils

Les packages

En python les mastodontes

L’intérêt de ces deux packages est qu’ils s’ouvrent sur les réseaux de neurones…

Pyro et Numpyro sont développés par Uber. Ils sont hyper puissants tout en gardanst une certaine simplicité d’utilisation. Numpyro est plus simple et repose sur Jax.

TensofFlow Probability est le package développé par Google. Les classes qui décrivent les distributions de probabilités ont été reprises par Pyro et Numpyro. C’est assez complexe !

Les outils

Les packages

En R

BRMS est le package de référence en R.

Les outils

Les packages

En Julia

Julia est un “nouveau langage” que je conseille vivement. Il s’apprend comme python ou R mais la syntaxe n’a pas hérité de dizaine d’années de modifications ce qui la rend simple et efficace.

Pour les curieux, Turing est un package vraiment complet. On peut écrire des modèles comme avec BRMS et la puissance de julia fait le reste… Si vous hésitez alors celui-là mérite tout votre intérêt