Intro

• JAGS (per ulteriori info ed esempi vedi articolo su R-bloggers)
• rjags

Modellino gaussiano analizzato con JAGS

• problema: dato un campione di osservazioni, inferire μ e σ della gaussiana;
  inferire anche il risultato di una ipotetica ulteriore osservazione (y);
• dati (100 valori): gaussian_sample.dat (generati con crea_dati.R);
• modello (vedi anche pagina di bugs): gauss_Nm_sigma_pred.bug;
• comandi di JAGS: comandi_jags.txt:
  • i risultati ('catena', prodotta dal comando ". coda *") sono salvati in
    • CODAchain1.txt;
    • CODAindex.txt
• script R per analizzare la catena risultante: gauss_Nm_sigma_pred.R
  (per una intro su R vedi per esempio qui)

Idem usando rjags

• problema, dati e modello sono esattamente come sopra;
• function per leggere i dati: read.jagsdata.R
• script per chiamare rjags e analizzare il risultato: do.R
  (l'analisi è soltanto abbozzata: vedi gauss_Nm_sigma_pred.R per fare altre cose).

Variante con script R che genera i dati, chiama rjags e analizza la catena

• modello esattamente come sopra, mentre il file dei dati non è usato;
• script per fare tutto: crea_do.R.

Esempio di fit

• modello: molla.bug;
• file R: molla.R
• Per info sui dati, tabella dei risultati e confronto con fit utilizzzando il metodo dei residui, vedi “analisi della molla” nella pagina di BUGS.

Esempi per il corso di dottorato (26° Ciclo, marzo 2011)

• Processi di Bernoulli
  • File di 'dati' in formato R: dati_Bernoulli.R
  • Inferenza di p e previsione del numero di successi in altre prove:
    • file con modello: inf_p_pred.bug;
    • script R (chiamate a rjags + analisi dei risultati): inf_p_pred.R;
    • idem, ma con modello incluso nello script R: inf_p_pred_tempfile.R.
  • Ancora inferenza su p, ma usando la sola informazione del numero di successi in n prove:
    • file con modello: inf_p_pred_sum.bug;
    • script R: inf_p_pred_sum.R.
  • Inferenza di n dato p:
    • file con modello: inf_n.bug;
    • script R: inf_n.R.
    [Attenzione: non funziona con JAGS 3.0 e 3.1! vedi qui]
  • Ancora inferenza di p, ma assumendo che l'esperimento sia stato eseguito prefissando il numero di successi:
    • file con modello: inf_p_pred_prove.bug;
    • script R: inf_p_pred_prove.R.
• Processi di Poisson
  • Inferenza di λ dai conteggi osservati e previsione del numero di conteggi di una ulteriore misura (stesso tempo e stessa intensità del processo):
    • file con modello: inf_lambda_pred.bug;
    • script R: inf_lambda_pred.R.
  • Inferenza dell'intensità del processo r:
    • file con modello: inf_r.bug;
    • script R: inf_r.R.
  • Inferenza dell'intensità del processo r, tendendo conto del background (eventualmente incerto):
    • file con modello: inf_r_bck.bug;
    • script R: inf_r_bck.R.
  • Simile al precedente, ma con background valutato da altra misura di conteggio:
    • file con modello: inf_r_bck_measured.bug;
    • script R: inf_r_bck_measured.R.
    (Per un articolo sullo stesso modello, risolto in modo analitico, vedi articolo su NIM di Loparco e Mazziotta.)
  • Suggerimento per ulteriori complicazioni: aggiungere qualche tipo di sistematica e altre incertezze, ad esempio immaginando di dover risalire dall'intensità r ad un flusso (ad esempio di raggi cosmici) per unità di tempo e di area, tenendo conto di:
    • efficienza del rivelatore, ricavata eventualmente da dati sperimentali con piccoli numeri (vedi inferenza processi di Bernoulli);
    • incertezze su tempo di misura e area del rivelatore.
• Campioni gaussiani con sistematiche
  • Tre campioni normali di tre grandezze fisiche omogenee, ottenuti con lo stesso strumento affetto da errore statistico ed errori sistematici di zero e di scala:
    • modello: norm-syst.bug;
    • script R: norm-syst.R.
  [Per altri semplici inferenze (+distr. predittive) da campioni gaussiani, vedi sopra.]
• Inferenza condizionata da un vincolo
  • Problema: vengono misurati in modo indipendente i tre angoli di un triangolo e, senza alcun vincolo sulla somma degli angoli interni si ottengono i seguenti valori, espressi in gradi, con incertezze standard: 58±2, 73±2 e 54±2. Rivalutare i valori degli angoli tenendo conto del vincolo.
  • Soluzione mediante distribuzioni multivariate ricondizionate: tr_cond.R.
  • Soluzione mediante rete bayesiana (assume che i valori siano ottenuti con likelihood gaussiane e prior uniformi fra 0 e 180 gradi).):
    • modello: triangolo.bug;
    • script R: triangolo.R.
    [Nota: l'esempio può girare anche sotto WinBUGS/OpenBUGS, rinominando il modello triangolo.bug.txt e usando il file dei dati triangolo.dat.txt.]