【问题标题】:Distance calculation on matrix using numpy使用 numpy 对矩阵进行距离计算
【发布时间】:2012-02-08 06:50:32
【问题描述】:

我正在尝试在 Python 中实现 K-means 算法(我知道有相关库,但我想学习如何自己实现它。)这是我遇到问题的函数:

def AssignPoints(points, centroids):
    """
    Takes two arguments:
    points is a numpy array such that points.shape = m , n where m is number of examples,
    and n is number of dimensions.

    centroids is numpy array such that centroids.shape = k , n where k is number of centroids.
    k < m should hold.

    Returns:
    numpy array A such that A.shape = (m,) and A[i] is index of the centroid which points[i] is assigned to.
    """

    m ,n = points.shape
    temp = []
    for i in xrange(n):
        temp.append(np.subtract.outer(points[:,i],centroids[:,i]))
    distances = np.hypot(*temp)
    return distances.argmin(axis=1)

这个函数的目的,给定 n 维空间中的 m 个点,以及 n 维空间中的 k 个质心,生成一个 (x1 x2 x3 x4 ... xm) 的 numpy 数组,其中 x1 是最接近的质心的索引第一点。这工作正常,直到我用 4 维示例进行了尝试。当我尝试放置 4 维示例时,出现此错误:

  File "/path/to/the/kmeans.py", line 28, in AssignPoints
    distances = np.hypot(*temp)
ValueError: invalid number of arguments

我该如何解决这个问题,或者如果我不能,你建议我如何计算我要在这里计算的内容?

我的答案

def AssignPoints(points, centroids):
    m ,n = points.shape
    temp = []
    for i in xrange(n):
        temp.append(np.subtract.outer(points[:,i],centroids[:,i]))
    for i in xrange(len(temp)):
        temp[i] = temp[i] ** 2
    distances = np.add.reduce(temp) ** 0.5
    return distances.argmin(axis=1)

【问题讨论】:

  • 您意识到len(temp) == nnp.hypot 只承认三个参数,不是吗?
  • 说实话,我对 np.hypot 的作用不是很清楚。当我第一次编写该函数时,我只是从另一个问题中复制了它。我认为无论尺寸如何,这都会起作用。你建议我应该如何编写这个函数?
  • np.hypot 计算三角形的斜边。它接受两个参数:leg1leg2,以及第三个可选参数,用于放置结果。当然,这个想法是一次性计算很多斜边,因此leg1leg2 应该是数组,包含你要计算的n 三角形的相应边。

标签: python numpy k-means


【解决方案1】:

试试这个:

np.sqrt(((points[np.newaxis] - centroids[:,np.newaxis]) ** 2).sum(axis=2)).argmin(axis=0)

或者:

diff = points[np.newaxis] - centroids[:,np.newaxis]
norm = np.sqrt((diff*diff).sum(axis=2))
closest = norm.argmin(axis=0)

不要问它在做什么:D

编辑:不,开个玩笑。中间的广播(points[np.newaxis] - centroids[:,np.newaxis])是从原始数组“制作”两个 3D 数组。结果是每个“平面”都包含所有点和一个质心之间的差异。我们就叫它diffs吧。

然后我们做通常的运算来计算欧几里得距离(差的平方的平方根):np.sqrt((diffs ** 2).sum(axis=2))。我们最终得到一个 (k, m) 矩阵,其中第 0 行包含到 centroids[0] 的距离等。因此,.argmin(axis=0) 为您提供了您想要的结果。

【讨论】:

  • 我认为这是最好的方法,正如@RicardoCárdenes 所展示的,您可以在一行中完成所有操作,但是为了任何必须阅读您的代码的人,请不要这样做。
  • OP 的解决方案为 2 个循环,我的为 1,你的为 0...所以你的是最好的。通常情况下,它也更难阅读
  • @Bago:实际上......我的意思不是只写一行;),而是为了避免循环。单线实际上只是一个副作用 O:)
【解决方案2】:

您需要在使用hypot的地方定义一个距离函数。通常在 K-means 中是 距离=总和((点质心)^2) 这是一些执行此操作的 matlab 代码……如果您不能,我可以移植它,但请试一试。就像你说的,学习的唯一途径。

function idx = findClosestCentroids(X, centroids)
%FINDCLOSESTCENTROIDS computes the centroid memberships for every example
%   idx = FINDCLOSESTCENTROIDS (X, centroids) returns the closest centroids
%   in idx for a dataset X where each row is a single example. idx = m x 1 
%   vector of centroid assignments (i.e. each entry in range [1..K])
%

% Set K
K = size(centroids, 1);

[numberOfExamples numberOfDimensions] = size(X);
% You need to return the following variables correctly.
idx = zeros(size(X,1), 1);


% Go over every example, find its closest centroid, and store
%               the index inside idx at the appropriate location.
%               Concretely, idx(i) should contain the index of the centroid
%               closest to example i. Hence, it should be a value in the 
%               range 1..K
%
for loop=1:numberOfExamples
    Distance = sum(bsxfun(@minus,X(loop,:),centroids).^2,2);
    [value index] = min(Distance);
    idx(loop) = index;
end;


end

更新

这应该返回距离,请注意上面的 matlab 代码只返回最近质心的距离(和索引)...您的函数返回所有距离,如下所示。

def FindDistance(X,centroids):
K=shape(centroids)[0]
examples, dimensions = shape(X)
distance = zeros((examples,K))
for ex in xrange(examples):
    distance[ex,:] = np.sum((X[ex,:]-centroids)**2,1)
return distance

【讨论】:

  • 所以我必须在一个循环中分别计算每个点的距离?
  • 你应该可以在你已经拥有的 for 循环中做到这一点,只需在 for 循环中定义 Distance[i] 而不是 temp...等一下,我会看看我有没有'没有这样做
  • @yasar11732 实际上...您应该可以一口气完成...但这取决于样本的大小(或者,如果您愿意,还取决于您的 RAM 数量可以致力于此)。让我们看看我能不能想出它...
  • 谢谢,我想出了一个稍微不同的方法。你可能想看看我对我的问题的回答。
猜你喜欢
  • 1970-01-01
  • 2021-03-05
  • 2014-10-02
  • 2016-06-27
  • 2018-08-05
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
相关资源
最近更新 更多