为什么 np.linalg.norm(..., axis=1) 比写出向量范数的公式慢？答案

【问题标题】：Why is np.linalg.norm(..., axis=1) slower than writing out the formula for vector norms?为什么 np.linalg.norm(..., axis=1) 比写出向量范数的公式慢？
【发布时间】：2021-03-05 00:04:21
【问题描述】：

为了将矩阵X 的行标准化为单位长度，我通常使用：

X /= np.linalg.norm(X, axis=1, keepdims=True)

尝试优化算法的这个操作，我很惊讶地发现在我的机器上写出规范化的速度大约快 40%：

X /= np.sqrt(X[:,0]**2+X[:,1]**2+X[:,2]**2)[:,np.newaxis]
X /= np.sqrt(sum(X[:,i]**2 for i in range(X.shape[1])))[:,np.newaxis]

怎么来的？ np.linalg.norm()的性能损失在哪里？

import numpy as np
X = np.random.randn(10000,3)

%timeit X/np.linalg.norm(X,axis=1, keepdims=True)
# 276 µs ± 4.55 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

%timeit X/np.sqrt(X[:,0]**2+X[:,1]**2+X[:,2]**2)[:,np.newaxis]
# 169 µs ± 1.38 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

%timeit X/np.sqrt(sum(X[:,i]**2 for i in range(X.shape[1])))[:,np.newaxis]
# 185 µs ± 4.17 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

我在支持 OpenBLAS 的 MacbookPro 2015 上观察到 (1) python3.6 + numpy v1.17.2 和 (2) python3.9 + numpy v1.19.3 的情况。

我不认为这是 this post 的重复，它解决了矩阵范数，而这个是关于向量的 L2 范数。

【问题讨论】：

查看源代码，它在后台处理了很多事情，这很可能就是原因所在。 “检查”这一点的一种快速方法是复制粘贴源，但删除所有不适用的杂物，然后再次运行测试
您可以尝试逐行分析。 stackoverflow.com/questions/3927628/…
@IanQuah this line 消耗了大部分时间。

标签： python numpy

【解决方案1】：

用于逐行 L2 范数的 source code 归结为以下代码行：

def norm(x, keepdims=False):
    x = np.asarray(x)
    s = x**2
    return np.sqrt(s.sum(axis=(1,), keepdims=keepdims))

简化的代码假定实值x 并利用np.add.reduce(s, ...) 等价于s.sum(...) 的事实。

因此，OP 问题与问为什么np.sum(x,axis=1) 比sum(x[:,i] for i in range(x.shape[1])) 慢：

%timeit X.sum(axis=1, keepdims=False)
# 131 µs ± 1.6 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
%timeit sum(X[:,i] for i in range(X.shape[1]))
# 36.7 µs ± 91.2 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

这个问题已经回答了here。简而言之，减少 (.sum(axis=1)) 伴随着开销成本，这些开销通常在浮点精度和速度方面得到回报（例如缓存机制、并行性），但在仅减少三列的特殊情况下不会.在这种情况下，与实际计算相比，开销相对较大。

如果X 有更多列，情况会发生变化。 numpy-boosted normalization 现在比使用 python for-loop 的归约快得多：

X = np.random.randn(10000,100)
%timeit X/np.linalg.norm(X,axis=1, keepdims=True)
# 3.36 ms ± 132 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
%timeit X/np.sqrt(sum(X[:,i]**2 for i in range(X.shape[1])))[:,np.newaxis]
# 5.92 ms ± 168 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

在这里找到另一个相关的 SO 线程：numpy ufuncs vs. for loop。

问题仍然是为什么 numpy 不明确地处理常见的缩减特殊情况（例如对具有低轴维度的矩阵的列或行的求和）。可能是因为这种优化的效果往往强烈依赖于目标机器，大大增加了代码复杂度。

【讨论】：