Appendice M — Introduzione al linguaggio Stan

Che cos’è Stan?

Stan è un linguaggio di programmazione progettato specificamente per l’analisi statistica bayesiana. Il nome rende omaggio a Stanislaw Ulam, uno dei pionieri dei metodi Monte Carlo. A differenza di software statistici come SPSS o Jamovi, Stan richiede di scrivere esplicitamente il modello statistico, ma in cambio offre una flessibilità straordinaria: qualunque modello bayesiano può essere implementato in Stan.

Non è necessario essere programmatori esperti per utilizzare Stan. Il linguaggio è stato progettato per essere leggibile e la sua sintassi rispecchia da vicino la notazione matematica utilizzata per descrivere i modelli statistici. Se si riesce a scrivere un modello in forma matematica, è possibile tradurlo in Stan.

M.1 Come funziona Stan

Quando eseguiamo un’analisi con Stan, il processo si svolge in tre fasi:

  1. Scrittura del modello: descriviamo il modello statistico in un file di testo (con estensione .stan).
  2. Compilazione: Stan traduce il nostro modello in codice ottimizzato per il computer. Questa operazione richiede alcuni secondi e viene eseguita una sola volta per ogni modello.
  3. Campionamento: Stan genera migliaia di campioni dalla distribuzione a posteriori utilizzando algoritmi MCMC.

M.2 Come accedere a Stan

Stan non viene utilizzato direttamente, ma tramite un’interfaccia che lo collega al nostro ambiente di lavoro. Le interfacce principali sono:

Interfaccia Linguaggio Note
cmdstanr R Raccomandata per R, usata in questo libro
brms R Sintassi semplificata, ideale per modelli standard
rstanarm R Modelli predefiniti con sintassi familiare
CmdStanPy Python Equivalente di cmdstanr per Python
Stan.jl Julia Per utenti Julia

In questo libro utilizziamo cmdstanr, che offre il miglior equilibrio tra controllo e semplicità. Per l’installazione, consultare l’Appendice K.

Esistono anche pacchetti complementari per l’analisi dei risultati:

  • bayesplot: grafici per visualizzare distribuzioni a posteriori e diagnostiche.
  • posterior: funzioni per manipolare i campioni MCMC.
  • loo: confronto tra modelli tramite validazione incrociata.

M.3 La struttura di un programma Stan

Un programma Stan è organizzato in blocchi, ciascuno con uno scopo preciso. Pensate ai blocchi come a sezioni di un modulo da compilare: ogni sezione richiede informazioni specifiche.

data {
  // Cosa osserviamo: i dati dell'esperimento
}

parameters {
  // Cosa vogliamo stimare: i parametri incogniti
}

model {
  // Come sono collegati: prior e verosimiglianza
}

generated quantities {
  // Cosa vogliamo calcolare in più (opzionale)
}

I blocchi devono apparire sempre in questo ordine, ma non tutti sono obbligatori. I tre blocchi fondamentali sono data, parameters e model.

M.3.1 Il blocco data: cosa abbiamo osservato

Nel blocco data dichiariamo tutte le informazioni che provengono dal nostro studio: il numero di partecipanti, le risposte osservate, le variabili indipendenti, e così via. Per ogni variabile dobbiamo specificare:

  • il tipo di dato (numero intero, numero reale, vettore, ecc.);
  • le dimensioni (quanti elementi contiene);
  • eventuali vincoli sui valori ammissibili.

Esempio commentato:

data {
  int<lower=1> N;              // Numero di partecipanti (almeno 1)
  int<lower=0> y;              // Risposte corrette osservate
  int<lower=0> n_trials;       // Numero totale di prove
  real<lower=0> alpha_prior;   // Parametro per il prior
  real<lower=0> beta_prior;    // Parametro per il prior
}

Perché specificare i vincoli? I vincoli (come lower=0) aiutano Stan a verificare che i dati siano sensati. Se per errore passiamo un numero negativo di partecipanti, Stan ci avviserà immediatamente invece di produrre risultati sbagliati.

M.3.2 Il blocco parameters: cosa vogliamo stimare

Nel blocco parameters dichiariamo i parametri del modello, ovvero le quantità incognite che vogliamo stimare a partire dai dati. Anche qui specifichiamo tipo e vincoli.

Esempio:

parameters {
  real<lower=0, upper=1> theta;  // Probabilità di successo
}

In questo esempio, theta rappresenta una probabilità, quindi deve essere compresa tra 0 e 1. I vincoli lower=0, upper=1 comunicano a Stan questa restrizione.

M.3.3 Il blocco model: come sono collegati dati e parametri

Il blocco model è il cuore del programma Stan. Qui definiamo:

  1. la distribuzione a priori per ogni parametro;
  2. la verosimiglianza, cioè come i dati dipendono dai parametri.

Esempio:

model {
  // Prior: cosa crediamo prima di vedere i dati
  theta ~ beta(alpha_prior, beta_prior);
  
  // Verosimiglianza: come i dati dipendono da theta
  y ~ binomial(n_trials, theta);
}

Come leggere questo codice: il simbolo ~ si legge “è distribuito come” oppure “segue una distribuzione”. Quindi:

  • theta ~ beta(2, 2) significa “\(\theta\) ha distribuzione Beta(2, 2)”;
  • y ~ binomial(n, theta) significa “\(y\) ha distribuzione Binomiale con parametri \(n\) e \(\theta\)”.

Questa sintassi corrisponde direttamente alla notazione matematica: \[ \theta \sim \text{Beta}(\alpha, \beta) \] \[ y \sim \text{Binomiale}(n, \theta) \]

Nota sui prior: se non specifichiamo un prior per un parametro, Stan assume un prior uniforme (non informativo). Tuttavia, è buona pratica dichiarare sempre esplicitamente i prior per rendere il modello più leggibile.

M.3.4 Il blocco generated quantities: calcoli aggiuntivi

Questo blocco opzionale permette di calcolare quantità derivate dai parametri stimati. Due usi comuni sono:

  1. posterior predictive check: generare dati simulati per verificare se il modello cattura le caratteristiche dei dati reali;
  2. log-verosimiglianza: calcolare valori utili per il confronto tra modelli.

Esempio:

generated quantities {
  int y_rep = binomial_rng(n_trials, theta);  // Dato simulato
  real log_lik = binomial_lpmf(y | n_trials, theta);  // Log-verosimiglianza
}

Le funzioni che terminano in _rng generano numeri casuali; quelle che terminano in _lpmf o _lpdf calcolano le log-probabilità.

M.4 I tipi di dati in Stan

Stan implementa un sistema di tipi di dati fortemente tipizzato. Questa rigidità, sebbene possa apparire scomoda, garantisce l’integrità del modello eliminando alla fonte numerose categorie di errori.

M.4.1 Tipi scalari: singoli valori

Tipo Descrizione Esempio
int Numero intero int N = 10;
real Numero reale (decimale) real sigma = 2.5;

M.4.2 Tipi composti: collezioni di valori

Tipo Descrizione Esempio
vector[N] Vettore colonna di N numeri reali vector[10] y;
row_vector[N] Vettore riga di N numeri reali row_vector[5] x;
matrix[R, C] Matrice R×C di numeri reali matrix[3, 4] A;
array[N] int Array di N numeri interi array[20] int counts;

Quando usare cosa:

  • usate int per conteggi (numero di partecipanti, numero di risposte corrette);
  • usate real per misure continue (tempi di reazione, punteggi su scale);
  • usate vector per sequenze di misure continue (punteggi di \(N\) partecipanti);
  • usate array[N] int per sequenze di conteggi.

M.4.3 Vincoli sui valori

I vincoli più comuni sono:

Vincolo Significato Uso tipico
lower=0 Valore minimo 0 Deviazioni standard, conteggi
upper=1 Valore massimo 1 Proporzioni
lower=0, upper=1 Tra 0 e 1 Probabilità
lower=1 Almeno 1 Numerosità campionarie

Esempio con vincoli:

data {
  int<lower=1> N;                    // Almeno 1 partecipante
  array[N] real<lower=0> rt;         // Tempi di reazione (positivi)
}

parameters {
  real<lower=0> mu;                  // Media (positiva per tempi)
  real<lower=0> sigma;               // Deviazione standard (positiva)
}

M.5 Passare i dati da R a Stan

Per eseguire un modello Stan da R, dobbiamo preparare i dati in una lista. I nomi degli elementi della lista devono corrispondere esattamente a quelli dichiarati nel blocco data.

Esempio in R:

# Supponiamo di avere un esperimento con 50 prove e 32 successi
data_list <- list(
  N = 1,              # Un partecipante
  y = 32,             # Successi osservati
  n_trials = 50,      # Prove totali
  alpha_prior = 2,    # Parametro prior
  beta_prior = 2      # Parametro prior
)

Errori comuni da evitare:

  • nomi scritti diversamente (ntrials invece di n_trials);
  • tipi sbagliati (passare "50" come stringa invece di 50 come numero);
  • dimensioni non corrispondenti.

M.6 Regole sintattiche essenziali

Stan ha alcune regole sintattiche che è importante conoscere.

M.6.1 Punto e virgola

Ogni istruzione deve terminare con un punto e virgola (;):

int N;                          // Corretto
theta ~ beta(2, 2);             // Corretto
real<lower=0> sigma             // ERRORE: manca il punto e virgola

M.6.2 Commenti

I commenti iniziano con // e vengono ignorati da Stan:

// Questo è un commento su una riga
int N;  // Commento a fine riga

/* Questo è un commento
   su più righe */

M.6.3 Indicizzazione

Stan numera gli elementi a partire da 1 (come R), non da 0 (come Python):

vector[3] v;
v[1]  // Primo elemento
v[3]  // Terzo e ultimo elemento

M.7 Un esempio completo

Mettiamo insieme tutti i concetti in un modello completo per stimare la probabilità di successo in un compito cognitivo:

// Modello Beta-Binomiale per stimare la probabilità di successo

data {
  int<lower=1> n_trials;           // Numero di prove
  int<lower=0, upper=n_trials> y;  // Successi osservati
  real<lower=0> alpha_prior;       // Prior Beta: parametro alpha
  real<lower=0> beta_prior;        // Prior Beta: parametro beta
}

parameters {
  real<lower=0, upper=1> theta;    // Probabilità di successo
}

model {
  // Distribuzione a priori
  theta ~ beta(alpha_prior, beta_prior);
  
  // Verosimiglianza
  y ~ binomial(n_trials, theta);
}

generated quantities {
  // Per verificare il modello: generiamo dati simulati
  int y_rep = binomial_rng(n_trials, theta);
  
  // Per confrontare modelli: calcoliamo la log-verosimiglianza
  real log_lik = binomial_lpmf(y | n_trials, theta);
}

Questo modello:

  1. riceve in input il numero di prove e i successi osservati;
  2. stima la probabilità di successo \(\theta\) con un prior Beta;
  3. genera dati simulati per verificare l’adeguatezza del modello;
  4. calcola la log-verosimiglianza per eventuali confronti tra modelli.

M.8 Blocchi opzionali avanzati

Oltre ai blocchi fondamentali, Stan offre blocchi aggiuntivi per modelli più complessi:

M.8.1 transformed data

Permette di pre-elaborare i dati prima dell’analisi. È utile per standardizzare variabili o calcolare quantità derivate dai dati.

transformed data {
  real y_mean = mean(y);           // Media dei dati
  vector[N] y_centered = y - y_mean;  // Dati centrati
}

M.8.2 transformed parameters

Permette di definire parametri derivati da altri parametri. Questi vengono salvati insieme ai parametri principali.

transformed parameters {
  real odds = theta / (1 - theta);  // Odds derivato da theta
}

M.8.3 functions

Permette di definire funzioni personalizzate da usare nel modello. Sono utili per modelli complessi con calcoli ripetuti.

M.9 Risorse per approfondire

Per una documentazione completa del linguaggio Stan: