二元信号探测的误检概率推导与其matlab验证

假设传输两个信号：$s_0[n]$和$s_1[n]$，传输过程中叠加了高斯白噪声$w[n]$，$w[n]\sim N(0,{\sigma }^2)$。在接收端，我们需要判断接收到的$x[n]$到底是哪个信号，就有了如下两个假设：
$$\begin{gathered} H_0:x[n]=s_0[n]+w[n] \\ {H}_1:x[n]=s_1[n]+w[n] \end{gathered}$$
假设发出两种信号的概率相等，当
$$\frac{\mathrm{Pr}(\vec{x}|H_1)}{\mathrm{Pr}(\vec{x}|H_0)}>\gamma =\frac{\mathrm{Pr}(H_0)}{\mathrm{Pr}(H_1)}=1$$
时，拒绝$H_0$，认为我们接收到了$s_1$信号。

为了进一步探讨该问题，将$\mathrm{Pr}(\vec{x}|H_i)$展开，因为接收到的信号$x[n]$每一位都满足$x[n]\sim N(s_i[n],{\sigma }^2)$，故整个信号$\vec{x}$有
$$\mathrm{Pr}(\vec{x}|H_i)=\frac{1}{(2\pi {\sigma }^2{)}^{\frac{N}{2}}}\exp [-\frac{1}{2{\sigma }^2}\sum _{n=0}^{N-1}(x[n]-s_i[n]{)}^2]$$
注意到对于两个假设$H_0$和$H_1$，$\mathrm{Pr}(\vec{x}|H_i)$的区别在于${\sum }_{n=0}^{N-1}(x[n]-s_i[n]{)}^2$，故设
$$D_i^2=\sum _{n=0}^{N-1}(x[n]-s_i[n]{)}^2=||\vec{x}-\overrightarrow{s_i}|{|}^2$$
不难看出，我们做出判断的依据就在于接收信号构成的向量$\vec{x}$与两已知信号向量$\overrightarrow{s_i}$之差的模长的平方。

再展开$D_i^2$，得
$$D_i^2=\sum _{n=0}^{N-1}{x}^2[n]-2\sum _{n=0}^{N-1}x[n]s_i[n]+\sum _{n=0}^{N-1}{s}_i^2[n]$$
发现$D_1^2$和$D_2^2$之间的差别在于$-2{\sum }_{n=0}^{N-1}x[n]s_i[n]+{\sum }_{n=0}^{N-1}{s}_i^2[n]$，故定义
$$T_i(x)=\sum _{n=0}^{N-1}x[n]s_i[n]-\frac{1}{2}\sum _{n=0}^{N-1}{s}_i^2[n]$$
其中，${\sum }_{n=0}^{N-1}{s}_i^2[n]$为信号$s_i$的能量，故$T_i$又可写为
$$T_i(x)=\sum _{n=0}^{N-1}x[n]s_i[n]-\frac{1}{2}{\epsilon }_i$$
可得出误检概率$P_e$的表达式为
$$\begin{array}{rl}{P}_e& =\mathrm{Pr}(H_1|H_0)\mathrm{Pr}(H_0)+\mathrm{Pr}(H_0|H_1)\mathrm{Pr}(H_1)\\ & =\frac{1}{2}[\mathrm{Pr}(H_1|H_0)+\mathrm{Pr}(H_0|H_1)]\\ & =\frac{1}{2}[\mathrm{Pr}(T_1(\vec{x})-T_0(\vec{x})>0|H_0)+\mathrm{Pr}(T_0(\vec{x})-T_1(\vec{x})>0|H_1)]\end{array}$$
故定义检验$T(\vec{x})=T_1(\vec{x})-T_0(\vec{x})$，即
$$T(\vec{x})=\sum _{n=0}^{N-1}x[n](s_1[n]-s_0[n])-\frac{1}{2}({\epsilon }_1-{\epsilon }_0)$$
通过计算可得$T$的期望与方差：
$$\begin{array}{rl}& E(T|H_0)=\sum _{n=0}^{N-1}{s}_0[n](s_1[n]-s_0[n])-\frac{1}{2}({\epsilon }_1-{\epsilon }_0)=-\frac{1}{2}||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2\\ & E(T|H_1)=\frac{1}{2}||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2\\ & Var(T|H_0)=Var(\sum _{n=0}^{N-1}x[n](s_1[n]-s_0[n])|H_0)={\sigma }^2||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2\\ & Var(T|H_1)={\sigma }^2||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2\end{array}$$
所以，$(T|H_0)$和$(T|H_1)$均满足正态分布，且
$$\begin{gathered} (T|H_0)\sim N(-\frac{1}{2}||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2,{\sigma }^2||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2) \\ (T|H_1)\sim N(\frac{1}{2}||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2,{\sigma }^2||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2) \end{gathered}$$
故误检概率可写为
$$\begin{array}{rl}{P}_e& =\mathrm{Pr}(T>0|H_0)\mathrm{Pr}(H_0)+\mathrm{Pr}(T<0|H_1)\mathrm{Pr}(H_1)\\ & =\frac{1}{2}[\mathrm{Pr}(T>0|H_0)+\mathrm{Pr}(T<0|H_1)]\\ & =\Phi (-\frac{\frac{1}{2}||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2}{\sqrt{{\sigma }^2||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2}})\\ & =\Phi (-\frac{1}{2}\sqrt{\frac{||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2}{{\sigma }^2}})\end{array}$$
可以看出，误检概率与两信号向量之差的模长成正比，与噪音的标准差成反比。

我们可以用matlab验证该推导结果，每次检验均生成两个随机信号$s_1,s_2$，分别叠加噪音计算$T$值，判断是否误检。

在这里只生成了长度为2的信号，一是因为长度太长每次检验的速度会增加，程序的总运行时间暴增；其次若信号长度太长，很难随机生成$||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2$较小的信号，使统计不准。
test.m

function res=test(sigma)
N=2;
s0=2*rand(1,N)-1;
s1=2*rand(1,N)-1;
w0=s0+sigma*randn(1,N);
w1=s1+sigma*randn(1,N);
T0=0;
T1=0;
E0=0;
E1=0;
for i=1:1:N
   E0=E0+s0(i)^2;
   E1=E1+s1(i)^2;
   T0=T0+w0(i)*(s1(i)-s0(i));
   T1=T1+w1(i)*(s1(i)-s0(i));
end
T0=T0-(E1-E0)/2;
T1=T1-(E1-E0)/2;
p=(s1-s0)*(s1-s0)';
%p=2*s1*s0'/(s1*s1'+s0*s0');
res=[p,(T0>0)+(T1<0)];

之后，统计两向量之差的模长与误检次数的关系，与推导结果比对，可得

在$||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2$偏大的部分，统计很不准确，这是因为对于随机生成的长度为2的信号，$||\overrightarrow{s_1}-\overrightarrow{s_0}|{|}^2$那么大的样本很少，统计结果就会非常不准确。可以考虑使用蒙特卡洛法消除部分偏差。

统计及生成gif代码：
statistic.m

clc;clear;
cnt=1;
for sigma=0.2:0.01:0.6
str='σ=';
N=10000000;
d=1000;
y=zeros(1,d);
tot=zeros(1,d);
data=zeros(N,2);
hb=0;
lb=400;
for i=1:1:N
    a=test(sigma);
    data(i,1)=a(1);
    data(i,2)=a(2);
    hb=max(hb,a(1));
    lb=min(lb,a(1));
end
hb=hb+0.01;
sep=(hb-lb)/d;
x=lb+sep/2:sep:hb-sep/2;
for i=1:1:N
    idx=floor((data(i,1)-lb)/sep)+1;
    y(idx)=y(idx)+data(i,2);
    tot(idx)=tot(idx)+2;
end
y=y./tot;
figure
plot(x,y,"linewidth",2);
hold on;
pd=makedist("Normal","mu",0,"sigma",1);
y1=cdf(pd,trans(x,sigma));
plot(x,y1,"linewidth",2);
xlabel("$||s_1-s_0||^2$","interpreter","latex","fontsize",16);
ylabel("prE","fontsize",16);
str1=[str num2str(sigma)];
title(str1,"fontsize",18);
print(1,'-dbmp',sprintf('img/%d',cnt));
cnt=cnt+1;
close;
end
%xlim([0.1 0.7]);
%ylim([-0.0025 0.03]);
%依次读取生成的所有图片
for j=1:41
    %获取当前图片
    A=imread(sprintf('img/%d.bmp',j));
    [I,map]=rgb2ind(A,256);
    %生成gif，并保存
    if(j==1)
        imwrite(I,map,'movefig.gif','DelayTime',0.1,'LoopCount',Inf)
    else
        imwrite(I,map,'movefig.gif','WriteMode','append','DelayTime',0.1)    
    end
end

trans.m

function re=trans(x,sigma)
re=-sqrt(x/sigma^2)/2;