Kullback-Leibler-Divergenz für zwei Proben

Ich habe versucht, eine numerische Schätzung der Kullback-Leibler-Divergenz für zwei Stichproben zu implementieren. Um die Implementierung zu debuggen, ziehen Sie die Stichproben aus zwei Normalverteilungen und N ( 1 , 2 ) .$\mathcal N (0,1)$$\mathcal N (1,2)$

Für eine einfache Schätzung habe ich zwei Histogramme erstellt und versucht, das Integral numerisch zu approximieren. Ich habe mich nicht mehr mit den Teilen des Histogramms befasst, bei denen die Fächer eines der Histogramme Null sind, sodass ich entweder durch Null oder den Logarithmus von Null dividiere. Wie gehe ich mit diesem Problem um?

Eine verwandte Frage kam mir in den Sinn: Wie kann man die KL-Divergenz zwischen zwei verschiedenen Gleichverteilungen genau berechnen? Muss ich das Integral auf die Vereinigung der Unterstützung beider Distributionen beschränken?

Die Kullback-Leibler-Divergenz ist definiert als Um dies aus empirischen Daten zu berechnen (abzuschätzen), benötigen wir möglicherweise einige Schätzungen der Dichtefunktionen p ( x ) , q

$$\DeclareMathOperator{\KL}{KL} \KL(P || Q) = \int_{-\infty}^\infty p(x) \log \frac{p(x)}{q(x)} \; dx$$ . Ein natürlicher Ausgangspunkt könnte also eine Dichteschätzung sein (und danach nur noch eine numerische Integration). Wie gut oder stabil eine solche Methode wäre, weiß ich nicht.$p(x), q(x)$

Aber zuerst Ihre zweite Frage, dann komme ich zur ersten zurück. Nehmen wir an, und q sind einheitliche Dichten für [ 0 , 1 ] bzw. [ 0 , 10 ] . Dann ist KL ( p | | q ) = log 10, während KL ( q | | p ) schwieriger zu definieren ist, aber der einzig vernünftige Wert ist ∞ , soweit ich sehen kann, da es die Integration von $p$$q$$[0,1]$$[0,10]$$\KL(p || q)=\log 10$$\KL(q || p)$$\infty$. Diese Ergebnisse sind aus der Interpretation, die ich inIntuition über die Kullback-Leibler (KL) -Divergenzgebe, vernünftig$\log(1/0)$was wir wählen können, um als log ∞ zu interpretieren$\log \infty$

Zurück zur Hauptfrage. Es wird sehr nichtparametrisch abgefragt, und es werden keine Annahmen über die Dichten gemacht. Wahrscheinlich sind einige Annahmen erforderlich. Unter der Annahme, dass die beiden Dichten als konkurrierende Modelle für dasselbe Phänomen vorgeschlagen werden, können wir wahrscheinlich annehmen, dass sie dasselbe dominierende Maß haben: Die KL-Divergenz zwischen einer kontinuierlichen und einer diskreten Wahrscheinlichkeitsverteilung wäre beispielsweise immer unendlich. Ein Artikel, der sich mit dieser Frage befasst, lautet wie folgt: https://pdfs.semanticscholar.org/1fbd/31b690e078ce938f73f14462fceadc2748bf.pdf Sie schlagen eine Methode vor, für die keine vorläufige Dichteschätzung erforderlich ist, und analysieren ihre Eigenschaften.

(Es gibt viele andere Papiere). Ich werde zurückkommen und einige Details aus diesem Papier veröffentlichen, die Ideen.

 EDIT               

Einige Ideen aus diesem Artikel, in denen es um die Abschätzung der KL-Divergenz mit iid-Proben aus absolut kontinuierlichen Verteilungen geht. Ich zeige ihren Vorschlag für eindimensionale Verteilungen, aber sie geben auch eine Lösung für Vektoren (unter Verwendung der Schätzung der Dichte des nächsten Nachbarn). Für Beweise lesen Sie das Papier!

Sie schlagen vor, eine Version der empirischen Verteilungsfunktion zu verwenden, die jedoch linear zwischen Stichprobenpunkten interpoliert wird, um eine kontinuierliche Version zu erhalten. Sie definieren

$$P_e(x) = \frac1{n}\sum_{i=1}^n U(x-x_i)$$ $U$$U(0)=0.5$$P_c$$c$ $$\hat{D}(P \| Q) = \frac1{n}\sum_{i=1}^n \log\left(\frac{\delta P_c(x_i)}{\delta Q_c(x_i)}\right)$$ $\delta P_c = P_c(x_i)-P_c(x_i-\epsilon)$$\epsilon$ eine Zahl ist, die kleiner als der kleinste Abstand der Proben ist.

R-Code für die Version der empirischen Verteilungsfunktion, die wir benötigen, ist

my.ecdf  <-  function(x)   {
    x   <-   sort(x)
    x.u <-   unique(x)
    n  <-  length(x) 
    x.rle  <-  rle(x)$lengths
    y  <-  (cumsum(x.rle)-0.5) / n
    FUN  <-  approxfun(x.u, y, method="linear", yleft=0, yright=1,
                           rule=2)
    FUN
}          

Beachten Sie, dass rleder Fall mit Duplikaten bearbeitet wirdx .

Dann ist die Schätzung der KL-Divergenz gegeben durch

KL_est  <-  function(x, y)   {
    dx  <-  diff(sort(unique(x)))
    dy  <-  diff(sort(unique(y)))
    ex  <-  min(dx) ; ey  <-  min(dy)
    e   <-  min(ex, ey)/2
    n   <-  length(x)    
    P  <-   my.ecdf(x) ; Q  <-  my.ecdf(y)
    KL  <-  sum( log( (P(x)-P(x-e))/(Q(x)-Q(x-e)))) / n
    KL              
}

Dann zeige ich eine kleine Simulation:

KL  <-  replicate(1000, {x  <-  rnorm(100)
                         y <- rt(100, df=5)
                         KL_est(x, y)})
hist(KL, prob=TRUE)

Dies ergibt das folgende Histogramm, das (eine Schätzung) der Stichprobenverteilung dieses Schätzers zeigt:

Zum Vergleich berechnen wir die KL-Divergenz in diesem Beispiel durch numerische Integration:

LR  <-  function(x) dnorm(x,log=TRUE)-dt(x,5,log=TRUE)
100*integrate(function(x) dnorm(x)*LR(x),lower=-Inf,upper=Inf)$value
[1] 3.337668

hmm ... der Unterschied ist groß genug, dass es hier viel zu untersuchen gibt!

Ich habe die Antwort von kjetil-b-halvorsen ein wenig erweitert und es tut mir leid, dass ich nicht kommentiert habe. Ich habe nicht den Ruf:

1. Ich habe das Gefühl, dass die analytische Berechnung sein sollte (ohne Multiplikation mit 100):

LR <- function(x) dnorm(x,log=TRUE)-dt(x,5,log=TRUE) integrate(function(x) dnorm(x)*LR(x),lower=-Inf,upper=Inf)$value

2. $\hat D(P||Q)$$\hat D(P||Q)-1 \to D(P||Q)$

Sobald diese beiden Korrekturen vorgenommen wurden, erscheinen die Ergebnisse realistischer.