python - 为什么 numpy 在数字化示例中比 matlab 慢得多？

Question

我正在比较 numpy 与 matlab 的性能，在一些情况下我观察到 numpy 明显更慢(索引、对数组的简单操作，例如绝对值、乘法、求和等)。让我们看一下下面的例子，它有点引人注目，涉及函数 digitize(我打算用它来同步时间戳):

import numpy as np
import time
scale=np.arange(1,1e+6+1)
y=np.arange(1,1e+6+1,10)
t1=time.time()
ind=np.digitize(scale,y)
t2=time.time()
print 'Time passed is %2.2f seconds' %(t2-t1)

结果是:

Time passed is 55.91 seconds

现在让我们使用等效函数 histc 尝试相同的示例 Matlab

scale=[1:1e+6];
y=[1:10:1e+6];
tic
[N,bin]=histc(scale,y);
t=toc;
display(['Time passed is ',num2str(t), ' seconds'])

结果是:

Time passed is 0.10237 seconds

快了 560 倍!

在学习使用 C++ 扩展 Python 时，我实现了自己的数字化版本(使用扩展库进行扩展):

import analysis # my C++ module implementing digitize
t1=time.time()
ind2=analysis.digitize(scale,y)
t2=time.time()
print 'Time passed is %2.2f seconds' %(t2-t1)
np.all(ind==ind2) #ok

结果是:

Time passed is 0.02 seconds

我的 digitize 版本假设输入都是单调的，这有一点作弊，这可以解释为什么它比 Matlab 更快。但是，对大小为 1e+6 的数组进行排序需要 0.16 秒(使用 numpy.sort)，因此与 Matlab 函数 历史。

所以问题是:

• 为什么 numpy.digitize 这么慢？这个函数不应该用编译和优化的代码编写吗？
• 为什么我自己的 digitize 版本比 numpy.digitize 快得多，但仍然比 Matlab 慢(假设输入已经排序，我非常有信心我使用可能最快的算法)？

我正在使用 Fedora 16，并且最近安装了 ATLAS 和 LAPACK 库(但性能发生了很大变化)。我应该重建 numpy 吗？我不确定我安装的 numpy 是否使用了适当的库来获得最大速度，也许 Matlab 使用了更好的库。

更新

根据目前的答案，我想强调的是，如果有人(在这种情况下像我一样)不关心直方图。我需要 hisc 的第二个输出，它是从输入值到所提供输入 bin 的索引的映射。这样的输出由 numpy 函数 digitize 和 searchsorted 提供。正如其中一个答案所说，searchsorted 比 digitize 快得多。但是，searchsorted 仍然比 Matlab 慢 2 倍:

t1=time.time()
ind3=np.searchsorted(y,scale,"right")
t2=time.time()
print 'Time passed is %2.2f seconds' %(t2-t1)

np.all(ind==ind3) #ok

结果是

Time passed is 0.21 seconds

那么现在的问题是:

1. 如果有一个等效函数 numpy.searchsorted 快 280 倍，那么使用 numpy.digitize 有什么意义?

2. 为什么 Matlab 函数 histc(它还提供 numpy.searchsorted 的输出)比 快 2 倍>numpy.searchsorted?

Answer 1

首先，让我们看看为什么 numpy.digitize 很慢。如果发现您的 bin 是单调的，则根据 bin 是非递减还是非递增调用这些函数之一(相关代码可在 numpy/lib/src/_compiled_base.c 中找到numpy git 仓库):

static npy_intp
incr_slot_(double x, double *bins, npy_intp lbins)
{
    npy_intp i;

    for ( i = 0; i < lbins; i ++ ) {
        if ( x < bins [i] ) {
            return i;
        }
    }
    return lbins;
}

static npy_intp
decr_slot_(double x, double * bins, npy_intp lbins)
{
    npy_intp i;

    for ( i = lbins - 1; i >= 0; i -- ) {
        if (x < bins [i]) {
            return i + 1;
        }
    }
    return 0;
}

如您所见，它正在进行线性搜索。线性搜索比二进制搜索慢得多，所以你可以回答为什么它慢。我会打开a ticket为此，在 numpy 跟踪器上。

其次，我认为 Matlab 实际上比您的 C++ 代码慢，因为 Matlab also assumes bins 是单调非递减的。