设有均值为零的随机信号向量
,其
自相关矩阵为

很明显,
是对称矩阵,且是非负定的(所有特征值都大于或等于0)。

现在,寻找一个线性变换
对
进行变换,即
,使得



上式的含义是:
y的各分量是不相关的,即
。通常将这个过程称为“空间解相关”、“空间白化”或“球化”。
称为空间解相关矩阵(空间白化矩阵、球化矩阵)。
![E[y_i y_j]=\delta_{ij}](https://chunqiu.blog.ustc.edu.cn/wp-content/plugins/latex/cache/tex_b3b63da6510f7d257532071b96c76da9.gif)

由
的性质可知,其存在特征值分解:




令
\begin{equation}\label{eq:B}B=\Sigma^{-1/2} Q^T\end{equation}
则有
因此,通过矩阵
线性变换后,
的各个分量变得不相关了。


对于
来说,特征值分解和奇异值分解是等价的,而奇异值分解的数值算法比特征值分解的数值算法具有更好的稳定性,因此一般都用奇异值分解来构造空间解相关矩阵
。


应该注意到,
“空间解相关”不能保证各分量信号之间的“独立性”,但它能够简化盲分离算法或改善分离算法的性能。
注:以上来自戴老师《盲信号处理》的课件。
最为熟知的例子是白噪声。元素
可以是一个时间序列在相继时间点
的值,且在噪声序列中没有时间上得相关性。术语“白”来自于白噪声的能谱在所有频率上是一个常数这一事实,就像含有各种颜色的白光谱一样。
白化的本质就是去相关加缩放。


式\eqref{eq:B}的解相关矩阵
肯定不是唯一的白化矩阵。容易看到,任何矩阵
(
为正交矩阵)也是白化矩阵。这是因为对
,下式成立:




一个重要的例子是矩阵
。这也是一个白化矩阵,因为它是用正交矩阵
左乘式\eqref{eq:B}的
得到的。这个矩阵称为
的逆均方根,并用
表示,因为它来自于均方根概念向矩阵的标准推广。





注:以上来自《Independent Component Analysis》
关于
白化代码的实现,其实很简单,下面给一个
function [z_w varargout] = myWhiten(z)
%--------------------------------------------------------------------------
% 语法:z_w = myWhiten(z);
% [z_w T] = myWhiten(z);
% 输入:z是一个mxn的矩阵,包含m维随机变量的各n个采样点。
% 输出:z_w是白化版本的z。T是mxm的白化变换矩阵。
%--------------------------------------------------------------------------
%% 计算样本协方差
R = cov(z'); % 1表示除以N来计算协方差
%% 白化z
[U D ~] = svd(R, 'econ'); % 用eig也行,[U, D] = eig(R);
%% 下面求白化矩阵
T = U * inv (sqrt(D)) * U'; % 称为协方差矩阵的逆均方根,inv计算不会太耗时间,因为D为对角阵。inv(sqrt(D))*U'也是一个可行白化矩阵
%% 乘以白化矩阵实现白化
z_w = T * z;
if (nargout == 2)
varargout{1} = T;
end
另外可以直接利用
FastICA的白化函数whitenv,它要配合其提供的PCA函数pcamat.m完成白化,其提供的例子中有提示如下
% EXAMPLE
% [E, D] = pcamat(vectors);
% [nv, wm, dwm] = whitenv(vectors, E, D);
之前写过关于FastICA工具箱的使用,再把使用
白化的例子给出来吧:
% 测试whitenv函数
clc
clear
close all
% 加载matlab自带的数据
load cities
stdr = std(ratings);
sr = ratings./repmat(stdr,329,1);
sr = sr';
figure
boxplot(sr','orientation','horizontal','labels',categories)
% 测试
firstEig = 1;
lastEig = 9;
s_interactive = 'off';
sr = remmean(sr); % 下面的pcamat和whitenv都没有去均值,这里先做取均值处理
[E, D] = pcamat(sr, firstEig, lastEig, s_interactive);
[nv, wm, dwm] = whitenv(sr, E, D);
figure
boxplot(nv','orientation','horizontal','labels',categories)
结果如下
