Intervalli di confidenza distribuzione asimmetrica
Voglio determinare gli intervalli di confidenza (al 68%) di una distribuzione, di cui non ho la forma analitica, ma un'istogramma dato da un alto numero di realizzazioni della variabile.
Se l'istogramma fosse gaussiano, stimerei la media e la deviazione standard dai dati, quindi l'intervallo di confidenza al 68% sarebbe centrato sulla media campionaria, con raggio pari alla deviazione standard.
Nel mio caso l'istogramma è evidentemente asimmetrico, tanto da non poterlo approssimare con una gaussiana. Come posso fare a determinare un analogo intervallo di confidenza?
Un'idea che ho avuto è questa:
- Eseguo un fit con un polinomio di quarto grado intorno al massimo dell'istogramma, per determinare il picco della distribuzione, cioè il valore più probabile.
- A questo punto conto tutte le realizzazioni della variabile che ricadono a sinistra del picco, quindi determino il punto a sinistra del picco tale che il 68% dei valori di sinistra ricada tra questo punto e il picco.
- Analogamente procedo a destra.
- In questo modo ho determinato i due estremi dell'intervallo di confidenza al 68%.
Concettualmente è corretto?
In pratica il metodo è un po' lungo però, perchè la determinazione del massimo è un po' arbitraria:
- Quanto posso scegliere larghi gli intervalli in cui suddivido l'istogramma?
- Quanti ne prendo per eseguire il fit di quarto grado?
Dovrei mettermi a fare prove per vedere se e quanto queste scelte influiscono nella determinazione del picco, e siccome non si tratta di studiare una distribuzione, bensì almeno alcune centinaia, non è una prospettiva allettante.
Esistono metodi più automatici per stimare gli intervalli di confidenza da una distribuzione campionaria? Per esempio pensavo a qualche formula che includesse i momenti della distribuzione fino a un certo ordine.
Se l'istogramma fosse gaussiano, stimerei la media e la deviazione standard dai dati, quindi l'intervallo di confidenza al 68% sarebbe centrato sulla media campionaria, con raggio pari alla deviazione standard.
Nel mio caso l'istogramma è evidentemente asimmetrico, tanto da non poterlo approssimare con una gaussiana. Come posso fare a determinare un analogo intervallo di confidenza?
Un'idea che ho avuto è questa:
- Eseguo un fit con un polinomio di quarto grado intorno al massimo dell'istogramma, per determinare il picco della distribuzione, cioè il valore più probabile.
- A questo punto conto tutte le realizzazioni della variabile che ricadono a sinistra del picco, quindi determino il punto a sinistra del picco tale che il 68% dei valori di sinistra ricada tra questo punto e il picco.
- Analogamente procedo a destra.
- In questo modo ho determinato i due estremi dell'intervallo di confidenza al 68%.
Concettualmente è corretto?
In pratica il metodo è un po' lungo però, perchè la determinazione del massimo è un po' arbitraria:
- Quanto posso scegliere larghi gli intervalli in cui suddivido l'istogramma?
- Quanti ne prendo per eseguire il fit di quarto grado?
Dovrei mettermi a fare prove per vedere se e quanto queste scelte influiscono nella determinazione del picco, e siccome non si tratta di studiare una distribuzione, bensì almeno alcune centinaia, non è una prospettiva allettante.
Esistono metodi più automatici per stimare gli intervalli di confidenza da una distribuzione campionaria? Per esempio pensavo a qualche formula che includesse i momenti della distribuzione fino a un certo ordine.
Risposte
"Sergio":
Se invece hai solo dei dati, se hai solo un alto numero di realizzazioni di una variabile aleatoria e nessuna ipotesi sulla sua legge di probabilità, non puoi parlare di intervalli di confidenza -- che, ripeto, si riferiscono a parametri che sono numeri, non variabili aleatorie. [1]
----------------------------
[1] Le cose cambiano nell'impostazione bayesiana, in cui i parametri sono visti come variabili aleatorie, ma anche in quel caso si muove sempre da un modello statistico.
Le distribuzioni di cui parlo sono proprio le posterior di parametri, che stimo attraverso le realizzazioni di una catena MonteCarlo.
Penso quindi che la risposta "calcola il quantile centrale" sia quella giusta. Il problema è come fare a calcolarlo? (a parte usando una funzione già pronta?) E come calcolare il valore più probabile?
In ogni caso, grazie per la risposta precedente.
Ripensandoci (a brevissima distanza di tempo) un'idea ce l'ho, ma che si può applicare solo se ho un modello della distribuzione:
Si tratta di trovare un intervallo sull'asse $x$ dell'istogramma, che:
- contenga il 68% dell'area della distribuzione;
- i valori della distribuzione agli estremi dell'intervallo (y) coincidono.
Fare questo senza un modello della distribuzione chiaramente diventa più complicato, ma almeno in principio è questo che dovrei fare?
Si tratta di trovare un intervallo sull'asse $x$ dell'istogramma, che:
- contenga il 68% dell'area della distribuzione;
- i valori della distribuzione agli estremi dell'intervallo (y) coincidono.
Fare questo senza un modello della distribuzione chiaramente diventa più complicato, ma almeno in principio è questo che dovrei fare?
"Sergio":
Semplice: metti i valori in ordine crescente, ne togli il 16% a sinistra, il 16% a destra, e ti rimane il 68% centrale.
Sai cosa non mi convince di questo metodo?
Il fatto è che in questo modo posso escludere dal mio intervallo a sinistra dei valori cui corrisponde una densità di probabilità più alta rispetto ad altri che includo a destra (pensa al caso in cui il picco della distribuzione è molto spostato verso sinistra).
Invece, come proponevo nel post precedente, prendo un intervallo che contiene tutti i valori al di sopra di una certa densità di probabilità, ed esclude tutti quelli al di sotto, e tale che complessivamente includa il 68% dei punti.
"Sergio":
Diciamo che sono un po' condizionato dal fatto che per me "catena MonteCarlo" vuol dire "statistica bayesiana". Forse anche per te, visto che parli di "posterior",
E' proprio la moda che intendo calcolare...
Ho tentato di creare un codice per determinare l'intervallo di credibilità al 68%, centrato in modo che i due estremi abbiano la stessa densità di probabilità.
Lo schema è il seguente:
1. Creo un istogramma, a passo costante, della catena di realizzazioni;
2. Inizializzo un valore di soglia a 0, quindi sommo i conteggi di tutti i canali che superano il valore di soglia;
3. Se la somma dei conteggi supera il 68% del totale (ciò ovviamente al primo step, visto che sarà il 100% del totale), aumento la soglia di 1;
4. Ripeto l'operazione 3, fino a quando la somma dei conteggi sarà minore o uguale al 68% del totale.
Alla fine, il primo canale oltre la soglia mi dà una stima del limite sinistro dell'intervallo che cerco, mentre l'ultimo mi dà una stima del limite destro dell'intervallo. Come incertezza attribuisco la larghezza di un canale.
Per testare il codice ho generato una catena fittizia, di $n=200000$ elementi, secondo una distribuzione normale centrata in 0 e di larghezza 1.
Primo problema: Qual è il numero ottimale di canali da usare?
Secondo problema: Perchè non stimo bene il limite destro?
Spero abbi la pazienza di leggere tutto quanto, e possibilmente suggerire miglioramenti della strategia.
Lo schema è il seguente:
1. Creo un istogramma, a passo costante, della catena di realizzazioni;
2. Inizializzo un valore di soglia a 0, quindi sommo i conteggi di tutti i canali che superano il valore di soglia;
3. Se la somma dei conteggi supera il 68% del totale (ciò ovviamente al primo step, visto che sarà il 100% del totale), aumento la soglia di 1;
4. Ripeto l'operazione 3, fino a quando la somma dei conteggi sarà minore o uguale al 68% del totale.
Alla fine, il primo canale oltre la soglia mi dà una stima del limite sinistro dell'intervallo che cerco, mentre l'ultimo mi dà una stima del limite destro dell'intervallo. Come incertezza attribuisco la larghezza di un canale.
Per testare il codice ho generato una catena fittizia, di $n=200000$ elementi, secondo una distribuzione normale centrata in 0 e di larghezza 1.
Primo problema: Qual è il numero ottimale di canali da usare?
Secondo problema: Perchè non stimo bene il limite destro?
Spero abbi la pazienza di leggere tutto quanto, e possibilmente suggerire miglioramenti della strategia.
Ho visto che forse non avevo fatto una scelta ottimale per includere o meno i canali agli estremi della zona di confidenza. Correggendo, ho ottenuto una stima migliore anche a destra ($0,96 pm 0,03$), ma stranamente sembra essere un po' più critica la scelta del numero di canali in cui suddivido l'istogramma.
Volendo aiutarsi analiticamente, io mi aspetterei delle distribuzioni a posteriori gaussiane, essendo le prior uniformi, e la funzione di verosimiglianza gaussiana.
Tuttavia c'è un motivo per cui vengono deformate:
Uno dei tre parametri è un angolo, il cui valore teorico è $90^o$, che cadrebbe al limite dell'intervallo possibile (tra $0$ e $90$). Siccome la posterior associata è per forza di cose asimmetrica, io ho ampliato il range oltre $90^o$, in modo da ottenere una distribuzione simmetrica, riuscendoci.
Tuttavia le altre variabili sono correlate al seno di quest'angolo, e ovviamente il trucco che ho usato non cambia il fatto che la distribuzione del seno non può essere gaussiana, essendo piccata in $1$, questo condiziona le distribuzioni delle altre variabili.
Non so come affrontare la cosa analiticamente però. Per il seno dell'angolo mi aspetterei una mezza gaussiana, in quanto prodotto di una gaussiana per un gradino. Per le altre variabili però non so come comportarmi, perchè oltre alla prior e alla verosimiglianza entrano in gioco le correlazioni tra variabili (volendo vale anche per l'angolo il discorso, ma in prima approssimazione ignorerei questo fatto).
Volendo aiutarsi analiticamente, io mi aspetterei delle distribuzioni a posteriori gaussiane, essendo le prior uniformi, e la funzione di verosimiglianza gaussiana.
Tuttavia c'è un motivo per cui vengono deformate:
Uno dei tre parametri è un angolo, il cui valore teorico è $90^o$, che cadrebbe al limite dell'intervallo possibile (tra $0$ e $90$). Siccome la posterior associata è per forza di cose asimmetrica, io ho ampliato il range oltre $90^o$, in modo da ottenere una distribuzione simmetrica, riuscendoci.
Tuttavia le altre variabili sono correlate al seno di quest'angolo, e ovviamente il trucco che ho usato non cambia il fatto che la distribuzione del seno non può essere gaussiana, essendo piccata in $1$, questo condiziona le distribuzioni delle altre variabili.
Non so come affrontare la cosa analiticamente però. Per il seno dell'angolo mi aspetterei una mezza gaussiana, in quanto prodotto di una gaussiana per un gradino. Per le altre variabili però non so come comportarmi, perchè oltre alla prior e alla verosimiglianza entrano in gioco le correlazioni tra variabili (volendo vale anche per l'angolo il discorso, ma in prima approssimazione ignorerei questo fatto).
Riporto una soluzione numerica al problema molto più elegante e agile della precedente, perchè non ha bisogno di utilizzare l'istogramma.
E' facile convincersi che cercare l'intervallo di credibilità al 68%, determinato come intervallo in cui la densità di probabilità è superiore a una certa soglia, è equivalente a cercare il più stretto tra gli intervalli di credibilità al 68%.
Dunque un possibile algoritmo è:
- ordinare i valori della catena in ordine crescente;
- considerare tutti le possibili successioni di elementi consecutivi pari al 68% del totale (siano $m$ elementi);
- misurare per ognuna la differenza tra ultimo e primo valore;
- individuare il minimo tra queste differenze;
- i corrispondenti estremi della successione sono gli estremi dell'intervallo di credibilità cercato.
Se si vuole stimare anche l'incertezza su questi valori, si può rieseguire il codice con $m+sqrt(m)$ ed $m-sqrt(m)$ elementi.
Testando il mio codice su 20 catene di 150000 elementi estratte secondo una distribuzione normale, ho trovato che gli estremi dell'intervallo di credibilità vengono abbastanza vicine a quelle vere, almeno al livello dei centesimi.
E' facile convincersi che cercare l'intervallo di credibilità al 68%, determinato come intervallo in cui la densità di probabilità è superiore a una certa soglia, è equivalente a cercare il più stretto tra gli intervalli di credibilità al 68%.
Dunque un possibile algoritmo è:
- ordinare i valori della catena in ordine crescente;
- considerare tutti le possibili successioni di elementi consecutivi pari al 68% del totale (siano $m$ elementi);
- misurare per ognuna la differenza tra ultimo e primo valore;
- individuare il minimo tra queste differenze;
- i corrispondenti estremi della successione sono gli estremi dell'intervallo di credibilità cercato.
Se si vuole stimare anche l'incertezza su questi valori, si può rieseguire il codice con $m+sqrt(m)$ ed $m-sqrt(m)$ elementi.
Testando il mio codice su 20 catene di 150000 elementi estratte secondo una distribuzione normale, ho trovato che gli estremi dell'intervallo di credibilità vengono abbastanza vicine a quelle vere, almeno al livello dei centesimi.
Accedi a tutti gli appunti
Tutor AI: studia meglio e in meno tempo