python - 在 Python 中向量化多维函数

Question

我经常潜伏在 Stack Overflow 上有一段时间了，每当我遇到编码问题时，我往往会从这里找到非常有用且清晰的信息。然而，我今天似乎找不到一个线索来解决我的具体问题。

今天早些时候，我学习了 Python 中的向量化函数，以加快计算速度。我目前正在尝试优化一个多月前编写的 python 程序。我的程序采用一个包含以下格式数据的文本文件:

<magnitude> <dmagnitude> <exposure_number>

然后，我将每一列分配给列表 mag , dmag ，和expnum .

我想要做的是创建 mag 的二维数组和dmag共享相同的值 expnum (具有相同的曝光数意味着 mag 和 dmag 指向相同的数据点)。

我对所有曝光数都这样做，最后，我取 mag 的中位数。和dmag ，以及 mag 的标准差对于每个基于曝光数的数组，并将它们全部组合成一个我可以绘制的数组。

目前，我有以下代码:

from numpy import loadtxt,array,asarray,append,std,median,empty,take

data = loadtxt(infile,usecols=(0,1,2))
mag = data1[:,2].tolist() 
dmag = data1[:,3].tolist() 
expnum = data1[:,4].tolist() 

#initialize variables
indexing = list() 
master_mag = list() 
master_dmag = list() 
sub_mag = list() 
sub_dmag = list() 

mag_std = array([]) 
mag_stdmed = array([]) 
mag_med = array([])  

while len(mag) > 0: 
    num=expnum[0] 
    for i in range(0,len(expnum)): 
        if expnum[i] == num: 
            sub_mag.append(mag[i]) 
            sub_dmag.append(dmag[i]) 
            indexing.append(i) 

    #add the sub lists to their master lists
    master_mag.append(sub_mag) 
    master_dmag.append(sub_dmag) 
    sub_mag=list() 
    sub_dmag=list()

    #remove from mag, dmag, and expnum the index referred to by indexing
    while len(indexing) > 0:    
        mag.pop(indexing[-1]) 
        dmag.pop(indexing[-1]) 
        expnum.pop(indexing[-1]) 
        indexing.pop() 

#make the master mag and dmag lists into numpy arrays 
master_mag=asarray(master_mag) 
master_dmag=asarray(master_dmag) 

#generate the mag and dmag median and mag std arrays 
for i in range(0,len(master_mag)): 
    mag_std=append(mag_std,std(master_mag[i])) 
    mag_med=append(mag_med,median(master_mag[i])) 
    mag_stdmed=append(mag_stdmed,median(master_dmag[i])) 

#create empty numpy arrays to be used for mag med vs. mag std 
#and mag med vs. dmag med 
med_std=empty([0,2]) 
med_dmed=empty([0,2]) 

#fill in those arrays 
for i in range(0,len(mag_std)): 
    med_std=append(med_std,[[mag_med[i],mag_std[i]]],axis=0) 
    med_dmed=append(med_dmed,[[mag_med[i],mag_stdmed[i]]],axis=0) 

#sort the median mag and dmag standard deviation arrays by median mag 
order_med_std=med_std[:,0].argsort() 
order_med_dmed=med_dmed[:,0].argsort() 

sorted_med_std=take(med_std,order_med_std,0) 
sorted_med_dmed=take(med_dmed,order_med_dmed,0) 

然后我准备绘制 sorted_med_dmed[:,0]与 sorted_med_dmed[:,1]和sorted_med_std[:,0]与 sorted_med_std[:,1]

这段代码可以工作，只是我觉得它太慢了(特别是当我需要处理超过 10,000 个数据点时)。我想尝试向量化这段代码以使其更快，但我不知道从哪里开始。

我需要一些帮助来弄清楚如何矢量化按曝光数匹配的组件。我想从一开始就创建一个多维数组，其格式为: array([[[mag],[dmag]],...])长度等于不同曝光次数的数量。有没有一种方法可以生成和更新这样的内联数组，而不必使用大量循环？

如果您需要进一步了解此代码的具体用途，请告诉我。

感谢您的宝贵时间。

Answer 1

解决此类问题的第一步应该始终是分析。我建议尝试 line_profiler ，因为它可以轻松地从视觉上找到热点。 (您也可以尝试 Python 的内置分析器，但我发现它的输出更难解析。)

这应该能让您了解哪些部分最能减慢您的代码速度。无需亲自尝试，我可以提供一些建议:

• 当 numpy 数组就足够时，尽量避免使用 Python 列表。如果您执行大量append，列表速度很快，但对于大多数其他操作，列表速度很慢，而且它们不支持矢量化。
• 与此相关，如果可以的话，请尽量避免调用 numpy.append。每次调用都涉及分配更多内存和复制，这在循环中可能会非常慢。

• 使用字典按键对数据进行分组。我发现 stdlib collections.defaultdict 对于这样的分组非常有用:

  groups = defaultdict(list)
  for a,b,key in data:
    groups[key].append((a,b))
  

• 使用 numpy 的自动向量化函数调用，而不是在循环中调用函数。例如，这段代码:

  #generate the mag and dmag median and mag std arrays 
  for i in range(0,len(master_mag)): 
    mag_std=append(mag_std,std(master_mag[i])) 
    mag_med=append(mag_med,median(master_mag[i])) 
    mag_stdmed=append(mag_stdmed,median(master_dmag[i]))
  

  写成这样会快得多:

  mag_std = numpy.std(master_mag, axis=0)
  mag_meg = numpy.median(master_mag, axis=0)
  mag_stdmed = numpy.median(master_dmag, axis=0)