阵列信源估计理论——AIC、MDL、HQ(EDC)准则估计

Poulen

于 2025-03-07 11:23:06 发布

阅读量1k

点赞数 20

分类专栏： DOA估计算法文章标签： MDL准则 EDC/HQ准则 AIC准则信源估计算法

本文链接：https://blog.csdn.net/Poulen/article/details/146090789

版权

DOA估计算法专栏收录该内容

8 篇文章

订阅专栏

1.简介

在空间谱估计技术中,估计信号源数是一个关键问题。空间谱估计技术中的大部分算法均需要知道入射信号数N。在实际应用场合,信号源数往往是一个未知数,往往需要先估计信号源的数目或假设已知,然后再估计信号源的方向。但当估计的信号源数目与真实的信号源数目不一致时,空间谱曲线中的峰值个数与实际源数不相同,而且往往也会对真实信号的估计产生严重的影响(如偏离真实信号方向等)。为什么信号源数的不同会对空间谱曲线带来如此大的影响?这主要是因为,大多数空间谱估计算法都是基于特征子空间的算法,也就是充分利用了信号子空间与噪声子空间的正交性,而当信号源数估计不准时,也就是信号子空间、噪声子空间估计不准,即两者之间不完全正交,就会造成估计信号源时的漏警或虚警，也就会造成在估计信号方向时的偏差。

根据特征空间的分析可知,在一定的条件下数据协方差矩阵的大特征值数对应于信号源数,而其他的小特征值是相等的(等于噪声功率)。这就说明可以直接根据数据协方差矩阵的大特征值来判断信号的源数,但在实际应用场合(包括数据仿真),由于快拍数据、信噪比等方面的限制,对实际得到的数据协方差矩阵进行特征分解后,不可能得到明显的大小特征值。如何才能从分别不明显的特征值中进行判断?一种经典的方法就是通过似然比来确定门限,但这种方法有一个明显的缺点就是门限的判断带有明显的主观性。

针对上述方法的缺点,很多学者提出了在信号数目估计方面较为有效的方法,包括信息论方法、平滑秩法、矩阵分解法、盖氏圆方法以及正则相关等方法。我将分为几篇短文依次介绍这几种常用的信号源数估计的方法，并对相应算法的性能进行综合测试。

2.信息论方法

信息论的方法是Wax M和Kailath T提出的,这些方法都是在Anderson Tw和Rissanen J提出的理论基础上发展起来的,如Akaike信息论(AIC)准则、最小描述长度(MDL)准则及有效检测（EDC）准则等方法。

信息论的方法有一个统一的表达形式

式中，表示似然函数，p(k)表示罚函数。通过对L(k)和p(k)的不同选择就可以得到不同的准则。下面介绍EDC准则的一般形式：

式中n为待估计信源个数，L为采样点，其中Λ(n)为似然函数，且

式(2)中C(L)需满足如下条件：

当C(L)满足上述条件时，准则EDC估计具有一致性。在式(2)中选择C(L)分别为1，lnL/2，ln(lnL)/2时，就可以得到AIC、MDL及HQ(EDC)等准则，即：

3.仿真实验：

下面通过仿真来比较AIC、MDL、HQ准则的估计性能。实验针对3个独立信源，线性均匀阵列阵元数为12，快拍数为128，信号入射方向为5°、30°、50°，阵元间距为半波长，独立实验500次。仿真结果见图1所示;

从仿真实验可以得出以下结论：

AIC准则不是一致性估计(显然C(L)=1时不满足上述的第2个条件),即在大快拍数的场合,它仍然有较大的误差概率;而MDL准则相对较好;HQ在低信噪比时具有相比于AIC、MDL更好的估计性能，但其在大信噪比时具有较大的误差概率。
MDL准则是一致性估计,也就是在高信噪比情况下该准则有较好的性能,但在小信噪比情况下该准则相比AIC有高的误差概率。

仿真源码：

%%  Author：Poulen
%%  shijian：2025.3.6
%%  该脚本分别用于验证和学习Akaike信息论（AIC）准则、最小描述长度（MLD）准则和有效检测（EDC or HQ）准则估计信源个数；
clear 
close all;
clc;

%产生信号
%仿真初始值设定
M=12;                              %阵元单元
c=3e8;                             %光速
f0=77e9;                           %初始频率
lambda=c/f0;                       %波长
slope=30e12;                       %调频斜率
time=60e-6;                        %60us
d=0:lambda/2:(M-1)*lambda/2;       %阵列天线
thita=[10,30,0];                  %波达方向
K = length(thita);

%产生原始信号
N=128;%信号长度
t=linspace(0,time,N);

A=zeros(M,K);%导向向量空间 M*K
S=zeros(K,N);%信号空间

SNR_Scale = -20:0.5:20;

result1 = zeros(1,length(SNR_Scale));
result2 = zeros(1,length(SNR_Scale));
result3 = zeros(1,length(SNR_Scale));

for idx =1: length(SNR_Scale)
    Moteclo_num = 500;
    Moteclo = Moteclo_num;
    count_1 = 0;
    count_2 = 0;
    count_3 = 0;
    f = f0;     %非相干信号的频率
    while(Moteclo)
        for i = 1:K
            A(:,i) = exp(-1j*2*pi/lambda*d(:)*sind(thita(i)));
            S(i,:) = exp(1j*2*pi*(f*t(:)));
            f = 2e9+f;           
        end
        S = A*S; %产生阵列接收数据
        %向数据添加白噪声
        S = S +(randn(size(S)).*std(S))/db2mag(SNR_Scale(idx));
        
        R = (1/N)*S*conj(S.');
        [V,D] = eig(R);                   %在matalb2021版本以上，特征值已经从小到大排列
%         Un = V(:,1:end-K);              %提取小特征对应的特征向量作为噪声子空间 M*(M-K)
        
        
         eigenvalue = flip(diag(D));    %提取特征值,从小到大排列
        
         n1 = AIC(N,M,eigenvalue);
         n2 = MLD(N,M,eigenvalue);
         n3 = EDC(N,M,eigenvalue);
         
         if(n1 == K)
             count_1 = count_1 + 1;
         end
         if(n2 == K)
             count_2 = count_2 + 1;
         end
         if(n3 == K)
             count_3 = count_3 + 1;
         end
         S=zeros(K,N);
         Moteclo = Moteclo -1;
    end
    result1(idx) = 100 * count_1 / Moteclo_num;
    result2(idx) = 100 * count_2 / Moteclo_num;
    result3(idx) = 100 * count_3 / Moteclo_num;
end
%% 绘图
figure;
plot(SNR_Scale,result1,'-.<');hold on;
plot(SNR_Scale,result2,'-.>');hold on;
plot(SNR_Scale,result3,'-.*');hold on;
xlabel('SNR/dB');
ylabel('成功概率/%');
legend('AIC准则','MLD准则','EDC准则');


 %% AIC算法实现
function [n] = AIC(N,M,eigenvalue)
    aicvalue = zeros(1,M);
    for i = 0 : M-1
        sumeigvalue = 0;
        multipeigvalue = 1;
        for j = i+1 : M
            sumeigvalue = sumeigvalue + eigenvalue(j);
            multipeigvalue = multipeigvalue * eigenvalue(j);
        end
        gama = (sumeigvalue * (1/(M-i))) / (multipeigvalue^(1/(M-i)));
        
        aicvalue(i+1) = 2*N*(M-i)*log(gama) + 2*i*(2*M-i);
    end
    
    [~,index] = min(aicvalue);
    n = index -1;
end

 %% MLD算法实现
function [n] = MLD(N,M,eigenvalue)
    aicvalue = zeros(1,M);
    for i = 0 : M-1
        sumeigvalue = 0;
        multipeigvalue = 1;
        for j = i+1 : M
            sumeigvalue = sumeigvalue + eigenvalue(j);
            multipeigvalue = multipeigvalue * eigenvalue(j);
        end
        gama = (sumeigvalue * (1/(M-i))) / (multipeigvalue^(1/(M-i)));
        
        aicvalue(i+1) = N*(M-i)*log(gama) + 0.5*i*(2*M-i)*log(N);
    end
    
    [~,index] = min(aicvalue);
    n = index -1;
end

 %% EDC算法实现
function [n] = EDC(N,M,eigenvalue)
    aicvalue = zeros(1,M);
    for i = 0 : M-1
        sumeigvalue = 0;
        multipeigvalue = 1;
        for j = i+1 : M
            sumeigvalue = sumeigvalue + eigenvalue(j);
            multipeigvalue = multipeigvalue * eigenvalue(j);
        end
        gama = (sumeigvalue * (1/(M-i))) / (multipeigvalue^(1/(M-i)));
        
        aicvalue(i+1) = N*(M-i)*log(gama) + 0.5*i*(2*M-i)*log(log(N));
    end
    
    [~,index] = min(aicvalue);
    n = index -1;
end