DBSCAN聚类算法

一、算法概述

　　DBSCAN（Density-Based Spatial Clustering of Applications with Noise，具有噪声的基于密度的聚类方法）是一种基于密度的空间聚类算法，相比其他的聚类方法，基于密度的聚类方法可以在有噪音的数据中发现各种形状和各种大小的簇。基于密度的聚类是寻找被低密度区域分离的高密度区域，这些高密度区域就是一个一个的簇，这里的密度指的是一个样本点的领域半径内包含的点的数目，基于密度的聚类方法可以用来过滤噪声孤立点数据，发现任意形状的簇（因为它将簇定义为密度相连的点的最大集合），与k-means算法的不同之处在于它不需要事先指定划分的簇的数目。

二、算法涉及的基本概念

E领域：给定对象半径为E内的区域称为该对象的E领域；
核心对象：如果给定对象E领域内的样本点数大于等于MinPts，则称该对象为核心对象；
直接密度可达：对于样本集合D，如果样本点q在p的E领域内，并且p为核心对象，那么对象q从对象p直接密度可达。
密度可达：对于样本集合D，给定一串样本点p1,p2….pn，p=p1,q=pn,假如对象pi从pi-1直接密度可达，那么对象q从对象p密度可达。
密度相连：存在样本集合D中的一点o，如果对象o到对象p和对象q都是密度可达的，那么p和q密度相联。

　　

 

　　可以发现，密度可达是直接密度可达的传递闭包，并且这种关系是非对称的，因为q不一定是核心对象，如果q为边界点，则q到p一定不密度可达。只有核心对象之间相互密度可达，然而密度相连是对称关系。

三、算法步骤

　　DBSCAN目的是找到密度相连对象的最大集合。我的理解是每次聚类都是寻找到某一个核心点，然后基于该核心点的领域半径内的点集A出发，寻找下一个核心点，其领域半径内的点集B合并到A里，再次从扩大的集合A里寻找新的核心点，直到这个集合A不再扩大为止，此时集合A即为密度相连点的最大集合，这个集合里的各点之间的密度相连都是对于第一次寻找到的核心点来说的。 

 

　　 当设置MinPts=4的时候，红点为高密度点，蓝点为异常点，黄点为边界点。红黄点串成一起成了一个簇。 

　　“其核心思想就是先发现密度较高的点，然后把相近的高密度点逐步都连成一片，进而生成各种簇。算法实现上就是，对每个数据点为圆心，以eps为半径画个圈（称为邻域），然后数有多少个点在这个圈内，这个数就是该点密度值。然后我们可以选取一个密度阈值MinPts，如圈内点数小于MinPts的圆心点为低密度的点，而大于或等于MinPts的圆心点高密度的点（称为核心点Core point）。如果有一个高密度的点在另一个高密度的点的圈内，我们就把这两个点连接起来，这样我们可以把好多点不断地串联出来。之后，如果有低密度的点也在高密度的点的圈内，把它也连到最近的高密度点上，称之为边界点。这样所有能连到一起的点就成一了个簇，而不在任何高密度点的圈内的低密度点就是异常点。”

四、代码实现

1、输入数据

2 2
3 1
3 4
3 14
5 3
8 3
8 6
9 8
10 4
10 7
10 10
10 14
11 13
12 8
12 15
14 7
14 9
14 15
15 8

2、代码

import scala.collection.mutable.ArrayBuffer
import scala.io.Source
import scala.util.control.Breaks._

object myDBSCAN {
  def main(args: Array[String]): Unit = {
    val minPts = 3 //密度阈值
    val ePs = 3 //领域半径
    val dim = 2 //数据集维度
    //处理输入数据
    val fileName = "input.txt"
    val lines = Source.fromFile(fileName).getLines()
    val points = lines.map(line => {
      //数据预处理
      val parts = line.split(" ").map(_.toDouble)
      var vector = Vector[Double]()
      for (i <- 0 to dim - 1)
        vector ++= Vector(parts(i))
      vector
    }).toArray
    println("数据集点的x坐标如下：")
    points.foreach(v => print(v(0) + ","))
    println()
    println("数据集点的y坐标如下：")
    points.foreach(v => print(v(1) + ","))
    println()
    val (cluster, types) = runDBSCAN(points, ePs, minPts)
    printResult(points, cluster, types)
  }

  def runDBSCAN(data: Array[Vector[Double]], ePs: Double, minPts: Int): (Array[Int], Array[Int]) = {
    val types = (for (i <- 0 to data.length - 1) yield -1).toArray //用于区分核心点1，边界点0，和噪音点-1(即cluster中值为0的点)
    val visited = (for (i <- 0 to data.length - 1) yield 0).toArray //用于判断该点是否处理过，0表示未处理过
    var number = 1 //用于标记类
    var xTempPoint = Vector(0.0, 0.0)
    var yTempPoint = Vector(0.0, 0.0)
    var distance = new Array[(Double, Int)](1)
    var distanceTemp = new Array[(Double, Int)](1)
    val neighPoints = new ArrayBuffer[Vector[Double]]()
    var neighPointsTemp = new Array[Vector[Double]](1)
    val cluster = new Array[Int](data.length) //用于标记每个数据点所属的类别
    var index = 0
    for (i <- 0 to data.length - 1) {
      //对每一个点进行处理
      if (visited(i) == 0) { //表示该点未被处理
        visited(i) = 1 //标记为处理过
        xTempPoint = data(i) //取到该点
        distance = data.map(x => (vectorDis(x, xTempPoint), data.indexOf(x))) //取得该点到其他所有点的距离Array{(distance,index)}
        neighPoints ++= distance.filter(x => x._1 <= ePs).map(v => data(v._2)) //找到半径ePs内的所有点(密度相连点集合)
        if (neighPoints.length > 1 && neighPoints.length < minPts) {
          breakable {
            for (j <- 0 to neighPoints.length - 1) {
              //此为非核心点，若其领域内有核心点，则该点为边界点
              var index = data.indexOf(neighPoints(j))
              if (types(index) == 1) {
                types(i) = 0 //边界点
                break
              }
            }
          }
        }
        else if (neighPoints.length >= minPts) {
          //核心点,此时neighPoints表示以该核心点出发的密度相连点的集合
          types(i) = 1
          cluster(i) = number
          while (neighPoints.isEmpty == false) { //对该核心点领域内的点迭代寻找核心点，直到所有核心点领域半径内的点组成的集合不再扩大（每次聚类 ）
            yTempPoint = neighPoints.head //取集合中第一个点
            index = data.indexOf(yTempPoint)
            if (cluster(index) == 0)
              cluster(index) = number //划分到与核心点一样的簇中
            if (visited(index) == 0) {
              //若该点未被处理，则标记已处理
              visited(index) = 1
              distanceTemp = data.map(x => (vectorDis(x, yTempPoint), data.indexOf(x))) //取得该点到其他所有点的距离Array{(distance,index)}
              neighPointsTemp = distanceTemp.filter(x => x._1 <= ePs).map(v => data(v._2)) //找到半径ePs内的所有点
              if (neighPointsTemp.length >= minPts) {
                types(index) = 1 //该点为核心点
                for (i <- 0 to neighPointsTemp.length - 1) {
                  //将其领域内未分类的对象划分到簇中,然后放入neighPoints
                  if (cluster(data.indexOf(neighPointsTemp(i))) == 0) {
                    cluster(data.indexOf(neighPointsTemp(i))) = number //只划分簇，没有访问到
                    if (!neighPoints.contains(neighPointsTemp(i)))
                      neighPoints += neighPointsTemp(i)
                  }
                }
              }
              else if (neighPointsTemp.length > 1 && neighPointsTemp.length < minPts) {
                breakable {
                  for (j <- 0 to neighPointsTemp.length - 1) {
                    //此为非核心点，若其领域内有核心点，则该点为边界点
                    var index1 = data.indexOf(neighPointsTemp(j))
                    if (types(index1) == 1) {
                      types(index) = 0 //边界点
                      break
                    }
                  }
                }
              }
            }
            neighPoints -= yTempPoint //将该点剔除
          } //end-while
          number += 1 //进行新的聚类
        }
        //清空neighPoints集合
        neighPoints.clear()
      }
    }
    (cluster, types)
  }

  def printResult(data: Array[Vector[Double]], cluster: Array[Int], types: Array[Int]) = {
    val result = data.map(v => (cluster(data.indexOf(v)), v)).groupBy(v => v._1) //Map[int,Array[(int,Vector[Double])]]
    //key代表簇号，value代表属于这一簇的元素数组
    result.foreach(v => {
      println("簇" + v._1 + "包含的元素如下:")
      v._2.foreach(v => println(v._2))
    })
    //val noise = cluster.zip(data).filter(v => v._1 == 0)
    //noise.foreach(v => types(data.indexOf(v._2)) = -1) //通过簇号0把噪音点在types中赋值-1,数据集中没有包含在任何簇中(即簇号为0)的数据点就构成异常点
    //val pointsTypes = data.map(v => (types(data.indexOf(v)), v)).groupBy(v => v._1) //Map[点类型int,Array[(点类型int,Vector[Double])]]
    //key代表点的类型号，value代表属于这一类型的元素数组
    //    pointsTypes.foreach(v => {
    //      if (v._1 == 1) println("核心点如下：")
    //      else if (v._1 == 0) println("边界点如下：")
    //      else println("噪音点如下：")
    //      v._2.foreach(v => println(v._2))
    //    })
    for (i <- 0 to data.length - 1) {
      if (cluster(i) == 0 && types(i) == (-1)) {
        println(data(i) + "是噪声点")
      } else if (cluster(i) > 0 && types(i) == 1) {
        println(data(i) + "是核心点")
      } else {
        println(data(i) + "是边界点")
      }
    }
  }

  //--------------------------自定义向量间的运算-----------------------------
  //--------------------------向量间的欧式距离-----------------------------
  def vectorDis(v1: Vector[Double], v2: Vector[Double]): Double = {
    var distance = 0.0
    for (i <- 0 to v1.length - 1) {
      distance += (v1(i) - v2(i)) * (v1(i) - v2(i))
    }
    distance = math.sqrt(distance)
    distance
  }
}

3、输出

数据集点的x坐标如下：
2.0,3.0,3.0,3.0,5.0,8.0,8.0,9.0,10.0,10.0,10.0,10.0,11.0,12.0,12.0,14.0,14.0,14.0,15.0,
数据集点的y坐标如下：
2.0,1.0,4.0,14.0,3.0,3.0,6.0,8.0,4.0,7.0,10.0,14.0,13.0,8.0,15.0,7.0,9.0,15.0,8.0,
簇2包含的元素如下:
Vector(10.0, 14.0)
Vector(11.0, 13.0)
Vector(12.0, 15.0)
Vector(14.0, 15.0)
簇1包含的元素如下:
Vector(2.0, 2.0)
Vector(3.0, 1.0)
Vector(3.0, 4.0)
Vector(5.0, 3.0)
Vector(8.0, 3.0)
Vector(8.0, 6.0)
Vector(9.0, 8.0)
Vector(10.0, 4.0)
Vector(10.0, 7.0)
Vector(10.0, 10.0)
Vector(12.0, 8.0)
Vector(14.0, 7.0)
Vector(14.0, 9.0)
Vector(15.0, 8.0)
簇0包含的元素如下:
Vector(3.0, 14.0)
噪音点如下：
Vector(3.0, 14.0)
核心点如下：
Vector(2.0, 2.0)
Vector(3.0, 1.0)
Vector(3.0, 4.0)
Vector(5.0, 3.0)
Vector(8.0, 3.0)
Vector(8.0, 6.0)
Vector(9.0, 8.0)
Vector(10.0, 4.0)
Vector(10.0, 7.0)
Vector(10.0, 10.0)
Vector(10.0, 14.0)
Vector(11.0, 13.0)
Vector(12.0, 8.0)
Vector(12.0, 15.0)
Vector(14.0, 7.0)
Vector(14.0, 9.0)
Vector(15.0, 8.0)
边界点如下：
Vector(14.0, 15.0)