Berechnung der Wahrscheinlichkeit, wenn

8

Ich versuche diese hintere Verteilung zu berechnen:

(θ|)=i=1npiyi(1pi)1yiallθ,pi|θi=1npiyi(1pi)1yi

Das Problem ist, dass der Zähler, der das Produkt einer Reihe von -Wahrscheinlichkeiten ist, zu klein ist. (Mein ist groß, ungefähr 1500).Bernoulli(pi,yi)n

Daher werden die posterioren Werte für all all zu 0 berechnet (ich berechne in R).θ

Zur Verdeutlichung hat jedes sein eigenes , zusammen bilden diese einen Vektor von Elementen für . Jedes hat seinen eigenen Element-Vektor von .yipipinn yθnpi

BEARBEITEN: Hinzufügen eines Wiedergabebeispiels (für den Zähler)

p <- sample(seq(0,1,by=0.01), 1500, replace=T)
y <- sample(c(0,1), 1500, replace=T)
dbern(y, p) # 1500-element vector, each element is < 1
prod(dbern(y, p)) # produce 0
exp(sum(log(dbern(y, p)))) # produce 0 since the sum is very negative
Heisenberg
quelle
Haben Sie stattdessen versucht, die Summe der Protokolle zu berechnen?
Ansari
1
Es gibt ähnliche Diskussion hier . Es gibt einige zusätzliche Erörterungen einiger Details solcher Berechnungen.
Glen_b -Reinstate Monica

Antworten:

7

Dies ist ein häufiges Problem bei der Berechnung von Wahrscheinlichkeiten für alle Arten von Modellen. Die üblichen Aufgaben bestehen darin, Protokolle zu bearbeiten und einen gemeinsamen Skalierungsfaktor zu verwenden, der die Werte in einen vernünftigeren Bereich bringt.

In diesem Fall würde ich vorschlagen:

Schritt 1: Wählen Sie eine ziemlich "typischen" , . Teilen Sie die Formel für Zähler und Nenner des allgemeinen Terms durch den Zähler für , um etwas zu erhalten, bei dem die Wahrscheinlichkeit eines Unterlaufs viel geringer ist.θθ0θ=θ0

Schritt 2: Arbeiten Sie an der Protokollskala. Dies bedeutet, dass der Zähler eine Exp der Summen der Unterschiede der Protokolle ist und der Nenner eine Summe der Exp der Summen der Unterschiede der Protokolle.

NB: Wenn eines Ihrer Ps 0 oder 1 ist, ziehen Sie diese separat heraus und führen Sie keine Protokolle dieser Begriffe. Sie sind einfach zu bewerten!

[Allgemeiner ausgedrückt kann diese Skalierung und das Arbeiten auf der Protokollskala so gesehen werden, dass sie eine Reihe von Protokollwahrscheinlichkeiten und dies tut: . Eine naheliegende Wahl für besteht darin, den größten Term 0 zu machen, wodurch wir : . Beachten Sie, dass Sie bei einem Zähler und einem Nenner für beide dasselbe , das dann abgebrochen wird. Oben entspricht dies der mit der höchsten Log-Wahrscheinlichkeit.]lilog(ieli)=c+log(ielic)clog(ieli)=maxi(li)+log(ielimaxi(li))cθ0

Die üblichen Begriffe im Zähler sind tendenziell moderater, und daher sind in vielen Situationen sowohl der Zähler als auch der Nenner relativ vernünftig.

Wenn der Nenner verschiedene Größen enthält, addieren Sie die kleineren, bevor Sie die größeren addieren.

Wenn nur wenige Begriffe stark dominieren, sollten Sie Ihre Aufmerksamkeit darauf richten, die Berechnung für diese relativ genau durchzuführen.

Glen_b -Reinstate Monica
quelle
Aber für alle Theta geht der Zähler immer auf 0. Wie dividiere ich dann den allgemeinen Term durch den Zähler? (Schritt 1)
Heisenberg
1
Schritt 1 ist Algebra, keine Computerberechnung. Der Zweck ist es, Ihnen in Schritt 2 etwas zu berechnen, das nicht unterläuft. Es sei denn, Sie sagen, es ist immer algebraisch Null. In diesem Fall tun Sie zweifellos etwas, was Sie nicht tun sollten.
Glen_b -State Monica
okay - ich werde es versuchen. Der Zähler ist nicht genau 0, nur sehr klein, was R nicht berechnen kann. Vielen Dank!
Heisenberg
3
Lieber Gott, du hast recht! Vielen herzlichen Dank. Jeder sagt immer wieder "benutze log.likelihood", aber nur du siehst das Problem wirklich.
Heisenberg
1

Versuchen Sie, die Eigenschaften der Verwendung der Logarithmen und der Summierung zu nutzen, anstatt das Produkt aus Dezimalzahlen zu verwenden. Verwenden Sie nach der Summierung einfach das Anti-Log, um es wieder in Ihre natürlichere Form zu bringen. Ich denke, so etwas sollte den Trick machen

exp(in(yilog(pi)+(1yi)log(1pi)))gexp(inyilog(pi)+(1yi)log(1pi))

Philchalmers
quelle
Der Zähler in Ihrem Vorschlag erzeugt immer noch eine 0, da die Summe innerhalb von exp () immer noch sehr negativ ist (<-1000). Mache ich etwas falsch Danke für Ihre Hilfe!
Heisenberg
Nun, wenn ein Wert in p tatsächlich 0 oder 1 ist, erzeugt das Protokoll automatisch -inf und protokolliert (1-p). Ansonsten denke ich, dass die Zahlen einfach zu klein werden, um wieder in die ursprüngliche Form gebracht zu werden.
Philchalmers
2
Beachten Sie, dass Sie jede Konstante den Begriffen in des obigen Ausdrucks hinzufügen und von diesen subtrahieren können, ohne das Ergebnis zu ändern. Die Einstellung von gleich dem Negativ des Maximalwerts von liefert die beste numerische Genauigkeitexp ( ) c log ( p ( θ | - ) )cexp()clog(p(θ|))
Wahrscheinlichkeitslogik