Home > voicebox > minspane.m

minspane

PURPOSE ^

MINSPANE calculate minimum spanning tree using euclidean distance [p,s]=X

SYNOPSIS ^

function [p,s]=minspane(x)

DESCRIPTION ^

MINSPANE calculate minimum spanning tree using euclidean distance [p,s]=X

 Inputs:  x(n,d)    d-dimensional data points, one per row

 Outputs: p(n-1,1)  indices of the parent of each node within the tree
                      (data point n is the root)
          s(n-1,1)  list of the edges in ascending order of euclidean
                      distance. Thus the shortest edge goes from node
                      s(1) to node p(s(1))

 The minimum spanning tree (or shortest spanning tree ) defines a set of
 n-1 links that interconnect n points (or nodes) with the minimum total
 length. We represent these links in the form of a tree with node n as
 the root. Each node (except node n) has a unique parent node that is
 given in the output vector p; it is possible for several nodes to share
 the same parent.
 It can be useful to know which of the links are the longest, so the output
 argument lists them in ascending order. s could be calculated directly by
 calculaing the length, l, of each link as follows:
       l=sqrt(sum((x(1:n-1,:) - x(p,:)).^2),2);
       [v,s]=sort(l);

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [p,s]=minspane(x)
0002 %MINSPANE calculate minimum spanning tree using euclidean distance [p,s]=X
0003 %
0004 % Inputs:  x(n,d)    d-dimensional data points, one per row
0005 %
0006 % Outputs: p(n-1,1)  indices of the parent of each node within the tree
0007 %                      (data point n is the root)
0008 %          s(n-1,1)  list of the edges in ascending order of euclidean
0009 %                      distance. Thus the shortest edge goes from node
0010 %                      s(1) to node p(s(1))
0011 %
0012 % The minimum spanning tree (or shortest spanning tree ) defines a set of
0013 % n-1 links that interconnect n points (or nodes) with the minimum total
0014 % length. We represent these links in the form of a tree with node n as
0015 % the root. Each node (except node n) has a unique parent node that is
0016 % given in the output vector p; it is possible for several nodes to share
0017 % the same parent.
0018 % It can be useful to know which of the links are the longest, so the output
0019 % argument lists them in ascending order. s could be calculated directly by
0020 % calculaing the length, l, of each link as follows:
0021 %       l=sqrt(sum((x(1:n-1,:) - x(p,:)).^2),2);
0022 %       [v,s]=sort(l);
0023 
0024 %      Copyright (C) Mike Brookes 2000-2009
0025 %      Version: $Id: minspane.m 1446 2012-02-22 09:15:35Z dmb $
0026 %
0027 %   VOICEBOX is a MATLAB toolbox for speech processing.
0028 %   Home page: http://www.ee.ic.ac.uk/hp/staff/dmb/voicebox/voicebox.html
0029 %
0030 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0031 %   This program is free software; you can redistribute it and/or modify
0032 %   it under the terms of the GNU General Public License as published by
0033 %   the Free Software Foundation; either version 2 of the License, or
0034 %   (at your option) any later version.
0035 %
0036 %   This program is distributed in the hope that it will be useful,
0037 %   but WITHOUT ANY WARRANTY; without even the implied warranty of
0038 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0039 %   GNU General Public License for more details.
0040 %
0041 %   You can obtain a copy of the GNU General Public License from
0042 %   http://www.gnu.org/copyleft/gpl.html or by writing to
0043 %   Free Software Foundation, Inc.,675 Mass Ave, Cambridge, MA 02139, USA.
0044 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0045 
0046 [np,nd]=size(x);
0047 
0048 % first do delauny tessalation to find feasible edges
0049 
0050 t=delaunayn(x); % nd+1 vertices per polytope
0051 nt=size(t,1);
0052 cn=choosenk(nd+1,2);
0053 nf=size(cn,1);
0054 ee=zeros(nt,nf,2);     % space for all the edges
0055 ee(:,:,1)=t(:,cn(:,1));
0056 ee(:,:,2)=t(:,cn(:,2));
0057 ee=reshape(ee,[],2);
0058 mk=ee(:,1)>ee(:,2);
0059 ee(mk,:)=ee(mk,2:-1:1); % make all edges in ascending order
0060 % ees=sparse(ee(:,1),ee(:,2),1);      % remove duplicates
0061 [er,ec]=find(sparse(ee(:,1),ee(:,2),1));      % remove duplicates
0062 ee=[er,ec];
0063 ne=size(ee,1);
0064 
0065 % now apply Kruskal shortest spanning tree algorithm
0066 
0067 sz=sum((x(ee(:,1),:)-x(ee(:,2),:)).^2,2);     % length^2 of each edge
0068 [vz,mz]=sort(sz);
0069 ee=ee(mz,:);        % sort edges into ascenging length order
0070 ts=ones(ne,1);      % size of each component
0071 tp=zeros(np,1);     % root node links
0072 ei=zeros(np-1,1);     % index of shortest spanning tree edges
0073 k=0;
0074 for i=1:ne
0075     i1=ee(i,1);
0076     j=tp(i1);
0077     while j         % find root node for ee(i,1)
0078         i1=j;
0079         j=tp(i1);
0080     end
0081     i2=ee(i,2);
0082     j=tp(i2);
0083     while j         % find root node for ee(i,2)
0084         i2=j;
0085         j=tp(i2);
0086     end
0087     if i1~=i2       % if they are different, merge them
0088         k=k+1;
0089         ei(k)=i;    % add to the shortest spanning tree
0090         if ts(i1)>ts(i2)
0091             tp(i2)=i1;          % make i2 a sub tree of i1
0092             ts(i1)=ts(i1)+ts(i2);
0093         else
0094             tp(i1)=i2;          % make i1 a sub tree of i2
0095             ts(i2)=ts(i1)+ts(i2);
0096         end
0097     end
0098 end
0099 ee=ee(ei,:);        % refine the edges to include only those in the minimum spanning tree
0100 
0101 % now arrange as a tree with point np as the head
0102 
0103 eet=sparse(ee(:,1),ee(:,2),1:np-1,np,np);
0104 p=zeros(np-1,1);              % points to parent node
0105 s=zeros(np-1,1);              % sorted index
0106 chn=np;                                 % start with the root node of the tree
0107 [rf,cf]=find(eet(:,chn));               % find any nodes that connect to it
0108 [rg,cg]=find(eet(chn,:));
0109 while ~isempty(rf) || ~isempty(rg)
0110     p(rf)=chn(cf);                      % set the parents
0111     rcf=rf(:)+np*(chn(cf(:))-1);
0112     s(eet(rcf))=rf(:);
0113     eet(rcf)=0;                         % delete the edges linking them to their parents
0114     p(cg)=chn(rg);
0115     rcg=chn(rg(:))+np*(cg(:)-1);
0116     s(eet(rcg))=cg(:);
0117     eet(rcg)=0;
0118     chn=[rf(:); cg(:)];                 % now search for their children
0119     [rf,cf]=find(eet(:,chn));
0120     [rg,cg]=find(eet(chn,:));
0121 end
0122

Generated on Tue 10-Oct-2017 08:30:10 by m2html © 2003