On combining estimated and analytic covariance matrices

Titolo: Combinazione di Matrici di Covarianza Stimata e Analitica

Data: 22 aprile 2026 (Preprint)
Autori: Alan Heavens, Lorne Whiteway, Elena Sellentin

1. Il Problema

L'analisi statistica dei dati cosmologici assume tradizionalmente una distribuzione di campionamento gaussiana per le statistiche di sintesi (es. spettro di potenza o funzione di correlazione). Tuttavia, questo approccio presenta due criticità fondamentali quando si utilizzano simulazioni per stimare la matrice di covarianza:

1. Incertezza sulla Covarianza: Quando la matrice di covarianza vera ($C$) è sconosciuta e viene stimata tramite un numero finito di simulazioni ($N_{sim}$), la distribuzione di verosimiglianza (likelihood) dei dati non è più gaussiana, ma segue una distribuzione Student-t multivariata. Questo deriva dalla marginalizzazione analitica sulla distribuzione a posteriori di Wishart inversa della vera covarianza.
2. Contributi Gaussiani Aggiuntivi: In molti scenari cosmologici reali (es. lensing debole, CMB), oltre alla covarianza stimata dalle simulazioni, esistono altre fonti di errore puramente gaussiane e analitiche (es. rumore strumentale, foreground, incertezze degli emulatori, o covarianza "super-sample").

Il problema centrale affrontato in questo lavoro è la convoluzione di una distribuzione Student-t (derivante dall'incertezza sulla covarianza stimata) con una distribuzione Gaussiana (derivante da errori analitici). Tale convoluzione non ha una forma chiusa semplice, rendendo difficile l'inferenza statistica. Le approcci attuali, come la correzione di Hartlap applicata a una likelihood gaussiana, falliscono nel riprodurre le "code pesanti" (heavy tails) della distribuzione reale, portando a stime di probabilità inaccurate, specialmente in presenza di tensioni tra dataset.

2. Metodologia

Gli autori propongono un approccio basato sull'approssimazione della convoluzione risultante con un'altra distribuzione Student-t multivariata, i cui parametri sono determinati tramite il matching dei momenti (covarianza e curtosi).

A. Derivazione dei Parametri Effettivi

Sia $X \sim t_\nu(0, \Sigma)$ la variabile Student-t (dove $\nu$ sono i gradi di libertà e $\Sigma$ la matrice di scala) e $Y \sim N(0, \Sigma_g)$ la variabile gaussiana indipendente. La somma $Z = X + Y$ viene approssimata come $Z \dot{\sim} t_{\nu^*}(0, \Sigma^*)$.

I parametri $\nu^*$ (gradi di libertà effettivi) e $\Sigma^*$ (matrice di scala effettiva) sono calcolati uguagliando:

1. La Covarianza Totale:
   $$\Sigma_Z = \Sigma_g + \left(\frac{\nu}{\nu-2}\right)\Sigma$$
   (valido per $\nu > 2$).
2. La Curtosi Multivariata di Mardia:
   La curtosi della distribuzione convoluta viene calcolata e uguagliata a quella di una Student-t con gradi di libertà $\nu^*$:
   $$\text{Kurt}(t_{\nu^*}) = \frac{p(p+2)(\nu^*-2)}{\nu^*-4}$$
   dove $p$ è la dimensionalità dei dati.

Risolvendo il sistema, si ottengono le formule per i nuovi parametri (Eq. 13 nel testo):
$$\nu^* = 4 + \frac{2p(p+2)}{\kappa}, \quad \Sigma^* = \left(\frac{\nu^*-2}{\nu^*}\right)\Sigma_Z$$
dove $\kappa$ è un termine che dipende dalla traccia e dal quadrato della traccia della matrice $N = \Sigma_Z^{-1}\Sigma$.

B. Campionamento Esatto (Alternativa)

Come alternativa all'approssimazione, gli autori descrivono un algoritmo di campionamento esatto basato sulla rappresentazione in miscela della distribuzione Student-t. La Student-t può essere vista come una miscela di Gaussiane con varianza scalata da una variabile casuale $\lambda$ distribuita Gamma. Questo permette di campionare la convoluzione esatta senza approssimazioni, sebbene sia computazionalmente più costoso.

3. Contributi Chiave

1. Formule Pratiche per l'Inferenza: Derivazione di espressioni analitiche compatte per i parametri della likelihood Student-t approssimata quando si combinano covarianze stimate e analitiche. Questo permette un'implementazione diretta ("drop-in replacement") nei pipeline di analisi cosmologica esistenti.
2. Generalizzazione del Framework: Estensione del lavoro di Sellentin & Heavens (2016) e Percival et al. (2022) per includere esplicitamente termini di rumore gaussiano aggiuntivi, un caso d'uso molto comune ma precedentemente privo di una soluzione analitica semplice.
3. Miglioramento della Copertura Statistica: Dimostrazione che l'approssimazione Student-t mantiene le code pesanti necessarie per una corretta copertura frequentista, a differenza delle approssimazioni gaussiane che sottostimano le probabilità nelle code.
4. Implementazione Software: Fornitura di un codice Python (Appendice D) per il calcolo della log-verosimiglianza basata sul matching dei momenti.

4. Risultati Sperimentali

Gli autori hanno testato il metodo attraverso esperimenti numerici e casi di studio cosmologici realistici:

• Dati Sintetici: In esperimenti di copertura (coverage tests), l'approssimazione Student-t ha mostrato una copertura P-P (Probability-Probability) molto più vicina all'ideale rispetto alle approssimazioni gaussiane (sia naive che con matching della covarianza), specialmente quando il contributo gaussiano non era dominante.
• Lensing Debole Tomografico: Applicazione a survey di shear cosmico (piccola e grande area) con covarianza super-sample analitica.
  • Le approssimazioni gaussiane (anche con covarianza totale corretta) hanno sottostimato significativamente gli errori e prodotto regioni di credibilità troppo strette.
  • L'approssimazione Student-t proposta ha fornito risultati quasi identici al campionamento esatto tramite modello gerarchico (BHM), confermando la sua accuratezza nel preservare la forma della distribuzione a posteriori.
• Robustezza: Il metodo funziona bene anche quando il numero di simulazioni è solo leggermente superiore al numero di punti dati, un regime in cui le correzioni gaussiane falliscono.

5. Significato e Implicazioni

Questo lavoro risolve un problema pratico critico nell'analisi dei dati cosmologici di nuova generazione (es. Euclid, LSST, CMB-S4).

• Accuratezza nelle Code: Fornisce un metodo per gestire correttamente l'incertezza sulle code della distribuzione, fondamentale per valutare la significatività statistica di eventuali tensioni tra diversi dataset cosmologici.
• Efficienza Computazionale: Offre una soluzione analitica che evita la necessità di costosi campionamenti gerarchici completi, permettendo di utilizzare likelihood standard con parametri modificati ($\nu^*$ e $\Sigma^*$).
• Versatilità: Il metodo è applicabile a qualsiasi scenario in cui si combinano incertezze stimate (simulazioni) e teoriche (rumore, emulatori), rendendolo uno strumento essenziale per l'inferenza bayesiana e frequentista in cosmologia.

In sintesi, il documento fornisce un ponte teorico e pratico tra l'uso di simulazioni per la covarianza e la necessità di includere errori analitici, garantendo che le inferenze cosmologiche rimangano statisticamente valide anche in regimi di alta precisione.