Ultimamente ho letto quasi solo articoli di MCMC. Sarà perché ho dovuto preparare un esame sull’argomento a Natale. C’è un articolo che azzera (!) il bias (e aumenta la varianza), uno che riduce di 10 o più volte la varianza (e vai con il bias!), uno che è più veloce di tutti (tranne sull’esempio difficile), uno che claima di essere fast con i gruppi di Lie ma non mi sembra mica.
Markov Chain Monte Carlo Methods, a survey with some frequent misunderstandings
Una piccola review basata su domande e risposte di Cross Validated. Ci ho trovato un paio di cose utili.
(1) Ho scoperto che quando campioni separatamente il priore e poi la sampling distribution tenendo il risultato solo se viene uguale (o circa uguale) ai dati, si chiama metodo ABC. E che se per efficienza usi una statistica dei dati, è meglio che abbia la stessa dimensionalità del parametro.
(2) L’importance resampling non è marginalmente esatto perché la somma dei pesi sul campione è casuale. All’inizio questa cosa non mi convinceva (che c’entra la somma dei pesi?), allora ho fatto i conti. Siano $p_1$ e $p_2$ densità di probabilità sullo stesso spazio. So campionare $p_2$ ma vorrei campionare $p_1$. Allora estraggo $\mathbf x = (x_1, \dots, x_n)$ secondo $p(\mathbf x) = \prod_{i=1}^n p_2(x_i)$ e definisco i pesi $w_i \equiv \frac {p_1(x_i)} {p_2(x_i)}$. Si ha che i valori attesi calcolati con i pesi sono uguali a quelli se avessi estratto con $p_1$, perché $$ E_{x \sim p_2} \left[ f(x) \frac {p_1(x)} {p_2(x)} \right] = E_{x \sim p_1} [f(x)]. $$
Però se estraggo $\tilde x$ da $\mathbf x$ con probabilità proporzionale a $w_i$, ovvero $$ p(\tilde x|\mathbf x) = \sum_{i=1}^n \delta(\tilde x – x_i) \frac {w_i} {\sum_{j=1}^n w_j}, $$ allora anche marginalizzando su $\mathbf x$ non viene $p(\tilde x) = p_1(\tilde x)$, bensì $$ p(\tilde x) = p_1(\tilde x) E_{x_k \sim p_2} \left[ \frac {n}{\frac {p_1(\tilde x)}{p_2(\tilde x)} + \sum_{k=1}^{n-1} \frac {p_1(x_k)}{p_2(x_k)}} \right], $$ che approssimando un paio di volte al primo ordine diventa $$ p(\tilde x) \approx p_1(\tilde x) \left( 1 + \frac 1n \left( 1 – \frac {p_1(\tilde x)}{p_2(\tilde x)} \right) \right). $$
In particolare $p(\tilde x) < p_1(\tilde x)$ se $p_2(\tilde x) < p_1(\tilde x)$, cioè il ripesamento non si “scorda” del tutto la distribuzione che ho usato per estrarre i campioni. Infatti nel caso $n = 1$ viene banalmente $p(\tilde x) = p_2(\tilde x)$, quindi ci deve essere una transizione continua da $p_2$ a $p_1$ aumentando $n$.
Expectation propagation as a way of life: A framework for Bayesian inference on partitioned data
Non l’ho letto tutto perché è lungo e dettagliato, comunque spiegano un metodo per fare un’inferenza bayesiana approssimata quando i dati sono divisi in vari sottoinsiemi separati (nel senso che la likelihood si fattorizza), l’inferenza completa sarebbe computazionalmente impegnativa, e le varie parti che hanno i dati non vogliono condividere i propri dati (ad esempio per ragioni di privacy).
Unbiased Markov chain Monte Carlo with couplings
Unbiased MCMC significa che ottengo una stima unbiased senza rimuovere i campioni iniziali, cioè senza bisogno che la catena converga alla distribuzione target. Non c’è trucco, non c’è magia, signori! È molto interessante ma come al solito c’è l’inghippo che quando tolgo il bias la varianza aumenta e nemmeno so quant’è. Il punto è che con una catena unbiased posso parallelizzare il conto: anziché girare una catena lunga, rimuovere i campioni iniziali e controllare la convergenza, giro tante catene corte unbiased e poi avendone tante posso stimare direttamente la varianza, e le posso mediare perché sono unbiased.
Il limite principale che intravedo è che posso solo fare buone stime di funzioni definite punto per punto. Ad esempio se devo fare una buona stima della varianza della distribuzione non saprei come procedere una volta che ho tante catenine senza convergenza. A questo proposito però danno qualche spunto per calcolare i quantili, anche se non lo esplorano.
Comunque mi ha genuinamente stupito che si possa azzerare il bias di un MCMC. Sono indeciso se considerare ciò come evidenza della provvidenza divina o come evidenza che i metodi frequentisti non hanno davvero senso.
Ensemble samplers with affine invariance
Sono arrivato a questo articolo passando da emcee. Un MCMC invariante per trasformazioni affini, costruito con tante catene che vengono aggiornate insieme. Anche questo a sentirlo sembra un po’ magico ma immagino che l’inghippo sia in come scegli lo stato iniziale, perché l’invarianza ovviamente si ha applicando la trasformazione anche allo stato iniziale.
Prediction diversity and selective attention in the wisdom of crowds
Questo articolo è… fuffa. Diciamo che è stato educativo imparare varie nuove fuffate.
La fuffata più bella è questa: che la diversità di un insieme di stime migliora l’accuratezza. Per un po’ lo dicono a parole e non capivo cosa volesse dire, poi tirano fuori la formula $$ \left( \frac 1N \sum_{i=1}^N x_i – X \right)^2 = \frac 1N \sum_{i=1}^N (x_i – X)^2 – \frac 1N \sum_{i=1}^N \left(x_i – \frac 1N \sum_{k=1}^N x_k \right)^2. $$ Insomma la diversità è la varianza campione. E l’idea sarebbe che visto che compare con il segno meno, fa calare l’errore della media… Comunque non sono scemi loro eh, lo tirano fuori per confutarlo.
Ha attirato la mia attenzione perché al primo anno avevo fatto un esperimento di questo tipo in cui gli studenti dovevano stimare la massa di un pennarello da lavagna (l’avevamo rubato in lab o in aula studio? non ricordo). Poi per calcolare la stima facevamo la media geometrica, perché la massa è una quantità positiva, e ricordo che la media geometrica veniva meglio di quella aritmetica, ed ero molto contento di questo risultato.
Qui usano la media, poi fanno dei rapporti e dei fit, e a quanto riportano il primo esperimento di questo tipo è stato fatto nel 1907 stimando la massa di un bue a una fiera agricola e lì usavano la mediana. Quindi immagino che ognuno usi la statistica che gli fa tornare le cose.
Encoding Musical Style with Transformer Autoencoders
Partendo da un paio di dataset di musica per pianoforte trascritta, in cui è anche segnato quali note sono la melodia, allenano un programma a generare musica nello stile del brano gli viene dato in input. Non ci ho capito più di tanto perché non sono pratico di apprendimento automatico, le cose che mi sono sembrate interessanti sono che il test di ascolto con le persone dà buoni risultati, e che riescono a interpolare tra due stili. Cioè: il programma codifica lo stile di qualcosa con un certo pacco di numeri, e facendo l’interpolazione lineare tra i numeri tirati fuori per due stili si ottiene effettivamente una via di mezzo sensata.
Topological density estimation
Dal nome e dall’abstract mi era sembrato più promettente di quello che è, però comunque utile. È un modo per scegliere la banda (la larghezza dei kernel) di una kernel density estimation in 1D. Usa un invariante topologico che in un certo senso conta le mode (non conta i massimi locali, non sarebbe invariante), e alla fine infatti funziona bene su distribuzioni molto multimodali. Ci sono un paio di scelte arbitrarie da fare e comunque è una normale KDE quindi non è invariante per omomorfismi.
Allega il codice all’articolo (bravo!) ma è in MATLAB (buuuh).
Fast Markov Chain Monte Carlo Algorithms via Lie Groups
Anche questo ha un titolo che sa di magia… Infatti è dello stesso autore del precedente. Cerca di costruire sistematicamente algoritmi Markov chain partendo da principi primi.
C’è una cosa che non mi torna: quando riesce a ricavare Metropolis dal suo macchinario di lemmi, non gli viene il requisito che il proposal sia simmetrico. Cioè, sembra che faccia Metropolis-Hastings però poi usa l’accettanza di Metropolis. Quel requisito serve ad avere il bilancio dettagliato; qui lui non richiede il bilancio dettagliato, però comunque immagino che usare bovinamente l’accettanza di Metropolis non faccia funzionare le cose, altrimenti si saprebbe che funziona lo stesso. Ma in generale qui lui non chiede mai nemmeno che la catena di Markov sia regolare, quindi forse il suo approccio è che il lettore deve arrangiarsi per avere la regolarità.
Guardando le figure con i test mi sembra che alla fine, tenendo conto del costo computazionale, sia meglio fare Metropolis-Hastings normale anziché i nuovi algoritmi che ricava in cui vengono proposti più stati a ogni passaggio. Considerando che usa proposal uniformi, suona in effetti sensato che proporre $d$ stati a ogni passaggio sia più o meno simile a proporre uno stato per $d$ passaggi. Immagino che si potrebbero trovare miglioramenti solo con proposal non banali che però dipendono dal problema specifico.
Gradient-based Adaptive Markov Chain Monte Carlo
Tirano fuori un algoritmo per adattare automaticamente il proposal di un Metropolis-Hastings, e lo testano con proposal gaussiani in cui tutta la matrice di covarianza è libera. Il criterio per adattare il proposal è massimizzare una combinazione di accettanza media ed entropia del proposal.
Fanno vari confronti e il loro algoritmo è più veloce di quelli comuni, in termini di campioni efficaci diviso tempo di calcolo. In particolare è più veloce di NUTS nella maggior parte dei test. Tuttavia proprio nel test con dimensionalità maggiore è peggio di NUTS (ma comunque molto meglio degli altri). Rispetto a NUTS è meno automatico, in particolare i numerelli che usano per gli iperparametri io non saprei proprio dire perché vanno scelti in quel modo, quindi se al tempo di calcolo ci aggiungo il tempo che serve a me per capire le cose direi che NUTS è ancora in vantaggio.
Il codice lo danno ma come al solito è in MATLAB!
Hamiltonian Monte Carlo Swindles
Applicano a HMC tecniche note per ridurre la varianza degli stimatori MCMC. Ancor prima di leggere in dettaglio, la domanda che mi è venuta in mente è: non è che poi c’è un bias che non so stimare? A tal proposito la Figura 1 mi sembra sospetta: nel grafico a destra, dove c’è la versione a varianza ridotta, i campioni stanno quasi sempre sopra il valore vero, e c’è un campione che esplode e va molto più lontano degli altri (sempre verso l’alto). Nel confronto con gli algoritmi non “migliorati” evidenziano le riduzioni della varianza, però non mi hanno convinto per bene che non ci siano problemi.
Comunque le due idee sono queste: (1) se la distribuzione è circa pari, produrre due catene anticorrelate e fare la media. La varianza della media di due variabili anticorrelate viene più piccola. (2) Costruire una distribuzione che approssima quella desiderata e su cui è più facile stimare i valori attesi, girare due catene correlate una sulla distribuzione target e l’altra sull’approssimante, e fare la differenza tra la catena target e la approssimante centrata sulla stima facile del valore atteso.