点云聚类学习 KMeans/DBSCAN

Overview

最近做的东西会处理一些Lidar的点云数据，虽然之前在看Autoware的时候有了解一些聚类的基本原理和实现，但还是稍微再学习一下聚类方法吧，这里就简单记录一下（大部分都是参考wiki），然后做个简单实现看一下两种方法具体有什么区别

KMeans

参考wiki，k-means clustering源于信号处理中的一种向量量化方法，现在则更多地作为一种聚类分析方法流行于数据挖掘领域。 其实简单来说在Lidar点云场景下，k-means就是把n个点云分成了k个聚类里去。其总体思路如下：

在具体的算法实践中，常见的k-means算法其实也是一种迭代的算法，还是参考wiki，其具体的算法流程为：

这个算法流程还是挺好理解的，其实就是根据距离把点分给最近的聚类中心，然后再更新聚类中心再分配再更新，最后达到收敛就是完成了。 那这里也就引出了另一个问题，其实我们是没有最初的初始聚类中心的，所以初始化步骤也是k-means聚类中的重要一步。通常使用的初始化方法有Forgy和随机划分(Random Partition)方法 。Forgy方法随机地从数据集中选择k个观测作为初始的均值点；而随机划分方法则随机地为每一观测指定聚类，然后运行“更新(Update)”步骤,即计算随机分配的各聚类的图心，作为初始的均值点。Forgy方法易于使得初始均值点散开，随机划分方法则把均值点都放到靠近数据集中心的地方。
最后来看一下简单实现，我这里实现的其实是一维的，也就是根据其中一个方向坐标进行的聚类，如果是3D的话其实就是初始化以及迭代的时候距离的更新需要改一下，还是很简单的，这里用到的就是随机划分的初始化方法

    // Initialize centroidsdouble min_val = *std::min_element(data.begin(), data.end());double max_val = *std::max_element(data.begin(), data.end());for (int i = 0; i < k; ++i){centroids[i] = min_val + (max_val - min_val) * i / (k - 1);}

然后直接迭代就行了

    for (int iter = 0; iter < max_iterations; ++iter){std::vector<std::vector<double>> clusters(k);// Assign points to the nearest centroidfor (const auto& point : data){int closest_centroid = 0;double min_distance = std::abs(point - centroids[0]);for (int j = 1; j < k; ++j){double distance = std::abs(point - centroids[j]);if (distance < min_distance){min_distance = distance;closest_centroid = j;}}clusters[closest_centroid].push_back(point);}// Update centroidsbool centroids_changed = false;for (int i = 0; i < k; ++i){if (!clusters[i].empty()){double new_centroid = std::accumulate(clusters[i].begin(), clusters[i].end(), 0.0) / clusters[i].size();if (std::abs(new_centroid - centroids[i]) > 1e-6){centroids[i] = new_centroid;centroids_changed = true;}}}if (!centroids_changed){break;}}

最后来看一下效果，可以看到整体效果还是不错的，把点云进一步地细分了，但是因为我这里用的是1D的聚类，所以在一些距离比较远当坐标有交集的时候不可避免地产生了误聚类的现象

DBSCAN

DBSCAN也是一种聚类方法，还是参考wiki，DBSCAN算法是以密度为本的：给定某空间里的一个点集合，这算法能把附近的点分成一组（有很多相邻点的点），并标记出位于低密度区域的局外点（最接近它的点也十分远）。 这个方法相比k-means算法复杂多了，首先所有的点被分成了三类：核心点可达点和局外点

如果 p 是核心点，则它与所有由它可达的点(包括核心点和非核心点)形成一个聚类，每个聚类拥有最少一个核心点，非核心点也可以是聚类的一部分，但它是在聚类的“边缘”位置，因为它不能达至更多的点。具体可以参考下面这个图

在这幅图里，minPts = 4，点 A 和其他红色点是核心点，因为它们的 ε-邻域（图中红色圆圈）里包含最少 4 个点（包括自己），由于它们之间相互相可达，它们形成了一个聚类。点 B 和点 C 不是核心点，但它们可由 A 经其他核心点可达，所以也属于同一个聚类。点 N 是局外点，它既不是核心点，又不由其他点可达。

所以DBSCAN算法的话，核心有两个参数，一个是选的范围的大小，另一个就是minPts这个值，它由一个任意未被访问的点开始，然后探索这个点的 ε-邻域，如果 ε-邻域里有足够的点，则建立一个新的聚类，否则这个点被标签为杂音。注意这个点之后可能被发现在其它点的 ε-邻域里，而该 ε-邻域可能有足够的点，届时这个点会被加入该聚类中。
如果一个点位于一个聚类的密集区域里，它的 ε-邻域里的点也属于该聚类，当这些新的点被加进聚类后，如果它(们)也在密集区域里，它(们)的 ε-邻域里的点也会被加进聚类里。这个过程将一直重复，直至不能再加进更多的点为止，这样，一个密度连结的聚类被完整地找出来。然后，一个未曾被访问的点将被探索，从而发现一个新的聚类或杂音。其伪代码流程如下：

当然DBSCAN由于相对复杂一些，所以我这里就没有自己去实现了，我是直接利用了PCL库来实现DBSCAN功能

std::vector<std::vector<Eigen::Vector3d>> LidarPreProcess::performDBSCANClustering(const std::vector<Eigen::Vector3d>& points, double clusterTolerance, int minClusterSize, int maxClusterSize)
{pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZ>);cloud->points.reserve(points.size());for (const auto& point : points){cloud->points.emplace_back(point.x(), point.y(), point.z());}cloud->width = cloud->points.size();cloud->height = 1;cloud->is_dense = true;pcl::search::KdTree<pcl::PointXYZ>::Ptr tree(new pcl::search::KdTree<pcl::PointXYZ>);tree->setInputCloud(cloud);std::vector<pcl::PointIndices> cluster_indices;pcl::EuclideanClusterExtraction<pcl::PointXYZ> ec;ec.setClusterTolerance(clusterTolerance);ec.setMinClusterSize(minClusterSize);ec.setMaxClusterSize(maxClusterSize);ec.setSearchMethod(tree);ec.setInputCloud(cloud);ec.extract(cluster_indices);std::vector<std::vector<Eigen::Vector3d>> clusters;for (const auto& indices : cluster_indices){std::vector<Eigen::Vector3d> cluster;for (int index : indices.indices){cluster.push_back(points[index]);}clusters.push_back(cluster);}return clusters;
}

这里就不放这个效果图了，大概的对比在下一节用文字描述一下

简单对比

其实把两个流程大概梳理一下，能发现第一个不同就是：k-means在使用的时候是需要指定聚类数的，而DBSCAN是没有这个需求的， 而且理论上，DBSCAN对那种非线性的点云处理效果应该是更好一些的，不过我这里因为点云在其中一个轴有更好地一种分布特性，用了DBSCAN之后反而会把一整根线分成多段(大概率也是我参数设置的不好)，在我的应用场景中，效果应该是单轴k-means > DBSCAN >>> k-means，所以最后实际的应用中，去结合点云的一些实际分布特性适当调整方法就可以了。

本文仅作个人简单学习记录，如有错误之处还烦请指正。