Numpy/Python与Matlab相比表现非常糟糕

Question

新手程序员在这里.我正在编写一个程序来分析点(单元格)的相对空间位置.程序从数组中获取边界和单元格类型,第1列中的x坐标,第2列中的y坐标和第3列中的单元格类型.然后,它会检查每个单元格的单元格类型以及与边界的适当距离.如果它通过,则计算它与阵列中每个其他单元格的距离,如果距离在指定的分析范围内,则将其添加到该距离的输出数组.

我的单元格标记程序是在wxpython中,所以我希望在python中开发这个程序并最终将它粘贴到GUI中.不幸的是,现在python在我的机器上运行核心循环需要大约20秒,而MATLAB可以执行~15循环/秒.由于我计划在大约30个案例中对几个探索性分析类型进行1000次循环(具有随机比较条件),这不是一个微不足道的差异.

我尝试运行一个分析器,数组调用是1/4的时间,几乎所有其余的都是未指定的循环时间.

这是主循环的python代码:

for basecell in range (0, cellnumber-1):
    if firstcelltype == np.array((cellrecord[basecell,2])):
        xloc=np.array((cellrecord[basecell,0]))
        yloc=np.array((cellrecord[basecell,1]))
        xedgedist=(xbound-xloc)
        yedgedist=(ybound-yloc)
        if xloc>excludedist and xedgedist>excludedist and yloc>excludedist and    yedgedist>excludedist:
            for comparecell in range (0, cellnumber-1):
                if secondcelltype==np.array((cellrecord[comparecell,2])):
                    xcomploc=np.array((cellrecord[comparecell,0]))
                    ycomploc=np.array((cellrecord[comparecell,1]))
                    dist=math.sqrt((xcomploc-xloc)**2+(ycomploc-yloc)**2)
                    dist=round(dist)
                    if dist>=1 and dist<=analysisdist:
                         arraytarget=round(dist*analysisdist/intervalnumber)
                         addone=np.array((spatialraw[arraytarget-1]))
                         addone=addone+1
                         targetcell=arraytarget-1
                         np.put(spatialraw,[targetcell,targetcell],addone)

这是主循环的matlab代码:

for basecell = 1:cellnumber;
    if firstcelltype==cellrecord(basecell,3);
         xloc=cellrecord(basecell,1);
         yloc=cellrecord(basecell,2);
         xedgedist=(xbound-xloc);
         yedgedist=(ybound-yloc);
         if (xloc>excludedist) && (yloc>excludedist) && (xedgedist>excludedist) && (yedgedist>excludedist);
             for comparecell = 1:cellnumber;
                 if secondcelltype==cellrecord(comparecell,3);
                     xcomploc=cellrecord(comparecell,1);
                     ycomploc=cellrecord(comparecell,2);
                     dist=sqrt((xcomploc-xloc)^2+(ycomploc-yloc)^2);
                     if (dist>=1) && (dist<=100.4999);
                         arraytarget=round(dist*analysisdist/intervalnumber);
                         spatialsum(1,arraytarget)=spatialsum(1,arraytarget)+1;
                    end
                end
            end            
        end
    end
end

谢谢!

Answer 1

以下是一些加速python代码的方法.

第一:当你只存储一个值时,不要制作np数组.您在代码中多次执行此操作.例如,

if firstcelltype == np.array((cellrecord[basecell,2])):

可以

 if firstcelltype == cellrecord[basecell,2]:

我会用一些时间陈述向你展示原因:

>>> timeit.Timer('x = 111.1').timeit()
0.045882196294822819
>>> t=timeit.Timer('x = np.array(111.1)','import numpy as np').timeit()
0.55774970267830071

这些调用之间的差异是一个数量级.

第二:以下代码:

arraytarget=round(dist*analysisdist/intervalnumber)
addone=np.array((spatialraw[arraytarget-1]))
addone=addone+1
targetcell=arraytarget-1
np.put(spatialraw,[targetcell,targetcell],addone)

可以替换为

arraytarget=round(dist*analysisdist/intervalnumber)-1
spatialraw[arraytarget] += 1

第三:你可以像菲利普analysisdist先前提到的那样摆脱平方.但是,既然你习惯analysisdist了arraytarget,你可能想要创建一个单独的变量,analysisdist2即analydist的平方,并将其用于比较.

第四:您正在寻找secondcelltype每次到达该点时匹配的单元格,而不是一次又一次地查找这些单元格.你可以定义一个数组:

comparecells = np.where(cellrecord[:,2]==secondcelltype)[0]

然后更换

for comparecell in range (0, cellnumber-1):
    if secondcelltype==np.array((cellrecord[comparecell,2])):

同

for comparecell in comparecells:

第五:使用psyco.它是一个JIT编译器.如果您使用的是最近的版本,Matlab有一个内置的JIT编译器.这应该会加快你的代码速度.

第六:如果在之前的所有步骤之后代码仍然不够快,那么您应该尝试对代码进行矢量化.这应该不会太难.基本上,numpy数组中的东西越多越好.这是我对矢量化的尝试:

basecells = np.where(cellrecord[:,2]==firstcelltype)[0]
xlocs = cellrecord[basecells, 0]
ylocs = cellrecord[basecells, 1]
xedgedists = xbound - xloc
yedgedists = ybound - yloc
whichcells = np.where((xlocs>excludedist) & (xedgedists>excludedist) & (ylocs>excludedist) & (yedgedists>excludedist))[0]
selectedcells = basecells[whichcells]
comparecells = np.where(cellrecord[:,2]==secondcelltype)[0]
xcomplocs = cellrecords[comparecells,0]
ycomplocs = cellrecords[comparecells,1]
analysisdist2 = analysisdist**2
for basecell in selectedcells:
    dists = np.round((xcomplocs-xlocs[basecell])**2 + (ycomplocs-ylocs[basecell])**2)
    whichcells = np.where((dists >= 1) & (dists <= analysisdist2))[0]
    arraytargets = np.round(dists[whichcells]*analysisdist/intervalnumber) - 1
    for target in arraytargets:
        spatialraw[target] += 1

您可以取出内部for循环,但必须小心,因为arraytargets的某些元素可能是相同的.此外,我实际上没有尝试所有的代码,因此可能存在错误或拼写错误.希望它能让你很好地了解如何做到这一点.哦,还有一件事.你创建analysisdist/intervalnumber一个单独的变量,以避免一遍又一遍地进行这种划分.