R应用在Python中的匹配速度

Question

我想有效地将复杂的函数应用于 Python 中的矩阵行（编辑：Python 3）。在 R 中，这是 apply 函数和它的近亲，它工作得很快。

在 Python 中，我知道这可以通过几种方式完成。列表理解、numpy.apply_along_axis、panas.dataframe.apply。

在我的编码中，这些 Python 方法非常慢。我应该使用另一种方法吗？或者也许我对这些 Python 方法的实现不正确？

这是一个例子。数学取自概率回归模型。要明确我的目标不是执行概率回归，我对一种有效的应用方法很感兴趣。

在 R 中：

> n = 100000
> p = 7
> x = matrix(rnorm(700000, 0 , 2), ncol = 7)
> beta = rep(1, p)

> start <- Sys.time()
> test <- apply(x, 1, function(t)(dnorm(sum(t*beta))*sum(t*beta)/pnorm(sum(t*beta))) )
> end <- Sys.time()
> print(end - start)
Time difference of 0.6112201 secs

在 Python 中通过理解：

import numpy as np
from scipy.stats import norm
import time

n = 100000
p = 7
x = norm.rvs(0, 2, n * p)
x = x.reshape( (n , p) )
beta = np.ones(p)

start = time.time()
test = [ 
norm.pdf(sum(x[i,]*beta))*sum(x[i,]*beta)/norm.cdf(sum(x[i,]*beta)) 
for i in range(100000) ]
end = time.time()
print (end - start)
23.316735982894897

在 Python 中通过 pandas.dataframe.apply：

frame = DataFrame(x)
f = lambda t: norm.pdf(sum(t))*sum(t)/norm.cdf(sum(t))
start = time.time()
test = frame.apply(f, axis = 1)
end = time.time()
print(end - start)
34.39404106140137

在this 问题中，最受好评的回答指出 apply_along_axis 不是为了速度。所以我不包括这种方法。

同样，我对快速执行这些计算很感兴趣。非常感谢您的帮助！

Answer 1

列表推导使用效率非常低的 Python 级别循环。您应该修改代码以利用 numpy 矢量化。如果您更改 time.time() 调用之间的内容

xbeta = np.sum(x * beta, axis=1)
test = norm.pdf(xbeta) * xbeta / norm.cdf(xbeta)

您会看到巨大的差异。对于我的机器，它在 0.02 秒内完成。 为了让您安心，我已经根据您的列表理解对其进行了测试，它们都给出了相同的结果。

xbeta 是你在计算中浪费了很多次的东西。通过沿第二轴求和，我们将其折叠为一维数组，即 100000 行。现在所有的计算都处理一维数组，所以我们让 numpy 处理剩下的事情。

Answer 2

这似乎是你想要的：

y = np.dot(x, beta)
test2 = norm.pdf(y) * y / norm.cdf(y)

# let's compare it to the expected output:
np.allclose(test2, np.array(test))
True

根据 ipython，这会计算“测试”列表中的所有 100000 个值（在数值容差范围内），但运行时间约为 11.5 毫秒：

%time y = np.dot(x, beta); test2 = norm.pdf(y) * y / norm.cdf(y)
CPU times: user 10 ms, sys: 0 ns, total: 10 ms
Wall time: 15.2 ms

这通过以下方式改进了您的版本：

1. 消除常见的sum(x[i,]*beta) 表达式
2. 将循环从 i=0 到 100000 向量化为点积。

在两者中，第二个更为重要。

我还注意到您上面的代码使用了sum，这是一个 python 内置函数，可以对任何迭代器求和。你几乎肯定应该使用np.sum，它是专门用于 numpy 数组的 numpy 矢量化版本。我不得不提一下这一点，因为一旦我以“批处理”形式重写了您的代码，sum 就隐含在点积中，因此np.sum 不会出现在最终版本中。