Discrepanza tra binomiale e beta in R?

Sep 14 2020

Sto ottenendo un risultato che non posso spiegare quando si utilizza la distribuzione beta.

Ho un risultato derivante da una distribuzione binomiale: 2 successi in 6 prove. Penso che lo stimatore di massima verosimiglianza per p sia 2/6 = 0,33?

dbinom(0:6, 6, 0.33)
[1] 0.090458382 0.267324771 0.329168562 0.216170399 0.079853991 0.015732428 0.001291468

Ma, quando uso la distribuzione beta, il punto più alto che ottengo è 0,25:

beta_df <- data.frame(PROB = seq(0, 1, 0.01), HEIGHT = dbeta(seq(0, 1, 0.01), 2, 4))
beta_df[which.max(beta_df$HEIGHT),] beta_df[which.max(beta_df$HEIGHT),]
   PROB   HEIGHT
26 0.25 2.109375

Non riesco a capire questo ... sto interpretando male i risultati o chiamando una di queste funzioni in modo errato? Grazie :)

Risposte

4 Tim Sep 14 2020 at 15:12

Perché ti aspetteresti di vedere risultati simili? Quelle sono distribuzioni diverse, usate per modellare cose completamente diverse. La prima è una distribuzione discreta, la seconda è una distribuzione continua. Rispondendo alla tua domanda, stai cercando la modalità delle distribuzioni. La modalità di distribuzione beta è$\frac{\alpha - 1}{\alpha + \beta - 2}$, così esattamente $0.25$ in caso dei valori forniti.

Per quanto riguarda il commento, nel caso binomiale stavi massimizzando solo la probabilità. Quando si utilizza il modello beta-binomiale bayesiano , quando lo si massimizza, si considera anche il precedente

$$ \hat p = \operatorname{arg\,max} \; \underbrace{p(X|\theta)}_\text{likelihood}\,\underbrace{p(\theta)}_\text{prior} $$

quindi la scelta del priore influenzerebbe il risultato. Quando si usa$\alpha=\beta=0$nel precedente, questo è un precedente di Haldane improprio che ha tutta la massa di probabilità sui valori$0$ e $1$(vedi l'immagine qui sotto presa in prestito da questo sito ).

Soprattutto quando la dimensione del campione è piccola , il precedente avrebbe un impatto sul risultato. In questo caso, trascinerà la massa di probabilità verso gli estremi. Per ottenere risultati paragonabili a MLE, potresti scegliere l'uniforme prima con$\alpha=\beta=1$.

2 LiKao Sep 14 2020 at 15:47

Stai utilizzando il parametro sbagliato per la distribuzione beta. Se hai un esperimento binomiale con$n$ successi e $m$ fallimenti, è necessario utilizzare il $beta(n+1,m+1)$distribuzione. Il motivo è che stai usando fondamentalmente un file$beta(1,1)$ (uniforme) prima, che devi aggiungere alla distribuzione (se usi un file $beta(a,b)$ prima invece si ottiene $beta(n+a,m+b)$).

Quindi se lo provi

beta_df <- data.frame(PROB = seq(0, 1, 0.01), HEIGHT = dbeta(seq(0, 1, 0.01), 3, 5))
beta_df[which.max(beta_df$HEIGHT),]

ottieni il risultato corretto, cioè $ 0,33 $ .

EDIT (più matematica):

Quindi la probabilità di $ n $ successi e $ m $ fallimenti è (fino a una costante moltiplicativa):

$ L (\ theta | \, m, n) \ sim \ theta ^ n (1- \ theta) ^ m $

Ma la densità beta è data come

$ f_ {a, b} (\ theta) \ sim \ theta ^ {(a-1)} (1- \ theta) ^ {(b-1)} $ .

Quindi se abbini $ a-1 = n $ e $ b-1 = m $ ottieni $ a = n + 1 $ e $ b = m + 1 $ , come dovresti.