作者:桂。
时间:2017-04-06 20:26:01
链接:
声明:欢迎被转载,不过记得注明出处哦~
前言
本文非负矩阵分解(Nonegative matrix factorization,NMF)系列第二篇,主要介绍最基本的NMF原理及代码实现,内容主要包括:
1)基于Euclidean距离的NMF推导及实现;
2)基于KL散度的NMF推导及实现;
3)NMF应用示例
开始之前,有两点需要补充:
- 前面分析用的是X=AS形式,感觉别扭,好多文章都是用V = WH,后续打算也采用这也表达方式;
- NMF其实是含有约束的优化问题,但乘法算法可以巧妙得让我们只需讨论:无约束优化问题。
一、基于Euclidean距离的NMF推导及实现
考虑无约束优化问题:
利用梯度下降:
其中:
如果直接梯度下降,对于无约束的优化问题,我们不能保证结果都是非负的,下面巧妙之处来了:将梯度下降法变为乘法算法。
令:
梯度下降法变换为乘法算法:
真是巧妙!一个复杂的约束性优化问题,就让一个简单的无约束给解决了。这样一来,如果原矩阵为非负,W、H初始值同样非负,结果自始至终都是非负,直至迭代到满足收敛条件。
收敛性证明可以参考:Lee D D, Seung H S. Algorithms for Non-negative Matrix Factorization[C]// NIPS. 2000:556--562.
给出对应的代码实现:
function [W, H] = nmf(V, K, MAXITER)%Euclidean distanceF = size(V,1);T = size(V,2); rand('seed',0)W = 1+rand(F, K);% W = W./repmat(sum(W),F,1);H = 1+rand(K, T); ONES = ones(F,T); for i=1:MAXITER H = H .* (W'*V)./(W'*W*H+eps) ; W = W .* (V*H')./(W*H*H'+eps);end
其实关键的就是循环里的两行。
二、基于KL散度的NMF推导及实现
整个思路与Euclidean distance下的求解思路如出一辙。
考虑无约束优化问题:
利用梯度下降算法:
其中:
根据梯度下降算法转化为乘法算法:
令:
梯度下降算法改写为乘法算法:
收敛性证明可以参考:Lee D D, Seung H S. Algorithms for Non-negative Matrix Factorization[C]// NIPS. 2000:556--562.
对应代码:
function [W, H] = nmf(V, K, MAXITER)%KL-divergenceF = size(V,1);T = size(V,2); rand('seed',0)W = 1+rand(F, K);% W = W./repmat(sum(W),F,1);H = 1+rand(K, T); ONES = ones(F,T);for i=1:MAXITER H = H .* (W'*( V./(W*H+eps))) ./ (W'*ONES); W = W .* ((V./(W*H+eps))*H') ./(ONES*H');end
三、NMF应用示例
对于一个混合语音,如鼓点和管乐器混合的单通道声音,可以利用非负矩阵进行分解,实现语音信号的分离。
思路:
语音的时频分析,得到的语谱图是一个二维数据矩阵,其中鼓点、管乐器的概率分布不同,利用NMF可以实现信号的分离。
对应代码(NMF调用上面任何一个都可以):
% Read in audio file[x0 fs] = audioread('test_audio.wav');x = x0(1:5*fs);%read 5s wavnw = 1024;ni = 256;X = fft(enframe(x,nw,ni)');V = abs(X); % NMFK = 2; % number of basis vectorsMAXITER = 200; % total number of iterations to run[W, H] = nmf(V, K, MAXITER); % get the mixture phasephi = angle(X);%ReconstructX1 = W(:,1)*H(1,:);X2 = W(:,2)*H(2,:);s1 = zeros(1,length(x));s2 = zeros(1,length(x));for i = 1:size(X1,2) nic = (1+(i-1)*ni):(nw+(i-1)*ni); s1(1,nic) = s1(1,nic)+real(ifft(X1(:,i).*exp(1j*phi(:,i))))'; s2(1,nic) = s2(1,nic)+real(ifft(X2(:,i).*exp(1j*phi(:,i))))';ends1 = s1/max(abs(s1));s2 = s2/max(abs(s2));
可以看出这里K是给定的,即NMF实现分解需要给出先验的类别数。这里给出语谱图,图片可以更加直观地观察分离效果,结果图:
看看时域的分离效果:
分离效果还是不错的。
参考:
- Lee D D, Seung H S. Algorithms for Non-negative Matrix Factorization[C]// NIPS. 2000:556--562.
- 张贤达:《矩阵分析与应用,第二版》.