KSG互信息估计器的原理详细推导(Kraskov, 2004)

KSG Estimators

KSG估计器是Kraskov在2004年提出的互信息估计器，其原文过于简略，我参考大量文献对细节进行了补充。原文位置：[https://arxiv.org/pdf/cond-mat/0305641.pdf]

KSG估计的基本方法是首先对互信息进行如下分解 (1)
$$I(X,Y)=H(X)+H(Y)-H(X,Y)$$
因此，只需要对$H(x)$等三项进行逐个估计即可。

KLE Esitmator

KLE估计器用于估计微分熵。

利用蒙特卡洛估计有(2)
$$H[x]=-\int \mu (x)log\mu (x)dx\approx -\frac{1}{N}log\hat{\mu }(x)$$
Then the key problem is estimating $\hat{\mu }$(x).

考虑 一个 $ϵ-ball$ 以 $x_i$为中心, 则 $(0,\frac{ϵ}{2})$的范围内有$k-1$个点，$(\frac{ϵ}{2},\frac{ϵ}{2}+dϵ)$中间有一个点（the kth neighbor），$(\frac{ϵ}{2}+dϵ)$之外有剩余的$N-k-1$​个点，因此首先在除了$x_i$的N-1个点中挑选k个点，其中一个点位于边缘上，这k个点中剩余选取k-1个点位于最中心范围内,则$p(ϵ)$是$x_i$周围$ϵ-ball$的概率质量(mass) (3)：
$$\begin{array}{rl}{p}_i(ϵ)=& \left(\genfrac{}{}{0ex}{}{N-1}{k}\right)\left(\genfrac{}{}{0ex}{}{k}{k-1}\right)(1-p_i(ϵ){)}^{N-k-1}{p}_i(ϵ{)}^{k-1}\frac{d{p}_i(ϵ{)}^1}{dϵ}\\ & =\frac{(N-1)!}{(N-k-1)!1!(k-1)!}(1-p_i(ϵ){)}^{N-k-1}{p}_i(ϵ{)}^{k-1}\frac{d{p}_i(ϵ{)}^1}{dϵ}\end{array}$$
令$u=p_i(ϵ)$，则$du=d{p}_i(ϵ)dϵ$
Beta函数定义为$beta(p,q)={\int }_0^1{x}^{p-1}(1-x{)}^{q-1}dx=\frac{(p-1)!(q-1)!}{(p+q-1)!=\frac{\frac{p+q}{pq}}{\left(\genfrac{}{}{0ex}{}{p+q}{p}\right)}}$
因此，验证积分为1： (4)
$$\begin{array}{rl}{\int }_0^{\infty }{p}_i(ϵ)dϵ& ={\int }_0^{1}k\left(\genfrac{}{}{0ex}{}{N-1}{k}\right){\mu }^{k-1}(1-\mu {)}^{N-k-1}d\mu =k(\genfrac{}{}{0ex}{}{N-1}{k})beta(k,N-k)\\ & =k\left(\genfrac{}{}{0ex}{}{N-1}{k}\right)\frac{N}{k(N-k)}\frac{1}{\left(\genfrac{}{}{0ex}{}{N}{k}\right)}=1\end{array}$$

Lemma 1.

$X\sim Beta(\alpha ,\beta ),\mathbb{E}[log(X)]=\psi (\alpha )-\psi (\alpha +\beta ),\psi ()=\frac{{\Gamma }^{\prime }()}{\Gamma ()}$

证明： (5)
$$\begin{array}{rl}E[log(X)]& ={\int }_0^{1}B(\alpha ,\beta {)}^{-1}log(x)x^{\alpha -1}(1-x{)}^{\beta -1}dx\\ & =B(\alpha ,\beta {)}^{-1}{\int }_0^1\frac{\partial }{\partial \alpha }{x}^{\alpha -1}(1-x{)}^{\beta -1}dx\end{array}$$
交换微分和积分的顺序（类似于莱布尼兹法则）: (6)
$$E[log(X)]=Beta(\alpha ,\beta {)}^{-1}\frac{\partial }{\partial \alpha }Beta(\alpha ,\beta )=\frac{\partial }{\partial \alpha }logB(\alpha ,\beta )$$
展开Beta函数： (7)
$$=\frac{\partial }{\partial \alpha }[log\Gamma (\alpha )+log\Gamma (\beta )-log\Gamma (\alpha +\beta )]=\psi (\alpha )-\psi (\alpha +\beta )$$

Proposition 2.

$E[log{p}_i(ϵ)]=\psi (k)-\psi (N)$

使用Lemma 1可得。

熵的估计

假设在$x_i$的邻域范围内，概率密度$\mu (x)$是均匀的，则有 (8)
$$p_i(ϵ)\approx c_d{ϵ}^d\mu (x_i)$$
其中，$c_d$是n维的单位球体体积，定义为${\int }_{\Vert s\Vert <1/2}ds$,对于某些范数。

如果是2范数，则n维球体积公式 (9)
$$V_n=\frac{{\pi }^{\frac{n}{2}}{R}^d}{\Gamma (1+\frac{n}{2})}$$
$R=1$时即为$c_d$，再乘以$R^d$。

因此，该区域概率密度$\mu (x_i)$乘以体积$c_d{ϵ}^d$，得到该区域概率质量$p_i(ϵ)$

根据（8）,有 (10)
$$\begin{array}{rl}E[-log\mu (X_i)]& \approx -E[log{p}_i(ϵ)]+log{c}_d+dE[logϵ]\\ & =-\varphi (k)+\varphi (N)+log{c}_d+dE[logϵ]\\ & \approx -\varphi (k)+\varphi (N)+log{c}_d+\frac{d}{N}\sum _{i=1}^{N}log{ϵ}_i\end{array}$$
其中，最后的${ϵ}_i$表示第$i$个点到他第kth最近邻的距离。

KSG1估计器

首先，$H(X)$和$H(Y)$可以直接使用（10）进行估计。

其次，使用(10)，联合熵的一个估计是 (11)：
$$\hat{H}(X,Y)=-\psi (k)+\psi (N)+log(c_{dz})+\frac{dz}{N}\sum log{ϵ}_i$$
其中，$dz=dx+dy$，即Z的维度是x和y维度之和。且对于球的体积，KSG采用的最大范数进行的度量，也就是最大范数下N维空间球的体积，因此区域看上去是n维的立方体。

其中可以注意到： (12)
$$c_{dz}=c_{dx}{c}_{dy}$$
我思考时类比于正方体的底面积乘以高，这一点在下面（23）式中使用。

然而，使用固定的k估计的$H(X)H(Y)H(X,Y)$直接对$I(X,Y)$进行估计，存在问题。原因是，（10）式的偏差主要来自于对（8）式假设的均匀概率密度的违反，因此这三项估计都因此会存在偏差。而对于一个固定的k来说，联合空间上的第k近邻的距离要远于边际空间，因此实际上同样是kth近邻，空间的大小不同，因此概率密度违反（8）式假设的大小也不同，因此会导致偏差非常不同从而无法抵消（且观看下图图示，以1近邻为例）。

为了避免这一点，由于（10）对任意k都成立，因此不必使用固定的k。KSG1的方法是首先找到Joint上kth近邻，然后在X marginnal space上确定其距离$ϵ(i)$，在该距离内的点的个数作为$n_x(i)$，则到joint上第kth点的距离正是$ϵ(i)$，而这个点在Xmarginal上则必然是第$n_x(i)+1$个点。对Y marginal space用同样大小的距离$ϵ(i)$做同等操作。具体的情况如下图所示，此时Joint上是1th近邻，而此时x方向上产生5近邻，而y方向上是3近邻。

此时对于$H(X)$来说，则有（考虑（2）式和（10）式中第一个约等号的蒙特卡洛近似，我未必写的非常非常详细，如果这里还需过于详细的介绍，则无需再看此文，而是重新补习统计相关知识）： (13)
$$\hat{H}(X)=-\frac{1}{N}\sum _{i=1}^N\psi (n_x(i)+1)+\psi (N)+log(c_{dz})+\frac{dx}{N}\sum log{ϵ}_i$$
对于Y，则利用同样的距离找到$n_y(i)+1$进行估计，导致Y的估计不准，因为y到其kth点的距离不是$\frac{ϵ}{2}$。

KSG2估计器

KSG2估计器更为复杂，其考虑在两个marginal的方向上采用不同的距离，此时图示joint的kth区域不再是超立方体区域，而是超长方体区域。首先根据joint上kth近邻的位置分别确定x marginal上的距离和y marginal上的距离（已经与KSG1不同，该方法确定x marginal上的距离后，直接将该距离应用于 y marginal上，因此是超正方体），此时点落在$x_i\pm {ϵ}_x(i)/2$和$y_i\pm {ϵ}_y(i)/2$之内。

此时存在两种情况，即${ϵ}_x(i)$和${ϵ}_y(i)$由同一点决定或由两点决定，具体见下图。

故考虑在将原先的超正方体概率替换为超长方体的概率，该超长方体的概率是两种情况的混合，则 (14)
$$P_k({ϵ}_x,{ϵ}_y)=P_k^b({ϵ}_x,{ϵ}_y)+P_k^c({ϵ}_x,{ϵ}_y)$$
其中第一种时，有k-1个点在最内，中间的微分区域（逻辑与（3）相同）有1个点，外部有N-1-k个点在最外，若不失一般性的假设${ϵ}_x<{ϵ}_y$，该种状况发生的概率表达为 (15)
$$\begin{array}{rl}{P}_k^b({ϵ}_x,{ϵ}_y)& =\left(\genfrac{}{}{0ex}{}{N-1}{k}\right)\left(\genfrac{}{}{0ex}{}{k}{k-1}\right)q_i^{k-1}(\frac{{\partial }^2{q}_i}{\partial {ϵ}_x\partial {ϵ}_y})(1-p_i({ϵ}_y){)}^{N-k-1}\\ & =k\left(\genfrac{}{}{0ex}{}{N-1}{k}\right)q_i^{k-1}(\frac{{\partial }^2{q}_i}{\partial {ϵ}_x\partial {ϵ}_y})(1-p_i({ϵ}_y){)}^{N-k-1}\end{array}$$
原文的表述将$q_i$吸进了偏导项中，其中$q_i$表示的是上述超长方体的区域，而$p_i$则仍然是由较大的${ϵ}_y$确定的超正方体区域。使用$p_i$是必须的，因为使用了最大范数，即最远的点的距离，保证了如果该正方形内的点都在矩形中，反之则不行。

对于第二种情况，此时有k-2个点在最内，而2个点在微分区域，N-k-1个点在外围区域，因此此种状况发生的概率为 (16)
$$\begin{array}{rl}{P}_k^c({ϵ}_x,{ϵ}_y)& =\left(\genfrac{}{}{0ex}{}{N-1}{k}\right)\left(\genfrac{}{}{0ex}{}{k}{k-2}\right)q_i^{k-2}(\frac{{\partial }^2{q}_i^2}{\partial {ϵ}_x\partial {ϵ}_y})(1-p_i({ϵ}_y){)}^{N-k-1}\\ & =\left(\genfrac{}{}{0ex}{}{N-1}{k}\right)\left(\genfrac{}{}{0ex}{}{k}{k-2}\right)2q_i{q}_i^{k-2}(\frac{{\partial }^2{q}_i}{\partial {ϵ}_x\partial {ϵ}_y})(1-p_i({ϵ}_y){)}^{N-k-1}\\ & =k(k-1)\left(\genfrac{}{}{0ex}{}{N-1}{k}\right)q_i^{k-1}(\frac{{\partial }^2{q}_i}{\partial {ϵ}_x\partial {ϵ}_y})(1-p_i({ϵ}_y){)}^{N-k-1}\end{array}$$
同时注意到，由于X和Y上是独立的，根据 (17)
$$u(\frac{\partial q}{\partial x\partial y})=(\frac{\partial q}{\partial x})(\frac{\partial q}{\partial y})$$
故（16）进一步拆解为 (18)
$$k(k-1)\left(\genfrac{}{}{0ex}{}{N-1}{k}\right)q_i^{k-2}(\frac{\partial q}{\partial x})(\frac{\partial q}{\partial y})(1-p_i({ϵ}_y){)}^{N-k-1}$$
验证$P_k({ϵ}_x,{ϵ}_y)=P_k^b({ϵ}_x,{ϵ}_y)+P_k^c({ϵ}_x,{ϵ}_y)$的积分为1. (19)
$$\begin{array}{rl}& {\int }_0^{\infty }{\int }_0^{{ϵ}_y}{P}_k({ϵ}_x,{ϵ}_y)d{ϵ}_{x}d{ϵ}_y\\ & =\left(\genfrac{}{}{0ex}{}{N-1}{k}\right){\int }_0^{\infty }(1-p_i({ϵ}_y){)}^{N-k-1}{\int }_0^{{ϵ}_y}\frac{\partial }{\partial x}(k{q}_i^{k-1}(\frac{\partial q_i}{\partial {ϵ}_y}))d{ϵ}_{x}d{ϵ}_y\\ & =k\left(\genfrac{}{}{0ex}{}{N-1}{k}\right){\int }_0^{\infty }(1-p_i({ϵ}_y){)}^{N-k-1}{p}_i({ϵ}_y{)}^{k-1}d{ϵ}_y\\ & =k\left(\genfrac{}{}{0ex}{}{N-1}{k}\right)Beta(k,N-k)=1\end{array}$$
下面证明$E[log{q}_i({ϵ}_x,{ϵ}_y)]=\psi (k)-\psi (N)-1/k$，首先在第二个积分号进行分部积分： (20)
E[logqi(ϵx,ϵy)]=∫0∞∫0ϵylog(qi)Pk(ϵx,ϵy)dϵxdϵy=(N−1k)∫0∞(1−pi(ϵy))N−k−1∫0ϵy∂∂x(kqik−1(∂qi∂ϵy))dϵxdϵy=(N−1k)∫0∞(1−pi(ϵy))N−k−1[log(qi)kqik−1(∂qi∂ϵy)∣0ϵy−∫0ϵy{kqik−2∂qi∂ϵx(∂qi∂ϵy)}dϵx]dϵy \begin{aligned} E[logq_i(\epsilon_x,\epsilon_y)]&=\int_0^{\infin}\int_0^{\epsilon_y}log(q_i)P_k(\epsilon_x,\epsilon_y)d\epsilon_xd\epsilon_y\\ &=\tbinom{N-1}{k}\int_0^{\infin}(1-p_i(\epsilon_y))^{N-k-1}\int_0^{\epsilon_y}\frac{\partial}{\partial x}({kq_i^{k-1}(\frac{\partial q_i}{\partial\epsilon_y}}))d\epsilon_xd\epsilon_y\\ &=\tbinom{N-1}{k}\int_0^{\infin}(1-p_i(\epsilon_y))^{N-k-1}[log(q_i)kq_i^{k-1}(\frac{\partial q_i}{\partial \epsilon_y})|_{0}^{\epsilon_y}\\&-\int_0^{\epsilon_y}\{kq_i^{k-2}\frac{\partial q_i}{\partial \epsilon_x}(\frac{\partial q_i}{\partial \epsilon_y})\}d\epsilon_x]d\epsilon_y\\ \end{aligned} E[logqi​(ϵx​,ϵy​)]​=∫0∞​∫0ϵy​​log(qi​)Pk​(ϵx​,ϵy​)dϵx​dϵy​=(kN−1​)∫0∞​(1−pi​(ϵy​))N−k−1∫0ϵy​​∂x∂​(kqik−1​(∂ϵy​∂qi​​))dϵx​dϵy​=(kN−1​)∫0∞​(1−pi​(ϵy​))N−k−1[log(qi​)kqik−1​(∂ϵy​∂qi​​)∣0ϵy​​−∫0ϵy​​{kqik−2​∂ϵx​∂qi​​(∂ϵy​∂qi​​)}dϵx​]dϵy​​
由最开始的Lemma1可知$X\sim Beta(\alpha ,\beta ),\mathbb{E}[log(X)]=\psi (\alpha )-\psi (\alpha +\beta ),\psi ()=\frac{{\Gamma }^{\prime }()}{\Gamma ()}$

可以观察到第一部分： (21)
$$\begin{array}{rl}& =\left(\genfrac{}{}{0ex}{}{N-1}{k}\right){\int }_0^{\infty }(1-p_i({ϵ}_y){)}^{N-k-1}[log(q_i)k{q}_i^{k-1}(\frac{\partial q_i}{\partial {ϵ}_y}){|}_0^{{ϵ}_y}d{ϵ}_{x}d{ϵ}_y\\ & =\left(\genfrac{}{}{0ex}{}{N-1}{k}\right){\int }_0^{\infty }(1-p_i({ϵ}_y){)}^{N-k-1}log(p_i({ϵ}_y))k{p}_i({ϵ}_y{)}^{k-1}(\frac{\partial p_i({ϵ}_y)}{\partial {ϵ}_y})d{ϵ}_y\\ & =\left(\genfrac{}{}{0ex}{}{N-1}{k}\right){\int }_0^{\infty }(1-p_i({ϵ}_y){)}^{N-k-1}log(p_i({ϵ}_y))k{p}_i({ϵ}_y{)}^{k-1}(\frac{\partial p_i({ϵ}_y)}{\partial {ϵ}_y})d{ϵ}_y\\ & =\left(\genfrac{}{}{0ex}{}{N-1}{k}\right){\int }_0^{\infty }(1-p_i({ϵ}_y){)}^{N-k-1}log(p_i({ϵ}_y))k{p}_i({ϵ}_y{)}^{k-1}d{p}_i({ϵ}_y)\\ & =\left(\genfrac{}{}{0ex}{}{N-1}{k}\right){\int }_0^1(1-u{)}^{N-k-1}log(u)k{u}^{k-1}du=\psi (k)-\psi (N)\end{array}$$
其中，注意到$q_i(0)=0$和当${ϵ}_x={ϵ}_y$时有$q({ϵ}_y,{ϵ}_y)$为一个超立方体，恰是$p_i({ϵ}_y)$，最后的部分由前面Prososition 2.的结论得出。

第二部分： (22)
−(N−1k)∫0∞(1−pi(ϵy))N−k−1(∫0ϵy{kqik−2∂qi∂ϵx(∂qi∂ϵy)}dϵx)dϵy=−(N−1k)∫0∞(1−pi(ϵy))N−k−1(∫0ϵy{kqik−2∂qi∂ϵx(∂qi∂ϵy)}dϵx)dϵy \begin{aligned} &-\tbinom{N-1}{k}\int_0^{\infin}(1-p_i(\epsilon_y))^{N-k-1}(\int_0^{\epsilon_y}\{kq_i^{k-2}\frac{\partial q_i}{\partial \epsilon_x}(\frac{\partial q_i}{\partial \epsilon_y})\}d\epsilon_x)d\epsilon_y\\ &=-\tbinom{N-1}{k}\int_0^{\infin}(1-p_i(\epsilon_y))^{N-k-1}(\int_0^{\epsilon_y}\{kq_i^{k-2}\frac{\partial q_i}{\partial \epsilon_x}(\frac{\partial q_i}{\partial \epsilon_y})\}d\epsilon_x)d\epsilon_y \end{aligned} ​−(kN−1​)∫0∞​(1−pi​(ϵy​))N−k−1(∫0ϵy​​{kqik−2​∂ϵx​∂qi​​(∂ϵy​∂qi​​)}dϵx​)dϵy​=−(kN−1​)∫0∞​(1−pi​(ϵy​))N−k−1(∫0ϵy​​{kqik−2​∂ϵx​∂qi​​(∂ϵy​∂qi​​)}dϵx​)dϵy​​
由局部均匀假设： (23)
$$\begin{array}{rl}& q_i({ϵ}_x,{ϵ}_y)\approx \mu (x_i,y_i)c_{dx}{c}_{dy}{ϵ}_x^{dx}{ϵ}_y^{dy}\\ \to & \frac{\partial }{\partial {ϵ}_y}{q}_i({ϵ}_x,{ϵ}_y)=\mu (x_i,y_i)c_{dx}{c}_{dy}{ϵ}_x^{dx}dy{ϵ}_y^{dy-1}=q_i({ϵ}_x,{ϵ}_y)\frac{dy}{{ϵ}_y}\end{array}$$
并且由刚才分部积分时提及的当${ϵ}_x={ϵ}_y$时有 $q({ϵ}_y,{ϵ}_y)$为一个超立方体，因此 (24)
$$\frac{\partial }{\partial {ϵ}_y}{q}_i({ϵ}_x,{ϵ}_y){|}_{{ϵ}_x={ϵ}_y}=p_i({ϵ}_y)\frac{dy}{{ϵ}_y}$$

则（22）继续简化有： (25)
=−(N−1k)∫0∞(1−pi(ϵy))N−k−1(∫0ϵy{kqik−2∂qi∂ϵx(∂qi∂ϵy)}dϵx)dϵy=−(N−1k)∫0∞(1−pi(ϵy))N−k−1(∫0ϵy{kqik−2∂qi∂ϵx(qidyϵy)}dϵx)dϵy=−(N−1k)∫0∞(1−pi(ϵy))N−k−1(dyϵy)(∫0ϵy{kqik−1∂qi∂ϵx}dϵx)dϵy=−(N−1k)∫0∞(1−pi(ϵy))N−k−1(dyϵy)pi(ϵy)k)dϵy=−(N−1k)∫0∞(1−pi(ϵy))N−k−1(dyϵy)pi(ϵy)k)dϵy \begin{aligned} &=-\tbinom{N-1}{k}\int_0^{\infin}(1-p_i(\epsilon_y))^{N-k-1}(\int_0^{\epsilon_y}\{kq_i^{k-2}\frac{\partial q_i}{\partial \epsilon_x}(\frac{\partial q_i}{\partial \epsilon_y})\}d\epsilon_x)d\epsilon_y\\ &=-\tbinom{N-1}{k}\int_0^{\infin}(1-p_i(\epsilon_y))^{N-k-1}(\int_0^{\epsilon_y}\{kq_i^{k-2}\frac{\partial q_i}{\partial \epsilon_x}(q_i\frac{dy}{\epsilon_y})\}d\epsilon_x)d\epsilon_y\\ &=-\tbinom{N-1}{k}\int_0^{\infin}(1-p_i(\epsilon_y))^{N-k-1}(\frac{dy}{\epsilon_y})(\int_0^{\epsilon_y}\{kq_i^{k-1}\frac{\partial q_i}{\partial \epsilon_x}\}d\epsilon_x)d\epsilon_y\\ &=-\tbinom{N-1}{k}\int_0^{\infin}(1-p_i(\epsilon_y))^{N-k-1}(\frac{dy}{\epsilon_y})p_i(\epsilon_y)^k)d\epsilon_y\\ &=-\tbinom{N-1}{k}\int_0^{\infin}(1-p_i(\epsilon_y))^{N-k-1}(\frac{dy}{\epsilon_y})p_i(\epsilon_y)^k)d\epsilon_y \end{aligned} ​=−(kN−1​)∫0∞​(1−pi​(ϵy​))N−k−1(∫0ϵy​​{kqik−2​∂ϵx​∂qi​​(∂ϵy​∂qi​​)}dϵx​)dϵy​=−(kN−1​)∫0∞​(1−pi​(ϵy​))N−k−1(∫0ϵy​​{kqik−2​∂ϵx​∂qi​​(qi​ϵy​dy​)}dϵx​)dϵy​=−(kN−1​)∫0∞​(1−pi​(ϵy​))N−k−1(ϵy​dy​)(∫0ϵy​​{kqik−1​∂ϵx​∂qi​​}dϵx​)dϵy​=−(kN−1​)∫0∞​(1−pi​(ϵy​))N−k−1(ϵy​dy​)pi​(ϵy​)k)dϵy​=−(kN−1​)∫0∞​(1−pi​(ϵy​))N−k−1(ϵy​dy​)pi​(ϵy​)k)dϵy​​
其中第三个等号注意到了$\frac{dy}{{ϵ}_y}$是常数，注意$dy$是$Y$的维度不是微分。第四个等号是将$k{q}_i^{k-1}$吸收到微分，再由积分所得。

进一步根据（24）有：(26)
$$\frac{d}{d{ϵ}_y}{p}_i({ϵ}_y)=\frac{d}{\partial {ϵ}_y}{q}_i({ϵ}_x,{ϵ}_y){|}_{{ϵ}_x={ϵ}_y}=q_i({ϵ}_x,{ϵ}_y){|}_{{ϵ}_x={ϵ}_y}\frac{dy}{d{ϵ}_y}=p_i({ϵ}_y)\frac{dy}{{ϵ}_y}$$
因此，（25）简化为(26)
$$\begin{array}{rl}& =-\left(\genfrac{}{}{0ex}{}{N-1}{k}\right){\int }_0^{\infty }(1-p_i({ϵ}_y){)}^{N-k-1}(\frac{dy}{{ϵ}_y})p_i({ϵ}_y{)}^k)d{ϵ}_y\\ & =-\left(\genfrac{}{}{0ex}{}{N-1}{k}\right){\int }_0^{\infty }(1-p_i({ϵ}_y){)}^{N-k-1}{p}_i({ϵ}_y{)}^{k-1}\frac{\partial }{\partial y}{p}_i({ϵ}_y)d{ϵ}_y\\ & =-\left(\genfrac{}{}{0ex}{}{N-1}{k}\right){\int }_0^{\infty }(1-p_i({ϵ}_y){)}^{N-k-1}{p}_i({ϵ}_y{)}^{k-1}d{p}_i=-\left(\genfrac{}{}{0ex}{}{N-1}{k}\right)Beta(k,N-k)\\ & =-\frac{1}{k}\end{array}$$
联立（21）和（26），得到(27)
$$E[log{q}_i({ϵ}_x,{ϵ}_y)]=\psi (k)-\psi (N)-\frac{1}{k}$$
进一步的，根据公式（23），自然可得(28)
$$\begin{array}{r}log\mu (x_i,y_i)\approx log{q}_i({ϵ}_x,{ϵ}_y)-log(c_{dx}{c}_{dy})-dxlog{ϵ}_x-dylog{ϵ}_y\end{array}$$
求期望有：(29)
$$\begin{gathered} E[log\mu (x_i,x_y)]= \\ \psi (k)-\psi (N)-\frac{1}{k}-log(c_{dx}{c}_{dy})-dx\frac{1}{N}\sum log{ϵ}_x-dy\frac{1}{N}\sum log{ϵ}_y \end{gathered}$$
再进一步利用（1）式可得KSG2估计器：(30)
$$I^{(2)}(X,Y)=\psi (k)+\psi (N)-\frac{1}{k}-\frac{1}{N}\sum \psi (n_x)-\frac{1}{N}\sum \psi (n_y)$$

注意，这里是$\psi (n_x)$，这是因为$n_x$此时定义为距离小于等于$\frac{{ϵ}_x}{2}$的点的个数。

至此，KSG估计器的形式推导完毕。