Подтвердить что ты не робот

Как работает уплотненная матрица расстояний? (Pdist)

scipy.spatial.distance.pdist возвращает сжатую матрицу расстояний. Из документация:

Возвращает матрицу сокращенных расстояний Y. Для каждого и (где) метрика dist (u = X [i], v = X [j]) вычисляется и сохраняется в записи ij.

Я думал, что ij означает i*j. Но я думаю, что я ошибаюсь. Рассмотрим

X = array([[1,2], [1,2], [3,4]])
dist_matrix = pdist(X)

то в документации указано, что dist(X[0], X[2]) должен быть dist_matrix[0*2]. Однако dist_matrix[0*2] равно 0, а не 2.8, как и должно быть.

Какую формулу я должен использовать для доступа к подобию двух векторов, заданных i и j?

4b9b3361

Ответ 1

Вы можете посмотреть на это следующим образом: предположим, что x есть m по n. Возможные пары строк m, выбранных по два за раз, составляют itertools.combinations(range(m), 2), например, для m=3:

>>> import itertools
>>> list(combinations(range(3),2))
[(0, 1), (0, 2), (1, 2)]

Итак, если d = pdist(x), k th кортеж в combinations(range(m), 2)) дает индексы строк x, связанных с d[k].

Пример:

>>> x = array([[0,10],[10,10],[20,20]])
>>> pdist(x)
array([ 10.        ,  22.36067977,  14.14213562])

Первый элемент dist(x[0], x[1]), второй - dist(x[0], x[2]), а третий - dist(x[1], x[2]).

Или вы можете просмотреть его как элементы в верхней треугольной части квадратной матрицы расстояния, нанизанные в 1D-массив.

например.

>>> squareform(pdist(x)) 
array([[  0.   ,  10.   ,  22.361],
       [ 10.   ,   0.   ,  14.142],
       [ 22.361,  14.142,   0.   ]])

>>> y = array([[0,10],[10,10],[20,20],[10,0]])
>>> squareform(pdist(y)) 
array([[  0.   ,  10.   ,  22.361,  14.142],
       [ 10.   ,   0.   ,  14.142,  10.   ],
       [ 22.361,  14.142,   0.   ,  22.361],
       [ 14.142,  10.   ,  22.361,   0.   ]])
>>> pdist(y)
array([ 10.   ,  22.361,  14.142,  14.142,  10.   ,  22.361])

Ответ 2

Вектор сжатой матрицы соответствует нижней треугольной области квадратной матрицы. Чтобы преобразовать точку в этой треугольной области, вам нужно вычислить количество точек слева в треугольнике и число, указанное в столбце.

Вы можете использовать следующую функцию для преобразования:

q = lambda i,j,n: n*j - j*(j+1)/2 + i - 1 - j

Check:

import numpy as np
from scipy.spatial.distance import pdist, squareform
x = np.random.uniform( size = 100 ).reshape( ( 50, 2 ) )
d = pdist( x )
ds = squareform( d )
for i in xrange( 1, 50 ):
    for j in xrange( i ):
        assert ds[ i, j ] == d[ q( i, j, 50 ) ]

Ответ 3

Матрица с конденсированным расстоянием до матрицы полного расстояния

Сжатая матрица расстояний, возвращаемая pdist, может быть преобразована в матрицу полного расстояния с помощью scipy.spatial.distance.squareform:

>>> import numpy as np
>>> from scipy.spatial.distance import pdist, squareform
>>> points = np.array([[0,1],[1,1],[3,5], [15, 5]])
>>> dist_condensed = pdist(points)
>>> dist_condensed
array([  1.        ,   5.        ,  15.5241747 ,   4.47213595,
        14.56021978,  12.        ])

Используйте squareform для преобразования в полную матрицу:

>>> dist = squareform(dist_condensed)
array([[  0.        ,   1.        ,   5.        ,  15.5241747 ],
       [  1.        ,   0.        ,   4.47213595,  14.56021978],
       [  5.        ,   4.47213595,   0.        ,  12.        ],
       [ 15.5241747 ,  14.56021978,  12.        ,   0.        ]])

Расстояние между точкой i, j сохраняется в dist [i, j]:

>>> dist[2, 0]
5.0
>>> np.linalg.norm(points[2] - points[0])
5.0

Индексы к конденсированному индексу

Можно преобразовать индексы, используемые для доступа к элементам квадратной матрицы, к индексу в конденсированной матрице:

def square_to_condensed(i, j, n):
    assert i != j, "no diagonal elements in condensed matrix"
    if i < j:
        i, j = j, i
    return n*j - j*(j+1)/2 + i - 1 - j

Пример:

>>> square_to_condensed(1, 2, len(points))
3
>>> dist_condensed[3]
4.4721359549995796
>>> dist[1,2]
4.4721359549995796

Сжатый индекс индексов

Также возможно другое направление без sqaureform, которое лучше с точки зрения времени выполнения и потребления памяти:

import math

def calc_row_idx(k, n):
    return int(math.ceil((1/2.) * (- (-8*k + 4 *n**2 -4*n - 7)**0.5 + 2*n -1) - 1))

def elem_in_i_rows(i, n):
    return i * (n - 1 - i) + (i*(i + 1))/2

def calc_col_idx(k, i, n):
    return int(n - elem_in_i_rows(i + 1, n) + k)

def condensed_to_square(k, n):
    i = calc_row_idx(k, n)
    j = calc_col_idx(k, i, n)
    return i, j

Пример:

>>> condensed_to_square(3, 4)
(1.0, 2.0)

Сравнение времени выполнения с квадратом

>>> import numpy as np
>>> from scipy.spatial.distance import pdist, squareform
>>> points = np.random.random((10**4,3))
>>> %timeit dist_condensed = pdist(points)
1 loops, best of 3: 555 ms per loop

Создание sqaureform оказывается очень медленным:

>>> dist_condensed = pdist(points)
>>> %timeit dist = squareform(dist_condensed)
1 loops, best of 3: 2.25 s per loop

Если мы ищем две точки с максимальным расстоянием, неудивительно, что поиск в полной матрице равен O (n), а в конденсированной форме только O (n/2):

>>> dist = squareform(dist_condensed)
>>> %timeit dist_condensed.argmax()
10 loops, best of 3: 45.2 ms per loop
>>> %timeit dist.argmax()
10 loops, best of 3: 93.3 ms per loop

Получение inideces для двух точек практически не требует времени в обоих случаях, но, конечно, для вычисления сконденсированного индекса есть некоторые накладные расходы:

>>> idx_flat = dist.argmax()
>>> idx_condensed = dist.argmax()
>>> %timeit list(np.unravel_index(idx_flat, dist.shape))
100000 loops, best of 3: 2.28 µs per loop
>>> %timeit condensed_to_square(idx_condensed, len(points))
100000 loops, best of 3: 14.7 µs per loop

Ответ 4

Если вы хотите получить доступ к элементу pdist, соответствующему (i, j) -тому элементу матрицы квадратного расстояния, математика выглядит следующим образом: Предположим i < j (иначе индексы переворота), если i == j, ответ равен 0.

X = random((N,m))
dist_matrix = pdist(X)

Тогда (i, j) -й элемент является dist_matrix [ind], где

ind = (N - array(range(1,i+1))).sum()  + (j - 1 - i).

Ответ 5

Это версия верхнего треугольника (я < j), которая должна быть интересна некоторым:

condensed_idx = lambda i,j,n: i*n + j - i*(i+1)/2 - i - 1

Это очень легко понять:

  • с i*n + j вы переходите к позиции в квадратной матрице;
  • с - i*(i+1)/2 вы удаляете нижний треугольник (включая диагональ) во всех строках до i;
  • с - i вы удаляете позиции в строке я перед диагональю;
  • с - 1 вы удаляете позиции в строке я по диагонали.

Check:

import scipy
from scipy.spatial.distance import pdist, squareform
condensed_idx = lambda i,j,n: i*n + j - i*(i+1)/2 - i - 1
n = 50
dim = 2
x = scipy.random.uniform(size = n*dim).reshape((n, dim))
d = pdist(x)
ds = squareform(d)
for i in xrange(1, n-1):
    for j in xrange(i+1, n):
        assert ds[i, j] == d[condensed_idx(i, j, n)]