当前位置：首页 > news >正文

基于最小方差无畸变响应滤波器组的谱相关密度估计（Matlab代码实现）

news 2026/4/29 2:25:25

💥💥💞💞欢迎来到本博客❤️❤️💥💥
🏆博主优势：🌞🌞🌞博客内容尽量做到思维缜密，逻辑清晰，为了方便读者。
⛳️座右铭：行百里者，半于九十。
📋📋📋本文内容如下：🎁🎁🎁

⛳️赠与读者

👨‍💻做科研，涉及到一个深在的思想系统，需要科研者逻辑缜密，踏实认真，但是不能只是努力，很多时候借力比努力更重要，然后还要有仰望星空的创新点和启发点。建议读者按目录次序逐一浏览，免得骤然跌入幽暗的迷宫找不到来时的路，它不足为你揭示全部问题的答案，但若能解答你胸中升起的一朵朵疑云，也未尝不会酿成晚霞斑斓的别一番景致，万一它给你带来了一场精神世界的苦雨，那就借机洗刷一下原来存放在那儿的“躺平”上的尘埃吧。
或许，雨过云收，神驰的天地更清朗.......🔎🔎🔎

💥第一部分——内容介绍

基于最小方差无畸变响应滤波器组的谱相关密度估计

摘要

现有循环平稳信号的谱相关密度估计大多采用加权重叠分段平均方法，该方法存在分辨率低、频率泄漏严重的问题，尤其在短数据观测场景下缺陷明显。本文将高分辨率的最小方差无畸变响应估计器应用于谱相关密度估计，基于匹配滤波器组框架，分别推导了传统加权重叠分段平均方法、最小方差无畸变响应方法以及归一化最小方差无畸变响应方法对应的谱相关密度估计器。匹配滤波器组为各类方法提供了统一的分析与对比基础，能够清晰揭示不同估计器之间的内在联系。通过数值仿真实验，本文从循环频率分辨率、估计偏差、方差以及均方误差等维度，对所提估计器进行定量验证与性能评估。结果表明，最小方差无畸变响应类方法在分辨率上显著优于传统方法，具备更优的工程应用潜力。

关键词

循环平稳信号；谱相关密度；最小方差无畸变响应；匹配滤波器组；高分辨率谱估计

一、引言

在通信、机械振动分析、气象、水文、声学等众多工程与科学领域中，大量信号由周期性调制、旋转、采样等过程生成，这类信号的统计特性随时间呈现周期性变化，无法用传统平稳模型准确描述，而循环平稳模型能够很好地适配其非平稳特征。谱相关密度作为循环平稳信号最核心的二阶统计描述量，能够同时反映信号的频谱结构与周期调制特性，是信号检测、特征提取、调制识别等任务的关键工具。

在实际应用中，通常需要从有限长度的观测数据中估计谱相关密度。目前主流的估计方法以加权重叠分段平均法为基础，该方法原理简单、实现便捷，但存在频率分辨率低、频谱泄漏严重的固有缺陷，在短数据、低信噪比等场景下性能急剧下降，难以满足高精度分析需求。

最小方差无畸变响应方法是一种高分辨率自适应谱估计方法，在阵列信号处理、波达方向估计等领域已得到广泛验证，其优势在于短数据条件下仍能保持高分辨率、低旁瓣的优良性能。将该方法引入谱相关密度估计领域，有望突破传统方法的性能瓶颈。

本文基于匹配滤波器组统一框架，构建基于加权重叠分段平均、最小方差无畸变响应、归一化最小方差无畸变响应的谱相关密度估计器，实现不同方法的直接对比，并通过数值实验系统评估各估计器在循环频率分辨率、偏差、方差等方面的性能，为循环平稳信号的高分辨率分析提供可靠方案。

二、谱相关密度基础理论

（一）循环平稳信号基本概念

如果一个随机信号的二阶统计量随时间呈现周期性变化，该信号即为循环平稳信号。与传统平稳信号不同，循环平稳信号具有额外的循环频率维度，能够刻画信号的周期调制、周期采样等内在非平稳特性。

（二）谱相关密度定义与意义

谱相关密度是循环自相关函数的频域表示，它以双频域形式同时描述信号在常规频率和循环频率上的功率分布。常规频率反映信号的频谱构成，循环频率反映信号的周期调制速率，二者结合能够完整揭示循环平稳信号的本质特征。

（三）循环解调与互谱等价关系

通过循环解调操作，可将原始信号变换为两组频率偏移后的新信号，此时谱相关密度等价于这两组新信号的互功率谱密度。这一关键转化使得传统功率谱、互功率谱估计方法可以直接用于谱相关密度估计，为本文方法的推导奠定基础。

（四）典型循环平稳信号

单频正弦信号叠加高斯白噪声是最典型的循环平稳信号，其谱相关密度具有明确的峰值位置，可作为标准测试信号用于验证各类估计器的有效性，也是本文数值实验所采用的信号模型。

三、基于匹配滤波器组的谱相关密度估计

（一）滤波器组谱估计原理

滤波器组谱估计的核心思想是将信号通过一组窄带滤波器，通过测量滤波器输出功率来反推信号的功率谱。不同谱估计方法的差异，主要体现在滤波器设计方式与等效带宽的选取规则上。匹配滤波器组框架为各类方法提供了统一实现形式，便于方法间的对比与改进。

（二）加权重叠分段平均方法

加权重叠分段平均方法采用固定的矩形窗滤波器，不依赖信号自身特征，属于非自适应方法。该方法通过对信号分段、加窗、重叠平均来平滑谱估计结果，实现简单、数值稳定性好，但分辨率低、频谱泄漏严重，短数据下性能较差。

（三）最小方差无畸变响应方法

最小方差无畸变响应方法是一种数据自适应滤波器设计方法，其目标是在目标频率上保持无畸变响应，同时最小化其他频率的输出功率。该方法根据信号自相关矩阵自适应调整滤波器权重，具备高分辨率、低旁瓣的优势，短数据场景下性能优势尤为突出。

在此基础上，通过不同的等效带宽归一化方式，可得到归一化最小方差无畸变响应方法，进一步优化分辨率性能。

（四）谱相关密度估计器构建

利用循环解调与互谱的等价关系，将信号进行循环解调处理后，分别代入加权重叠分段平均、最小方差无畸变响应、归一化最小方差无畸变响应互谱估计框架，即可得到对应的谱相关密度估计器，实现从传统方法到高分辨率自适应方法的完整覆盖。

四、数值仿真与性能评估

（一）实验设置

本文采用单频谐波叠加白噪声的标准循环平稳信号作为测试对象，设置固定采样率、信号频率与数据长度，构建短数据仿真环境。自相关矩阵采用前后向平滑估计，以提升小样本下的估计精度。实验从循环频率分辨率概率、估计偏差、方差、均方误差四个维度，全面评估三种方法的性能。

（二）循环频率分辨率分析

循环频率分辨率是衡量谱相关密度估计器的核心指标，代表方法区分邻近循环频率分量的能力。实验结果显示，所有方法均存在阈值效应，当循环频率间隔小于阈值时无法分辨分量；在相同滤波器长度或分段长度下，分辨率由高到低依次为归一化最小方差无畸变响应、最小方差无畸变响应、加权重叠分段平均。即使增加分段长度，传统方法的分辨率仍远低于高分辨率自适应方法。

（三）估计偏差与方差分析

偏差反映估计结果与理论值的偏离程度，方差反映估计结果的波动稳定性。实验表明，所有估计器均存在一定偏差，其中加权重叠分段平均方法偏差最大；最小方差无畸变响应方法与传统方法偏差趋势相近，但整体更优；归一化最小方差无畸变响应方法在短滤波器长度下偏差与方差均最小，但随着长度增加性能快速退化。

（四）均方误差综合对比

均方误差综合反映偏差与方差的整体性能。结果显示，最小方差无畸变响应可视为传统方法的自适应改进版本，中等长度下二者均方误差接近；归一化最小方差无畸变响应在短长度下均方误差最优，长长度下性能显著下降。整体存在明显的性能权衡：滤波器或分段长度增加，分辨率提升，但均方误差增大；长度减小，均方误差降低，但分辨率下降。

五、结论与展望

（一）研究结论

最小方差无畸变响应及其归一化方法可有效应用于谱相关密度估计，在循环频率分辨率上显著优于传统加权重叠分段平均方法。
最小方差无畸变响应方法在偏差、方差、稳定性之间取得良好平衡，综合性能优于传统方法，适合通用工程场景。
归一化最小方差无畸变响应方法在短滤波器长度下分辨率与均方误差最优，适合密集循环频率分离的高精度需求，但长滤波器下性能退化明显。
匹配滤波器组框架能够有效统一各类谱相关密度估计器，简化方法推导与实现，便于方法间的对比分析。

（二）工程应用建议

通用高分辨率分析场景，优先选用最小方差无畸变响应方法。
密集循环频率分离、超分辨需求场景，选用短滤波器长度的归一化最小方差无畸变响应方法。
低计算资源、实时性要求高的场景，可选用传统加权重叠分段平均方法。
最优混合方案：先用归一化最小方差无畸变响应定位循环频率峰值位置，再用最小方差无畸变响应或传统方法计算对应位置的谱相关密度幅度，兼顾高分辨率与低计算量。

（三）未来研究方向

优化最小方差无畸变响应类方法的矩阵求逆运算，降低计算复杂度，提升实时性。
开展理论分析，建立偏差、方差、分辨率的闭式表达，完善算法理论支撑。
将方法扩展到多分量、时变循环频率、非高斯噪声等复杂实际场景。
结合通信、雷达、机械故障诊断等实际应用，开展实测数据验证与工程化落地。

📚第二部分——运行结果

部分代码：

fh = figure(15);

set(fh, 'color', 'white'); % sets the color to white

fig1 = semilogy(snr,mse(:,7),'o-',snr,mse(:,8),'^-',snr,mse(:,9),'s-');

set(fig1, 'LineWidth', 1.5, 'MarkerSize', 8.0);

legend('NMVDR N/2', 'NMVDR N/4', 'NMVDR N/8');

title('MSE of the NMVDR SCD Estimator','FontSize',16,'FontWeight', 'bold');

ylabel('MSE','FontSize',16,'FontWeight', 'bold');

xlabel('Signal to Noise Ratio (dB)','FontSize',16,'FontWeight', 'bold');

grid on;

set(gca, 'Box', 'off','TickDir', 'out', 'FontSize',16 ); % here gca means get current axis

%print -deps -painters img/scd_mse03

fh = figure(16);

set(fh, 'color', 'white'); % sets the color to white

fig1 = semilogy(snr,mse(:,10),'o-',snr,mse(:,11),'^-',snr,mse(:,12),'s-');

set(fig1, 'LineWidth', 1.5, 'MarkerSize', 8.0);

legend('APES N/2', 'APES N/4', 'APES N/8');

title('MSE of the APES SCD Estimator','FontSize',16,'FontWeight', 'bold');

ylabel('MSE','FontSize',16,'FontWeight', 'bold');

xlabel('Signal to Noise Ratio (dB)','FontSize',16,'FontWeight', 'bold');

grid on;

set(gca, 'Box', 'off','TickDir', 'out', 'FontSize',16 ); % here gca means get current axis

%print -deps -painters img/scd_mse04

fh = figure(18);

set(fh, 'color', 'white'); % sets the color to white

fig1 = semilogy(snr,mse(:,1),'o-',snr,mse(:,2),'o--',snr,mse(:,3),'o-.',snr,mse(:,4),'s-',snr,mse(:,5),'s--',snr,mse(:,6),'s-.',snr,mse(:,7),'^-',snr,mse(:,8),'^--',snr,mse(:,9),'^-.',snr,variance(:,10),'h-',snr,variance(:,11),'h--',snr,variance(:,12),'h-.');

set(fig1, 'LineWidth', 1.5, 'MarkerSize', 8.0);

legend( 'WOSA N', 'WOSA N/2', 'WOSA N/4','MVDR N/2', 'MVDR N/4', 'MVDR N/8','NMVDR N/2', 'NMVDR N/4', 'NMVDR N/8','APES N/2', 'APES N/4', 'APES N/8');

title('MSE Comparison of SCD Estimators','FontSize',16,'FontWeight', 'bold');

ylabel('MSE','FontSize',16,'FontWeight', 'bold');

xlabel('Signal to Noise Ratio (dB)','FontSize',16,'FontWeight', 'bold');

grid on;

set(gca, 'Box', 'off','TickDir', 'out', 'FontSize',16 ); % here gca means get current axis

%print -deps -painters img/scd_mse06

fh = figure(19);

%set(0,'DefaultAxesColorOrder',[0 0 0]); % makes the graph black and white

set(fh, 'color', 'white'); % sets the color to white

fig1 = semilogy(snr,mse(:,1),'o-',snr,mse(:,2),'o--',snr,mse(:,3),'o-.',snr,mse(:,4),'s-',snr,mse(:,5),'s--',snr,mse(:,6),'s-.',snr,mse(:,7),'^-',snr,mse(:,8),'^--',snr,mse(:,9),'^-.');

set(fig1, 'LineWidth', 1.5, 'MarkerSize', 8.0);

legend( 'WOSA N', 'WOSA N/2', 'WOSA N/4','MVDR N/2', 'MVDR N/4', 'MVDR N/8','NMVDR N/2', 'NMVDR N/4', 'NMVDR N/8');

title('MSE Comparison of SCD Estimators','FontSize',16,'FontWeight', 'bold');

ylabel('MSE','FontSize',16,'FontWeight', 'bold');

xlabel('Signal to Noise Ratio (dB)','FontSize',16,'FontWeight', 'bold');