Быстрый (<n ^ 2) алгоритм кластеризации
У меня есть 1 миллион 5-мерных точек, которые мне нужно сгруппировать в k кластеров с k < 1 миллион. В каждом кластере две точки не должны быть слишком далеко друг от друга (например, они могут быть ограничивающими сферами с заданным радиусом). Это означает, что, вероятно, должно быть много кластеров размером 1.
Но! Мне нужно, чтобы время работы было значительно ниже n ^ 2. n log n или так должно быть хорошо. Причина, по которой я делаю эту кластеризацию, заключается в том, чтобы избежать вычисления матрицы расстояния всех n точек (которая занимает n ^ 2 раза или много часов), вместо этого я хочу просто вычислить расстояния между кластерами.
Я попробовал алгоритм k-means pycluster, но быстро понял, что он слишком медленный. Я также пробовал следующий жадный подход:
-
Разделите пространство на 20 частей в каждом измерении. (так что всего 20 ^ 5 штук). Я буду хранить кластеры в этих сетчатых коробках, согласно их центроидам.
-
Для каждой точки извлекайте сетчатые ячейки, находящиеся в пределах r (радиус максимальной ограничивающей сферы). Если есть достаточно кластер, добавьте его в этот кластер, иначе создайте новый кластер.
Однако, похоже, это дает мне больше кластеров, чем я хочу. Я также реализовал аналогичные ему методы дважды, и они дают очень разные ответы.
Существуют ли какие-либо стандартные подходы к кластеризации быстрее, чем n ^ 2 раза? Вероятностные алгоритмы в порядке.
Ответы
Ответ 1
Рассмотрим алгоритм приближенного ближайшего соседа (ANN) или чувствительного к местоположению хэширования (LSH). Они не решают проблему кластеризации напрямую, но они смогут сказать вам, какие точки "близки" друг к другу. Изменяя параметры, вы можете определить, что они находятся близко друг к другу. И это быстро.
Точнее, LSH может обеспечить хэш-функцию, h
, такую, что для двух точек x
и y
и метрики расстояния d
,
d(x,y) <= R1 => P(h(x) = h(y)) >= P1
d(x,y) >= R2 => P(h(x) = h(y)) <= P2
где R1 < R2
и P1 > P2
. Так что да, это вероятностно. Вы можете обработать полученные данные, чтобы получить истинные кластеры.
Вот информация о LSH, в том числе руководство E2LSH. ANN похож по духу; Дэвид Маунт имеет информацию здесь или попробуйте FLANN (имеет привязки Matlab и Python).
Ответ 2
Вам может понравиться попробовать мой исследовательский проект K-tree. Он хорошо масштабируется с большими входами относительно k-средних и формирует иерархию кластеров. Компромисс заключается в том, что он создает кластеры с более высоким искажением. Он имеет среднее время выполнения O (n log n) и наихудший случай O (n ** 2), что происходит, только если у вас есть какая-то странная топология. Более подробная информация о анализе сложности приведена в моей тесте Мастера. Я использовал его с очень объемными текстовыми данными и не имел проблем.
Иногда в дереве могут происходить плохие расщепления, где все данные идут в одну сторону (кластер). Строка в SVN имеет дело с этим иначе, чем текущая версия. Он случайным образом разбивает данные, если есть плохие расщепления. Предыдущий метод может заставить дерево стать слишком глубоким, если есть плохие расщепления.
Ответ 3
Поместите данные в индекс, например R * -tree (Wikipedia), тогда вы можете запустить множество алгоритмов кластеризации на основе плотности ( такие как DBSCAN (Wikipedia) или ОПТИКА (Википедия)) в O(n log n)
.
Кластеризация на основе плотности (Wikipedia), кажется, именно то, что вы хотите ( "не слишком далеко друг от друга" )
Ответ 4
У людей создается впечатление, что k-средства медленны, но медлительность действительно является только проблемой для алгоритма EM (Lloyd's). Методы стохастического градиента для k-средних на порядки выше EM (см. Www.eecs.tufts.edu/~dsculley/papers/fastkmeans.pdf).
Реализация здесь: http://code.google.com/p/sofia-ml/wiki/SofiaKMeans
Ответ 5
У меня есть модуль Perl, который делает именно то, что вы хотите Algorithm:: ClusterPoints.
Во-первых, он использует алгоритм, который вы описали в своем сообщении, чтобы разделить точки в многомерных секторах, а затем использует грубую силу для поиска кластеров между точками в смежных секторах.
Сложность варьируется от O (N) до O (N ** 2) в очень деградированных случаях.
Обновление
@Denis: нет, это намного хуже:
Для размеров d
размер сектора (или небольшого гиперкуба) s
определяется так, что его диагональ l
- это минимальное расстояние c
, разрешенное между двумя точками в разных кластерах.
l = c
l = sqrt(d * s * s)
s = sqrt(c * c / d) = c / sqrt(d)
Затем вы должны рассмотреть все сектора, которые касаются гиперсферы диаметром r = 2c + l
, центрированным в сводном секторе.
Грубо, мы должны рассматривать ceil(r/s)
строки секторов во всех направлениях, а это означает n = pow(2 * ceil(r/s) + 1, d)
.
Например, для d=5
и c=1
получаем l=2.236
, s=0.447
, r=3.236
и n=pow(9, 5)=59049
На самом деле мы должны проверять меньшие соседние сектора, так как здесь мы рассматриваем те, которые касаются гиперкуба размера (2r+1)/s
, и нам нужно только проверить тех, кто касается описанной гиперсферы.
Учитывая биективный характер отношения "находятся в одном кластере", мы можем также уменьшить число секторов, которые необходимо проверить.
В частности, Algorithm:: ClusterPoints для случая, когда d=5
проверяет 3903 сектора.
Ответ 6
Ниже представлен небольшой тестовый стенд, чтобы узнать, как быстро
scipy.spatial.cKDTree
находится на ваших данных,
и получить приблизительное представление о том, как рассеиваются расстояния между соседними точками.
Хороший способ запуска K-кластера для различных K
состоит в том, чтобы построить MST ближайших пар и удалить самый длинный K-1; видеть
Уэйн, Жадные алгоритмы.
Визуализация кластеров будет забавной - проецировать на 2d с помощью PCA?
(Любопытно, это ваш K 10, 100, 1000?)
Добавлено 17 декабря: реальное время автономной работы: 100000 x 5 10 секунд, 500000 x 5 60 секунд
#!/usr/bin/env python
# time scipy.spatial.cKDTree build, query
from __future__ import division
import random
import sys
import time
import numpy as np
from scipy.spatial import cKDTree as KDTree
# http://docs.scipy.org/doc/scipy/reference/spatial.html
# $scipy/spatial/kdtree.py is slow but clean, 0.9 has cython
__date__ = "2010-12-17 dec denis"
def clumpiness( X, nbin=10 ):
""" how clumpy is X ? histogramdd av, max """
# effect on kdtree time ? not much
N, dim = X.shape
histo = np.histogramdd( X, nbin )[0] .astype(int) # 10^dim
n0 = histo.size - histo.astype(bool).sum() # uniform: 1/e^lambda
print "clumpiness: %d of %d^%d data bins are empty av %.2g max %d" % (
n0, nbin, dim, histo.mean(), histo.max())
#...............................................................................
N = 100000
nask = 0 # 0: ask all N
dim = 5
rnormal = .9
# KDtree params --
nnear = 2 # k=nnear+1, self
leafsize = 10
eps = 1 # approximate nearest, dist <= (1 + eps) * true nearest
seed = 1
exec "\n".join( sys.argv[1:] ) # run this.py N= ...
np.random.seed(seed)
np.set_printoptions( 2, threshold=200, suppress=True ) # .2f
nask = nask or N
print "\nkdtree: dim=%d N=%d nask=%d nnear=%d rnormal=%.2g leafsize=%d eps=%.2g" % (
dim, N, nask, nnear, rnormal, leafsize, eps)
if rnormal > 0: # normal point cloud, .9 => many near 1 1 1 axis
cov = rnormal * np.ones((dim,dim)) + (1 - rnormal) * np.eye(dim)
data = np.abs( np.random.multivariate_normal( np.zeros(dim), cov, N )) % 1
# % 1: wrap to unit cube
else:
data = np.random.uniform( size=(N,dim) )
clumpiness(data)
ask = data if nask == N else random.sample( data, sample )
t = time.time()
#...............................................................................
datatree = KDTree( data, leafsize=leafsize ) # build the tree
print "%.1f sec to build KDtree of %d points" % (time.time() - t, N)
t = time.time()
distances, ix = datatree.query( ask, k=nnear+1, eps=eps )
print "%.1f sec to query %d points" % (time.time() - t, nask)
distances = distances[:,1:] # [:,0] is all 0, point to itself
avdist = distances.mean( axis=0 )
maxdist = distances.max( axis=0 )
print "distances to %d nearest: av" % nnear, avdist, "max", maxdist
# kdtree: dim=5 N=100000 nask=100000 nnear=2 rnormal=0.9 leafsize=10 eps=1
# clumpiness: 42847 of 10^5 data bins are empty av 1 max 21
# 0.4 sec to build KDtree of 100000 points
# 10.1 sec to query 100000 points
# distances to 2 nearest: av [ 0.07 0.08] max [ 0.15 0.18]
# kdtree: dim=5 N=500000 nask=500000 nnear=2 rnormal=0.9 leafsize=10 eps=1
# clumpiness: 2562 of 10^5 data bins are empty av 5 max 80
# 2.5 sec to build KDtree of 500000 points
# 60.1 sec to query 500000 points
# distances to 2 nearest: av [ 0.05 0.06] max [ 0.13 0.13]
# run: 17 Dec 2010 15:23 mac 10.4.11 ppc