【问题标题】:Division by zero in numpy (sub)arrays在numpy(子)数组中除以零
【发布时间】:2014-01-14 09:46:08
【问题描述】:

我有三个数组,它们用数学函数处理以获得最终结果数组。一些数组包含 NaN,一些包含 0。但是除以零逻辑上会引发警告,使用 NaN 的计算会给出 NaN。所以我想对涉及零的数组的某些部分进行某些操作:

r=numpy.array([3,3,3])
k=numpy.array([numpy.nan,0,numpy.nan])
n=numpy.array([numpy.nan,0,0])
1.0*n*numpy.exp(r*(1-(n/k)))

例如在 k == 0 的情况下,我想得到结果 0。在所有其他情况下,我想计算上面的函数。那么有什么方法可以对数组的某些部分(通过索引)进行这样的计算以获得最终的单个结果数组?

【问题讨论】:

    标签: python arrays numpy divide-by-zero


    【解决方案1】:
    import numpy
    r=numpy.array([3,3,3])
    k=numpy.array([numpy.nan,0,numpy.nan])
    n=numpy.array([numpy.nan,0,0])
    indxZeros=numpy.where(k==0)
    indxNonZeros=numpy.where(k!=0)
    d=numpy.empty(k.shape)
    d[indxZeros]=0
    d[indxNonZeros]=n[indxNonZeros]/k[indxNonZeros]
    print d
    

    【讨论】:

    • 请使用 edit 按钮将此作为编辑添加到您的原始答案中。
    • 而不是使用np.where 来计算数字索引数组,使用布尔索引通常更快,即使用indxZeros = k == 0; indxNonZeros = ~idxZeros 也可以,可能更快。
    • 感谢tobias_k,我对在这里回答问题还是很陌生。我通常只是浏览已经问过的问题来回答我的问题。谢谢你的建议,下次我会记住这一点
    【解决方案2】:

    是否遵循您的需要?

    >>> rv = 1.0*n*numpy.exp(r*(1-(n/k)))
    >>> rv[k==0] = 0
    >>> rv
    array([ nan,   0.,  nan])
    

    【讨论】:

    • 基本上是的,但是这种方法首先一次性处理然后“清理”仍然会引发警告。我也想避免警告,例如排除之前已经有 k==0 的情况
    【解决方案3】:

    所以,你可能认为解决这个问题的方法是使用numpy.where,但是如下:

    numpy.where(k==0, 0, 1.0*n*numpy.exp(r*(1-(n/k))))
    

    仍然会给出警告,因为表达式实际上是针对 k 为零的情况进行评估的,即使这些结果没有被使用。

    如果这真的让您感到困扰,您可以使用 numexpr 作为这个表达式,它实际上会在 where 语句上分支,而不是评估 k==0 的情况:

    import numexpr
    numexpr.evaluate('where(k==0, 0, 1.0*n*exp(r*(1-(n/k))))')
    

    另一种方式,根据您要求的索引,会导致一些易读性损失

    result = numpy.zeros_like(k)
    good = k != 0
    result[good] = 1.0*n[good]*numpy.exp(r[good]*(1-(n[good]/k[good])))
    

    这可以通过定义gaussian 函数来绕过:

    def gaussian(r, k, n):
        return 1.0*n*numpy.exp(r*(1-(n/k)))
    
    result = numpy.zeros_like(k)
    good = k != 0
    result[good] = gaussian(r[good], k[good], n[good])
    

    【讨论】:

      猜你喜欢
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      • 2012-06-26
      • 1970-01-01
      • 1970-01-01
      • 1970-01-01
      相关资源
      最近更新 更多