如何在不使用 for 循环的情况下将矢量化函数映射到 numpy 数组？答案

【问题标题】：How can I map a vectorized function to a numpy array without using a for loop?如何在不使用 for 循环的情况下将矢量化函数映射到 numpy 数组？
【发布时间】：2018-10-14 04:19:04
【问题描述】：

这就是我已经拥有的：

import numpy as np
import matplotlib.pyplot as plt

def monteCarloPi(n):
    np.random.seed()                 #seed the random number generator
    y = np.random.rand(n)*2 - 1      #n random samples on (-1,1)
    x = np.linspace(-1,1,n)          #x axis to plot against

    square = np.array([x,y])         #collecting axes as a single object

    mask1 = ((x**2 + y**2) < 1)      #filters

    hits = np.sum(mask1)             #calculating approximation
    ratio = hits/n
    pi_approx = ratio * 4

    return pi_approx

这是我希望做的事情：

x = np.arange(100,1000)
y = monteCarloPi(x)
plt.scatter(x,y)

但是，当我运行上述代码块时，出现以下错误：

---------------------------------------------------------------------
TypeError                           Traceback (most recent call last)
<ipython-input-52-bf4dcedaa309> in <module>()
      1 x = np.arange(100,1000)
----> 2 y = monteCarloPi(x)
      3 plt.scatter(x,y)

    <ipython-input-51-8d5b36e22d4b> in monteCarloPi(n)
      1 def monteCarloPi(n):
      2     np.random.seed()                 #seed the random number generator
----> 3     y = np.random.rand(n)*2 - 1      #n random samples on (-1,1)
      4     x = np.linspace(-1,1,n)          #x axis to plot against
      5 

mtrand.pyx in mtrand.RandomState.rand()

mtrand.pyx in mtrand.RandomState.random_sample()

mtrand.pyx in mtrand.cont0_array()

TypeError: only integer scalar arrays can be converted to a scalar index

根据我对numpy 中广播如何工作的理解，这应该工作。我可以只使用for 循环，但随着样本数量的增加，它会变得非常慢。

暂停

【问题讨论】：

random.rand 只接受标量值、数字、它应该生成的数组的形状。您正在向它传递一个包含许多值的数组x。以下带有linspace 的行也会有这个n 的问题。我不明白 broadcasting 与此有什么关系。
好吧...所以你要告诉我的是，对许多 n 值运行此模拟的唯一方法是使用 for 循环？
我想我认为这与广播有关的原因是该函数输出一个标量；我想要发生的是放入一个数组并取出一个新数组。
快速编译的 numpy 代码适用于矩形数组 - 每行具有相同的列数，对于更高的维度以此类推。广播让我们可以对具有各种维度组合的数组进行数学运算。严格意义上的向量化意味着编写可以对具有不同维数的数组进行操作的代码。它仍然循环，但它将这些循环向下移动到已编译的代码块中。
在你的函数内部，n 的不同值会产生不同大小的数组，y、x、square。 sum 将它们压缩回一个标量。它不能同时处理多个n，因为这将涉及y 等不同大小的数组。 y for n=100 不能与 y for n=200 组合成一个二维数组。 x 等也一样。在 numpy 中，“参差不齐的数组”并不快。

标签： python numpy array-broadcasting numpy-ufunc

【解决方案1】：

这是一种选择，其中基于最大样本大小，然后在start>0 时进行子采样（不包括错误处理）。

import numpy as np
import matplotlib.pyplot as plt

def monteCarloPi(n,start=0,stride=1):
    np.random.seed()                 # seed the random number generator
    y = np.random.rand(n)*2 - 1      # n random samples on (-1,1)
    x = np.linspace(-1,1,n)          # x axis to plot against
    mask = ( x**2 + y**2 ) < 1       # masking
    samples = {}
    inds = arange(n)
    for k in range(n-start,n+1,stride):
        sub_inds = np.random.choice(inds,k,replace=False)
        sub_mask = mask[sub_inds]
        sub_hits = np.sum(sub_mask)
        ratio    = sub_hits/n
        pi_approx = ratio * 4
        samples[k]=pi_approx
    return pi_approx

这仍然需要一个 for 循环，但它可以在方法内部快速处理，因为您是从一个大的随机样本中进行二次抽样。重新创建您的原始呼叫（从n=100 运行到n=1000 [注意，我将在此处上升到n=1000]）：

estimates = monteCarloPi(1000,start=900)
plt.plot(estimates.keys(),estimates.values())

您当然可以传递原始的x=arange(100,1001)，但随后需要在方法中进行错误检查（以确保传递了数组或列表），然后n 将等于最后一个元素x (n=x[-1])，最后，循环将在 x (for k in x:) 的元素上完成。

【讨论】：