Gibbs Sampling using Anti-correlation Gaussian Data Augmentation, with Applications to L1-ball-type Models

Riepilogo Tecnico: Campionamento Gibbs utilizzando l'Arricchimento dei Dati con Gaussiana Anti-correlata

Enunciato del Problema
Il lavoro affronta le sfide computazionali nell'inferenza bayesiana per modelli che impiegano priori di tipo sfera L1. Questi priori generalizzano i modelli "spike-and-slab" trasformando una distribuzione precursore continua (ad esempio, una Gaussiana o un processo Gaussiano) sul bordo di una sfera L1, inducendo zeri esatti con probabilità positiva. Questo quadro supporta la sparsità strutturata, dove le occorrenze di zeri possono essere dipendenti (ad esempio, gruppi temporali o spaziali) e gli elementi non nulli possono esibire correlazione (ad esempio, regolarità).

Gli algoritmi esistenti per tali modelli affrontano ostacoli significativi:

1. Mancanza di Coniugazione: I campionatori Gibbs standard faticano perché la trasformazione di soft-thresholding rompe la coniugazione e le varianze Gaussiane latenti sono spesso non fisse.
2. Lenta Miscelazione: L'aggiornamento dei parametri elemento per elemento in presenza di correlazioni porta a una lenta miscelazione della catena di Markov.
3. Costo Computazionale: I metodi basati sul gradiente come Hamiltonian Monte Carlo (HMC) o il No-U-Turn Sampler (NUTS) possono gestire la non differenziabilità della funzione di soft-thresholding, ma soffrono di elevati costi computazionali per iterazione a causa delle multiple valutazioni del gradiente.

Metodologia
Gli autori propongono un Campionatore Gibbs a Blocchi che utilizza una tecnica di arricchimento dei dati innovativa denominata "Gaussiana Anti-correlata".

1. Meccanismo Principale: Il metodo introduce variabili latenti, $r$ e $t$, estratte da distribuzioni Gaussiane specificamente progettate per annullare i termini di correlazione quadratica ($\theta'M\theta$ e $\beta'H\beta$) nell'esponente della densità a posteriori.

   • Nello specifico, se $M$ è la matrice di precisione per $\theta$, una variabile latente $r$ viene campionata da $N((dI - M)\theta, dI - M)$, dove $d > \lambda_p(M)$.
   • Questa trasformazione rende la distribuzione a posteriori condizionata dei parametri $(\beta, \theta)$ condizionalmente indipendente tra le dimensioni date le variabili latenti.

2. Schema di Campionamento:

   • Aggiornamento a Blocchi: Invece di aggiornare un elemento alla volta, l'algoritmo aggiorna l'intero vettore $\beta$ (e di conseguenza $\theta$) in un singolo blocco.
   • Struttura di Mixture: La distribuzione a posteriori condizionata per ogni elemento $\beta_j$ diventa una mixture a tre componenti (corrispondenti a $\theta_j > 0$, $\theta_j = 0$, o $\theta_j < 0$). L'algoritmo campiona un indicatore discreto $b_j$ e successivamente campiona $\beta_j$ da una corrispondente distribuzione Gaussiana troncata.
   • Campionamento Efficiente delle Variabili Latenti: Per contesti di regressione in cui la covarianza coinvolge un prodotto $X'\Omega X$, gli autori derivano un algoritmo non iterativo per campionare la Gaussiana anti-correlata. Questo metodo utilizza una decomposizione ai valori singolari (SVD) di $X$ pre-calcolata ed evita decomposizioni matriciali di complessità $O(p^3)$ durante le iterazioni MCMC, raggiungendo una complessità di $O(\max(n, p))$ per passo.

3. Estensioni: Il quadro è esteso a modelli Gaussiani latenti generali (ad esempio, quelli con verosimiglianze non Gaussiane) e al campionamento da Gaussiane multivariate troncate con vincoli separabili.

Contributi Chiave

• Nuovo Arricchimento dei Dati: Introduzione della variabile latente "Gaussiana anti-correlata", che annulla i termini di interazione quadratica per abilitare aggiornamenti a blocchi.
• Efficienza Algoritmica: Sviluppo di una procedura di campionamento rapida e non iterativa per le variabili latenti nei contesti di regressione, evitando costose inversioni matriciali ad ogni passo.
• Garanzia Teorica: Dimostrazione dell'ergodicità geometrica per il campionatore proposto nei modelli lineari, garantendo una rapida convergenza della catena di Markov.
• Ampia Applicabilità: Dimostrazione che il metodo si applica oltre le priori di tipo sfera L1, estendendosi a modelli Gaussiani latenti generali e Gaussiane troncate.

Risultati
Gli autori validano il metodo attraverso esperimenti numerici:

• Regressione Lineare: In simulazioni con predittori ad alta dimensionalità ($p$ fino a 5000) e disegni correlati, il campionatore Gibbs a blocchi anti-correlato supera significativamente NUTS, i campionatori a slice per componente e i campionatori a slice basati su componenti principali empiriche in termini di Dimensione del Campione Effettivo per unità di tempo (ESS/tempo). Sebbene NUTS misceli più velocemente per iterazione, il suo elevato costo computazionale per iterazione risulta in una minore efficienza complessiva.
• Selezione delle Variabili: Il metodo ottiene tassi di falsi positivi (FPR) e tassi di falsi negativi (FNR) competitivi rispetto a NUTS, specialmente in contesti ad alta correlazione dove altri campionatori non convergono o miscelano male.
• Processo Gaussiano con Soft-Thresholding (Lisciatura delle Immagini): Applicato a dati fMRI, il campionatore proposto ha completato 10.000 iterazioni in circa 270 minuti con un elevato ESS/ora. Al contrario, un'implementazione in Stan di NUTS ha richiesto 68 ore per lo stesso numero di iterazioni, con un ESS/ora significativamente inferiore. Le stime a posteriori risultanti erano lisce e sparse, mentre il soft-thresholding diretto dei dati grezzi produceva risultati rumorosi e non lisci.

Significato e Affermazioni
Il lavoro afferma che il metodo proposto offre un miglioramento sostanziale nell'efficienza computazionale per la stima a posteriori nei modelli di sparsità strutturata. Sfruttando l'arricchimento dei dati per raggiungere l'indipendenza condizionata, gli autori abilitano aggiornamenti a blocchi che riducono drasticamente il costo computazionale per iterazione mantenendo una rapida miscelazione.

Gli autori collocano questo lavoro come soluzione al compromesso tra campionatori marginali (miscelazione rapida, costo per iterazione lento) e campionatori con arricchimento dei dati (costo per iterazione rapido, potenzialmente miscelazione più lenta). Sostengono che il loro specifico arricchimento raggiunge il miglior equilibrio, fornendo un'alternativa a basso costo e ad alta efficienza rispetto ai metodi basati sul gradiente come NUTS, specialmente in contesti ad alta dimensionalità o quando si affrontano strutture latenti complesse come i processi Gaussiani. La garanzia di ergodicità geometrica fornisce inoltre una base teorica per l'affidabilità del campionatore nei modelli lineari.