optics matlab实现

关于optics算法的一些基本概念，在此一一忽略。

先求得所有节点的核心距离，用cd矩阵表示；

然后对每个节点进行处理，这个时候不需要考虑该节点是不是核心对象，按顺序取节点，如果该拓展点是核心对象，处理该节点，然后更新其他所有节点的可达距离。这个时候该核心对象的邻居节点求得的可达距离肯定是比较小的，所以优先处理这些邻居节点，最后得到的排序结果也是紧挨着的；处理完所有的邻居节点，这个时候，那就是得到了相对来说离这团数据很近的一个节点，其实，这个相当于还是盲取，这个本身就是要在剩下的节点中随机去一个数据，所以，都关系不大啦，然后得到该节点的邻居直到处理完所有的节点。
如果当前节点不是核心对象，那么还是任意取值。反正就是要把所有的节点处理完，所以影响不大。
如果处理的节点不是某一团数据的一个边缘值，那么总能处理到该团数据的核心值，这样也能改变其他未被处理过的对象的可达距离，他们最终还是一类；

最后得到RD矩阵，就是可达距离，而输出序列就是order矩阵。根据order矩阵得到输出序列，然后rd矩阵对应可达距离，这样即可得到如下图：

关于上图的解释：

每一个凹槽对应于一个聚类结果，因为核心对象的可达距离肯定很大，对吧，而在一个聚类簇中，其他对象到该核心对象的距离肯定很小，所以，每一个凹槽的开始就是一个新的聚类的开始。

代码如下：

% -------------------------------------------------------------------------
% Function: 
% [RD,CD,order]=optics(x,k)
% -------------------------------------------------------------------------
% Aim: 
% Ordering objects of a data set to obtain the clustering structure 
% -------------------------------------------------------------------------
% Input: 
% x - data set (m,n); m-objects, n-variables
% k - number of objects in a neighborhood of the selected object
% (minimal number of objects considered as a cluster)
% -------------------------------------------------------------------------
% Output: 
% RD - vector with reachability distances (m,1)
% CD - vector with core distances (m,1)
% order - vector specifying the order of objects (1,m)
% -------------------------------------------------------------------------
% Example of use:
% x=[randn(30,2)*.4;randn(40,2)*.5+ones(40,1)*[4 4]];
% [RD,CD,order]=optics(x,4)
% -------------------------------------------------------------------------
% 

function [RD,CD,order]=optics(x,k)

[m,n]=size(x);
CD=zeros(1,m);%核心距离
RD=ones(1,m)*10^10;%可达距离

% Calculate Core Distances 得到每个点的核心距离 计算每个点和其他节点的距离，然后进行排序，得到前
% k个点 故而得到核心距离 
for i=1:m	
    D=sort(dist(x(i,:),x))
    CD(i)=D(k);  
end
order=[];
seeds=[1:m];
ind=1;
while ~isempty(seeds)
    ob=seeds(ind);      
    seeds(ind)=[] ;
    order=[order ob];
    mm=max([ones(1,length(seeds))*CD(ob);dist(x(ob,:),x(seeds,:))]);
    ii=(RD(seeds))>mm;
    RD(seeds(ii))=mm(ii);
    [i1 ind]=min(RD(seeds));
end   

%RD(1)=max(RD(2:m))+.1*max(RD(2:m));
kk=RD>k;
RD(kk)=10;
plot(RD(order))
function [D]=dist(i,x)

% function: [D]=dist(i,x)
%
% Aim: 
% Calculates the Euclidean distances between the i-th object and all objects in x	 
% Input: 
% i - an object (1,n)
% x - data matrix (m,n); m-objects, n-variables	    
%                                                                 
% Output: 
% D - Euclidean distance (m,1)

[m,n]=size(x);
D=(sum((((ones(m,1)*i)-x).^2)'));
D=sqrt(D);
if n==1
   D=abs((ones(m,1)*i-x))';
end

　测试数据如下：

x=

2 3
2 4
1 4
1 3
2 2
3 2
8 7
8 6
7 7
7 6
8 5
100 2
8 20
8 19
7 18
7 17
8 21

k=4