python/numpy 的时间和内存高效随机抽样

Question

我需要计算单位正方形内2个随机均匀分布点之间距离的期望值。这个问题可以解析求解，答案大约是 0.521405。我想要完成的是使用随机数生成器和 Python 使用和不使用 numpy 来得出这个数字。使用以下代码

dist=list()
random.seed()
for each in range(100000000):
    x1=random.random()
    x2=random.random()
    y1=random.random()
    y2=random.random()
    dist.append(math.sqrt((x1-x2)**2+(y1-y2)**2))
print(round(sum(dist)/float(len(dist)),6))

它在我的机器上 125 秒 迭代了 1 亿次，但只有 4 个小数位是正确的。现在，我用 numpy 创建了以下代码

dist=list()
start_time = time.time()
np.random.seed()
for each in range(100000000):
    x = np.array((np.random.random(),np.random.random()))
    y = np.array((np.random.random(),np.random.random()))
    dist.append(np.linalg.norm(x-y))
print(round(np.mean(dist),6))

并且需要 1111 秒来迭代相同的 1 亿 次！

由于结果只有 4 个小数位是正确的，所以我尝试将迭代次数增加到 10 亿 使用以前的版本，没有 numpy。 我想，那由于每个浮点数最多为 64 位（我使用的是 64 位 Python），因此该列表大约需要 8 GB。但是，该程序用完了 26GB 内存，并且当列表有 7.9 亿项

时出现异常错误

所以我在以下方面征求您的意见：

有没有办法利用各种 numpy 优化并实际使其运行得更快？
有没有办法让程序更节省内存？我意识到列表是一种比我需要的更复杂的数据结构
我是否正确地假设为了得到正确的 6 位十进制数字我需要接近 10^12 的迭代次数？（因为 N 次测量的标准误差以 1/sqrt(N) 递减）

提前致谢！

Answer 1

为了回答前两个子问题，这里有一种方法一次性生成所有随机数，然后利用 np.einsum 替换大部分 np.linalg.norm 的工作（平方微分和沿行求和），保持代码的其余部分不变。实现看起来像这样 -

N = 100000 # Edit this to change number of random points
x,y = np.random.random((2,N,2))
diffs = x-y
out = round(np.mean(np.sqrt(np.einsum('ij,ij->i',diffs,diffs))),6)

Answer 2

我建议这个变体 没有 numpy，假设你使用的是 Python 3:

random.seed()
dist_count = 100000000
dist_sum = 0

for _ in range(dist_count):
    dx = random.random() - random.random()  # x1 - x2
    dy = random.random() - random.random()  # y1 - y2

    dist += math.sqrt( dx*dx + dy*dy )

dist_average = dist_sum / dist_count
print(round(dist_average, 6))

首先，如果我们只是计算并加总距离，为什么要将所有距离存储在一个列表中？直接将每个随机距离添加到一个整数变量会更快。此外，我们已经知道我们创建了多少个随机距离，因为这是我们在 for 循环的范围内指定的，所以我们不需要 len(dist) 之类的东西或单独的计数器。

此外，我们不必为每个坐标指定名称，我们可以即时计算 dx 和 dy 的差异。这也有助于我们进行下一步。

将相同的值与其自身相乘比将其乘以 2 的幂 (further reading...) 快一个数量级。因此，我们当然将 a**2 替换为 a*a。现在我们直接存储上面的差异就有用了。

最后我们将距离总和除以计数并显示一次结果。

Answer 3

至于内存问题，您对'float takes 64 bits'的假设不太正确。每个浮动对象（实际上，python 中的盒装对象）将占用 24 个字节。您可以使用 sys.getsizeof(0.0) 检查它。因此，您的程序需要的 space 应该是您估计的 3 倍，这大约是您实际经历过的。

Answer 4

好吧，hypot 在我的笔记本电脑上要快一点，13.5 秒对 einsum 的 15.2 秒内存方面，它可能不适用于 10 亿，当前版本需要 ~6G

代码：

import numpy as np

np.random.seed(12345)

N = 100000000 # Edit this to change number of random points
x,y = np.random.random((2,N,2))
diffs = x - y

out = round(np.mean( np.hypot(diffs[:,0], diffs[:,1]) ),6)

#out = round(np.mean(np.sqrt(np.einsum('ij,ij->i',diffs,diffs))),6)

print(out)

python/numpy 的时间和内存高效随机抽样

Time and memory-efficient random sampling with python/numpy

python

random

numpy

list

out-of-memory