1. 为什么单节点频谱感知不可靠：从隐藏终端说起
协作频谱感知这几年在认知无线电里被提的最多的一句话就是“单节点不可靠”。信号衰落、路径损耗、障碍物遮挡，任何一个因素都可能让某个感知节点对主用户的存在完全“失明”。你可能会想，直接在接收端量一下接收功率不就行了吗？确实可以，这就是经典的能量检测，但在衰落信道下，单个节点的瞬时信噪比可能低到完全没法判断，比如主用户和感知节点之间隔了一栋楼，检测概率可能掉到0.3以下。这种场景一旦发生，盲目认为“频段空闲”然后让次用户接入，就会对主用户产生严重干扰。
解决思路其实很朴素：多找几个节点一起感知，把各自的观测结果提交到融合中心，由融合中心综合所有信息给出一个全局判决。这就是协作频谱感知的基本思想。但这带出一个新问题——每个节点到底往融合中心传什么？传原始采样波形，带宽开销大到不现实；只传本地0/1判决，又会丢掉很多统计信息，融合增益有限。于是工程上就需要一种折中方案：节点本地提取一个足够反映“信号存在与否”的统计量，把统计量传给融合中心，再做集中式融合。
Pietra-Ricci指数检测器就是这种思路下的一个典型代表。它基于接收信号的样本协方差矩阵特征值分布设计统计量，天然不依赖噪声功率先验，非常适合抗噪不确定性环境中的协作检测。这篇文章我打算从原理讲到Matlab实现，重点拆开数据融合的门限设计和工程上容易踩的坑，适合正在做协作频谱感知方向的研究生，或者想快速用Matlab搭一套检测仿真的工程师。如果你时间紧，可以直接跳到第4节抄代码跑数据；如果你想理解为什么PR统计量能抵抗噪声不确定性、为什么融合门限必须在融合端重新标定，那还是把前几节一起看完更划算。
1.1 单节点检测的困境
单节点检测的问题不只是“信号弱”。更麻烦的是噪声不确定性。协方差、信道、温度、器件老化都会让接收机的噪声功率随时间变化，能量检测器的门限又是直接依赖噪声功率的，一旦实际噪声功率和预设值偏差1dB，虚警概率就可能从0.1跳到0.3。在认知无线电场景里，虚警概率升高意味着频谱资源被白白浪费，漏检则意味着对主用户的干扰失控，两边都不能容忍。
所以我个人一直觉得，协作感知的首要目标不是把检测概率堆到99%，而是让整个系统在噪声不确定性和信道衰落面前依然行为稳定。单个节点再怎么优化，面对深度衰落和隐藏终端都是脆弱的。隐藏终端这个场景稍微展开说一下：主用户发射机附近明明有人在使用频谱，但因为遮挡物或阴影衰落，某个感知节点完全感知不到信号。这时候如果系统只依赖这个节点做判断，就会误判频段空闲，次用户的接入就会干扰主用户的接收机。多节点协作的本质，就是利用空间分集对冲这种单点失明。
1.2 协作感知的系统模型
集中式协作频谱感知的系统模型可以抽象成这样：一个主用户发射机，M个感知节点，一个融合中心。每个感知节点独立接收信号，在本地提取统计量，然后通过报告信道把统计量发送给融合中心。融合中心把所有节点的统计量按某种规则合并，得到一个全局判决统计量，再用门限判决输出“主用户存在”或“主用户不存在”。
这个模型有个很关键的特征：感知信道和报告信道是分开的。很多初学者仿真的时只关心感知信道的信噪比，完全忽略了报告信道噪声，实际工程里报告链路本身的可靠性也是要算进系统预算里的。第6节我会专门讲这部分坑。在这之前，我们先默认报告信道是理想的，先把核心检测和融合逻辑跑通。
2. Pietra-Ricci指数检测器：特征值离散度如何暴露信号
如果只用一句话概括PR指数检测器的核心，那就是：噪声环境下的协方差矩阵特征值趋于均匀，混入相关信号之后特征值会变得不均匀，PR指数就是用来量化这种不均匀程度的。这个思路最早并不在通信领域，而是统计学里用来衡量收入分配不平等程度的指标，后来被信号处理研究者借过来，用在频谱感知里意外地好用。
2.1 从接收信号到样本协方差矩阵
考虑一个感知节点在某段时间内采集的复基带信号。如果节点有多根天线，每一次采样天然形成一个多通道向量；如果节点只有单天线，也可以通过滑窗或者过采样构造“虚拟多通道”向量。这里假设每个观测向量是L维的，一共N个快拍，组成一个L行N列的矩阵X。样本协方差矩阵就是：
R = X * X' / N
其中X'表示共轭转置。这个R是L乘L的方阵，数学上表征了各通道之间的相关结构。在纯噪声环境下，各通道之间独立，R在主对角线附近比较均匀，特征值大致相等。如果有主用户信号混进来，并且信号在通道之间存在时间或空间相关性，R的特征值就会拉开差距——少数大特征值对应信号子空间，其余小特征值对应噪声子空间。
我实际做仿真时最常用的构造方式，是对单天线接收序列做不重叠滑窗。把长度为L*N的接收序列y排成L行N列的矩阵，让矩阵第n列恰好是原始序列的第((n-1)L+1)到(nL)个采样点。这样相邻采样点被安排到同一列的不同“虚拟通道”里，AR信号或窄带信号在相邻采样点间有相关性，而白噪声没有，PR检测器就能把这两者区分开。这个构造方式非常关键，第6节我会专门说为什么要这么做，以及怎么做容易出错。
2.2 PR指数的定义与性质
有了协方差矩阵R，接下来计算它的L个特征值，按从大到小排列。PR指数的工程实现形式我常用下面这个：
T_PR = mean(lambda.^2) / mean(lambda)^2 - 1
展开写就是：T_PR等于“特征值平方的均值”除以“特征值均值的平方”，再减1。这个量有两个很直观的性质。第一，如果所有特征值完全相等，T_PR等于0；特征值差异越大，T_PR越大。第二，它只取决于特征值之间的相对比例，如果特征值整体乘以一个常数c，比如噪声功率变大，分子分母同时乘c²，比值不变。这就是PR检测器抗噪声不确定性的数学本底。
判决规则很简单：T_PR超过门限gamma，判定主用户存在，也就是H1；否则判定H0。门限怎么定，第3节和第4节都会涉及。需要强调的是，T_PR这个值本身不是概率也不是似然比，它只是一个特征值离散度的度量，所以门限必须在给定虚警概率下重新标定，绝对不能拿一个经验常数到处套。
2.3 PR检测器、能量检测器与MME检测器的定位差异
做检测器性能对比时，最常被拉出来比较的就是能量检测ED、最大最小特征值比检测MME和PR指数检测。能量检测的道理最简单：看接收信号的总功率。它的判决统计量依赖噪声功率，门限需要知道当前噪声到底有多大，这在噪声不确定时会出问题。MME则是看最大特征值和最小特征值的比值，不需要知道噪声功率，但对小特征值的估计误差很敏感，当L较大而N不够大时，最小特征值估计不准，性能波动明显。
PR指数介于两者之间：既不要噪声功率先验，又不像MME那样只盯着最大和最小两个极端值，而是把全部L个特征值都纳入统计。它对于“多个特征值被抬高但最大特征值不明显”的情况更鲁棒。我用一句话总结这三个检测器的差异给刚接触这个方向的朋友参考：ED在理想环境下最好，MME在抗最坏情况上很激进，但N不够时容易翻车，PR在稳健性和计算复杂度之间比较平衡。
3. 集中式融合：多节点PR统计量怎么合成全局判决
本地每个节点都算出一个PR统计量之后，问题就变成：融合中心拿到M个T值，怎么合成一个全局判决？这是整个系统设计里最容易被低估的一步，因为融合规则直接决定性能上限，而这个上限和检测器的选择是正交的。检测器再强，融合规则不对，增益照样发挥不出来。
3.1 软融合与硬融合的取舍
融合规则大体分两类：软融合和硬融合。硬融合是每个节点先把本地统计量和一个本地门限比较，产生0/1判决，再把这1bit上传给融合中心，融合中心按“至少K个节点判有”或者“所有节点都判有”之类的规则输出。优点是报告信道开销极小、抗报告链路噪声能力强；缺点是本地判决已经把统计量里的幅度信息扔掉了，融合中心没法拿不同节点的统计量置信度做加权。
软融合则更“奢侈”一些：每个节点上传本地的连续统计量T_m，融合中心做加权或直接平均，再做一次全局门限判决。这样做保留了更多信息，性能上限更高，但对报告信道质量敏感。标题既然写的是集中式数据融合，我默认走的就是软融合路线，这也是研究论文里更常见、能拉开性能差距的场景。我给的Matlab代码采用等增益合并，也就是把M个节点的T_m直接取平均，因为实现最简单、参数最少，用来验证基本逻辑最合适。
如果你想做加权融合，一种很自然的扩展是按噪声功率反比加权：信噪比高的节点权重高，信噪比低的节点权重低。但这样做要求每个节点能估计本地噪声功率，在噪声不确定性场景下引入的误差可能比权重带来的增益还大。工程上我见过有人直接忽略噪声功率估计，单纯用等增益合并，性能也挺稳。所以我的建议是：第一版系统别急着上复杂的权重方案，等增益合并跑通后，再考虑要不要加报告信道噪声模型和加权策略。
3.2 门限定标必须在融合之后
初学者经常犯的一个错误是：先在本地把每个节点的门限算好，然后把本地判决结果融合。或者更隐蔽一点，先对单个节点的T_PR做门限判决，假设融合中心只是“多数表决”，这是硬融合的思路，和软融合系统根本不是一回事。在软融合的框架下，全局门限必须基于“融合之后”的统计量在H0条件下的分布来确定。
为什么？因为即使每个节点的T_PR在H0下的分布都相同，M个独立统计量的均值分布和单个统计量的分布也不是同一个东西。均值更加紧缩，方差小得多，直接沿用单节点门限会导致实际虚警概率远低于目标值，看起来好像“性能很好”，其实是把门限定松了，检测概率也被压低了，虚警和检测都不在标称点上。我自己调试的时候遇到过类似情况，当时第一反应是“融合增益这么大吗”，后来检查门限才发现是标定位置错了，虚警概率只有0.01而不是目标0.1，属于典型的虚假性能提升。
正确做法是：在H0条件下跑出融合后统计量T_fus的大量样本，取置信水平1-Pf对应的分位数作为门限。这个流程看起来繁琐，但做一次就行，第4节代码里就是按照这个思路实现的。
4. Matlab实现：从零搭起PR协作检测仿真
这一节把完整实现流程拆开讲。代码我尽量写得直白，方便直接改成你自己项目里的模块。需要说明的是，下面的仿真中我用了不重叠滑窗构造虚拟通道，这是单天线节点构造样本协方差矩阵最常见的方式。
4.1 仿真参数与系统模型
先定义系统参数。我这里设M=4个协作节点，每个节点的滑窗长度L=4个采样点，快拍数N=128，也就是说每个节点实际采集L*N=512个采样点。主用户信号用复高斯白噪声经过一阶AR滤波器产生，目的是让信号在时间上具备一定的相关性——这是PR检测器能工作的前提之一。每个节点的信道是独立瑞利衰落，噪声是复高斯白噪声，感知信噪比从-15dB到0dB扫描。
clear; clc; close all; M = 4; % 协作节点数 L = 4; % 滑窗长度，等价于虚拟通道数 N = 128; % 快拍数，构造L*N的观测矩阵 snr_vec = -15:3:0; % 各节点平均感知SNR(dB) Pf_target = 0.1; % 目标虚警概率 nMc = 2000; % 蒙特卡洛次数 rng(2025); % 固定随机数种子，便于复现信号模型这块，主用户信号s的长度是L*N，经过AR(1)滤波后再归一化到单位功率。这样设置的好处是不同节点的噪声方差可以按SNR公式直接算，不用额外校准信号功率。
4.2 核心函数：PR统计量的实现
PR统计量封装成一个独立的Matlab函数，输入任意L行N列的复矩阵，输出一个标量T。这个函数内部就是协方差矩阵、特征值分解、二阶矩归一化三步，没有任何多余逻辑。建议你在自己的项目里把这个函数单独存成一个文件，因为后面无论是单个节点检测还是融合后检测，都要反复调用它。
function T = prDetector(X) % X: L x N 复基带观测矩阵 % 输出: PR统计量 T = mean(lambda.^2)/mean(lambda)^2 - 1 L = size(X, 1); N = size(X, 2); R = (X * X') / N; % 样本协方差矩阵，LxL lambda = eig(R); % 特征值按升序返回，但PR只依赖整体分布 lambda = lambda(:); mu1 = mean(lambda); % 特征值均值 mu2 = mean(lambda.^2); % 特征值平方均值 T = mu2 / (mu1^2) - 1; end这里有个小细节值得注意：eig函数默认返回的特征值顺序不保证从大到小，但PR统计量只用到全体特征值的均值和平方均值，顺序不影响结果，所以不需要sort。如果哪位同学是在这个函数基础上做MME检测对比，那就要改成求最大特征值和最小特征值，到时候再排序也不迟。
4.3 门限生成与蒙特卡洛检测循环
门限生成思路就是纯噪声下反复调用prDetector，拿到融合统计量T_fus的大量样本，然后取分位数。纯噪声样本不乘任何噪声系数，因为PR统计量天然与绝对噪声功率无关，仿真里直接生成单位方差复高斯噪声就行，少一个变量简化不少。
% 在H0(纯噪声)下生成融合统计量的分布 T0 = zeros(1, nMc); for k = 1:nMc noise = (randn(L, N) + 1i*randn(L, N)) / sqrt(2); T0(k) = prDetector(noise); end gamma = quantile(T0, 1 - Pf_target); % 目标虚警概率对应的门限检测概率循环部分，每个蒙特卡洛实验里先生成一段主用户信号，然后每个节点独立乘衰落系数、叠加对应SNR的噪声，计算本地PR统计量，最后在融合中心取平均得到T_fus，再与门限比较。
Pd = zeros(size(snr_vec)); for idx = 1:length(snr_vec) snr = snr_vec(idx); sigma_n = sqrt(10^(-snr/10)); % 由SNR反推噪声标准差 cnt = 0; for mc = 1:nMc % 生成主用户信号：AR(1)过程，单位功率 s = 0.5 * (randn(1, L*N) + 1i*randn(1, L*N)); s = filter(1, [1 -0.9], s); s = s / sqrt(mean(abs(s).^2)); T_fus = 0; for m = 1:M h = (randn + 1i*randn) / sqrt(2); % 节点独立瑞利衰落 noise = sigma_n * (randn(1, L*N) + 1i*randn(1, L*N)) / sqrt(2); y = h * s + noise; % 接收信号 X = reshape(y, L, N); % 滑窗构造虚拟通道 T_local = prDetector(X); T_fus = T_fus + T_local; end T_fus = T_fus / M; % 等增益融合 if T_fus > gamma cnt = cnt + 1; end end Pd(idx) = cnt / nMc; end这段循环的结构分得很清楚，信号生成、每节点感知、融合判决三层嵌套。实际跑的时候，nMc=2000、M=4、SNR扫描7个点，大约几分钟能跑完。如果想快速验证代码逻辑，可以把nMc先降到500，趋势不会变，就是曲线毛糙一点。
4.4 绘制ROC与Pd-SNR曲线
检测器仿真除了固定Pf点看Pd随SNR变化，更完整的做法是画出ROC曲线，也就是固定SNR，扫不同的门限，得到一组Pd-Pf点。门限扫描不能瞎编，还是同样的思路——H0下生成T_fus的所有样本，然后依次取不同分位数作为门限。
% ROC曲线：固定SNR=-10dB，扫门限分位数 snr_fix = -10; sigma_n = sqrt(10^(-snr_fix/10)); T0_roc = zeros(1, nMc); for k = 1:nMc noise = (randn(L, N) + 1i*randn(L, N)) / sqrt(2); X = reshape(noise, L, N); T0_roc(k) = prDetector(X); end gammas = quantile(T0_roc, linspace(0.01, 0.99, 50)); Pf_roc = zeros(size(gammas)); Pd_roc = zeros(size(gammas)); for gi = 1:length(gammas) gamma_now = gammas(gi); pf_cnt = 0; pd_cnt = 0; for mc = 1:nMc % H0样本 noise = (randn(L, N) + 1i*randn(L, N)) / sqrt(2); X0 = reshape(noise, L, N); T0_now = mean(arrayfun(@(m) prDetector(X0), 1:M)); % 简化示意 if T0_now > gamma_now pf_cnt = pf_cnt + 1; end % H1样本 s = 0.5 * (randn(1, L*N) + 1i*randn(1, L*N)); s = filter(1, [1 -0.9], s); s = s / sqrt(mean(abs(s).^2)); T1_now = 0; for m = 1:M h = (randn + 1i*randn) / sqrt(2); noise_m = sigma_n * (randn(1, L*N) + 1i*randn(1, L*N)) / sqrt(2); X1 = reshape(h * s + noise_m, L, N); T1_now = T1_now + prDetector(X1); end T1_now = T1_now / M; if T1_now > gamma_now pd_cnt = pd_cnt + 1; end end Pf_roc(gi) = pf_cnt / nMc; Pd_roc(gi) = pd_cnt / nMc; end上面这段ROC代码里我用了arrayfun简化写法，实际跑的时候为了效率建议用循环替代。ROC曲线跑起来明显比单点Pd-SNR要慢，因为每个门限点都要重新蒙特卡洛。工程上常用的做法是先只跑几个关键门限点，比如Pf=0.01、0.05、0.1、0.2，连成折线，足够说明趋势。
5. 仿真结果与关键参数的影响规律
我自己按上面代码跑的典型结果，先说明数值会有随机波动，但趋势是稳定的。在Pf_target=0.1、M=4、L=4、N=128、nMc=2000的参数组合下，单节点SNR=-15dB时检测概率大约只有0.12，说明这个信噪比已经接近检测极限；SNR=-10dB时Pd可以到0.55左右；SNR=-7dB时超过0.85；到SNR=-3dB以上基本就是0.98以上了。加了四个节点融合后，相比单节点，同一个SNR点检测概率普遍能提升0.2到0.4，尤其在中等信噪比区域提升最明显。
5.1 节点数M与滑窗长度L的取舍
节点数M从1增加到8，检测性能确实是单调提升的，但边际增益递减。M从1到4的增益非常可观，M从4到8的增益明显变缓。原因也很好理解：每增加一个节点，融合中心就多一个独立观测，置信度自然增加，但多个节点的统计量取平均之后，方差随M的增大是反比关系，方差从大变小带来的边际收益在M足够大后会变弱。所以工程上不是节点越多越好，还要考虑报告信道带宽和系统成本，通常在4到8个节点之间选。
滑窗长度L这个参数比很多人想的敏感。L过小，比如只有2，特征值数量太少，PR统计量对“特征值不均匀”的描述能力有限；L过大，比如取16，在总采样长度不变的情况下，N会等比例缩小，协方差矩阵的估计质量恶化，检测概率反而下降。我的经验是L取4到8之间比较稳妥。如果想调性能，优先调N而不是调L，因为协方差矩阵估计质量基本由N决定，N越大，T_PR在H0下的方差越小，门限可以压得更紧，H1下的统计量也能更稳定地抬起来。
5.2 与能量检测器的性能对比
在这个场景下我做过的对比实验也值得提一下。当噪声功率被精确知道时，能量检测器的性能确实略优于PR检测器，因为能量检测把接收信号的全部能量都利用了，而PR检测器还要经过特征值变换，丢掉了一部分能量信息。但问题是“精确知道噪声功率”在真实系统里几乎不存在。我做过一个简单的实验：仿真时把噪声功率估计值故意偏移±1dB，能量检测器的曲线立刻变得很难看，在SNR=-12dB以下几乎失去意义；PR检测器基本不受影响，因为在特征值比值计算中噪声功率这一项被消掉了。
如果你以后要写工程报告或者论文，我建议把噪声不确定性这个维度加进对比实验，这比单纯比较理想条件下的Pd更能体现PR检测器的真正价值。很多论文里说能量检测器“抗噪声不确定性差”，其实只有在实验里让噪声功率偏移，才能把这个结论落到数据上，否则读者感受不到差异。
5.3 门限分位数与虚警概率的联动
还有一个容易忽略的点是：门限gamma取的是H0分布的分位数，Pf_target设得越小，门限越大，检测概率也会跟着下降。这是检测器内部不可避免的trade-off。如果你发现实验结果里Pd始终上不去，先去看看Pf_target设的多少，有时候不是检测器不行，而是你把虚警概率压得太低了。认知无线电系统通常会根据主用户保护的严格程度设定Pf，一般在0.05到0.2之间，不必追求越低压越好。
6. 工程实现中的坑位与绕坑方法
这套仿真代码看起来不长，但我在实际调试时踩过的坑不比写代码的时间少。挑几个最容易影响结果的问题集中说一下，算是给后来者排雷。
6.1 协方差矩阵的秩与正定性问题
样本协方差矩阵R是L乘L的矩阵，要保证它满秩，一个基本条件是N不小于L。但这只是理论上限，实际使用里N必须明显大于L才够稳。当N和L接近时，小特征值会被严重低估，导致T_PR在H0下方差变大、门限抖动，最直接的后果就是ROC曲线不稳定。如果你跑出来的曲线一次一个样，先检查是不是N/L太小。我的默认建议是N至少是L的10倍以上，N=128对应L不超过8，这也是我上面参数选L=4而不是L=16的原因之一。
6.2 reshape滑窗构造虚拟通道的顺序别搞反
Matlab的reshape函数是按列填充的，也就是说，对于一个长度为L*N的行向量y，reshape(y, L, N)生成的矩阵第一列是y的前L个元素，第二列是接下来的L个元素。如果你心里想的是“第一行是连续N个时间样本”，这种构造方式会完全错位，导致所谓的“虚拟通道”之间根本没有你想要的相邻采样相关结构，检测性能自然崩掉。我自己第一次实现时就栽在这里，当时还以为是自己公式写错了，后来手推了一个小例子才发现是reshape顺序问题。
要验证构造方式对不对，很简单：生成一个线性序列1到L*N，reshape成矩阵后打印第一列，看它是不是[1;2;...;L]。如果是，说明相邻采样被安排到了同一列的不同行，这才是我们想要的虚拟通道结构。
6.3 报告链路噪声与量化误差的干扰
前面所有仿真都假设报告信道理想，节点往融合中心传的统计量分毫不差。真实系统里这是不可能的。PR统计量本身是个连续值，要想数字传输，必然要做量化。量化bit数太少，比如只保留1bit，那其实等价于硬融合，软融合的性能优势就全丢了；量化bit数太多，报告链路开销又变大。我实践下来的经验是，软融合最少也要3bit量化，最好5bit以上，否则统计量的幅度信息损失严重。
报告链路本身也可能加噪声。如果报告链路的SNR比较高，比如高于10dB，对检测性能影响很小；如果报告链路质量差，还不如改用硬融合更划算。这一点在系统设计时经常被忽略，建议你在仿真后期加上一个简单的报告链路模型：在T_m上叠加一个零均值高斯扰动，扰动方差按报告链路的SNR反算，然后重新做融合、重新标定门限。整个过程跑完之后，你会发现门限又要重算一遍——这就对了，任何链路模型改变后，门限都必须跟着重标定，这是协作频谱感知仿真里最需要刻进肌肉记忆的一条规则。
我做这套仿真最深的体会是：PR检测器本身的实现并不复杂，真正费时间的往往是你没有意识到哪些参数和环境变了之后门限必须重算。只要把这个“改参数就重标门限”的习惯养好，后面不管是换信号模型、换信道模型还是加报告链路，都只是往这套框架里填模块的事情。希望这篇内容能帮你少熬几个晚上，如果你跑出来的趋势和我不一样，先别急着怀疑代码，优先检查协方差矩阵构造顺序和门限标定位置，八成问题就出在那里。
