本节内容摘自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 [ ∂ ln p ( 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
−
1
n=0,1,\ldots,N-1
n=0,1,…,N−1。试说明正则条件
E
[
∂
ln
p
(
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
=
0
N
−
1
p
(
x
[
n
]
;
θ
)
,
p({\bf x};\theta)=\prod_{n=0}^{N-1}p(x[n];\theta),
p(x;θ)=n=0∏N−1p(x[n];θ),显然,如果我们能够证明
E
[
∂
ln
p
(
x
;
θ
)
∂
θ
]
≠
0
,
f
o
r
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;\theta
p(y;θ的图像如下图所示。看图的时候,可以考虑
y
y
y是固定值。显然应该有
y
=
x
[
n
]
>
0
y=x[n]>0
y=x[n]>0。这时,如果
θ
<
y
\theta<y
θ<y,函数值为0;
θ
>
y
\theta>y
θ>y,则有函数值为
1
/
θ
1/\theta
1/θ。因此,有
E
[
∂
ln
p
(
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估计。
我们很自然会想到用样本均值
θ ^ = 1 N ∑ n = 0 N − 1 x [ 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])=\theta E(x[n])=θ。对估计求均值,得到
E ( θ ^ ) = θ . {\rm E}(\hat \theta)=\theta. E(θ^)=θ.因此 θ ^ \hat \theta θ^为无偏估计,其方差为
v a r ( θ ^ ) = 1 N v a r ( x [ n ] ) = β 2 12 N = θ 2 3 N . (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
=
0
N
−
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
β
N
0
≤
x
[
n
]
≤
β
,
n
=
0
,
1
,
…
,
N
−
1
0
o
t
h
e
r
w
i
s
e
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
β
N
max
x
[
n
]
<
β
,
min
x
[
n
]
>
0
0
o
t
h
e
r
w
i
s
e
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
β
N
u
(
β
−
max
x
[
n
]
)
u
(
min
x
[
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
β
N
u
(
β
−
max
x
[
n
]
)
,
h
(
x
)
=
u
(
min
x
[
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
)
=
max
x
[
n
]
T({\bf x})=\max x[n]
T(x)=maxx[n]为充分统计量,并且可以证明其为完备的(证明略)。
(2) 找到
T
(
x
)
T({\bf x})
T(x)的函数,得到无偏估计
上面我们找到完备的充分统计量为
T
(
x
)
=
max
x
[
n
]
T({\bf x})=\max x[n]
T(x)=maxx[n]。我们首先来确定它的均值。显然,
T
T
T为顺序统计量(order statistics)。我们先来看累积分布函数
P
r
{
T
≤
ξ
}
=
Pr
{
x
[
0
]
≤
ξ
,
x
[
1
]
≤
ξ
,
…
,
x
[
N
−
1
]
≤
ξ
}
=
∏
n
=
0
N
−
1
Pr
{
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为
p
T
(
ξ
)
=
d
P
r
{
T
≤
ξ
}
d
ξ
=
N
Pr
{
x
[
n
]
≤
ξ
}
N
−
1
d
P
r
{
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]≤ξ}.注意到
d
P
r
{
x
[
n
]
≤
ξ
}
d
ξ
\frac{d{\rm Pr}\{x[n]\le \xi\}}{d\xi}
dξdPr{x[n]≤ξ}为
x
[
n
]
x[n]
x[n]的概率密度函数
p
x
[
n
]
(
ξ
;
θ
)
p_{x[n]}(\xi; \theta)
px[n](ξ;θ),且
p
x
[
n
]
(
ξ
;
θ
)
=
{
1
β
0
<
ξ
<
β
0
o
t
h
e
r
w
i
s
e
.
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<ξ<βξ>β.最后,我们得到
p
T
(
ξ
)
=
{
0
ξ
<
0
N
(
ξ
β
)
N
−
1
1
β
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
)
=
∫
−
∞
∞
ξ
p
T
(
ξ
)
d
ξ
=
∫
0
β
ξ
N
(
ξ
β
)
N
−
1
1
β
d
ξ
=
N
N
+
1
β
=
2
N
N
+
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
+
1
2
N
T
\hat \theta=\frac{N+1}{2N}T
θ^=2NN+1T,因此最终无偏估计为
θ
^
=
N
+
1
2
N
max
x
[
n
]
.
\hat \theta=\frac{N+1}{2N}\max x[n].
θ^=2NN+1maxx[n].
有点出人意料的是,样本均值不是均匀分布噪声的MVU估计!
下面我们来看看二者方差情况。
(3)最小估计方差
我们可以得到最小估计方差为
v
a
r
(
θ
^
)
=
(
N
+
1
2
N
)
v
a
r
(
T
)
,
{\rm var}(\hat \theta)=\left(\frac{N+1}{2N}\right){\rm var}(T),
var(θ^)=(2NN+1)var(T),其中
v
a
r
(
T
)
=
∫
0
β
ξ
2
N
ξ
N
−
1
β
N
d
ξ
−
(
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.最终,我们可以得到最小方差为
v
a
r
(
θ
^
)
=
β
2
4
N
(
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)重写如下
v
a
r
(
θ
^
)
=
1
N
v
a
r
(
x
[
n
]
)
=
β
2
12
N
(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
≥
2
N\ge 2
N≥2,则(5.10)更小。