Bag-Of-Words中K-Means聚类的效率优化

转自：http://www.zhizhihu.com/html/y2010/2257.html

 

        最初的Bag of words，也叫做“词袋”，在信息检索中，Bag of words model假定对于一个文本，忽略其词序和语法，句法，将其仅仅看做是一个词集合，或者说是词的一个组合，文本中每个词的出现都是独立的，不依赖于其他词 是否出现，或者说当这篇文章的作者在任意一个位置选择一个词汇都不受前面句子的影响而独立选择的。

        穿越到2011年10月24号写了篇博客《一个用BoW|Pyramid BoW+SVM进行图像分类的Matlab Demo》，希望有帮助。

       

        现在Computer Vision中的Bag of words来表示图像的特征描述也是很流行的。大体思想是这样的，假设有5类图像，每一类中有10幅图像，这样首先对每一幅图像划分成patch（可以是刚性分割也可以是像SIFT基于关键点检测的），这样，每一个图像就由很多个patch表示，每一个patch用一个特征向量来表示，咱就假设用Sift表示的，一幅图像可能会有成百上千个patch，每一个patch特征向量的维数128。

        接下来就要进行构建Bag of words模型了，假设Dictionary词典的Size为100，即有100个词。那么咱们可以用K-means算法对所有的patch进行聚类，k=100，我们知道，等k-means收敛时，我们也得到了每一个cluster最后的质心，那么这100个质心（维数128）就是词典里的100个词了，词典构建完毕。

        词典构建完了怎么用呢？是这样的，先初始化一个100个bin的初始值为0的直方图h。每一幅图像不是有很多patch么？我们就再次计算这些patch和和每一个质心的距离，看看每一个patch离哪一个质心最近，那么直方图h中相对应的bin就加1，然后计算完这幅图像所有的patches之后，就得到了一个bin=100的直方图，然后进行归一化，用这个100维的向量来表示这幅图像。对所有图像计算完成之后，就可以进行分类聚类训练预测之类的了。

        那么，这里影响效率的一个方面是构建词典时的K-means聚类，我在用的时候遇到了两个问题：1、内存溢出。这是由于一般的K-means函数的输入是待聚类的完整的矩阵，在这里就是所有patches的特征向量f合成的一个大矩阵，由于这个矩阵太大，内存不顶了。我内存为4G。2、效率低。因为需要计算每一个patch和每一个质心的欧拉距离，还有比较大小，那么要是循环下来这个效率是很低的。

        为了解决这个问题，我采用一下策略，不使用整一个数据矩阵X作为输入的k-means，而是自己写循环，每次处理一幅图像的所有patches，对于效率的问题，因为matlab强大的矩阵处理能力，可以有效避免耗时费力的自己编写的循环迭代。

        我自己的代码下载：

        代码下载链接：PG_BOW_DEMO.zip,如果无法下载，当然可以再发邮件问我要，pagelee.sd at gmail.com。

        Demo中的图像是我自己研究中用到的一些Action的图像，我都采集的简单的一共6类，每一类60幅，40训练20测试。请注意图像的版权问题，自己研究即可，不能商用。

分类器用的是libsvm，最好自己mex重新编译一下。

        Demo中的版权我也已经注明。

        具体我用的代码形式：

</span>   %% for n=1:niters % Save old centres to check for termination if n>1 e2=max(max(abs(centres - old_centres))) ; elsee2=0; end  old_centres = centres; tempc=centres; num_points=ones(1,ncentres);   for f = 1:nimages %省略了一些读取特征的细节 data = features.data;%每个图像的所有patch的特征矩阵 [ndata, data_dim] = size(data);   id = eye(ncentres); d2 = distance(data,centres);%计算两个矩阵的欧式距离，这样计算以及下面几句代码都大大提高计算效率 % Assign each point to nearest centre [minvals, index] = min(d2', [], 1); post = id(index,:); % matrix, if word i is in cluster j, post(i,j)=1, else 0;   num_points = num_points + sum(post, 1);   for j = 1:ncentres tempc(j,:) =  tempc(j,:)+sum(data(find(post(:,j)),:), 1); end  end  for j = 1:ncentres if num_points(j)>0 centres(j,:) =  tempc(j,:)/num_points(j); endendif n > 1 % Test for termination   newError=max(max(abs(centres - old_centres))); if newError<minError minError=newError; end  if newError < 0.01 %某个聚类要收敛的阈值 fprintf('Saving texton dictionary\n'); save (['/','centres'],'centres');      % save the settings of descriptor in opts.globaldatapath break; end  fprintf('The %d th interation finished \n',n); endend

下面是求欧式距离的distance函数：

</span>   function d = distance(a,b) % DISTANCE - computes Euclidean distance matrix % % E = distance(A,B) % %    A - (MxD) matrix %    B - (NxD) matrix % % Returns: %    E - (MxN) Euclidean distances between vectors in A and B % % % Description : %    This fully vectorized (VERY FAST!) m-file computes the %    Euclidean distance between two vectors by: % %                 ||A-B|| = sqrt ( ||A||^2 + ||B||^2 - 2*A.B ) % % Example : %    A = rand(100,400); B = rand(200,400); %    d = distance(A,B);   % Author   : Roland Bunschoten %            University of Amsterdam %            Intelligent Autonomous Systems (IAS) group %            Kruislaan 403  1098 SJ Amsterdam %            tel.(+31)20-5257524 %            bunschot@wins.uva.nl % Last Rev : Oct 29 16:35:48 MET DST 1999 % Tested   : PC Matlab v5.2 and Solaris Matlab v5.3 % Thanx    : Nikos Vlassis   % Copyright notice: You are free to modify, extend and distribute %    this code granted that the author of the original code is %    mentioned as the original author of the code.   if (nargin ~= 2) b=a; end  if (size(a,2) ~= size(b,2)) error('A and B should be of same dimensionality'); end  %aa=sum(a.*a,1); bb=sum(b.*b,1); ab=a'*b; %d = sqrt(abs(repmat(aa',[1 size(bb,2)]) + repmat(bb,[size(aa,2) 1]) - 2*ab)); aa=sum(a.*a,2); bb=sum(b.*b,2); ab=a*b'; d = sqrt(abs(repmat(aa,[1 size(bb,1)]) + repmat(bb',[size(aa,1) 1]) - 2*ab));

另外从水木上看到一个人写的短小精悍的kmeans，大家参考下吧：

function label = litekmeans(X, k) n = size(X,2); last = 0; label = ceil(k*rand(1,n));  % random initialization while any(label ~= last) [~,~,label] = unique(label);   % remove empty clusters E = sparse(1:n,label,1,n,k,n);  % transform label into indicator matrix center = X*(E*spdiags(1./sum(E,1)',0,k,k));    % compute center of each cluster last = label; [~,label] = max(bsxfun(@minus,center'*X,0.5*sum(center.^2,1)'),[],1); % assign samples to the nearest centers end

 

另外在最新的一篇博文《一个用BoW|Pyramid BoW+SVM进行图像分类的Matlab Demo》详细讲了怎么结合SVM分类的问题和Matlab Demo，很直观，希望有用。

下面再次说一下基本的理解：

1、个人感觉Bag of Words/Bag of Features原理都差不多，个人还是当一个概念来理解的，可能bag of features涵盖的内容更加具体，特征更有代表性吧。

2、聚类的初始点可以自己选取，也可以随机选取，其实对于聚类的结果虽然有影响，但是最后用于分类的时候影响没那么明显。

-------------

刚性分割 其实术语是叫 dense sift 你搜一下；

代码也不是纯碎我写的 主要还是集成了 http://www.di.ens.fr/willow/pdfs/cvpr06b.pdf 的代码 稍微做了自己分为方便的改动。