DBSCAN算法

注：调用该算法时，需要将SetOfPoints设为全局变量并在用到该变量时需要显式的声明。% 算法描述%   基于密度的聚类算法% 参数%   Eps: 一个阈值,表示以某个样本点为圆心所形成的圆的半径%   MinPts: 一个阈值,表示以某个样本点为圆心以Eps为半径所形成的圆内至少要包含的点的数目% 返回值%   void% 注示：%   -1 表示未分类标号%   0 表示噪音，噪音点即从任何点出发都不会密度可达的点，处于团边界上的点不一定为噪音点%   正整数表示簇标号function void = DBSCAN( Eps, MinPts )     global SetOfPoints;    % 未分类标号和噪音标号    global UNCLASSIFIED;    UNCLASSIFIED = -1;    global NOISE;    NOISE = 0;    normal_classid = 1;      % 数据集长度和属性数目    setLength = size(SetOfPoints, 1);    % 预处理    for i = 1 : setLength        % 将所有数据点初始化为未分类状态        SetOfPoints(i, 3) = UNCLASSIFIED;        % 在每一条样本数据后追加1个序号属性        SetOfPoints(i, 4) = i;    end        % 开始分类    for i = 1 : setLength        if SetOfPoints(i, 3) == UNCLASSIFIED            if ExpandCluster(SetOfPoints(i, :), normal_classid, Eps, MinPts)                normal_classid = normal_classid + 1;            end        end    end    visualization(SetOfPoints);end% 算法描述%   寻找从当前样本点p密度相连的最大的数据点集，直到没有新的数据点密度可达为止% 参数描述%   Point: 待处理的当前样本点%   ClId: 当前簇的标号%   Eps: 一个阈值,表示以某个样本点为圆心所形成的圆的半径%   MinPts: 一个阈值,表示以某个样本点为圆心以Eps为半径所形成的圆内至少要包含的点的数目% 返回值%   1: 当前点p是核心点%   0: 当前点为噪音点function res = ExpandCluster(Point, ClId, Eps, MinPts)    global NOISE;    global UNCLASSIFIED;    seeds = regionQuery(Point, Eps);    seedsLen = size(seeds, 1);        % 点Point可能是噪音点    if seedsLen < MinPts        changeClId(Point, NOISE);        res = 0;        return;    else        res = 1;        % 将种子集合中的点的簇标号更新为ClId        changeClIds(seeds, ClId);        seeds = deleteElementByValue(seeds, Point(4), 4);        while seedsLen ~= 0            % 取种子集合中的第一个点            currentP = seeds(1,:);            result = regionQuery(currentP, Eps);            resultLen = size(result, 1);            % 若currentP是核心点            if resultLen >= MinPts                for i = 1 : resultLen                    resultP = result(i, :);                    if resultP(3) == NOISE || resultP(3) == UNCLASSIFIED                        % 若resultP未分类，则resultP可能是核心点，故将其加入种子集合中                        if resultP(3) == UNCLASSIFIED                            seeds(size(seeds, 1) + 1, :) = resultP;                        end                        % 无论点resultP之前的簇标号是噪音还是未分类，它们当前都是从currentP这个核心点密度可达的，所以肯定不是噪音点                        changeClId(resultP, ClId);                    end                end            end            % 从种子集合中删除currentP            seeds = deleteElementByValue(seeds, currentP(4), 4);            seedsLen = size(seeds, 1);        end        return;    endend% 算法描述%   寻找点p在Eps邻域内的所有点% 参数%   p: 当前待处理的点%   Eps: 一个阈值,表示以某个样本点为圆心所形成的圆的半径% 返回值%   seeds: 从点p出发直接密度可达的所有点function seeds = regionQuery(p, Eps)    global SetOfPoints;    seeds = zeros(0, 4);    % 种子集合长度    seedsLen = 0;    setLength = size(SetOfPoints, 1);        % 搜索所有在点p的Eps邻域内的点    for i = 1 : setLength        if abs(p(1) - SetOfPoints(i, 1)) <= Eps && abs(p(2) - SetOfPoints(i, 2)) <= Eps            if sqrt((p(1) - SetOfPoints(i, 1)) * (p(1) - SetOfPoints(i, 1)) + (p(2) - SetOfPoints(i, 2)) * (p(2) - SetOfPoints(i, 2))) <= Eps                % 将第i个点加入种子集合中                seedsLen = seedsLen + 1;                seeds(seedsLen, :) = SetOfPoints(i, :);            end        end    endend% 算法描述%   从一个集合中删除给定列等于给定值的元素% 参数%   set: 目标集合%   elementValue: 给定值%   columnIndex: 给定列索引% 返回值%   newSet: 新的集合function newSet = deleteElementByValue(set, elementValue, columnIndex)    newSet = [];    setLen = size(set, 1);    if setLen == 0        return;    end    for i = 1 : setLen        if set(i, columnIndex) == elementValue            set(i, :) = [];            break;        end    end    newSet = set;end% 算法描述%   修改单个样本点的簇标号% 参数%   Point: 待修改簇标号的点%   ClId: 簇标号% 返回值%   voidfunction void = changeClId(Point, ClId)    global SetOfPoints;    SetOfPoints(Point(4), 3) = ClId;end% 算法描述%   修改多个样本点的簇标号% 参数%   PointSet: 待修改簇标号的点集%   ClId: 簇标号% 返回值%   voidfunction void = changeClIds(PointSet, ClId)    pointSetLen = size(PointSet, 1);    % 将PointSet中的点的簇标号修改为ClId    for i = 1 : pointSetLen        changeClId(PointSet(i,:), ClId);    endend% 算法描述%   可视化聚类结果，即为不同的簇中的点标注不同的颜色和形状% 参数描述%   表示数据集的全局变量名称% 返回值%   voidfunction void = visualization( rings )    ringsLen = size(rings, 1);    colorMatrix = ['.y','.m','.c','.r',',g','.b'];        for i = 1 : ringsLen        if rings(i, 3) == 0            % 如果是噪音点则设置为黑色的点            plot(rings(i, 1), rings(i, 2), '.k');        else            % 如果存在多个簇，则循环的为每个簇设置形状及颜色            colorIndex = mod(rings(i, 3), size(colorMatrix, 2)) + 1;            switch colorIndex                case 1                    plot(rings(i, 1), rings(i, 2), 'b+');                   case 2                    plot(rings(i, 1), rings(i, 2), 'm*');                 case 3                    plot(rings(i, 1), rings(i, 2), 'ro');                 case 4                    plot(rings(i, 1), rings(i, 2), 'cs');                 case 5                    plot(rings(i, 1), rings(i, 2), 'gd');                 case 6                    plot(rings(i, 1), rings(i, 2), 'yp');             end        end        hold all;    endend