参数估计有两种形式:点估计与区间估计。本文选择几种常用的点估计方法作一些讨论。
用于估计未知参数的统计量称为点估计(量)。参数 θ \theta θ 的估计量常用 θ ^ = θ ^ ( x 1 , x 2 , … , x n ) \hat{\theta} = \hat{\theta}(x_{1},x_{2}, \dots, x_{n}) θ^=θ^(x1,x2,…,xn) 表示,参数 θ \theta θ 的可能取值范围称为参数空间,记为 Θ = { θ } \Theta = \{\theta\} Θ={ θ}。
文章目录
最大似然估计
最大似然估计,即对似然函数最大化,其关键是从样本 x x x 和含有位置参数 θ \theta θ 的分布 p ( x , θ ) p(x,\theta) p(x,θ) 获得似然函数。设 x = ( x 1 , x 2 , … , x n ) x=(x_{1},x_{2},\dots,x_{n}) x=(x1,x2,…,xn) 是来自含有未知参数的某分布 p ( x , θ ) p(x,\theta) p(x,θ) 的一个样本,那么其联合分布为:
p ( x , θ ) = ∏ i = 1 n p ( x i , θ ) p(x,\theta) = \prod_{i=1}^{n}p(x_{i},\theta) p(x,θ)=i=1∏np(xi,θ) 其中 p ( x i , θ ) p(x_{i},\theta) p(xi,θ) 在连续场合是指密度函数在 x i x_{i} xi 处的值,在离散场合为分布列中的一个概率 P θ ( X = x i ) P_{\theta}(X=x_{i}) Pθ(X=xi) 。对样本分布 p ( x , θ ) p(x,\theta) p(x,θ) 我们知道:
- 样本如何产生?先有 θ \theta θ 后有 x x x,即先有一个给定的 θ \theta θ 的值 θ 0 \theta_{0} θ0,然后由分布 p ( x , θ 0 ) p(x,\theta_{0}) p(x,θ0) 经过随机抽样产生样本观察值 x x x。
- 如今我们有了 x x x 如何追溯参数 θ 0 \theta_{0} θ0 呢?当给定样本观察值 x x x 时样本分布 p ( x , θ ) p(x,\theta) p(x,θ) 仅是 θ \theta θ 的函数,可记为 L ( θ , x ) L(\theta,x) L(θ,x) 或 L ( θ ) L(\theta) L(θ),并称其为似然函数。对于不同的 θ 1 , θ 2 ∈ Θ \theta_{1},\theta_{2}\in\Theta θ1,θ2∈Θ,可使得样本观察值 x x x 出现的机会不同。若 L ( θ 1 ) > L ( θ 2 ) L(\theta_{1}) > L(\theta_{2}) L(θ1)>L(θ2),表明 θ 1 \theta_{1} θ1 会使 x x x 出现的机会比 θ 2 \theta_{2} θ2 更大些,即 θ 1 \theta_{1} θ1 比 θ 2 \theta_{2} θ2 更像真值 θ 0 \theta_{0} θ0。也就是说 L ( θ ) L(\theta) L(θ) 成为了度量 θ \theta θ 更像真值的程度,其值越大越像。按此思路,在参数空间 Θ \Theta Θ 中使 L ( θ ) L(\theta) L(θ) 最大的 θ ^ \hat{\theta} θ^ 就是最像 θ 0 \theta_{0} θ0 的真值,这个 θ ^ \hat{\theta} θ^ 就是 θ \theta θ 的最大似然估计。
这里给出两个实例。
1.伯努利分布实例
假设 P ( X = 1 ) = p , P ( X = 0 ) = 1 − p P(X=1)=p,P(X=0)=1-p P(X=1)=p,P(X=0)=1−p 综合起来就有
P ( X ) = p X ( 1 − p ) 1 − X P(X)=p^{X}(1-p)^{1-X} P(X)=pX(1−p)1−X
此时如果有一组数据 D D D 是从这个随机变量中采样得到的,那么就有
m a x p log P ( D ) = max p log ∏ i N P ( D i ) = max p ∑ i N log P ( D i ) = max p ∑ i N [ D i log p + ( 1 − D i ) log ( 1 − p ) ] \begin{aligned} \ max_{p}\log P(D)&= \max_{p}\log\prod_{i}^{N}P(D_{i}) \\ &=\max_{p}\sum_{i}^{N}\log P(D_{i}) \\ &=\max_{p}\sum_{i}^{N}[D_{i}\log p+(1-D_{i})\log(1-p)] \end{aligned} maxplogP(D)=pmaxlogi∏NP(Di)=pmaxi∑NlogP(Di)=pmaxi∑N[Dilogp+(1−Di)log(1−p)]
对上式求导,则有
∇ p max p log P ( D ) = ∑ i N [ D i 1 p + ( 1 − D i ) 1 p − 1 ] \nabla_{p}\max_{p}\log P(D)=\sum_{i}^{N}[D_{i}\frac{1}{p}+(1-D_{i})\frac{1}{p-1}] ∇ppmaxlogP(D)=i∑N[Dip1+(1−Di)p−11]
求极值,令导数为 0 0 0,就有
∑ i N [ D i 1 p + ( 1 − D i ) 1 p − 1 ] = 0 ∑ i N [ D i ( p − 1 ) + ( 1 − D i ) p ] = 0 ∑ i N ( p − D i ) = 0 p = 1 N ∑ i N D i \begin{aligned} & \sum_{i}^{N}[D_{i}\frac{1}{p}+(1-D_{i})\frac{1}{p-1}]=0 \\ & \sum_{i}^{N}[D_{i}(p-1)+(1-D_{i})p]=0 \\ & \sum_{i}^{N}(p-D_{i})=0 \\ & p=\frac{1}{N}\sum_{i}^{N}D_{i} \end{aligned} i∑N[Dip1+(1−Di)p−11]=0i∑N[Di(p−1)+(1−Di)p]=0i∑N(p−Di)=0p=N1i∑NDi
即全部采样的平均值。
2.高斯分布实例
令 p ( x ) = 1 2 π σ 2 e − ( x − μ ) 2 2 σ 2 p(x)=\frac{1}{\sqrt{2\pi\sigma^{2}}}e^{-\frac{(x-\mu)^{2}}{2\sigma^{2}}} p(x)=2πσ2 1e−2σ2(x−μ)2,采用同样的方法有
max p log P ( D ) = max p log ∏ i N P ( D i ) = max p ∑ i N log P ( D i ) = max p ∑ i N [ − 1 2 log ( 2 π σ 2 ) − ( D i − μ ) 2 2 σ 2 ] = max [ − N 2 log ( 2 π σ 2 ) − 1 2 σ 2 ∑ i N ( D i − μ ) 2 ] \begin{aligned} \max_{p}\log P(D) &= \max_{p}\log\prod_{i}^{N}P(D_{i}) \\ &= \max_{p}\sum_{i}^{N}\log P(D_{i}) \\ &= \max_{p}\sum_{i}^{N}[-\frac{1}{2}\log(2\pi\sigma^{2})-\frac{(D_{i}-\mu)^{2}}{2\sigma^{2}}] \\ &= \max[-\frac{N}{2}\log(2\pi\sigma^{2})-\frac{1}{2\sigma^{2}}\sum_{i}^{N}(D_{i}-\mu)^{2}] \end{aligned} pmaxlogP(D)=pmaxlogi∏NP(Di)=pmaxi∑NlogP(Di)=pmaxi∑N[−21log(2πσ2)−2σ2(Di−μ)2]=max[−2Nlog(2πσ2)−2σ21i∑N(Di−μ)2]
此处包含两个参数,分别估计。
首先对 μ \mu μ 求导,有
∂ max μ log P ( D ) ∂ μ = − 1 σ 2 ∑ i N ( μ − D i ) \frac{\partial\max_{\mu}\log P(D)}{\partial \mu} = -\frac{1}{\sigma^{2}}\sum_{i}^{N}(\mu-D_{i}) ∂μ∂maxμlogP(D)=−σ21i∑N(μ−Di)
令导数为 0 0 0,有
− 1 σ 2 ∑ i N ( μ − D i ) = 0 , μ = 1 N ∑ i N D i -\frac{1}{\sigma^{2}}\sum_{i}^{N}(\mu-D_{i})=0,\quad \mu=\frac{1}{N}\sum_{i}^{N}D_{i} −σ21i∑N(μ−Di)=0,μ=N1i∑NDi
注意很容易看出这个结果与最小二乘计算完全相同,实质上最小二乘法可以与极大似然法中假定误差遵循正态分布的特殊情况相对应。
其次对 σ 2 \sigma^{2} σ2 求导,有
∂ max σ 2 log P ( D ) ∂ σ 2 = − N 2 σ 2 + 1 2 σ 4 ∑ i N ( D i − μ ) 2 \frac{\partial\max_{\sigma^{2}}\log P(D)}{\partial\sigma^{2}} = -\frac{N}{2\sigma^{2}}+\frac{1}{2\sigma^{4}}\sum_{i}^{N}(D_{i}-\mu)^{2} ∂σ2∂maxσ2logP(D)=−2σ2N+2σ41i∑N(Di−μ)2
令导数为 0 0 0,有
− N 2 σ 2 + 1 4 σ 4 ∑ i N ( D i − μ ) 2 = 0 -\frac{N}{2\sigma^{2}}+\frac{1}{4\sigma^{4}}\sum_{i}^{N}(D_{i}-\mu)^{2}=0 −2σ2N+4σ41i∑N(Di−μ)2=0
σ 2 = 1 N ∑ i N ( D i − μ ) 2 \sigma^{2} = \frac{1}{N}\sum_{i}^{N}(D_{i}-\mu)^{2} σ2=N1i∑N(Di−μ)2
可见最终计算结果与期望方差计算方式完全一致。注意最大似然估计并不一定具有无偏性。
对似然函数添加或剔除一个与参数 θ \theta θ 无关的量 c ( x ) > 0 c(x)>0 c(x)>0,不影响寻求最大似然估计的最终结果,故 c ( x ) L ( θ ) c(x)L(\theta) c(x)L(θ) 仍然是 θ \theta θ 的似然函数。例如,对于正态分布而言:
L ( μ , σ 2 ) = ∏ i = 1 n 1 2 π σ 2 e − ( x i − μ ) 2 2 σ 2 ∝ ( σ 2 ) − n 2 exp { − 1 2 σ 2 ∑ i = 1 n ( x i − μ ) 2 } L(\mu,\sigma^{2}) = \prod_{i=1}^{n}\frac{1}{\sqrt{2\pi}\sigma^{2}}e^{-\frac{(x_{i}-\mu)^{2}}{2\sigma^{2}}} \propto (\sigma^{2})^{-\frac{n}{2}}\exp\left\{-\frac{1}{2\sigma^{2}}\sum_{i=1}^{n}(x_{i}-\mu)^{2}\right\} L(μ,σ2)=i=1∏n2π σ21e−2σ2(xi−μ)2∝(σ2)−2nexp{ −2σ21i=1∑n(xi−μ)2}
不变原理: 设 X ∼ p ( x , θ ) , θ ∈ Θ X\sim p(x,\theta), \theta\in\Theta X∼p(x,θ),θ∈Θ,若 θ \theta θ 的最大似然估计为 θ ^ \hat{\theta} θ^ ,则对任意函数 γ = g ( θ ) \gamma=g(\theta) γ=g(θ), γ \gamma γ 的最大似然估计为 γ ^ = g ( θ ^ ) \hat{\gamma}=g(\hat{\theta}) γ^=g(θ^)。
贝叶斯估计
统计学中有两个主要学派:频率学派(又称经典学派)和贝叶斯学派。前述最大似然估计属于经典统计学范畴。频率学派利用总体信息和样本信息进行统计推断,贝叶斯学派与之的区别在于还用到了先验信息。
贝叶斯学派最基本的观点是:任一未知量 θ \theta θ 都可以看做随机变量,可用一个概率分布区描述,这个分布称为先验分布 (记为 π ( θ ) \pi(\theta) π(θ))。因为任一未知量都有不确定性,而在表述不确定性地程度时,概率与概率分布是最好的语言。依赖于参数 θ \theta θ 的密度函数在经典统计学中记为 p ( x , θ ) p(x,\theta) p(x,θ),它表示参数空间 Θ \Theta Θ 中不同的 θ \theta θ 对应不同的分布。在贝叶斯统计中应记为 p ( x ∣ θ ) p(x|\theta) p(x∣θ) ,表示随机变量 θ \theta θ 给定某个值时, X X X 的条件密度函数。
从贝叶斯观点看,样本 x x x 的产生要分两步进行:首先,设想从先验分布 π ( θ ) \pi(\theta) π(θ) 中产生一个样本 θ ′ \theta' θ′ ,这一步人是看不到的,所以是“设想”;再从 p ( x ∣ θ ′ ) p(x|\theta') p(x∣θ′) 中产生一个样本 x = ( x 1 , x 2 , x 3 , … , x n ) x=(x_{1},x_{2},x_{3},\dots,x_{n}) x=(x1,x2,x3,…,xn) 。这时样本 x x x 的联合条件密度函数为:
p ( x ∣ θ ′ ) = ∏ i = 1 n p ( x i ∣ θ ′ ) p(x|\theta')=\prod_{i=1}^{n}p(x_{i}|\theta') p(x∣θ′)=i=1∏np(xi∣θ′) 这个联合分布综合了总体信息和样本信息,又称为似然函数。它与极大似然估计中的似然函数没有什么区别。 θ ′ \theta' θ′ 仍然是未知的,它是按照先验分布 π ( θ ) \pi(\theta) π(θ) 产生的,为了把先验信息综合进去,不能只考虑 θ ′ \theta' θ′,对 θ \theta θ 的其它值发生的可能性也要加以考虑,故要用 π ( θ ) \pi(\theta) π(θ) 进行综合。这样一来,样本 x x x 和参数 θ \theta θ 的联合分布为:
h ( x , θ ) = p ( x ∣ θ ) π ( θ ) h(x,\theta)=p(x|\theta)\pi(\theta) h(x,θ)=p(x∣θ)π(θ) 这个联合分布综合了总体信息、样本信息和先验信息。
我们的核心目标是对 θ \theta θ 进行估计,若把 h ( x , θ ) h(x,\theta) h(x,θ) 作如下分解:
h ( x , θ ) = π ( θ ∣ x ) m ( x ) h(x,\theta) = \pi(\theta|x)m(x) h(x,θ)=π(θ∣x)m(x) 其中 m ( x ) m(x) m(x) 是 X X X 的边际密度函数:
m ( x ) = ∫ Θ h ( x , θ ) d θ = ∫ Θ p ( x ∣ θ ) π ( θ ) d θ m(x) = \int_{\Theta}h(x,\theta)\mathrm{d}\theta = \int_{\Theta}p(x|\theta)\pi(\theta)\mathrm{d}\theta m(x)=∫Θh(x,θ)dθ=∫Θp(x∣θ)π(θ)dθ 它与 θ \theta θ 无关。因此,能用来对 θ \theta θ 进行估计的只有条件分布 π ( θ ∣ x ) \pi(\theta|x) π(θ∣x),它的计算公式是:
π ( θ ∣ x ) = h ( x , θ ) m ( x ) = p ( x ∣ θ ) π ( θ ) m ( x ) = p ( x ∣ θ ) π ( θ ) ∫ Θ p ( x ∣ θ ) π ( θ ) d θ \pi(\theta|x)=\frac{h(x,\theta)}{m(x)} = \frac{p(x|\theta)\pi(\theta)}{m(x)} = \frac{p(x|\theta)\pi(\theta)}{\int_{\Theta}p(x|\theta)\pi(\theta)\mathrm{d}\theta} π(θ∣x)=m(x)h(x,θ)=m(x)p(x∣θ)π(θ)=∫Θp(x∣θ)π(θ)dθp(x∣θ)π(θ) 这就是贝叶斯公式的密度函数形式。 这个条件分布称为 θ \theta θ 的后验分布,它集中了总体信息、样本信息和先验信息中有关 θ \theta θ 的一切信息。也可以说是总体和样本对先验分布 π ( θ ) \pi(\theta) π(θ) 作调整的结果,比先验分布更接近 θ \theta θ 的实际情况。上述公式是在 x x x 和 θ \theta θ 都是连续随机变量场合下的贝叶斯公式。其它场合下的贝叶斯公式如下:
- x x x 离散, θ \theta θ 连续: π ( θ ∣ x j ) = p ( x j ∣ θ ) π ( θ ) ∫ Θ p ( x j ∣ θ ) π ( θ ) d θ \pi(\theta|x_{j})=\frac{p(x_{j}|\theta)\pi(\theta)}{\int_{\Theta}p(x_{j}|\theta)\pi(\theta)\mathrm{d}\theta} π(θ∣xj)=∫Θp(xj∣θ)π(θ)dθp(xj∣θ)π(θ)
- x x x 连续, θ \theta θ 离散: π ( θ i ∣ x ) = p ( x ∣ θ i ) π ( θ i ) ∑ i p ( x ∣ θ i ) π ( θ i ) \pi(\theta_{i}|x) =\frac{p(x|\theta_{i})\pi(\theta_{i})}{\sum_{i}p(x|\theta_{i})\pi(\theta_{i})} π(θi∣x)=∑ip(x∣θi)π(θi)p(x∣θi)π(θi)
- x x x 离散, θ \theta θ 离散: π ( θ i ∣ x j ) = p ( x j ∣ θ i ) π ( θ i ) ∑ i p ( x j ∣ θ i ) π ( θ i ) \pi(\theta_{i}|x_{j}) =\frac{p(x_{j}|\theta_{i})\pi(\theta_{i})}{\sum_{i}p(x_{j}|\theta_{i})\pi(\theta_{i})} π(θi∣xj)=∑ip(xj∣θi)π(θi)p(xj∣θi)π(θi)
先验分布的确定十分关键,其原则有二:一是要根据先验信息;二是要使用方便,即在数学上处理方便。先验分布的确定有一些比较成熟的方法,如共轭先验分布法,此处不做详细讨论。
回到我们的核心目标,寻求参数 θ \theta θ 的估计 θ ^ \hat{\theta} θ^ 只需要从后验分布 π ( θ ∣ x ) \pi(\theta| x) π(θ∣x) 中合理提取信息即可。常用的提取方式是用后验均方误差准则,即选择这样的统计量
θ ^ = θ ^ ( x 1 , x 2 , … , x n ) \hat{\theta} = \hat{\theta}(x_{1},x_{2},\dots,x_{n}) θ^=θ^(x1,x2,…,xn) 使得后验均方误差达到最小,即
min M S E ( θ ^ ∣ x ) = min E θ ∣ x ( θ ^ − θ ) 2 \min\mathrm{MSE}(\hat{\theta} | x) =\min E^{\theta|x}(\hat{\theta}-\theta)^{2} minMSE(θ^∣x)=minEθ∣x(θ^−θ)2 这样的估计 θ ^ \hat{\theta} θ^ 称为 θ \theta θ 的贝叶斯估计,其中 E θ ∣ x E^{\theta|x} Eθ∣x 表示用后验分布 π ( θ ∣ x ) \pi(\theta|x) π(θ∣x) 求期望。求解上式并不困难,
KaTeX parse error: No such environment: split at position 7: \begin{̲s̲p̲l̲i̲t̲}̲ E^{\theta|x}(\… 这是关于 θ ^ \hat{\theta} θ^ 的二次三项式,二次项系数为正,必有最小值:
θ ^ = ∫ Θ θ π ( θ ∣ x ) d θ = E ( θ ∣ x ) \hat{\theta} = \int_{\Theta}\theta\pi(\theta|x)\mathrm{d}\theta=E(\theta|x) θ^=∫Θθπ(θ∣x)dθ=E(θ∣x) 也就是说,在均方误差准则下, θ \theta θ 的贝叶斯估计 θ ^ \hat{\theta} θ^ 就是 θ \theta θ 的后验期望 E ( θ ∣ x ) E(\theta|x) E(θ∣x)。
类似的可证,在已知后验分布为 π ( θ ∣ x ) \pi(\theta|x) π(θ∣x) 的情况下,参数函数 g ( θ ) g(\theta) g(θ) 在均方误差下的贝叶斯估计为 $\hat{g}(\theta)=E[g(\theta)|x] $。
贝叶斯公式中, m ( x ) m(x) m(x) 为样本的边际分布,它不依赖于 θ \theta θ ,在后验分布计算中仅起到一个正则化因子的作用,加入把 m ( x ) m(x) m(x) 省略,贝叶斯公式可改写为如下形式:
π ( θ ∣ x ) ∝ p ( x ∣ θ ) π ( θ ) \pi(\theta|x) \propto p(x|\theta)\pi(\theta) π(θ∣x)∝p(x∣θ)π(θ) 上式右边虽然不是 θ \theta θ 的密度函数或分布列,但在需要时利用正则化立即可以恢复密度函数或分布列原型。这时,可把上式右端称为后验分布的核,加入其中还有不含 θ \theta θ 的因子,仍可剔去,使核更为精炼。