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

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

基本框架

假设观察值是$y_i$，它服从正态分布
$$y_i{\sim }_{iid}N({\theta }_i,{\sigma }^2)$$

我们引入下面的分为三层的先验：
$$\begin{gathered} {\theta }_i|{\lambda }_i\sim N(0,{\lambda }_i^2) \\ {\lambda }_i|\tau \sim C^{+}(0,\tau ) \\ \tau |\sigma \sim C^{+}(0,\sigma ) \end{gathered}$$

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

half-Cauchy分布

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

我们先验证一下归一性：
$${\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$$

当然这个分布的期望也是不存在的
$${\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$$

为什么它叫马蹄估计量

考虑${\lambda }_i$的边缘先验分布，
$$\begin{gathered} 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\tau \end{gathered}$$

定义${\kappa }_i=1/(1+{\lambda }_i^2)$，这个量在Bayesian shrinkage中非常重要，我们在下一个小标题介绍它的意义，但我们可以先分析它的先验分布。现在我们只想做一点定性分析，了解一下${\kappa }_i$的先验的形状，所以简单起见假设$\sigma =\tau =1$，于是
$$\begin{gathered} p({\lambda }_i)=\frac{2}{\pi ({\lambda }_i^2+1)},{\lambda }_i>0 \\ p(k_i)=p({\lambda }_i({\kappa }_i))|{\lambda }_i^{\prime }|=\frac{1}{\pi }{x}^{-1/2}(1-x{)}^{-1/2},x\in (0,1) \end{gathered}$$

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

后验均值、shrinkage与$\kappa$

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

定义$\kappa =1/(1+{\lambda }^2)$，则${\lambda }^2=\frac{1-\kappa }{\kappa }$
p ( k ) = p ( λ ( κ ) ) ∣ λ ′ ∣ = 1 π κ − 1 / 2 ( 1 − κ ) − 1 / 2 , κ ∈ ( 0 , 1 ) p ( θ ∣ κ ) = 1 2 π λ 2 e − θ 2 2 λ 2 = κ 2 π ( 1 − κ ) e − κ θ 2 2 ( 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}{\s