兩個樣本的 Kullback-Leibler 散度
我試圖對兩個樣本的 Kullback-Leibler 散度進行數值估計。要調試實現，請從兩個正態分佈中抽取樣本和.
對於一個簡單的估計，我生成了兩個直方圖並嘗試在數值上逼近積分。我一直在處理直方圖的那些部分，其中一個直方圖的 bin 為零，這樣我要么最終除以零，要么對數為零。我該如何處理這個問題？
我想到了一個相關的問題：如何準確計算兩個不同均勻分佈之間的 KL-Divergence？我是否必須將積分限制為兩個發行版支持的聯合？
Kullback-Leibler 散度定義為
所以為了從經驗數據中計算（估計）這個，我們可能需要對密度函數的一些估計. 所以一個自然的起點可以是通過密度估計（然後，只是數值積分）。我不知道這種方法有多好或多穩定。 但首先是你的第二個問題，然後我會回到第一個問題。讓我們說和是均勻的密度和分別。然後儘管更難定義，但給出它的唯一合理值是，據我所知，因為它涉及整合我們可以選擇將其解釋為. 從我在 Kullback-Leibler (KL) Divergence的直覺中給出的解釋來看，這個結果是合理的
回到主要問題。它以非常非參數的方式提出，並且沒有對密度進行任何假設。可能需要一些假設。但是假設這兩種密度被提議作為同一現象的競爭模型，我們可能會假設它們具有相同的主導度量：例如，連續和離散概率分佈之間的 KL 散度總是無窮大。一篇解決這個問題的論文如下： https ://pdfs.semanticscholar.org/1fbd/31b690e078ce938f73f14462fceadc2748bf.pdf 他們提出了一種不需要初步密度估計的方法，並分析了它的特性。
（還有很多其他論文）。我會回來發布那篇論文的一些細節，這些想法。
EDIT該論文的一些想法是關於使用來自絕對連續分佈的 iid 樣本估計 KL 散度。我展示了他們對一維分佈的建議，但他們也給出了向量的解決方案（使用最近鄰密度估計）。如需校樣，請閱讀論文！
他們建議使用經驗分佈函數的版本，但在樣本點之間進行線性插值以獲得連續版本。他們定義
在哪裡是 Heavyside 階躍函數，但定義為. 然後該函數線性插值（並水平擴展超出範圍）是(為連續）。然後他們建議通過以下方式估計 Kullback-Leibler 散度
在哪裡和是一個小於樣本最小間距的數字。 我們需要的經驗分佈函數版本的 R 代碼是
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 }請注意，
rle它用於處理x.然後KL散度的估計由下式給出
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 }然後我展示了一個小模擬：
KL <- replicate(1000, {x <- rnorm(100) y <- rt(100, df=5) KL_est(x, y)}) hist(KL, prob=TRUE)它給出了以下直方圖，顯示了該估計器的採樣分佈的（估計）：
為了比較，我們在這個例子中通過數值積分來計算 KL 散度：
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嗯…差異足夠大，以至於這裡有很多需要調查的地方！
