稀疏数据分析：马蹄估计量及其理论性质

基本框架
- half-Cauchy分布
- 为什么它叫马蹄估计量
- 后验均值、shrinkage与κ\kappaκ
- θ\thetaθ的边缘先验分布的阶
马蹄估计的一致性
收敛速率
- Clarke-Barron(1990)引理
- 马蹄估计的最优收敛速率

这是对The horseshoe estimator for sparse signal这篇论文的回顾，这篇论文在Modern Bayesian statistics与Bayesian Machine Learning领域比较重要，它提出了一种新的先验结构——horseshoe prior，基于这种先验得到的均值的后验估计在shrinkage上具有类似LASSO的性质，也就是保留数据中稀疏的信号、让噪声缩减为0。

基本框架

假设观察值是yiy_iyi，它服从正态分布
yi∼iidN(θi,σ2)y_i \sim_{iid} N(\theta_i,\sigma^2)yi∼iidN(θi,σ2)

我们引入下面的分为三层的先验：
θi∣λi∼N(0,λi2)λi∣τ∼C+(0,τ)τ∣σ∼C+(0,σ)\theta_i|\lambda_i \sim N(0,\lambda_i^2) \\ \lambda_i|\tau \sim C^+(0,\tau) \\ \tau|\sigma \sim C^+(0,\sigma)θi∣λi∼N(0,λi2)λi∣τ∼C+(0,τ)τ∣σ∼C+(0,σ)

其中C+(0,a)C^+(0,a)C+(0,a)是尺度参数为aaa的half-Cauchy分布。假设σ\sigmaσ的先验是Jeffrey先验，即密度函数与1/σ1/\sigma1/σ成正比。

half-Cauchy分布

如果X∼C+(0,a),a>0X \sim C^+(0,a),a>0X∼C+(0,a),a>0，则称XXX服从尺度参数为aaa的half-Cauchy分布，它的密度函数是
f(x)=2aπ(x2+a2)f(x)=\frac{2a}{\pi(x^2+a^2)}f(x)=π(x2+a2)2a

我们先验证一下归一性：
∫0∞f(x)dx=∫0∞2aπ(x2+a2)dx=2πarctan⁡(x/a)∣0∞=1\int_0^{\infty}f(x)dx = \int_0^{\infty}\frac{2a}{\pi(x^2+a^2)}dx = \frac{2}{\pi}\arctan(x/a)|_0^{\infty}=1∫0∞f(x)dx=∫0∞π(x2+a2)2adx=π2arctan(x/a)∣0∞=1

当然这个分布的期望也是不存在的
∫0+∞xf(x)dx=∫0∞2axπ(x2+a2)dx=aπln⁡(x2+a2)∣0+∞=+∞\int_0^{+\infty}xf(x)dx = \int_0^{\infty} \frac{2ax}{\pi(x^2+a^2)}dx=\frac{a}{\pi}\ln (x^2+a^2)|_0^{+\infty}=+\infty∫0+∞xf(x)dx=∫0∞π(x2+a2)2axdx=πaln(x2+a2)∣0+∞=+∞

为什么它叫马蹄估计量

考虑λi\lambda_iλi的边缘先验分布，
p(λi,τ,σ)∝τπ(λi2+τ2)σπ(τ2+σ2)1σ∝τ(λi2+τ2)(τ2+σ2)p(λi)∝∫0∞∫0∞τ(λi2+τ2)(τ2+σ2)dσdτp(\lambda_i,\tau,\sigma) \propto \frac{\tau}{\pi(\lambda_i^2+\tau^2)}\frac{\sigma}{\pi(\tau^2+\sigma^2)}\frac{1}{\sigma}\propto \frac{\tau}{(\lambda_i^2+\tau^2)(\tau^2+\sigma^2)} \\ p(\lambda_i) \propto \int_0^{\infty}\int_0^{\infty}\frac{\tau}{(\lambda_i^2+\tau^2)(\tau^2+\sigma^2)}d\sigma d\taup(λi,τ,σ)∝π(λi2+τ2)τπ(τ2+σ2)σσ1∝(λi2+τ2)(τ2+σ2)τp(λi)∝∫0∞∫0∞(λi2+τ2)(τ2+σ2)τdσdτ

定义κi=1/(1+λi2)\kappa_i=1/(1+\lambda_i^2)κi=1/(1+λi2)，这个量在Bayesian shrinkage中非常重要，我们在下一个小标题介绍它的意义，但我们可以先分析它的先验分布。现在我们只想做一点定性分析，了解一下κi\kappa_iκi的先验的形状，所以简单起见假设σ=τ=1\sigma=\tau=1σ=τ=1，于是
p(λi)=2π(λi2+1),λi>0p(ki)=p(λi(κi))∣λi′∣=1πx−1/2(1−x)−1/2,x∈(0,1)p(\lambda_i)=\frac{2}{\pi(\lambda_i^2+1)},\lambda_i>0 \\ p(k_i)=p(\lambda_i(\kappa_i))|\lambda_i'|=\frac{1}{\pi}x^{-1/2}(1-x)^{-1/2},x \in (0,1)p(λi)=π(λi2+1)2,λi>0p(ki)=p(λi(κi))∣λi′∣=π1x−1/2(1−x)−1/2,x∈(0,1)

因此ki∼Beta(1/2,1/2)k_i \sim Beta(1/2,1/2)ki∼Beta(1/2,1/2)，懒得自己画图我就扒了百度百科的图，看α=β=0.5\alpha=\beta=0.5α=β=0.5（粉红色）那条，那就是我们κi\kappa_iκi的先验分布，是不是非常像一个马蹄铁的形状，所以这种先验结构被称为马蹄先验，基于这种先验的贝叶斯方法被称为马蹄估计。

后验均值、shrinkage与κ\kappaκ

现在来填一个小坑，κ\kappaκ为什么重要？我们可以做一点简单的推导来理解κ\kappaκ的含义，考虑非常简单的情况，固定τ=σ=1\tau=\sigma=1τ=σ=1，先验可以被简化为两层θ∣λ∼N(0,λ2)λ∼C+(0,1)\theta|\lambda \sim N(0,\lambda^2) \\ \lambda \sim C^+(0,1)θ∣λ∼N(0,λ2)λ∼C+(0,1)

定义κ=1/(1+λ2)\kappa=1/(1+\lambda^2)κ=1/(1+λ2)，则λ2=1−κκ\lambda^2=\frac{1-\kappa}{\kappa}λ2=κ1−κ
p(k)=p(λ(κ))∣λ′∣=1πκ−1/2(1−κ)−1/2,κ∈(0,1)p(θ∣κ)=12πλ2e−θ22λ2=κ2π(1−κ)e−κθ22(1−κ)p(k)=p(\lambda(\kappa))|\lambda'|=\frac{1}{\pi}\kappa^{-1/2}(1-\kappa)^{-1/2},\kappa \in (0,1) \\ p(\theta|\kappa) = \frac{1}{\sqrt{2\pi\lambda^2}}e^{-\frac{\theta^2}{2\lambda^2}} =\sqrt{\frac{\kappa}{2\pi(1-\kappa)}} e^{-\frac{\kappa \theta^2}{2(1-\kappa)}}p(k)=p(λ(κ))∣λ′∣=π1κ−1/2(1−κ)−1/2,κ∈(0,1)p(θ∣κ)=2πλ21e−2λ2θ2=2π(1−κ)κe−2(1−κ)κθ2

我们尝试写出后验均值的表达式，如果y∼N(θ,1)y \sim N(\theta,1)y∼N(θ,1)
E[θ∣y]=∫−∞+∞θp(θ∣y)dθ=∫−∞+∞θ∫01p(θ,κ∣y)dκdθ=Fubini∫01dκ∫−∞+∞θp(θ,κ∣y)dθE[\theta|y] = \int_{-\infty}^{+\infty} \theta p(\theta|y)d\theta = \int_{-\infty}^{+\infty} \theta \int_0^1 p(\theta,\kappa|y)d\kappa d\theta \\ =_{Fubini} \int_0^1 d\kappa \int_{-\infty}^{+\infty} \theta p(\theta,\kappa|y)d\thetaE[θ∣y]=∫−∞+∞θp(θ∣y)dθ=∫−∞+∞θ∫01p(θ,κ∣y)dκdθ=Fubini∫01dκ∫−∞+∞θp(θ,κ∣y)dθ

下面我们用Hierarchical技巧：
p(θ,κ∣y)=p(θ∣κ,y)p(κ∣y)p(\theta,\kappa|y)=p(\theta|\kappa,y)p(\kappa|y)p(θ,κ∣y)=p(θ∣κ,y)p(κ∣y)

于是
∫01dκ∫−∞+∞θp(θ,κ∣y)dθ=∫01p(κ∣y)dκ∫−∞+∞θp(θ∣κ,y)dθ=∫01E[θ∣κ,y]p(κ∣y)dκ\int_0^1 d\kappa \int_{-\infty}^{+\infty} \theta p(\theta,\kappa|y)d\theta = \int_0^1 p(\kappa|y)d\kappa \int_{-\infty}^{+\infty} \theta p(\theta|\kappa,y)d\theta \\ = \int_0^1 E[\theta|\kappa,y]p(\kappa|y)d\kappa ∫01dκ∫−∞+∞θp(θ,κ∣y)dθ=∫01p(κ∣y)dκ∫−∞+∞θp(θ∣κ,y)dθ=∫01E[θ∣κ,y]p(κ∣y)dκ

因为θ∣κ∼N(0,1−κκ),y∣θ∼N(θ,1)\theta|\kappa \sim N(0,\frac{1-\kappa}{\kappa}),y|\theta \sim N(\theta,1)θ∣κ∼N(0,κ1−κ),y∣θ∼N(θ,1)，这是一个经典的共轭分布族，直接用正态分布共轭分布族的结论：
E[θ∣κ,y]=(1−κ)yE[\theta|\kappa,y]=(1-\kappa)yE[θ∣κ,y]=(1−κ)y

所以
∫01E[θ∣κ,y]p(κ∣y)dκ=∫01(1−κ)yp(κ∣y)dκ=(1−E[κ∣y])y\int_0^1 E[\theta|\kappa,y]p(\kappa|y)d\kappa = \int_0^1 (1-\kappa)yp(\kappa|y)d\kappa = (1-E[\kappa|y])y∫01E[θ∣κ,y]p(κ∣y)dκ=∫01(1−κ)yp(κ∣y)dκ=(1−E[κ∣y])y

综上
E[θ∣y]=(1−E[κ∣y])yE[\theta|y]=(1-E[\kappa|y])yE[θ∣y]=(1−E[κ∣y])y

这是θ\thetaθ的后验均值，也是它贝叶斯估计，我们把yyy理解为信号的一个观测，θ\thetaθ代表信号，N(0,1)N(0,1)N(0,1)是噪声，于是这个估计的作用是从观测中去掉噪声还原信号，不难发现
y−E[θ∣y]y=E[κ∣y]\frac{y-E[\theta|y]}{y} = E[\kappa|y]yy−E[θ∣y]=E[κ∣y]

也就是说我们基于马蹄估计还原的信号比观测更小，于是马蹄估计是一种shrinkage estimation，relative amount of shrinkage等于E[κ∣y]E[\kappa|y]E[κ∣y]，我们希望信号的E[κ∣y]E[\kappa|y]E[κ∣y]非常小，噪声的E[κ∣y]E[\kappa|y]E[κ∣y]非常大，这样我们就可以用这个方法把信号从噪声观测中还原出来，下文会分析马蹄估计理论上具有这个性质。

θ\thetaθ的边缘先验分布的阶

θ\thetaθ的边缘先验分布为
p(θ)=∫0∞12πλe−θ22λ22π(1+λ2)dλp(\theta)=\int_0^{\infty} \frac{1}{\sqrt{2\pi}\lambda}e^{-\frac{\theta^2}{2\lambda^2}} \frac{2}{\pi(1+\lambda^2)}d\lambdap(θ)=∫0∞2πλ1e−2λ2θ2π(1+λ2)2dλ

它关于θ−1\theta^{-1}θ−1的阶是我们比较关注的量，因为现代贝叶斯统计认为θ\thetaθ的先验在000处的值越大，并且在尾部关于θ−1\theta^{-1}θ−1的阶越小，用这种先验来做稀疏数据分析的效果就越好。

定理1 关于这个边缘先验分布p(θ)p(\theta)p(θ)，我们有下面两个结论：

lim⁡θ→0p(θ)=∞\lim_{\theta \to 0}p(\theta) = \inftylimθ→0p(θ)=∞
K2log⁡(1+4θ2)<p(θ)<Klog⁡(1+2θ2),K=(2π3)−1/2\frac{K}{2}\log (1+\frac{4}{\theta^2})<p(\theta)<K\log(1+\frac{2}{\theta^2}),K=(2\pi^3)^{-1/2}2Klog(1+θ24)<p(θ)<Klog(1+θ22),K=(2π3)−1/2

证明
第一条，如果θ→0\theta \to 0θ→0，讨论积分
2K∫0∞1λ1(1+λ2)dλ=2K[ln⁡∣x∣−0.5ln⁡(1+x2)]∣0+∞=∞2K\int_0^{\infty} \frac{1}{\lambda}\frac{1}{(1+\lambda^2)}d\lambda=2K[\ln |x|-0.5\ln(1+x^2)]|_0^{+\infty}=\infty2K∫0∞λ1(1+λ2)1dλ=2K[ln∣x∣−0.5ln(1+x2)]∣0+∞=∞

所以lim⁡θ→0p(θ)=∞\lim_{\theta \to 0}p(\theta) = \inftylimθ→0p(θ)=∞。

第二条，做换元u=1/λ2u=1/\lambda^2u=1/λ2，
p(θ)=K∫0+∞11+ue−θ2u2dup(\theta)=K\int_0^{+\infty}\frac{1}{1+u}e^{-\frac{\theta^2u}{2}}dup(θ)=K∫0+∞1+u1e−2θ2udu

再做换元z=1+uz=1+uz=1+u，
p(θ)=Keθ22∫1∞1ze−zθ22dz=Keθ22E1(θ22)p(\theta) = Ke^{\frac{\theta^2}{2}}\int_1^{\infty}\frac{1}{z}e^{-\frac{z\theta^2}{2}}dz = Ke^{\frac{\theta^2}{2}}E_1(\frac{\theta^2}{2})p(θ)=Ke2θ2∫1∞z1e−2zθ2dz=Ke2θ2E1(2θ2)

这里的E1E_1E1是一个特殊函数，它被称为exponential integral function，关于它有一个不等式
0.5e−tlog⁡(1+2t)<E1(t)<e−tlog⁡(1+1t),∀t>00.5e^{-t}\log(1+\frac{2}{t})<E_1(t) < e^{-t}\log(1+\frac{1}{t}),\forall t>00.5e−tlog(1+t2)<E1(t)<e−tlog(1+t1),∀t>0

所以
K2log⁡(1+4θ2)<p(θ)<Klog⁡(1+2θ2),K=(2π3)−1/2\frac{K}{2}\log (1+\frac{4}{\theta^2})<p(\theta)<K\log(1+\frac{2}{\theta^2}),K=(2\pi^3)^{-1/2}2Klog(1+θ24)<p(θ)<Klog(1+θ22),K=(2π3)−1/2

当θ\thetaθ足够大时，上界与下界关于θ−1\theta^{-1}θ−1的阶都是2，因为
log⁡(1+2θ2)=2θ−2+o(θ−2)log⁡(1+4θ2)=4θ−2+o(θ−2)\log(1+\frac{2}{\theta^2}) = 2\theta^{-2}+o(\theta^{-2}) \\ \log(1+\frac{4}{\theta^2}) = 4\theta^{-2}+o(\theta^{-2})log(1+θ22)=2θ−2+o(θ−2)log(1+θ24)=4θ−2+o(θ−2)

马蹄估计的一致性

定理2 假设p(∣y−θ∣)p(|y-\theta|)p(∣y−θ∣)是给定样本数据为yyy的似然函数，这个式子表示分布的参数是θ\thetaθ，并且它是一个位置参数；假设θ\thetaθ的边缘先验分布为p(θ)p(\theta)p(θ)，它是一个混合高斯分布，即θ∣λ∼N(0,λ2)\theta|\lambda \sim N(0,\lambda^2)θ∣λ∼N(0,λ2)，其中λ\lambdaλ的先验为p(λ),λ>0p(\lambda),\lambda>0p(λ),λ>0，假设似然函数与边缘先验分布p(θ)p(\theta)p(θ)使yyy的边缘密度m(y)m(y)m(y)有界，定义下面三个量
m∗(y)=∫p(∣y−θ∣)p∗(θ)dθp∗(θ)=∫λ>0p(θ∣λ)p∗(λ)dλp∗(λ)=λ2p(λ)m^*(y)=\int p(|y-\theta|)p^*(\theta)d\theta \\ p^*(\theta) = \int_{\lambda>0}p(\theta|\lambda)p^*(\lambda)d\lambda \\ p^*(\lambda)=\lambda^2p(\lambda)m∗(y)=∫p(∣y−θ∣)p∗(θ)dθp∗(θ)=∫λ>0p(θ∣λ)p∗(λ)dλp∗(λ)=λ2p(λ)

则
E[θ∣y]=1m(y)ddym∗(y)E[\theta|y]=\frac{1}{m(y)}\frac{d}{dy}m^*(y)E[θ∣y]=m(y)1dydm∗(y)

如果p(∣y−θ∣)p(|y-\theta|)p(∣y−θ∣)是正态的，那么
E[θ∣y]=y+ddylog⁡m(y)E[\theta|y]=y+\frac{d}{dy}\log m(y)E[θ∣y]=y+dydlogm(y)

评注
这个定理不但看上去比较长，而且也是一个非常重要的定理。在原文中，作者们还证明了如果y∼N(θ,1),θ∣λ∼N(0,λ2),λ∼C+(0,τ)y \sim N(\theta,1),\theta|\lambda \sim N(0,\lambda^2),\lambda \sim C^+(0,\tau)y∼N(θ,1),θ∣λ∼N(0,λ2),λ∼C+(0,τ)，τ\tauτ为已知的常数，则∣y−E[θ∣y]∣|y-E[\theta|y]|∣y−E[θ∣y]∣有界，并且
lim⁡∣y∣→∞ddylog⁡m(y)→0\lim_{|y| \to \infty} \frac{d}{dy}\log m(y) \to 0∣y∣→∞limdydlogm(y)→0

这个结果说明了马蹄估计的一致性。这个结果的证明并不复杂，它的要点在于如何计算ddylog⁡m(y)\frac{d}{dy}\log m(y)dydlogm(y)这个值，其中
m(y)=K∫0∞e−y22(1+τ2λ2)11+λ2τ211+λ2dλm(y) = K\int_0^{\infty} e^{-\frac{y^2}{2(1+\tau^2\lambda^2)}}\frac{1}{\sqrt{1+\lambda^2\tau^2}}\frac{1}{1+\lambda^2}d\lambdam(y)=K∫0∞e−2(1+τ2λ2)y21+λ2τ211+λ21dλ

原文的思路是做换元，z=11+τ2λ2z=\frac{1}{1+\tau^2\lambda^2}z=1+τ2λ21，然后把积分凑成Φ1\Phi_1Φ1函数，这个函数被称为Humbert series，
Φ1(a,b,c;x,y)=Γ(c)Γ(a)Γ(c−a)∫01ta−1(1−t)c−a−1(1−xt)−beytdt\Phi_1(a,b,c;x,y)=\frac{\Gamma(c)}{\Gamma(a)\Gamma(c-a)}\int_0^1t^{a-1}(1-t)^{c-a-1}(1-xt)^{-b}e^{yt}dtΦ1(a,b,c;x,y)=Γ(a)Γ(c−a)Γ(c)∫01ta−1(1−t)c−a−1(1−xt)−beytdt

根据Gordy (1998)等式，Humbert series可以用Kummer第一类函数M(a,b;x)M(a,b;x)M(a,b;x)表示，
Φ1(a,b,c;x,y)=ex∑n=0∞(a)n(b)n(c)nynn!M(c−a,c+n,−x)\Phi_1(a,b,c;x,y)=e^{x}\sum_{n=0}^{\infty}\frac{(a)_n(b)_n}{(c)_n}\frac{y^n}{n!}M(c-a,c+n,-x)Φ1(a,b,c;x,y)=exn=0∑∞(c)n(a)n(b)nn!ynM(c−a,c+n,−x)

Kummer第一类函数M(a,b;x)M(a,b;x)M(a,b;x)满足下面的等式(Slater 1960 Chapter 4)，
M(a,b;x)={Γ(a)Γ(b)exxa−b(1+O(x−1)),x>0Γ(a)Γ(b−a)(−x)−a(1+O(x−1)),x<0M(a,b;x)=\begin{cases} \frac{\Gamma(a)}{\Gamma(b)}e^xx^{a-b}(1+O(x^{-1})),x>0 \\ \frac{\Gamma(a)}{\Gamma(b-a)}(-x)^{-a}(1+O(x^{-1})),x<0 \end{cases}M(a,b;x)={Γ(b)Γ(a)exxa−b(1+O(x−1)),x>0Γ(b−a)Γ(a)(−x)−a(1+O(x−1)),x<0

基于这几个特殊函数进行计算即可验证一致性。

证明
有两个有用的等式
ddyp(y−θ)=−ddθp(y−θ)λ2ddθN(θ∣0,λ2)=−θN(θ∣0,λ2)\frac{d}{dy}p(y-\theta)=-\frac{d}{d\theta}p(y-\theta) \\ \lambda^2 \frac{d}{d\theta}N(\theta|0,\lambda^2)=-\theta N(\theta|0,\lambda^2)dydp(y−θ)=−dθdp(y−θ)λ2dθdN(θ∣0,λ2)=−θN(θ∣0,λ2)

第一个等式就是Schwartz定理，因为ppp是概率密度，已经是一阶导数了，而二阶导数与求导次序无关，所以第一个等式成立；第二个等式就是对正态分布求导。

所以对于
E[θ∣y]=∫θp(θ∣y)dθE[\theta|y]=\int \theta p(\theta|y)d\thetaE[θ∣y]=∫θp(θ∣y)dθ

我们用前文使用过的Hierarchical技巧，
p(θ∣y)=1m(y)p(y∣θ)p(θ∣λ)p(λ)=1m(y)p(y−θ)N(θ∣0,λ2)p(λ)p(\theta|y) = \frac{1}{m(y)}p(y|\theta)p(\theta|\lambda)p(\lambda) \\ = \frac{1}{m(y)}p(y-\theta)N(\theta|0,\lambda^2)p(\lambda)p(θ∣y)=m(y)1p(y∣θ)p(θ∣λ)p(λ)=m(y)1p(y−θ)N(θ∣0,λ2)p(λ)

于是
∫θp(θ∣y)dθ=∫1m(y)θp(y−θ)N(θ∣0,λ2)p(λ)dθ=1m(y)∫ddyp(y−θ)N(θ∣0,λ2)p∗(λ)dθ=1m(y)ddym∗(y)\int \theta p(\theta|y)d\theta = \int \frac{1}{m(y)}\theta p(y-\theta)N(\theta|0,\lambda^2)p(\lambda)d\theta \\ = \frac{1}{m(y)}\int \frac{d}{dy}p(y-\theta)N(\theta|0,\lambda^2)p^*(\lambda)d\theta = \frac{1}{m(y)}\frac{d}{dy}m^*(y)∫θp(θ∣y)dθ=∫m(y)1θp(y−θ)N(θ∣0,λ2)p(λ)dθ=m(y)1∫dydp(y−θ)N(θ∣0,λ2)p∗(λ)dθ=m(y)1dydm∗(y)

收敛速率

作为理论研究者，我们除了关心统计模型的渐近性质（比如一致性）外，我们还想了解模型的非渐近性质(non-asymptotic properties)，比如收敛速率，下面给出的是马蹄估计的后验分布趋近于真实分布的速率。

Clarke-Barron(1990)引理

用θ0\theta_0θ0表示参数的真实值，p(y∣θ)p(y|\theta)p(y∣θ)表示样本的概率密度，用L(p1,p2)L(p_1,p_2)L(p1,p2)表示从p2p_2p2到p1p_1p1的Kullback-Leibler divergence，
L(p1,p2)=Ep1[log⁡(p1/p2)]L(p_1,p_2) = E_{p_1}[\log(p_1/p_2)]L(p1,p2)=Ep1[log(p1/p2)]

∀ϵ>0\forall \epsilon>0∀ϵ>0，定义θ0\theta_0θ0的Kullback-Leibler信息邻域为
Aϵ={θ:L(pθ0,pθ)≤ϵ}A_{\epsilon} = \{\theta:L(p_{\theta_0},p_{\theta}) \le \epsilon\}Aϵ={θ:L(pθ0,pθ)≤ϵ}

用μn\mu_nμn表示后验分布，nnn为样本数，μ\muμ表示先验概率分布，密度函数的贝叶斯估计为
p^n=∫pθμn(dθ)\hat p_n = \int p_{\theta} \mu_n(d\theta)p^n=∫pθμn(dθ)

定义贝叶斯估计的Cesaro average risk为
Rn=1n∑j=1nL(pθ0,p^j)R_n = \frac{1}{n}\sum_{j=1}^n L(p_{\theta_0},\hat p_j)Rn=n1j=1∑nL(pθ0,p^j)

Clarke-Barron(1990)引理
如果先验满足∀ϵ>0\forall \epsilon>0∀ϵ>0，μ(Aϵ)>0\mu(A_{\epsilon})>0μ(Aϵ)>0，则
Rn≤ϵ−1nμ(Aϵ)R_n \le \epsilon-\frac{1}{n}\mu(A_{\epsilon})Rn≤ϵ−n1μ(Aϵ)

马蹄估计的最优收敛速率

定理4 假设先验满足∀ϵ>0\forall \epsilon>0∀ϵ>0，μ(Aϵ)>0\mu(A_{\epsilon})>0μ(Aϵ)>0

如果θ\thetaθ的边缘先验分布连续有上界且支撑集包含θ0\theta_0θ0的至少一个邻域，则Rn=O(n−1log⁡n)R_n=O(n^{-1}\log n)Rn=O(n−1logn)
马蹄估计量的Cesaro average risk具有下面的性质: 如果θ0=0\theta_0=0θ0=0，Rn=O(n−1(log⁡n−blog⁡log⁡n))R_n=O(n^{-1}(\log n - b \log \log n))Rn=O(n−1(logn−bloglogn))，其中bbb是一个常数；如果θ0≠0\theta_0 \ne 0θ0=0，Rn=O(n−1log⁡n)R_n = O(n^{-1}\log n)Rn=O(n−1logn)

证明
第一条以及第二条θ0≠0\theta_0 \ne 0θ0=0，取ϵn=1/n\epsilon_n=1/nϵn=1/n，如果∀ϵ>0\forall \epsilon>0∀ϵ>0，μ(Aϵ)>0\mu(A_{\epsilon})>0μ(Aϵ)>0，则
μ(Aϵn)=∫Aϵnp(θ)dθ\mu(A_{\epsilon_{n}})=\int_{A_{\epsilon_n}} p(\theta)d\thetaμ(Aϵn)=∫Aϵnp(θ)dθ

这里的AϵnA_{\epsilon_n}Aϵn是Kullback-Leibler信息邻域
Aϵn={θ:L(pθ0,pθ)≤1n}={θ:Epθ0[log⁡(pθ0/pθ)]≤1n}⊂{θ:∣θ−θ0∣≤1n}A_{\epsilon_n} = \{\theta:L(p_{\theta_0},p_{\theta}) \le \frac{1}{n}\}=\{\theta:E_{p_{\theta_0}}[\log(p_{\theta_0}/p_{\theta})] \le \frac{1}{n}\} \\ \subset \{\theta:|\theta-\theta_0| \le \frac{1}{\sqrt{n}}\}Aϵn={θ:L(pθ0,pθ)≤n1}={θ:Epθ0[log(pθ0/pθ)]≤n1}⊂{θ:∣θ−θ0∣≤n1}

因为p(θ)p(\theta)p(θ)有上界，∃C>0\exists C>0∃C>0，p(θ)<Cp(\theta)<Cp(θ)<C
μ(Aϵn)≤∫θ0−1/nθ0+1/nCdθ=2C/n=O(1/n)\mu(A_{\epsilon_{n}}) \le \int_{\theta_0-1/\sqrt{n}}^{\theta_0+1/\sqrt{n}}Cd\theta = 2C/\sqrt{n} = O(1/\sqrt{n})μ(Aϵn)≤∫θ0−1/nθ0+1/nCdθ=2C/n=O(1/n)

根据Clarke-Barron(1990)引理，
Rn≤1n−log⁡(Cn−1/2)n=O(log⁡n/n)R_n \le \frac{1}{n}-\frac{\log(Cn^{-1/2})}{n }=O(\log n/n)Rn≤n1−nlog(Cn−1/2)=O(logn/n)

第二条θ0=0\theta_0 = 0θ0=0，根据定理一，
p(θ)≥Klog⁡(1+4θ−2)p(\theta) \ge K\log(1+4\theta^{-2})p(θ)≥Klog(1+4θ−2)

所以
μ(Aϵ)≥K∫0ϵlog⁡(1+4θ−2)dθ\mu(A_{\epsilon}) \ge K \int_0^{\sqrt{\epsilon}}\log(1+4\theta^{-2})d\thetaμ(Aϵ)≥K∫0ϵlog(1+4θ−2)dθ

做换元u=1/θ2u=1/\theta^2u=1/θ2，
μ(Aϵ)≥K∫4/ϵ∞log⁡(1+u)u3/2du=Kϵ1/2log⁡(1+4ϵ−1)+2K∫4/ϵ∞1u1/2(1+u)du\mu(A_{\epsilon}) \ge K\int_{4/\epsilon}^{\infty} \frac{\log(1+u)}{u^{3/2}}du \\ = K\epsilon^{1/2}\log(1+4\epsilon^{-1})+2K\int_{4/\epsilon}^{\infty} \frac{1}{u^{1/2}(1+u)}duμ(Aϵ)≥K∫4/ϵ∞u3/2log(1+u)du=Kϵ1/2log(1+4ϵ−1)+2K∫4/ϵ∞u1/2(1+u)1du

取ϵ=1/n\epsilon=1/nϵ=1/n，根据Clarke-Barron(1990)引理，可以得到Rn=O(n−1(log⁡n−blog⁡log⁡n))R_n=O(n^{-1}(\log n - b \log \log n))Rn=O(n−1(logn−bloglogn))。

稀疏数据分析：马蹄估计量及其理论性质相关推荐

数据分析和数据挖掘的理论研究必要性
2019独角兽企业重金招聘Python工程师标准>>> 数据分析,并不抽象,传统的数据分析,包括很多.例如信号处理中的DCT,滤波,IDCT变换.由于确定了滤波窗口的特性,使得对一个 ...
数据分析面试题——统计理论
1. 扑克牌54张,平均分成2份,求这2份都有2张A的概率. 解: 排列公式 : A n m = n ( n − 1 ) ( n − 1 ) . . . ( n − m + 1 ) = n ! ( n ...
数据分析的理论与实践
学习重在与理论在实践中的运用,数据分析更为如此.而我们所知的数据分析包括很多的理论知识,简单的有统计学中的参数与非参,复杂的有机器学习中的神经网络和支持向量机.在学习的过程中我们不可能在一开始就把这些 ...
UA MATH567 高维统计专题1 稀疏信号及其恢复1 L0-norm minimization
UA MATH567 高维统计专题1 稀疏信号及其恢复1 L0-norm minimization L0L^0L0-norm L0L_0L0-norm minimization Exhaustive ...
【Python】时间序列数据分析与预测之Python工具汇总
本文中总结了十多种时间序列数据分析和预测工具和python库,在我们处理时间序列项目时,可以翻开本文,根据需要选择合适的工具,将会事半功倍! 在处理时间序列项目时,数据科学家或 ML 工程师通常会使用 ...
信号检测与估计理论_论文解读 | 利用脑功能连接实现疲劳驾驶检测
©PaperWeekly 原创 · 作者|张玮玮学校|东北大学硕士生研究方向|脑电情绪识别论文标题:Driving Fatigue Recognition with Functional Con ...
(附源码)Python云顶之弈数据分析系统毕业设计451545
目录摘要 1 1 绪论 1 1.1研究背景 1 1.2开发意义 1 1.3系统开发技术的特色 1 1.4论文结构与章节安排 1 2 基于Python云顶之弈数据分析系统系统分析 3 2.1 可行性 ...
时间序列数据分析与预测之Python工具汇总
‍ ‍ 大家好,我是辰哥‍ ‍ 本文中硬核总结了十多种时间序列数据分析和预测工具和python库,在我们处理时间序列项目时,可以翻开本文,根据需要选择合适的工具,将会事半功倍! 在处理时间序列项目时, ...
数据分析思维学习1：业务指标、数据获取、数据仓库、数据治理、数据分析方法
数据分析入门实战最近上了一个关于数据分析思维的网课,内容比较偏业务上的实际应用,关注点更多在于数据思维.以下是我在学习过程中的笔记整理. 文章目录数据分析入门实战第一章数据分析与数据分析师数据 ...

稀疏数据分析：马蹄估计量及其理论性质

稀疏数据分析：马蹄估计量及其理论性质

基本框架

half-Cauchy分布

为什么它叫马蹄估计量

后验均值、shrinkage与κ\kappaκ

θ\thetaθ的边缘先验分布的阶

马蹄估计的一致性

收敛速率

Clarke-Barron(1990)引理

马蹄估计的最优收敛速率

稀疏数据分析：马蹄估计量及其理论性质相关推荐

最新文章

热门文章