了凡春秋USTC
谱聚类
http://chunqiu.blog.ustc.edu.cn/?p=505
最近忙着写文章,好久不写博客了。最近看到一个聚类方法--谱聚类,号称现代聚类方法,看到它简洁的公式推导、实现代码,不禁要尝试一把。关于它的理论,google一搜有很多博客讲,这里就不赘述了,反正最后还是归结为一个SVD分解问题,参考网址如下
看到下面的图,乐了
上网搜索,看到了几个相关的工具箱,下面一一介绍:
1.这个是加利福尼亚大学一个实验室为他们的论文《Parallel Spectral Clustering in Disributed Systems》而写的测试代码,这篇文章发在IEEE Transactions on Pattern Analysis and Machine Intelligence上,机器学习方面的*期刊。不过看了看作者,怎么感觉都是中国人啊,下载地址见http://alumni.cs.ucsb.edu/~wychen/sc.html#FAQ
压缩包内有Readme,有工具箱的使用介绍,不过由于是为文章效果而写的,所以感觉并不能称得上工具箱吧,我自己编了个测试例子试了试,效果、速度确实不错:
原始待聚类数据
聚类后
代码如下:
% 谱聚类算法测试
clc
clear
close all
% 构造原始数据
ang = 0:pi/200:pi/2*3;
x1 = sin(ang) + randn(size(ang))*0.1;
y1 = cos(ang) + randn(size(ang))*0.1;
x2 = 2.3*sin(ang+pi) + randn(size(ang))*0.1;
y2 = 2.3*cos(ang) + randn(size(ang))*0.1;
x3 = 6*sin(ang+pi/2*3+pi/6) + randn(size(ang))*0.1;
y3 = 5*cos(ang+pi/2*3) + randn(size(ang))*0.1;
X = [x1, x2, x3; y1, y2, y3]';
figure
plot(X(:, 1), X(:, 2), '.')
axis equal
% 谱聚类
num_neighbors = 20;
block_size = 5;
gen_nn_distance(X, num_neighbors, block_size, 0); % 计算相似矩阵,保存到NN_sym_distance.mat文件中
filename = [num2str(num_neighbors), '_NN_sym_distance'];
load(filename)
num_clusters = 3;
sigma = 10;
cluster_labels = sc(A, sigma, num_clusters); % 谱聚类
figure
plot(X(cluster_labels == 1, 1), X(cluster_labels == 1, 2), 'r.', ...
X(cluster_labels == 2, 1), X(cluster_labels == 2, 2), 'b.', ...
X(cluster_labels == 3, 1), X(cluster_labels == 3, 2), 'g.')
axis equal
2.有个matlab exchange上的工具箱,是德国斯图加特大学一个博士整齐其博士论文的研究而写的GUI,下载链接如下(压缩包内有其博士论文)
同样,Readme里有详细使用说明,只要运行MAIN.m文件即可弹出GUI。运行可能会出个警告
影响不大,Readme里是这么说的:
recommended: export_fig configured to save to pdf files (Mathworks File Exchange ID #23629)
如果你不用保存图片,不用理会。
GUI如下
对其给的twomoons的例子做测试,结果如下
额,是不是有点自己打脸的赶脚!自己给的例子,效果也不够完美。。。估计调整一些参数能得到理想的效果吧!
Matlab exchange上还有个谱聚类的文件,提供了谱聚类各种算法的简单实现,对于学习谱聚类具体实现挺有帮助!
3.很多介绍谱聚类的博客都有提到算法的matlab实现,因为实现代码太简单了,如下
% 谱聚类简单实现
function idx = spectral_clustering(W, k)
D = diag(sum(W));
L = D-W;
opt = struct('issym', true, 'isreal', true);
[V, ~] = eigs(L, D, k, 'SM', opt);
idx = kmeans(V, k);
end
谱聚类复杂的理论最终归结为代码就这几句,当然这是考虑最简化的情况,如果考虑细节,就得靠上面提到的工具箱了。对于这个代码,我做了测试,效果还行(使用第一个代码里的计算相似矩阵的函数),不过多次测试的效果说明,不同数据对效果是有影响的,有时就得不到好的分离效果。
% 谱聚类算法测试
clc
clear
close all
% 构造原始数据
ang = 0:pi/200:pi/2*3;
x1 = sin(ang) + randn(size(ang))*0.1;
y1 = cos(ang) + randn(size(ang))*0.1;
x2 = 2.3*sin(ang+pi) + randn(size(ang))*0.1;
y2 = 2.3*cos(ang) + randn(size(ang))*0.1;
x3 = 6*sin(ang+pi/2*3+pi/6) + randn(size(ang))*0.1;
y3 = 5*cos(ang+pi/2*3) + randn(size(ang))*0.1;
X = [x1, x2, x3; y1, y2, y3]';
figure
plot(X(:, 1), X(:, 2), '.')
axis equal
% 谱聚类
num_neighbors = 20;
block_size = 5;
gen_nn_distance(X, num_neighbors, block_size, 0); % 计算相似矩阵,保存到NN_sym_distance.mat文件中
filename = [num2str(num_neighbors), '_NN_sym_distance'];
load(filename)
idx = spectral_clustering(A, 3);
figure
plot(X(idx == 1, 1), X(idx == 1, 2), 'r.', ...
X(idx == 2, 1), X(idx == 2, 2), 'b.', ...
X(idx == 3, 1), X(idx == 3, 2), 'g.')
axis equal
如果对那个博士论文的例子做测试,效果如下
貌似这个图确实很难分,这个效果不如那个德国博士的效果。
% 谱聚类算法测试
clc
clear
close all
% 构造原始数据
M = csvread('twomoon-2d-50s.test.csv');
X = M(:, 2:3);
figure
plot(X(:, 1), X(:, 2), '.')
axis equal
% 谱聚类
num_neighbors = 20;
block_size = 5;
gen_nn_distance(X, num_neighbors, block_size, 0); % 计算相似矩阵,保存到NN_sym_distance.mat文件中
filename = [num2str(num_neighbors), '_NN_sym_distance'];
load(filename)
idx = spectral_clustering(A, 2);
figure
plot(X(idx == 1, 1), X(idx == 1, 2), 'r.', ...
X(idx == 2, 1), X(idx == 2, 2), 'b.')
axis equal
如果用第一个工具箱的方法再测试这个例子
效果也就那样!
% 谱聚类算法测试
clc
clear
close all
% 构造原始数据
M = csvread('twomoon-2d-50s.test.csv');
X = M(:, 2:3);
figure
plot(X(:, 1), X(:, 2), '.')
axis equal
% 谱聚类
num_neighbors = 20;
block_size = 5;
gen_nn_distance(X, num_neighbors, block_size, 0); % 计算相似矩阵,保存到NN_sym_distance.mat文件中
filename = [num2str(num_neighbors), '_NN_sym_distance'];
load(filename)
num_clusters = 2;
sigma = 10;
cluster_labels = sc(A, sigma, num_clusters); % 谱聚类
figure
plot(X(cluster_labels == 1, 1), X(cluster_labels == 1, 2), 'r.', ...
X(cluster_labels == 2, 1), X(cluster_labels == 2, 2), 'b.')
axis equal