最大信噪比盲源分离算法研究

最大信噪比盲源分离算法研究
1 程序说明
该算法参考文为“基于最大信噪比的盲源分离算法”（《计算机仿真》，2006），程序的Matlab代码如下：
function[ys,w]=SNR_Max(x)
[n,T]=size(x);
x=x-mean(x')'*ones(1,T); %去均值
%---------------白化处理-------------------------
[F,D]=eig((x*x')/T);
v=F*(D^(-0.5))*F';
x=v*x;
p=80;
a=ones(1,p)/p;
x=x';
S=filter(a,1,x); %滑动平均处理
U=cov(S-x,1); %计算协方差
V=cov(x,1);
[W d]=eig(V,U); %计算特征值
ys=(x*W)'; %源信号估计
mplot(x);
figure(2);
mplot(ys);
上述程序中引用了一个画图子程序mplot(s)，能画小于6个信号的图形，该程序代码如下：
function mplot(s)
[dim,sample]=size(s);
if dim>sample
s=s';
[dim,sample]=size(s);
end
if dim>6
error('dim can not more than six');
end
for i=1:dim
subplot(dim,1,i);
plot(s(i,:));
axis('tight');
set(gca,'XTick',[]);
end
程序调用说明：输入混合信号x为一个n*T阶数据矩阵，输出的分离信号ys是n*T阶
估计数据矩阵，W 为n*n 阶分离矩阵，调用格式为[ys,W]=SNR_Max(x)。

2 算法说明
基于最大信噪比的盲源分离算法，以盲源分离效果越好时信噪比越大这一特点，建立信噪比目标函数，把求优过程转化成广义特征值求解，用求出的广义特征值构成特征向量矩阵----分离矩阵，该算法是全局最优的盲源分离算法，具有低的计算复杂度。

设T N n s n s n s )](),....,([)(1=为N 维源信号向量，)(n x 为N 维混合信号向量，A 为N*N 阶瞬时线性混合矩阵，信号的混合模型可表示为：
)()(n As n x = （2-1）
盲源分离就是仅有观测信号)(n x 和源信号)(n s 的概率分布先验知识来恢复出)(n s 。

即寻找一个N*N 阶的分离矩阵W ，使其输出
)()()()(n Gs n WAs n Wx n y === （2-2）
为)(n s 的一个估计，称)(n y 为估计信号或分离信号。

这里的G 为全局变换矩阵。

下面介绍以信噪比最大化为准则，建立目标函数以学习的方式确定分离矩阵W 的算法。

2.1 信噪比目标函数
根据盲分离式（2-2），把源信号s 与估计信号y 的误差e=s-y 作为噪声信号，建立信噪比函数为：
T
T
T T y s y s s s e e s s SNR )
()(log 10log 10-∙-∙=∙∙= （2-3） 由于源信号s 是未知的，考虑到y(n)含有噪声，因此用估计信号y(n)的滑动平均~
y 代替源信号s ，式（2-3）变为
T
T
T T
y y y y y
y e e s s SNR )
()(log 10log
10~
~~~
-∙-∙=∙∙= （2-4） 式中：∑=-=-=p
j i i p j j n y p n y 0
~
1,....,2,1,0),(1)(，滑动平均长度P 可根据信号的噪声特点
选取（P 可为小于100的整数，程序中设定p=80），为了简化计算，将式（2-4）分子中的~
y 用y 代替，因此得到最大信噪比目标函数为 T
T y y y y y y SNR y F )
)((log
10)(~
~--∙== （2-5）
式中：;;~
~
x W y Wx y == W 为分离矩阵；~
x 为混合信号经滑动平均处理后的信号，即
∑=-=-=p
j i p i j n x p n x 0
~
1,.....,2,1,0),(1)( （2-6）
方程（2-5）可写成
U
V W C W WCW
W x x x x W W Wxx y y y y y y x W F T
T
T
T T
T T
T log
10log
10))((log
10)
)((log 10),(~
~
~
~
~
==--=--∙= （2-7）
式中：T
T
x x x x C xx C ))((,~
~
~
--==为相关矩阵；T
T
W C W U WCW V ~
;==。

2.2 分离算法的推导
以分离矩阵W 对式（2-7）两边求梯度，得
~
22C U
W C V W W F -=∂∂ （2-8） 由于目标函数F(W,x)的极点值为式（2-8）的零点，因此得
~C W U
V
WC = （2-9）
通过求解（2-9）就可以得分离矩阵∧W ，∧W 为矩阵1
~
-∙C
C 的特征向量。

（可用Matlab 中eig()
函数求得）。

当然分离的源信号向量为x W y ∧
=，其中y 的每一行代表一个分离信号，或称y 是源信号s 的估计。

只要分理出的信号不相关时，通过∧
W 得到的信号是统计独立的，因此上述算法可分离独立的源信号。

3仿真结果与分析
选择两路信号，一路为超高斯分布信号，一路为亚高斯分布的正弦信号为源信号，随机产生的线性混合矩阵A={0.1509 0.3784}{0.6979 0.8600}，取数据长度L=1000。

这两个信号分别按如下方式产生
⎥⎥⎦
⎤
⎢⎢⎣⎡∆+-=⎥⎦⎤⎢⎣⎡=)100sin(3]9
11)23,mod([
215t t S S S π 其中t 取离散值1,2，…，1000，4
102-⨯=∆，信号S1的峭度大于零，属于超高斯信号，S2的峭度小于零，属于亚高斯信号，这两路信号的波形如图1所示。

（a ）原始信号—超高斯信号
-1
-0.8
-0.6-0.4-0.200.2
0.40.60.81
（b ）原始信号---亚高斯信号
图1 原始信号
两路原始信号经过混合矩阵混合后的信号如图2所示。

图2 混合信号
经过最大信噪比盲源分离算法分离后的波形如图3所示。

图3 基于最大信噪比盲源分离算法得到的分离信号
由图3可以看出，经过最大信噪比盲源分离算法的操作，得到了非常好的分离效果。

具体的Matlab实现程序如下：
clear;
t=1:1000;
E=0.0002;
S1=((mod(t,23)-11)/9).^5+3;
plot(S1)
A=[0.1509 0.6979 ;0.3784 0.8600]
S2=sin(100*pi*E*t);
figure(2)
plot(S2)
S=[S1; S2];
x=A*S;
[ys,W]=SNR_Max(x)。