本节内容摘自Steven M. Kay,《Fundamentals of Statistical Signal Processing: Estimation Theory》。
我们来个例子,如何用完备统计量找到MVU估计。
【例5.8】均匀噪声的均值
- 问题描述
我们有数据
x[n]=w[n],n=0,1,…,N−1(1)\tag{1} x[n]=w[n],\quad n=0,1,\ldots,N-1 x[n]=w[n],n=0,1,…,N−1(1)这里w[n]∼U(0,β)w[n]\sim {\mathcal U}(0,\beta)w[n]∼U(0,β),β>0\beta>0β>0。我们希望能够得到均值θ=β/2\theta=\beta/2θ=β/2的MVU估计。 - CRLB为何不适用?
从【定理3.1】我们知道,要使用CRLB,需要概率密度函数满足正则条件
E[∂lnp(x;θ)∂θ]=0.(2)\tag{2} {\rm E}\left[\frac{\partial \ln p({\bf x};\theta)}{\partial \theta}\right]=0. E[∂θ∂lnp(x;θ)]=0.(2)现在我们来证明,对于均匀分布噪声均值来说,正则条件(2)不满足。
下面我们来看看习题3.1。
顺便说一句,这本书的答案可见(https://download.csdn.net/download/qianlongchen/9433628?utm_source=iteye)。虽然是手写不太清楚,而且VIP才能下载,不过还是很解决问题啦。
【Problem 3.1】如果x[n]∼U(0,θ)x[n]\sim {\mathcal U}(0,\theta)x[n]∼U(0,θ),n=0,1,…,N−1n=0,1,\ldots,N-1n=0,1,…,N−1。试说明正则条件
E[∂lnp(x;θ)∂θ]=0
{\rm E}\left[\frac{\partial \ln p({\bf x};\theta)}{\partial \theta}\right]=0
E[∂θ∂lnp(x;θ)]=0不成立,因而CRLB不适用。
【解答】
p(x[n];θ)=1θ(u(x[n])−u(x[n]−θ))
p(x[n];\theta)=\frac{1}{\theta}\left(u(x[n])-u(x[n]-\theta) \right)
p(x[n];θ)=θ1(u(x[n])−u(x[n]−θ))这里的u(x)u(x)u(x)为阶跃函数。因此,我们可以得到
p(x;θ)=∏n=0N−1p(x[n];θ),
p({\bf x};\theta)=\prod_{n=0}^{N-1}p(x[n];\theta),
p(x;θ)=n=0∏N−1p(x[n];θ),显然,如果我们能够证明
E[∂lnp(x;θ)∂θ]≠0, for n=0,…,N−1
{\rm E}\left[\frac{\partial \ln p({\bf x};\theta)}{\partial \theta}\right]\ne 0,{\rm \ for \ }n=0,\ldots,N-1
E[∂θ∂lnp(x;θ)]=0, for n=0,…,N−1则获证。我们令y=x[n]y=x[n]y=x[n],有
p(y;θ)=1θ(u(y)−u(y−θ)),
p(y;\theta)=\frac{1}{\theta}\left(u(y)-u(y-\theta) \right),
p(y;θ)=θ1(u(y)−u(y−θ)),则可以得到p(y;θp(y;\thetap(y;θ的图像如下图所示。看图的时候,可以考虑yyy是固定值。显然应该有y=x[n]>0y=x[n]>0y=x[n]>0。这时,如果θ<y\theta<yθ<y,函数值为0;θ>y\theta>yθ>y,则有函数值为1/θ1/\theta1/θ。因此,有
E[∂lnp(y;θ)∂θ]=E[∂ln(1/θ)∂θ]=−1θ≠0.
{\rm E}\left[\frac{\partial \ln p({y};\theta)}{\partial \theta}\right]={\rm E}\left[\frac{\partial \ln (1/\theta)}{\partial \theta}\right]=-\frac{1}{\theta}\ne 0.
E[∂θ∂lnp(y;θ)]=E[∂θ∂ln(1/θ)]=−θ1=0.
- 下面我们来找到一个无偏估计,并且用充分统计量的完备性,来考察这个无偏估计是不是MVU估计。
我们很自然会想到用样本均值
θ^=1N∑n=0N−1x[n] \hat \theta=\frac{1}{N}\sum_{n=0}^{N-1}x[n] θ^=N1n=0∑N−1x[n]作为估计。这里x[n]∼U(0,2θ)x[n]\sim {\mathcal U}(0,2\theta)x[n]∼U(0,2θ),因此E(x[n])=θ{\rm E}(x[n])=\thetaE(x[n])=θ。对估计求均值,得到
E(θ^)=θ. {\rm E}(\hat \theta)=\theta. E(θ^)=θ.因此θ^\hat \thetaθ^为无偏估计,其方差为
var(θ^)=1Nvar(x[n])=β212N=θ23N.(5.9)\tag{5.9} {\rm var}(\hat \theta)=\frac{1}{N}{\rm var}(x[n])=\frac{\beta^2}{12N}=\frac{\theta^2}{3N}. var(θ^)=N1var(x[n])=12Nβ2=3Nθ2.(5.9)
(1) 找到完备的充分统计量
我们可以用阶跃函数来表示均匀分布的PDF,即
p(x[n];θ)=1β[u(x[n])−u(x[n]−β)],
p(x[n];\theta)=\frac{1}{\beta}[u(x[n])-u(x[n]-\beta)],
p(x[n];θ)=β1[u(x[n])−u(x[n]−β)],则可以得到对所有数据的PDF为
p(x;θ)=1βN∏n=0N−1[u(x[n])−u(x[n]−β)].
p({\bf x};\theta)=\frac{1}{\beta^N}\prod_{n=0}^{N-1}[u(x[n])-u(x[n]-\beta)].
p(x;θ)=βN1n=0∏N−1[u(x[n])−u(x[n]−β)].与Problem 3.1类似,我们可以得到
p(x;θ)={1βN0≤x[n]≤β,n=0,1,…,N−10otherwise
p({\bf x};\theta)=\left\{\begin{matrix} \frac{1}{\beta^N} & 0\le x[n]\le \beta,n=0,1,\ldots,N-1 \\0& {\rm otherwise}\end{matrix}\right.
p(x;θ)={βN100≤x[n]≤β,n=0,1,…,N−1otherwise进一步,可以得到p(x;θ)={1βNmaxx[n]<β,minx[n]>00otherwise
p({\bf x};\theta)=\left\{\begin{matrix} \frac{1}{\beta^N} & \max{x[n]}<\beta,\min x[n]>0 \\0& {\rm otherwise}\end{matrix}\right.
p(x;θ)={βN10maxx[n]<β,minx[n]>0otherwise因此
p(x;θ)=1βNu(β−maxx[n])u(minx[n])=g(T(x),θ)h(x),
p({\bf x};\theta)=\frac{1}{\beta^N}u(\beta-\max x[n])u(\min x[n])=g(T({\bf x}),\theta)h({\bf x}),
p(x;θ)=βN1u(β−maxx[n])u(minx[n])=g(T(x),θ)h(x),其中
g(T(x),θ)=1βNu(β−maxx[n]),h(x)=u(minx[n]).
g(T({\bf x}),\theta)=\frac{1}{\beta^N}u(\beta-\max x[n]),h({\bf x})=u(\min x[n]).
g(T(x),θ)=βN1u(β−maxx[n]),h(x)=u(minx[n]).
通过Neyman-Fisher因式分解定理,T(x)=maxx[n]T({\bf x})=\max x[n]T(x)=maxx[n]为充分统计量,并且可以证明其为完备的(证明略)。
(2) 找到T(x)T({\bf x})T(x)的函数,得到无偏估计
上面我们找到完备的充分统计量为T(x)=maxx[n]T({\bf x})=\max x[n]T(x)=maxx[n]。我们首先来确定它的均值。显然,TTT为顺序统计量(order statistics)。我们先来看累积分布函数
Pr{T≤ξ}=Pr{x[0]≤ξ,x[1]≤ξ,…,x[N−1]≤ξ}=∏n=0N−1Pr{x[n]≤ξ}=Pr{x[n]≤ξ}N.
\begin{aligned}
{\rm Pr}\{T\le \xi\}&=\Pr\{x[0]\le \xi,x[1]\le\xi, \ldots,x[N-1]\le \xi\}\\
&=\prod_{n=0}^{N-1}\Pr\{x[n]\le \xi\}\\
&=\Pr\{x[n]\le \xi\}^N.
\end{aligned}
Pr{T≤ξ}=Pr{x[0]≤ξ,x[1]≤ξ,…,x[N−1]≤ξ}=n=0∏N−1Pr{x[n]≤ξ}=Pr{x[n]≤ξ}N.因此得到PDF为
pT(ξ)=dPr{T≤ξ}dξ=NPr{x[n]≤ξ}N−1dPr{x[n]≤ξ}dξ.
\begin{aligned}
p_{T}(\xi)&=\frac{d{\rm Pr}\{T\le \xi\}}{d\xi}\\
&=N\Pr\{x[n]\le \xi\}^{N-1}\frac{d{\rm Pr}\{x[n]\le \xi\}}{d\xi}.
\end{aligned}
pT(ξ)=dξdPr{T≤ξ}=NPr{x[n]≤ξ}N−1dξdPr{x[n]≤ξ}.注意到dPr{x[n]≤ξ}dξ\frac{d{\rm Pr}\{x[n]\le \xi\}}{d\xi}dξdPr{x[n]≤ξ}为x[n]x[n]x[n]的概率密度函数px[n](ξ;θ)p_{x[n]}(\xi; \theta)px[n](ξ;θ),且
px[n](ξ;θ)={1β0<ξ<β0otherwise.
p_{x[n]}(\xi; \theta)=\left\{\begin{matrix}\frac{1}{\beta}&0<\xi<\beta\\0&{\rm otherwise}.\end{matrix} \right.
px[n](ξ;θ)={β100<ξ<βotherwise.积分后得到
Pr{x[n]≤ξ}={0ξ<0ξβ0<ξ<β1ξ>β.
\Pr\{x[n]\le\xi\}=\left\{\begin{matrix}0&\xi<0\\\frac{\xi}{\beta}&0<\xi<\beta\\1&{\xi>\beta}.\end{matrix} \right.
Pr{x[n]≤ξ}=⎩⎨⎧0βξ1ξ<00<ξ<βξ>β.最后,我们得到
pT(ξ)={0ξ<0N(ξβ)N−11β0<ξ<β0ξ>β.
p_T(\xi)=\left\{\begin{matrix}0&\xi<0\\
N\left(\frac{\xi}{\beta}\right)^{N-1}\frac{1}{\beta} & 0<\xi<\beta\\
0&{\xi>\beta}.\end{matrix} \right.
pT(ξ)=⎩⎪⎨⎪⎧0N(βξ)N−1β10ξ<00<ξ<βξ>β.因此,可以求得
E(T)=∫−∞∞ξpT(ξ)dξ=∫0βξN(ξβ)N−11βdξ=NN+1β=2NN+1θ
\begin{aligned}
{\rm E}(T)&=\int_{-\infty}^{\infty}\xi p_T({\xi})d\xi\\
&=\int_{0}^{\beta}\xi N\left(\frac{\xi}{\beta}\right)^{N-1}\frac{1}{\beta}d\xi\\
&=\frac{N}{N+1}\beta\\
&=\frac{2N}{N+1}\theta
\end{aligned}
E(T)=∫−∞∞ξpT(ξ)dξ=∫0βξN(βξ)N−1β1dξ=N+1Nβ=N+12Nθ为了得到无偏估计,我们令θ^=N+12NT\hat \theta=\frac{N+1}{2N}Tθ^=2NN+1T,因此最终无偏估计为
θ^=N+12Nmaxx[n].
\hat \theta=\frac{N+1}{2N}\max x[n].
θ^=2NN+1maxx[n].
有点出人意料的是,样本均值不是均匀分布噪声的MVU估计!
下面我们来看看二者方差情况。
(3)最小估计方差
我们可以得到最小估计方差为
var(θ^)=(N+12N)var(T),
{\rm var}(\hat \theta)=\left(\frac{N+1}{2N}\right){\rm var}(T),
var(θ^)=(2NN+1)var(T),其中
var(T)=∫0βξ2NξN−1βNdξ−(NβN+1)2=Nβ2(N+1)2(N+2).
\begin{aligned}
{\rm var}(T)&=\int_{0}^{\beta}\xi^2 \frac{N\xi^{N-1}}{\beta^N}d\xi-\left(\frac{N\beta}{N+1}\right)^2\\
&=\frac{N\beta^2}{(N+1)^2(N+2)}.
\end{aligned}
var(T)=∫0βξ2βNNξN−1dξ−(N+1Nβ)2=(N+1)2(N+2)Nβ2.最终,我们可以得到最小方差为
var(θ^)=β24N(N+2).(5.10)\tag{5.10}
{\rm var}(\hat \theta)=\frac{\beta^2}{4N(N+2)}.
var(θ^)=4N(N+2)β2.(5.10)我们把样本均值的方差(5.9)重写如下
var(θ^)=1Nvar(x[n])=β212N(5.9)\tag{5.9}
{\rm var}(\hat \theta)=\frac{1}{N}{\rm var}(x[n])=\frac{\beta^2}{12N}
var(θ^)=N1var(x[n])=12Nβ2(5.9)显然如果N≥2N\ge 2N≥2,则(5.10)更小。