Sunforger

Sunforger

指数分布族初探

论文读到指数型分布族相关内容,遂做一个整理。不尽准确,望大佬们不吝指教。

指数分布族 Exponential Family 也称为指数族、指数族分布,是统计中最重要的参数分布族

学习指数型分布族,应当与指数分布 Exponential Distribution 区分开来。两者并不是同一个东西。

所谓 “族”,英语里中称为 family,通常是具有相似特征的一类。指数型分布族是一组分布,其概率密度函数和概率分布函数随着分布参数的变化而变化。

常见的指数型分布族包括:
正态分布、卡方分布、二项分布、多项分布、Poisson 分布、Pascal 分布、 β\beta分布、 Γ\Gamma分布、对数正态分布等等。具体可见维基百科词条知乎专栏

指数分布族#

指数型分布族 exponential family 的概率密度函数有以下表现形式:

fX(x;θ)=h(x)exp{η(θ),T(x)A(θ)}f_\mathbf{X}(x;\theta) = h(x)\exp\{\langle\eta(\theta), T(x)\rangle-A(\theta)\}

其中,θ\theta 是唯一的参数,满足上式的所有 θ\theta 构成参数空间 Θ\Theta,对应的参数分布族 {Pθ:θΘ}\{P_\theta:\theta\in\Theta\} 是指数族。必须注意,这里的参数 θ\theta 不局限于狭隘的实数,也可以是 nn 维向量 θRn\theta\in \mathbb{R}^n

参数 θ\theta 改变,分布 XX 的形态(概率密度函数、概率分布函数以及对应的图像)也将改变。随机变量 xx 服从于分布 XX。函数 T(x),h(x),η(θ),A(θ)T(x), h(x), \eta(\theta), A(\theta) 都是已知函数。函数 h(x)h(x) 为非负函数。

h(x)h(x) 通常被称为是基础度量 base measure.

T(x)T(x) 是充分统计量 sufficient statistic.

A(θ)A(\theta) 是累积量生成函数 cumulant generating function 或称为对数配分函数(配分函数的对数) log-partition function 。显然 A(θ)A(\theta) 是实值函数,返回一个实数。

其中,η(θ)\eta(\theta)T(x)T(x) 可以是实数,也可以是向量。

由定义式子可以看出,根据指数函数的性质 exp{}=e{}>0\exp\{\cdot\} = e^{\{\cdot\}} > 0 是非负的。所以指数族的支撑集只与 h(x)h(x) 有关。也就是说,只与 xx 有关,而与未知参数 θ\theta 无关。我们利用这一点,可以排除非指数型分布族(比如均匀分布)。

这里需要简单补充一下所谓支撑集的概念。简单来说,对于实值函数 ff 而言, ff 的支撑集定义为:

supp(f)={xX:f(x)0}\text{supp}(f)=\{x\in X:f(x)\neq0\}

支撑集是函数 ff 原定义域 XX 的一个子集。了解更多,可以参考 维基百科词条CSDN 博客。在概率密度函数中,由于概率非负,可以定义随机变量的支撑集为(见 知乎专栏):

supp(X)={xR:fX(x)>0}\text{supp}(X) =\{x\in R : f_X(x)\gt 0\}

几种等价形式#

基于指数的运算法则,通过等价变换,给出指数族的两种等价形式:

fX(x;θ)=h(x)g(θ)exp{η(θ),T(x)}f_\mathbf{X}(x;\theta) = h(x)g(\theta)\exp\{\langle\eta(\theta), T(x)\rangle\}

fX(x;θ)=exp{η(θ),T(x)A(θ)+B(x)}f_\mathbf{X}(x;\theta) = \exp\{\langle\eta(\theta), T(x)\rangle-A(\theta)+B(x)\}

对应的替换关系是:A(θ)=lng(θ)-A(\theta) = \ln g(\theta)B(x)=lnh(x)B(x)=\ln h(x)

特别地,我们取 Z(θ)=1g(θ)Z(\theta) = \dfrac{1}{g(\theta)}, 可以得到另一中非常常见的指数族的表达式如下。其中 Z(θ)Z(\theta) 是这个分布的配分函数 partition function。

fX(x;θ)=1Z(θ)h(x)exp{η(θ),T(x)}f_\mathbf{X}(x;\theta) = \frac{1}{Z(\theta)}h(x)\exp\{\langle\eta(\theta), T(x)\rangle\}

规范形式#

上述定义式中,η(θ)\eta(\theta) 是关于参数 θ\theta 的一个函数。在指数族中,我们要求 η()\eta(\cdot) 是一个双射函数(也就是一一对应函数)。双射意味着函数单调可导,且存在反函数。

我们借助双射函数的特性,就可以简化指数族的形式。令 θ^=η(θ)\hat\theta = \eta(\theta),当然,这个变化是可逆的 θ=η1(θ^)\theta = \eta^{-1}(\hat\theta)。 于是,我们得到:fX(x;θ^)=h(x)exp{θ^,T(x)A(θ^)}f_\mathbf{X}(x;\hat\theta) = h(x)\exp\{\langle\hat\theta, T(x)\rangle-A^\prime(\hat\theta)\}

等价替换一下符号,我们得到指数族规范形式 Canonical Form 如下:

fX(x;θ)=h(x)exp{θ,T(x)A(θ)}f_\mathbf{X}(x;\theta) = h(x)\exp\{\langle\theta, T(x)\rangle-A(\theta)\}

我们通常把更新过后的这个参数 θ\theta 称为指数族的规范参数(或自然参数)。

自然形式#

虽然在定义上各有各的说法,但一般认为指数族的自然形式 Natural Form 和规范形式两者等同或几乎等同。例如斯坦福大学材料伯克利大学材料麻省理工课件博客知乎专栏 1 以及 知乎专栏 2

维基百科提供了另一种理解,这里不做介绍。

自然参数空间#

介绍自然参数空间 Natural Parameter Space 之前,我们首先介绍对数配分函数 log-partition function A(θ)A(\theta)

A(θ)=log(Xh(x)exp{θ,T(x)} dx) A(\theta) = \log\left(\int_X h(x)\exp\{\langle\theta, T(x)\rangle\}~{\rm d}x\right)

所谓配分函数,可以理解为归一化常数的一种特殊形式。

对数配分函数 A(θ)A(\theta) 使得 fX(x;θ)f_\mathbf{X}(x;\theta) 得以归一化,也就是说,保证了 fX(x;θ)f_\mathbf{X}(x;\theta) 是一个概率密度函数。理解这个归一化,可以参考上文 几种等价形式 小节中的这个表达式

fX(x;θ)=1Z(θ)h(x)exp{η(θ),T(x)}f_\mathbf{X}(x;\theta) = \frac{1}{Z(\theta)}h(x)\exp\{\langle\eta(\theta), T(x)\rangle\}

其中Z(θ)=Xh(x)exp{θ,T(x)} dxZ(\theta) = \int_X h(x)\exp\{\langle\theta, T(x)\rangle\}~{\rm d}x是与 xx 无关的一个函数。然后把两边式子同时积分,得到

XfX(x;θ)=1Z(θ)Xh(x)exp{η(θ),T(x)}=1\int_Xf_\mathbf{X}(x;\theta) = \frac{1}{Z(\theta)}\int_X h(x)\exp\{\langle\eta(\theta), T(x)\rangle\} = 1

所谓自然参数空间,就是使配分函数有限 (<\lt \infty) 时的参数 θ\theta 的集合,即:

N={θ:Xh(x)exp{θ,T(x)} dx<}={θ:Z(θ)<}\mathcal N = \left\{\theta:\int_X h(x)\exp\{\langle\theta, T(x)\rangle\}~{\rm d}x \lt \infty\right\} = \left\{\theta:Z(\theta) \lt \infty\right\}

自然参数空间有这一些特殊的性质。首先,自然参数空间 N\mathcal N 是一个凸集 Convex Set,对数配分函数 A(θ)A(\theta) 是一个凸函数 Convex Function。证明如下:

考虑两个不同的参数 θ1N, θ2N\theta_1\in\mathcal N,~\theta_2\in\mathcal N,给定 0<λ<10\lt\lambda\lt 1 证明 θ=λθ1+(1λ)θ2\theta=\lambda\theta_1+(1-\lambda)\theta_2 也在自然参数空间 N\mathcal N 内(即证明 θN\theta\in\mathcal N 也成立)

Z(θ)=exp{A(θ)}=exp{A(λθ1+(1λ)θ2)}=Xh(x)exp{(λθ1+(1λ)θ2),T(x)} dx=X(h(x)λexp{λθ1,T(x)})(h(x)1λexp{(1λ)θ2,T(x)}) dx(Xh(x)exp{1λλθ1,T(x)} dx)λ(Xh(x)exp{11λ(1λ)θ2,T(x)} dx)1λ=Z(θ1)λZ(θ2)1λ\begin{aligned} Z(\theta) &= \exp\{A(\theta)\} = \exp\{A(\lambda\theta_1+(1-\lambda)\theta_2)\}\\ &=\int_X h(x)\exp\{\langle(\lambda\theta_1+(1-\lambda)\theta_2), T(x)\rangle\}~{\rm d}x \\ & = \int_X \left(h(x)^{\lambda}\exp\{\langle\lambda\theta_1, T(x)\rangle \}\right)\left(h(x)^{1-\lambda}\exp\{\langle(1-\lambda)\theta_2, T(x)\rangle\}\right)~{\rm d}x \\ &\leq \left(\int_X h(x)\exp\{\frac1\lambda\langle\lambda\theta_1, T(x)\rangle\} ~{\rm d}x \right)^\lambda \left(\int_X h(x)\exp\{\frac1{1-\lambda}\langle(1-\lambda)\theta_2, T(x)\rangle\} ~{\rm d}x \right)^{1-\lambda} \\ &=Z(\theta_1)^\lambda \cdot Z(\theta_2)^{1-\lambda} \end{aligned}

上式中的 \leq 来自于赫尔德不等式 Hölder's inequality 其定义可以参考 Wolfram MathWorld知乎专栏。 值得一提的是,著名的数学软件 Mathematica 就是 Wolfram Research 公司开发的。

因为 θ1,θ2N\theta_1,\theta_2\in \mathcal NZ(θ1),Z(θ2)<Z(\theta_1),Z(\theta_2)\lt\infty 成立。所以 Z(θ)=Z(θ1)λZ(θ2)1λ<Z(\theta) = Z(\theta_1)^\lambda \cdot Z(\theta_2)^{1-\lambda} \lt \infty 也成立,根据定义,即可得 θN\theta\in\mathcal N。于是可以证明自然参数空间 N\mathcal N 是一个凸集。

将上式取对数,得到:

A(θ)=A(λθ1+(1λ)θ2)λA(θ1)+(1λ)A(θ2)A(\theta) = A(\lambda\theta_1+(1-\lambda)\theta_2) \leq \lambda A(\theta_1) + (1-\lambda)A(\theta_2)

于是可以证明,对数配分函数 A(θ)A(\theta) 是一个凸函数。θ1θ2\theta_1\neq\theta_2 时,Hölder's inequality 无法取到等号,A(θ)A(\theta) 是严格的凸函数

关于凸集、凸函数的定义,可以参看凸优化教程 知乎专栏 或凸优化经典教材 cvxbook by Stephen Boyd

指数分布族实例#

回顾指数族的规范形式,下面我们证明几种常见的分布,属于指数族。

fX(x;θ)=h(x)exp{θ,T(x)A(θ)}f_\mathbf{X}(x;\theta) = h(x)\exp\{\langle\theta, T(x)\rangle-A(\theta)\}

伯努利分布(两点分布)#

伯努利分布的概率质量函数(伯努利分布是离散的,所以是概率质量函数)为:

p(x;λ)=λx(1λ)(1x)p(x;\lambda) = \lambda^x\cdot (1-\lambda)^{(1-x)}

其中,λ\lambda 是这个伯努利分布的参数(事件发生的概率),x=0x =0 (事件不发生),x=1x =1 (事件发生)。不存在其他的 xx 取值。

我们将式子改写:

p(x;λ)=λx(1λ)(1x)=exp{log(λ1λ)x+log(1λ)}\begin{aligned} p(x;\lambda) &= \lambda^x\cdot (1-\lambda)^{(1-x)}\\ &\color{red}=\exp\left\{\log\left(\frac{\lambda}{1-\lambda}\right)x+\log(1-\lambda) \right\} \end{aligned}

我们取

θ=λ1λ,T(x)=x,A(θ)=log(1λ)=log(1+eθ),h(x)=1\theta = \frac{\lambda}{1-\lambda}, \quad T(x)=x,\quad A(\theta) = -\log(1-\lambda) = \log(1+e^\theta),\quad h(x) = 1

可证伯努利分布属于单参数指数族。

泊松分布#

泊松分布的概率质量函数如下:

p(x;λ)=λxeλx!=1x!exp{xlogλλ}\begin{aligned} p(x;\lambda) &= \frac{\lambda^xe^{-\lambda}}{x!} \\ &\color{red} = \frac{1}{x!}\exp\{x\log\lambda-\lambda\} \end{aligned}

θ=logλ,T(x)=x,A(θ)=λ=eθ,h(x)=1x!\theta = \log\lambda,\quad T(x) = x,\quad A(\theta)=\lambda=e^\theta,\quad h(x)=\frac{1}{x!}

可证泊松分布属于单参数指数族。

高斯分布(正态分布)#

高斯分布的概率密度函数如下

p(x;μ,σ2)=12πσexp{12σ2(xμ)2}=12πexp{μσ2x12σ2x212σ2μ2logσ}\begin{aligned} p(x;\mu,\sigma^2) &= \frac{1}{\sqrt{2\pi}\sigma}\exp\left\{ -\frac{1}{2\sigma^2}(x-\mu)^2 \right\}\\ & \color{red}=\frac{1}{\sqrt{2\pi}}\exp\left\{ \frac{\mu}{\sigma^2}x-\frac{1}{2\sigma^2}x^2-\frac{1}{2\sigma^2}\mu^2-\log\sigma \right\} \end{aligned}

θ=[μ/σ21/2σ2],T(x)=[xx2],A(θ)=μ22σ2+logσ=θ124θ212log(2θ2),h(x)=12π\theta = \begin{bmatrix}\mu / \sigma^2 \\ \\ -1/2\sigma^2\end{bmatrix},\quad T(x) = \begin{bmatrix} x \\ \\ x^2\end{bmatrix},\quad A(\theta) = \frac{\mu^2}{2\sigma^2}+\log\sigma=-\frac{\theta_1^2}{4\theta_2}-\frac12\log(-2\theta_2),\quad h(x)=\frac{1}{\sqrt{2\pi}}

可证高斯分布属于多参数指数族。

指数族的性质#

充分统计量#

关于充分统计量的理解,在本文之外,还可以参看 知乎专栏博客。这些材料对理解内容也会有很大帮助。本文的笔记内容,也部分来自于这些材料。

X1,,XnX_1,\cdots,X_nXX 的一组样本。在观测前,样本X1,,XnX_1,\cdots,X_n 是随机变量,在观测后,样本 X1,,XnX_1,\cdots,X_n 是具体的值。

数理统计的角度,我们希望通过样本,来推断原分布。充分统计量 Sufficient Statistic 作为统计量 Statistic 是定义在样本空间上的可测函数,记作 T(X1,,X2)T(X_1,\cdots, X_2) 很多地方也直接写作 T(X)T(X)。作为统计量,它缩减了原来随机变量包含的信息。

比如说,求样本均值时,样本值的顺序,是我们不关心的信息。

对于一组样本,样本本身会存在一个联合概率密度函数,可以记作 f(x)f(x) 。如果这个分布本身不存在参数(或参数已知),那么这个函数本质上就刻画了这一组样本包含的所有信息。

上述的联合概率密度函数,若存在未知的参数 θ\theta,则记作 f(x;θ)f(x;\theta)fθ(x)f_\theta(x)。给定统计量 TT 的值 T=tT=t,如果对应的条件分布 Fθ(XT=t)F_\theta(X|T=t) 是一个与未知参数 θ\theta 无关的分布(也就是说,是一个确定的分布)那么这个统计量 TT 就是一个充分统计量 Sufficient Statistic

充分统计量,保留了关于参数 θ\theta 的全部有用信息,消除了无用信息。

在充分统计量的基础上,我们更近一步,介绍极小充分统计量 Minimum Sufficient Statistic。在直觉上,我们肯定希望充分统计量的形式越简单越好,极小充分统计量的定义就是这么来的。

如果 T=T(X)T^\star = T^\star(X) 是一个充分统计量,对于任意的充分统计量 T=T(X)T=T(X),存在可测函数 φ\varphi,使得 T=φ(T)T^\star = \varphi(T) 那么, TT^\star 是极小充分统计量。

这个定义的逻辑在于:如果 TT^\star 是充分统计量,那么 TT 一定是充分统计量。

导数与期望#

学习期望时,我们知道,求解期望是在算一个积分。但指数分布族的特殊性质能将期望与导数联系起来。而求导一般比积分简单,因此我们会更喜欢导数。

我们对累积量生成函数 Cumulant Generating Function A(θ)A(\theta) 求一阶导,可以得到充分统计量 TT 的期望。

A(θ)θT=θT{logXh(x)exp{θ,T(x)} dx}=XT(x)h(x)exp{θ,T(x)} dxXh(x)exp{θ,T(x)} dx=1Z(θ)XT(x)h(x)exp{θ,T(x)} dx=XT(x)h(x)exp{θ,T(x)A(θ)} dx=XT(x)fX(x;θ) dx=E[T(X)]\begin{align*} \frac{\partial A(\theta)}{\partial \theta^T} &= \frac{\partial}{\partial \theta^T}\left\{ \log \int_Xh(x)\exp\{\langle\theta, T(x)\rangle\}~{\rm d}x\right\} \\ &= \frac{\int_X T(x)h(x)\exp\{\langle\theta, T(x)\rangle\}~{\rm d}x}{\int_Xh(x)\exp\{\langle\theta, T(x)\rangle\}~{\rm d}x} \\ &=\frac{1}{Z(\theta)} \int_X T(x)h(x)\exp\{\langle\theta, T(x)\rangle\}~{\rm d}x\\ &=\int_X T(x)h(x)\exp\{\langle\theta,T(x)\rangle - A(\theta)\}~{\rm d}x \\ &=\int_X T(x)f_\mathbf{X}(x;\theta)~{\rm d}x \\ &=\mathbb E[T(X)] \end{align*}

这里的公式比较绕,有几个点需要注意:

  1. 为什么对于 θT\theta^T 求导?可以简单理解成为了保证应用求导的链式法则时, θ,T(x)\langle\theta,T(x) \rangle 求导出来的东西是 T(x)T(x) 而不是 T(x)TT(x)^T
  2. 为什么求导符号和积分符号可以换?这里满足勒贝格控制收敛定理 Dominated Convergence Theorem
  3. 为什么式子里又多出了一个 A(θ)A(\theta) ?我们发现分母部分正好可以提出一个配分函数 Z(θ)Z(\theta),这个量是与积分变量 xx 无关的,因此我们可以按照指数的运算法则,把它移进去。这一步可以参考上文的 几种等价形式 那个小节。
  4. 最后一步怎么变成期望的?因为 fX(x;θ)f_\mathbf{X}(x;\theta) 是概率分布。

导数与方差#

对累积量生成函数 A(θ)A(\theta) 求二阶导,可以得到充分统计量 TT 的方差。

θ(A(θ)θT)=θXT(x)h(x)exp{θ,T(x)A(θ)} dx=XT(x)h(x)exp{θ,T(x)A(θ)}(T(x)TθA(θ)) dx=XT(x)(T(x)θTA(θ))Th(x)exp{θ,T(x)A(θ)} dx=XT(x)(T(x)E[T(X)])Th(x)exp{θ,T(x)A(θ)} dx=XT(x)T(x)Th(x)exp{θ,T(x)A(θ)} dxE[T(X)]TXT(x)h(x)exp{θ,T(x)A(θ)} dx=E[T(X)T(X)T]E[T(X)]E[T(X)]T=Var[T(X)]\begin{align*} \frac{\partial}{\partial \theta}(\frac{\partial A(\theta)}{\partial \theta^T}) &= \frac{\partial}{\partial \theta} \int_X T(x)h(x)\exp\{\langle\theta,T(x)\rangle - A(\theta)\}~{\rm d}x \\ &= \int_X T(x)h(x)\exp\{\langle\theta,T(x)\rangle - A(\theta)\}\left(T(x)^T - \frac{\partial}{\partial \theta}A(\theta)\right)~{\rm d}x \\ &= \int_X T(x)\left(T(x) - \frac{\partial}{\partial \theta^T}A(\theta)\right)^T h(x)\exp\{\langle\theta,T(x)\rangle - A(\theta)\}~{\rm d}x \\ &= \int_X T(x)\left(T(x) - \mathbb E[T(X)]\right)^T h(x)\exp\{\langle\theta,T(x)\rangle - A(\theta)\}~{\rm d}x \\ &= \int_X T(x)T(x)^T h(x)\exp\{\langle\theta,T(x)\rangle - A(\theta)\}~{\rm d}x \\ &\quad- \mathbb E[T(X)]^T \int_X T(x) h(x)\exp\{\langle\theta,T(x)\rangle - A(\theta)\}~{\rm d}x \\ &= \mathbb E[T(X)T(X)^T]-\mathbb E[T(X)]\cdot\mathbb E[T(X)]^T \\ &= Var[T(X)] \end{align*}

与上一节类似,这里也用到了导数积分互换,具体细节可以参看勒贝格控制收敛定理。

关于矩阵、向量的求导和转置,可以参看 博客园博客 。博客以及文中的引用链接给出了详细的解释。

参数化#

所谓参数化 parameterization 意味着用参数来表示。

如果指数族的参数 θ\theta 的元素是线性无关的,充分统计量 T(x)T(x) 的元素也是线性无关的,那么我们可以称这个指数族为最小指数族 minimal exponential family

似乎没有 minimal exponential family 的对应中文翻译,所以粗暴字面翻译最小指数族。但感觉翻译成最简指数族可能更加贴切一点。原因如下:
对于那些非 minimal 的指数族 (non-minimal) 的指数族,我们可以通过某种合适的参数替换或参数变换,得到一个最小指数族。

最小指数族对数配分函数 A(θ)A(\theta) 是严格的凸函数,满足 Fenchel's inequality。在介绍 Fenchel 不等式之前,首先引入凸共轭。

参考 维基百科
凸共轭 Convex Conjugate(也称 Fenchel Conjugate)
对于原空间 XX 上的扩充实值函数 extended real-valued function

f:XR  {,+}f: X\rightarrow\mathbb R~\cup~\{-\infty, +\infty\}

它在对偶空间 dual space XX^* 上的共轭函数 conjugate function 记作

f=XR  {,+}f^*=X^*\rightarrow\mathbb R~\cup~\{-\infty, +\infty\}

我们定义对偶空间中的点 xXx^*\in X^* 与原空间中的点 xXx\in X 的对应关系是:

f(x)=sup{x,xf(x)}f^*(x^*)=\sup\{\langle x^*,x\rangle-f(x)\}

其中,sup\supsupremum\rm supremum,即最小上界(上确界)。还有 inf(infimum)\inf (\rm infimum) 是最大下界(下确界)与 max(maximum)\max \rm(maximum)min(minimum)\min (\rm minimum) 的区别在于:
CSDN 博客知乎专栏

  1. 实值函数的必有上确界 / 下确界(一定能取到)。但不一定有最大值或最小值(可能最大 / 最小值点在定义上会取不到)。比如 f(x)=sinxxf(x)=\frac{\sin x}{x}
  2. 如果最大值 / 最小值能取到,就是上确界 / 下确界。

对于 A(θ)A(\theta) 来说,其凸共轭 A(θ)A^*(\theta^*)

A(θ)=sup{θ,θA(θ)}A^*(\theta^*) = \sup\{\langle\theta^*,\theta\rangle-A(\theta)\}

我们定义 μ=E[T(X)]\mu = \mathbb E[T(X)],于是 θT(θ,θA(θ))=θμ\dfrac{\partial}{\partial\theta^T}\left(\langle\theta^*,\theta\rangle-A(\theta)\right) = \theta^*-\mu

因此,当 θ=μ\theta^*=\mu 时,导数值为零,取到上确界。对应凸共轭为 A(μ)=μ,θA(θ)A^*(\mu)=\langle\mu,\theta\rangle-A(\theta),稍作变形我们得到

A(μ)+A(θ)=μ,θA^*(\mu) +A(\theta) = \langle\mu,\theta\rangle

Fenchel 不等式 Fenchel's inequality
另一方面,根据 Fenchel 不等式,任意 xX, xXx\in X,~x^*\in X^*

f(x)+f(x)x,xf(x)+f^*(x^*)\geq\langle x^*, x\rangle

由于 μA(θ)\mu\in\partial A(\theta) 上式取到等号。

均值表示法 指数族可以采用标准参数化 canonical parameterization 来表示,也可以以均值参数化 mean parameterization 来表示。因为 θ\theta 与均值 μ\mu 是一一对应的。即,既可以看做是 θ\theta 的函数,也可以看做是均值 μ\mu 的函数。

统计推断#

最大似然估计求总体均值#

首先回顾一下最大似然估计 Maximum Likelihood Estimation 的理念。

有一个未知的分布,我们有一系列的样本观测值。于是我们要拿着这些样本观测值,去反猜最有可能的分布。这就出现了两个问题:

  1. 模型确定吗?一般为了简化问题,会告知模型。实际问题中,如果没有告知模型,可能需要逐个模型去尝试。
  2. 参数确定吗?参数是不确定的。如果模型已知,那么一般操作是以这一组样本观测值去拟合模型,然后反推参数。

以最大似然估计求总体均值 μ\mu 。步骤:

  1. 已知 n 次重复采样的独立同分布样本观测值构成集合 D=(x1,x2,,xN)\mathcal D=(x_1,x_2,\cdots,x_N)
  2. 写出似然函数。方法是直接把这些样本值带入概率密度函数,并将结果相乘。
L(θD)=i=1Nf(xi;θ)=i=1Nh(xi)exp{η(θ),T(xi)A(θ)}L(\theta|\mathcal D) =\prod_{i=1}^N f(x_i;\theta) = \prod_{i=1}^N h(x_i)\exp\{\langle\eta(\theta), T(x_i)\rangle-A(\theta)\}
  1. 对似然函数取对数,并求导,得到 score function
l(θD)=logL(θD)=log(i=1Nh(xi))+θT(i=1NT(xi))NA(θ)θl(θD)=i=1NT(xi)NθA(θ)\begin{align*} &l(\theta|\mathcal D) = \log L(\theta|\mathcal D) = \log\left(\prod_{i=1}^N h(x_i)\right) + \theta^T\left( \sum_{i=1}^N T(x_i) \right) - NA(\theta) \\ & \nabla_\theta l(\theta|\mathcal D)= \sum_{i=1}^N T(x_i) - N\nabla_\theta A(\theta) \end{align*}
  1. 令导数为 0,解似然方程。
θl(θD)=0θA(θ^)=1Ni=1NT(xi)\nabla_\theta l(\theta|\mathcal D) = 0 \quad \Longrightarrow \quad \nabla_\theta A(\hat\theta) = \frac1N\sum_{i=1}^N T(x_i)

最大似然估计,本质是让似然函数取到极大值。但也有特殊情况:

  1. 如果对数似然函数单调,导致导数零点不存在
  2. 或由于样本太少,导致导数零点虽然存在但取不到的情况发生

一般会取端点值。

我们定义总体均值 μ=E[T(X)]\mu = \mathbb E[T(X)],结合上式,我们得到

μ^MLE=E[T(X)] = θA(θ^)=1Ni=1NT(xi)\hat\mu_{MLE} = \mathbb E[T(X)] ~{\color{red} = }~ \nabla_\theta A(\hat\theta) = \frac1N\sum_{i=1}^N T(x_i)

这个等式(红色等号)之所以能够成立,我们在上面的 导数与期望 小节中已经证明。

μ^MLE\hat\mu_{MLE} 是无偏的。因为

E[μ^MLE]=1Ni=1NE[T(Xi)]=1NNμ=μ\mathbb E [\hat\mu_{MLE}] = \frac1N\sum_{i=1}^N\mathbb E[T(X_i)] = \frac1N N\mu = \mu

μ^MLE\hat\mu_{MLE} 是有效的。可以证明 μ^MLE\hat\mu_{MLE} 是最小方差无偏估计 uniformly minimum-variance unbiased estimator (UMVUE)

上面我们说到,对数似然函数的一阶导数也称为 score function 记作

S(X;θ)=θl(X;θ)=θlogL(X;θ)=θlogi=1Nf(xi;θ)=i=1NT(xi)NθA(θ)\begin{align*} S(X;\theta) &= \nabla_\theta l(X;\theta) = \nabla_\theta\log L(X;\theta) = \nabla_\theta\log \prod_{i=1}^N f(x_i;\theta) \\ & = \sum_{i=1}^N T(x_i) - N\nabla_\theta A(\theta) \end{align*}

其中, XX 是样本序列 {X1,X2,,Xn}\{X_1,X_2,\cdots, X_n\},对应的样本观测值为 {x1,x2,,xn}\{x_1, x_2,\cdots,x_n\}

参考 知乎问答,我们引入费舍尔信息 Fisher Information。费舍尔信息是 score function 的二阶矩 second moment。

I(θ)=E[S2(X;θ)]I(\theta) = \mathbb E[S^2(X;\theta)]

Fisher Information 是用于衡量参数估计的精确度的。
N 次观测得到的 Fisher Information 是单次观测得到的 Fisher Information 的 N 倍
后文中我们将以单次观测的 Fisher Information 为例

score function 是关于 θ\theta 的函数,显然这个 fisher information 矩阵也是关于 θ\theta 的。 参考 维基百科网络博客 易证明:

E[S(X;θ)]=XS(X;θ)f(x;θ) dx=Xθf(x;θ)f(x;θ)f(x;θ) dx=θXf(x;θ) dx=θ1=0\begin{align*} \mathbb E[S(X;\theta)] & = \int_X S(X;\theta) f(x;\theta) ~{\rm d}x = \int_X\frac{\frac{\partial}{\partial \theta} f(x;\theta)}{f(x;\theta)}f(x;\theta) ~{\rm d}x\\ &=\color{red} \frac{\partial}{\partial \theta}\int_X f(x;\theta)~{\rm d}x = \frac{\partial}{\partial \theta} 1 = 0 \end{align*}

此处根据勒贝格收敛定理,导数与积分发生互换。
离散情况将积分号替换成求和即可。可能会产生一个 N 倍关系。

因此 I(θ)=E[S2(X;θ)]E2[S(X;θ)]=Var[S(X;θ)]I(\theta) = \mathbb E[S^2(X;\theta)] - \mathbb E^2[S(X;\theta)] = Var[S(X;\theta)]。即 fisher information 是 score function 的方差。

由于此处 S(X;θ)S(X;\theta) 是二阶可导的,因此可以证明:

E[S2(X;θ)]=E[2θ2logL(X;θ)]\begin{align*} \mathbb E[S^2(X;\theta)] =-\mathbb E\left[\frac{\partial^2}{\partial\theta^2}\log L(X;\theta)\right] \end{align*}

证明过程类似,由于

E[2θ2logL(X;θ)]=X2θ2logL(X;θ)f(x;θ) dx=XθS(X;θ)f(x;θ) dx=Xθ(θf(x;θ)f(x;θ))f(x;θ) dx=X(2θ2f(x;θ)f(x;θ)(θf(x;θ)f(x;θ))2)f(x;θ) dx=0X(θlogL(X;θ))2f(x;θ) dx=XS2(X;θ)f(x;θ) dx=E[S2(X;θ)]\begin{align*} \mathbb E\left[\frac{\partial^2}{\partial\theta^2}\log L(X;\theta)\right] &= \int_X \frac{\partial^2}{\partial\theta^2}\log L(X;\theta)f(x;\theta)~{\rm d}x = \int_X\frac{\partial}{\partial\theta}S(X;\theta)f(x;\theta)~{\rm d}x \\ &=\int_X\frac{\partial}{\partial\theta}\left(\frac{\frac{\partial}{\partial\theta}f(x;\theta)}{f(x;\theta)}\right)f(x;\theta)~{\rm d}x\\ &=\int_X\left(\frac{\frac{\partial^2}{\partial\theta^2}f(x;\theta)}{f(x;\theta)}-\left(\frac{\frac{\partial}{\partial\theta}f(x;\theta)}{f(x;\theta)}\right)^2\right)f(x;\theta)~{\rm d}x \\ &={\color{red} 0}-\int_X\left(\frac{\partial}{\partial\theta}\log L(X;\theta)\right)^2 f(x;\theta)~{\rm d}x \\ &=-\int_X S^2(X;\theta)f(x;\theta)~{\rm d}x\\ &=-\mathbb E[S^2(X;\theta)] \end{align*}

红色部分积分为 0 是因为积分号与二阶导互换后,求导得 0
离散情况将积分号替换成求和即可。可能会产生一个 N 倍关系。

我们于此总结一下 Fisher Information 的几个等价替换式

I(θ)=E[S2(X;θ)]=E[2θ2logL(X;θ)]=E[θS(X;θ)]=Var[S(X;θ)]I(\theta) = \mathbb E[S^2(X;\theta)] = -\mathbb E\left[\frac{\partial^2}{\partial\theta^2}\log L(X;\theta)\right] = -\mathbb E\left[\frac{\partial}{\partial\theta}S(X;\theta)\right] = Var[S(X;\theta)]

另一方面,我们有 L(θ)=fX(x;θ)=h(x)exp{θ,T(x)A(θ)}L(\theta) = f_X(x;\theta) = h(x)\exp\{\langle\theta, T(x)\rangle-A(\theta)\}
取对数然后求二阶导后,我们得到:

2θ2logL(X;θ)=2θ2A(θ)\frac{\partial^2}{\partial\theta^2}\log L(X;\theta) = -\frac{\partial^2}{\partial\theta^2} A(\theta)

所以,可以得到:

I(θ)=E[2θ2logL(X;θ)]=E[2θ2A(θ)]=Var[T(X)]I(\theta) = -\mathbb E\left[\frac{\partial^2}{\partial\theta^2}\log L(X;\theta)\right] = -\mathbb E\left[-\frac{\partial^2}{\partial\theta^2} A(\theta) \right] =Var[T(X)]

我们发现,自然参数 θ\theta 的 Fisher Information 正好是充分统计量的方差 Var[T(X)]Var[T(X)]

另一方面,

【未完待续】

Loading...
Ownership of this post data is guaranteed by blockchain and smart contracts to the creator alone.