本文主要内容摘自下面文献:
[1] Fadel F. Digham (2007). “On the Energy Detection of Unknown Signals Over Fading Channels.” IEEE TRANSACTIONS ON COMMUNICATION 55(1): 4.
[2] H.Urkowitz (1967). “Energy Detection of Unknown Deterministic Signals.” Proceedings of the IEEE 55(4): 9.【参见博文《合作频谱感知:未知确定信号的能量检测》】
[3] Zhi Quan, S. C., and Ali H. Sayed (2008). “Optimal Linear Cooperation for Spectrum Sensing in Cognitive Radio Networks.” IEEE JOURNAL OF SELECTED TOPICS IN SIGNAL PROCESSING, 2(1): 13.
[4] 博文《卡方分布》
1、系统模型
Urkowitz (1967)一文中,分别讨论了低通和带通两种情况。这里主要集中讨论带通情况。设发送的低通信号为
S
L
P
(
t
)
=
s
c
(
t
)
+
j
s
s
(
t
)
S_{\rm LP}(t)=s_c(t)+js_s(t)
SLP(t)=sc(t)+jss(t),信道为
h
=
α
e
j
θ
h=\alpha e^{j\theta}
h=αejθ,则带通接收信号为
r
(
t
)
=
{
R
{
[
h
s
L
P
(
t
)
+
n
L
P
(
t
)
]
e
j
2
π
f
c
t
}
:
H
1
R
{
n
L
P
(
t
)
e
j
2
π
f
c
t
}
:
H
0
(1)
\tag{1} r(t)={\Large\{ }\begin{aligned} {\mathcal R}\{[hs_{\rm LP}(t)+n_{\rm LP}(t)]e^{j2\pi f_ct}\}:\mathcal H_1\\ {\mathcal R}\{n_{\rm LP}(t)e^{j2\pi f_ct}\}:\mathcal H_0 \end{aligned}
r(t)={R{[hsLP(t)+nLP(t)]ej2πfct}:H1R{nLP(t)ej2πfct}:H0(1)其中
n
L
P
=
n
c
(
t
)
+
j
n
s
(
t
)
n_{\rm LP}=n_c(t)+jn_s(t)
nLP=nc(t)+jns(t)为等效低通噪声。此时定义信噪比为
γ
=
α
2
E
s
N
0
\gamma =\frac{\alpha^2 E_s}{N_0}
γ=N0α2Es。
根据Urkowitz (1967),我们可以得到判决变量为
r
=
2
N
0
∫
0
T
r
2
(
t
)
d
t
.
(2)
\tag{2} r=\frac{2}{N_0}\int_{0}^{T}r^2(t)dt.
r=N02∫0Tr2(t)dt.(2)下面我们分别针对
H
0
\mathcal H_0
H0和
H
1
\mathcal H_1
H1进行推导。
-
H
0
\mathcal H_0
H0
带通噪声可以表示为
n ( t ) = n c ( t ) cos ( 2 π f c t ) − n s ( t ) sin ( 2 π f c t ) , n(t)=n_c(t)\cos(2\pi f_ct)-n_s(t)\sin(2\pi f_ct), n(t)=nc(t)cos(2πfct)−ns(t)sin(2πfct),因此,判决变量为
r = 2 N 0 ∫ 0 T n 2 ( t ) d t = 1 N 0 ∫ 0 T [ n c 2 ( t ) + n s 2 ( t ) ] d t \begin{aligned} r&=\frac{2}{N_0}\int_{0}^{T}n^2(t)dt\\ &=\frac{1}{N_0}\int_{0}^{T}\left[ n_c^2(t)+n_s^2(t)\right]dt\\ \end{aligned} r=N02∫0Tn2(t)dt=N01∫0T[nc2(t)+ns2(t)]dt进一步,我们把 n c ( t ) n_c(t) nc(t)和 n s ( t ) n_s(t) ns(t)分别用有限个抽样样值表示,注意这时的抽样速率为 W W W而非 2 W 2W 2W,若 N = 2 W T N=2WT N=2WT,则有 n c ( t ) = ∑ i = 1 N / 2 n c i s i n c ( W t − i ) n s ( t ) = ∑ i = 1 N / 2 n s i s i n c ( W t − i ) , \begin{aligned} n_c(t)&=\sum_{i=1}^{N/2}n_{ci}{\rm sinc}\left(Wt-i \right)\\ n_s(t)&=\sum_{i=1}^{N/2}n_{si}{\rm sinc}\left(Wt-i \right), \end{aligned} nc(t)ns(t)=i=1∑N/2ncisinc(Wt−i)=i=1∑N/2nsisinc(Wt−i),利用sinc函数的性质,我们可以得到
r = 1 N 0 W ∑ i = 1 N / 2 [ n c i 2 + n s i 2 ] = ∑ i = 1 N / 2 [ ( n c i N 0 W ) 2 + ( n s i N 0 W ) 2 ] r=\frac{1}{N_0W}\sum_{i=1}^{N/2}\left[ n_{ci}^2+n_{si}^2\right]=\sum_{i=1}^{N/2}\left[ \left(\frac{n_{ci}}{\sqrt{N_0W}}\right)^2+\left(\frac{n_{si}}{\sqrt{N_0W}}\right)^2\right] r=N0W1i=1∑N/2[nci2+nsi2]=i=1∑N/2[(N0Wnci)2+(N0Wnsi)2]满足中心 χ 2 \chi^2 χ2分布,自由度为 2 T W 2TW 2TW。
-
H
1
\mathcal H_1
H1
低通等效信号为
s L P ( t ) = s c ( t ) + j s s ( t ) , s_{\rm LP}(t)=s_c(t)+js_s(t), sLP(t)=sc(t)+jss(t),等效低通信道为 h L P = α c + j α s h_{\rm LP}=\alpha_c+j\alpha_s hLP=αc+jαs,则经过衰落信道后的低通等效接收信号(不含噪声)为
s L P ′ ( t ) = α c s c ( t ) + j α c s s ( t ) + j α s s c ( t ) − α s s s ( t ) , s'_{\rm LP}(t)=\alpha_cs_c(t)+j\alpha_c s_s(t)+j\alpha_s s_c(t)-\alpha_ss_s(t), sLP′(t)=αcsc(t)+jαcss(t)+jαssc(t)−αsss(t),则带通信号可以表示为
s ′ ( t ) = R [ s L P ′ ( t ) e j 2 π f c t ] = [ α c s c ( t ) − α s s s ( t ) ] cos ( 2 π f c t ) − [ α c s s ( t ) + α s s c ( t ) ] sin ( 2 π f c t ) \begin{aligned} s'(t)&={\mathcal R}\left[ s'_{\rm LP}(t)e^{j2\pi f_ct} \right]\\ &=\left[ \alpha_cs_c(t)-\alpha_ss_s(t)\right]\cos(2\pi f_ct)-\left[\alpha_c s_s(t)+\alpha_s s_c(t)\right]\sin(2\pi f_ct)\\ \end{aligned} s′(t)=R[sLP′(t)ej2πfct]=[αcsc(t)−αsss(t)]cos(2πfct)−[αcss(t)+αssc(t)]sin(2πfct)由此,可以得到接收信号为
r ( t ) = [ α c s c ( t ) − α s s s ( t ) + n c ( t ) ] cos ( 2 π f c t ) − [ α c s s ( t ) + α s s c ( t ) + n s ( t ) ] sin ( 2 π f c t ) = r c ( t ) cos ( 2 π f c t ) − r s ( t ) sin ( 2 π f c t ) \begin{aligned} r(t)&=\left[ \alpha_cs_c(t)-\alpha_ss_s(t)+n_c(t)\right]\cos(2\pi f_ct)-\left[\alpha_c s_s(t)+\alpha_s s_c(t)+n_s(t)\right]\sin(2\pi f_ct)\\ &=r_c(t)\cos(2\pi f_ct)-r_s(t)\sin(2\pi f_ct) \end{aligned} r(t)=[αcsc(t)−αsss(t)+nc(t)]cos(2πfct)−[αcss(t)+αssc(t)+ns(t)]sin(2πfct)=rc(t)cos(2πfct)−rs(t)sin(2πfct)同理,可以用抽样样值表示
r c ( t ) = ∑ i = 1 N / 2 [ α c s c i − α s s s i + n c i ] s i n c ( W t − i ) r s ( t ) = ∑ i = 1 N / 2 [ α c s s i + α s s c i + n s i ] s i n c ( W t − i ) , \begin{aligned} r_c(t)&=\sum_{i=1}^{N/2}[\alpha_c s_{ci}-\alpha_s s_{si}+n_{ci}]{\rm sinc}\left(Wt-i \right)\\ r_s(t)&=\sum_{i=1}^{N/2}[\alpha_c s_{si}+\alpha_s s_{ci}+n_{si}]{\rm sinc}\left(Wt-i \right), \end{aligned} rc(t)rs(t)=i=1∑N/2[αcsci−αsssi+nci]sinc(Wt−i)=i=1∑N/2[αcssi+αssci+nsi]sinc(Wt−i),由此,我们可以得到
r = 2 N 0 ∫ 0 T r 2 ( t ) d t = 1 N 0 ∫ 0 T [ r c 2 ( t ) + r s 2 ( t ) ] d t , (2) \tag{2} \begin{aligned} r=&\frac{2}{N_0}\int_{0}^{T}r^2(t)dt\\ &=\frac{1}{N_0}\int_{0}^{T}\left[ r_c^{2}(t)+r_s^{2}(t)\right]dt\\ \end{aligned}, r=N02∫0Tr2(t)dt=N01∫0T[rc2(t)+rs2(t)]dt,(2)即
H 1 : r = 1 N 0 W [ ∑ i = 1 N / 2 ( α c s c i − α s s s i + n c i ) 2 + ∑ i = 1 N / 2 ( α c s s i + α s s c i + n s i ) 2 ] = ∑ i = 1 N / 2 ( α c s c i − α s s s i + n c i N 0 W ) 2 + ∑ i = 1 N / 2 ( α c s s i + α s s c i + n s i N 0 W ) 2 \begin{aligned} \mathcal H_1:\qquad r&=\frac{1}{N_0W}\left[ \sum_{i=1}^{N/2}(\alpha_c s_{ci}-\alpha_s s_{si}+n_{ci})^2+\sum_{i=1}^{N/2}(\alpha_c s_{si}+\alpha_s s_{ci}+n_{si})^2\right]\\ &=\sum_{i=1}^{N/2}\left(\frac{\alpha_c s_{ci}-\alpha_s s_{si}+n_{ci}}{\sqrt{N_0W}}\right)^2+\sum_{i=1}^{N/2}\left(\frac{\alpha_c s_{si}+\alpha_s s_{ci}+n_{si}}{\sqrt{N_0W}}\right)^2 \end{aligned} H1:r=N0W1⎣⎡i=1∑N/2(αcsci−αsssi+nci)2+i=1∑N/2(αcssi+αssci+nsi)2⎦⎤=i=1∑N/2(N0Wαcsci−αsssi+nci)2+i=1∑N/2(N0Wαcssi+αssci+nsi)2满足非中心 χ 2 \chi^2 χ2分布,方差 σ 2 = 1 \sigma^2=1 σ2=1,非中心参数为 μ = 2 γ \mu=2\gamma μ=2γ, γ = α 2 E s N 0 \gamma=\frac{\alpha^2 E_s}{N_0} γ=N0α2Es。
我们来推导下非中心参数:
μ
=
∑
i
=
1
N
/
2
(
α
c
2
s
c
i
2
N
0
W
+
α
s
2
s
s
i
2
N
0
W
)
+
∑
i
=
1
N
/
2
(
α
c
2
s
s
i
2
N
0
W
+
α
s
2
s
c
i
2
N
0
W
)
=
α
c
2
N
0
∑
i
=
1
N
/
2
s
c
i
2
W
+
α
s
2
N
0
∑
i
=
1
N
/
2
s
s
i
2
W
+
α
c
2
N
0
∑
i
=
1
N
/
2
s
s
i
2
W
+
α
s
2
N
0
∑
i
=
1
N
/
2
s
c
i
2
W
=
2
α
2
E
s
N
0
\begin{aligned} \mu&=\sum_{i=1}^{N/2}\left(\frac{\alpha_c^2 s^2_{ci}}{{N_0W}}+\frac{\alpha_s^2 s^2_{si}}{{N_0W}}\right)+\sum_{i=1}^{N/2}\left(\frac{\alpha_c^2 s^2_{si}}{{N_0W}}+\frac{\alpha_s^2 s^2_{ci}}{{N_0W}}\right)\\ &=\frac{\alpha_c^2}{N_0} \sum_{i=1}^{N/2}\frac{s^2_{ci}}{{W}}+\frac{\alpha_s^2}{N_0} \sum_{i=1}^{N/2}\frac{s^2_{si}}{{W}}+\frac{\alpha_c^2}{N_0} \sum_{i=1}^{N/2}\frac{s^2_{si}}{{W}}+\frac{\alpha_s^2}{N_0} \sum_{i=1}^{N/2}\frac{s^2_{ci}}{{W}}\\ &=\frac{2\alpha^2 E_s}{N_0} \end{aligned}
μ=i=1∑N/2(N0Wαc2sci2+N0Wαs2ssi2)+i=1∑N/2(N0Wαc2ssi2+N0Wαs2sci2)=N0αc2i=1∑N/2Wsci2+N0αs2i=1∑N/2Wssi2+N0αc2i=1∑N/2Wssi2+N0αs2i=1∑N/2Wsci2=N02α2Es
因此,
r
r
r的概率密度函数可以表示为
p
R
(
r
)
=
{
r
N
2
−
1
e
−
r
2
σ
2
σ
N
2
N
2
Γ
(
N
2
)
,
H
0
1
2
σ
2
(
r
μ
)
N
−
2
4
e
−
μ
+
r
2
σ
2
I
N
2
−
1
(
μ
r
σ
2
)
,
H
1
,
p_R(r)=\left\{\begin{aligned} &\frac{r^{\frac{N}{2}-1}e^{-\frac{r}{2\sigma^2}}}{\sigma^N2^{\frac{N}{2}}\Gamma(\frac{N}{2})},&\quad {\mathcal H}_0\\ &\frac{1}{2\sigma^2}(\frac{r}{\mu})^{\frac{N-2}{4}}e^{-\frac{\mu+r}{2\sigma^2}}I_{\frac{N}{2}-1}(\frac{\sqrt{\mu r}}{\sigma^2}),&\quad {\mathcal H}_1, \end{aligned}\right.
pR(r)=⎩⎪⎪⎪⎨⎪⎪⎪⎧σN22NΓ(2N)r2N−1e−2σ2r,2σ21(μr)4N−2e−2σ2μ+rI2N−1(σ2μr),H0H1,这里的
I
v
(
⋅
)
I_v(\cdot)
Iv(⋅)为
v
v
v阶第一类修正贝塞尔函数。
2、AWGN信道上的检测概率与虚警概率
为了计算检测与虚警概率,有很多近似方法。下面我们主要讨论高斯分布来近似。我们知道,随着样本数
N
N
N的增大,根据中心极限定理,可以认为
r
r
r的分布接近高斯分布。
我们知道检测变量
r
r
r满足如下分布:
r
∼
{
χ
N
2
,
H
0
χ
N
2
(
μ
)
,
H
1
r\sim \left\{\begin{aligned} \chi_N^2,&\quad \mathcal H_0\\ \chi_N^2(\mu),&\quad \mathcal H_1 \end{aligned} \right.
r∼{χN2,χN2(μ),H0H1这里
μ
=
2
E
s
α
2
N
0
=
2
γ
\mu=\frac{2E_s\alpha^2}{N_0}=2\gamma
μ=N02Esα2=2γ,显然
γ
\gamma
γ为信噪比。
用高斯分布近似卡方分布:
若 m = ∑ i = 1 N r i 2 m=\sum_{i=1}^{N}r_i^2 m=∑i=1Nri2,其中 r i ∼ N ( 0 , σ 2 ) r_i\sim {\mathcal N}(0,\sigma^2) ri∼N(0,σ2),显然 m σ 2 \frac{m}{\sigma^2} σ2m为自由度等于 N N N的卡方分布。当 N N N足够大时,我们可以把卡方分布近似成高斯分布。由于 ( r i σ 2 ) 2 ∼ N ( 0 , 1 ) (\frac{r_i}{\sigma^2})^2\sim \mathcal N(0,1) (σ2ri)2∼N(0,1),因此 m σ 2 \frac{m}{\sigma^2} σ2m的均值为 N N N,方差为 2 N 2N 2N。将 m m m近似为高斯分布,有
m ∼ N ( N σ 2 , 2 N σ 4 ) . m\sim {\mathcal N}(N\sigma^2,2N\sigma^4). m∼N(Nσ2,2Nσ4).
若 m = ∑ i = 1 N r i 2 m=\sum_{i=1}^{N}r_i^2 m=∑i=1Nri2,其中 r i ∼ N ( μ , σ 2 ) r_i\sim {\mathcal N}(\mu,\sigma^2) ri∼N(μ,σ2),显然 m σ 2 \frac{m}{\sigma^2} σ2m为自由度等于 N N N的非中心卡方分布,其均值为 N + μ N+\mu N+μ,方差为 2 ( N + 2 μ ) 2(N+2\mu) 2(N+2μ),这里 μ \mu μ为非中心参数。同样当 N N N足够大时,我们将 m m m近似为高斯分布,有
m ∼ N ( ( N + μ ) σ 2 , 2 ( N + 2 μ ) σ 4 ) ) . m\sim {\mathcal N}\left((N+\mu)\sigma^2,2(N+2\mu)\sigma^4)\right). m∼N((N+μ)σ2,2(N+2μ)σ4)).
当
N
N
N足够大时,我们可以把卡方分布近似成均值为
m
r
=
{
N
σ
2
;
H
0
(
N
+
μ
)
σ
2
;
H
1
m_r=\left\{ \begin{aligned} N\sigma^2;&\quad H_0\\ (N+\mu)\sigma^2;&\quad H_1 \end{aligned}\right.
mr={Nσ2;(N+μ)σ2;H0H1方差为
σ
r
2
=
{
2
N
σ
4
;
H
0
2
(
N
+
2
μ
)
σ
4
;
H
1
\sigma^2_r=\left\{ \begin{aligned} 2N\sigma^4;&\quad H_0\\ 2(N+2\mu)\sigma^4;&\quad H_1 \end{aligned}\right.
σr2={2Nσ4;2(N+2μ)σ4;H0H1
的高斯分布,即
r
∼
N
(
m
r
,
σ
r
2
)
r\sim {\mathcal N}(m_r,\sigma^2_r)
r∼N(mr,σr2)。
若判决门限设为
η
\eta
η,则检测概率为
P
d
=
P
r
(
r
>
η
∣
H
1
)
=
Q
[
η
−
E
(
r
∣
H
1
)
V
a
r
(
r
∣
H
1
)
]
,
P_d={\rm Pr}(r>\eta|H_1)=Q\left[\frac{\eta-{\rm E}(r|H_1)}{\sqrt{{\rm Var}(r|H_1)}} \right],
Pd=Pr(r>η∣H1)=Q[Var(r∣H1)η−E(r∣H1)],虚警概率为
P
f
=
P
r
(
r
>
η
∣
H
0
)
=
Q
[
η
−
E
(
r
∣
H
0
)
V
a
r
(
r
∣
H
0
)
]
.
P_f={\rm Pr}(r>\eta|H_0)=Q\left[\frac{\eta-{\rm E}(r|H_0)}{\sqrt{{\rm Var}(r|H_0)}} \right].
Pf=Pr(r>η∣H0)=Q[Var(r∣H0)η−E(r∣H0)].
显然,不同门限值大小,会影响
P
d
P_d
Pd与
P
f
P_f
Pf。此外,我们还有漏检概率为
P
m
=
1
−
P
d
.
P_m=1-P_d.
Pm=1−Pd.
图1给出了漏检概率
P
m
P_m
Pm随虚检概率(
P
f
P_f
Pf)变化的曲线,这里没有考虑衰落,即
α
=
1
\alpha=1
α=1。显然,
P
m
P_m
Pm随着
P
f
P_f
Pf的增加而减少。我们还可以发现,随着信噪比
γ
\gamma
γ的增大,曲线下面积在减少,因此检测性能变好。
图1 漏检概率(
P
m
P_m
Pm)与虚检概率(
P
f
P_f
Pf)关系曲线