爲什麼scipy.stats.nanmean會給出numpy.nansum的不同結果？

>>> import numpy as np 
>>> from scipy import stats 
>>> a = np.r_[1., 2., np.nan, 4., 5.] 
>>> stats.nanmean(a) 
2.9999999999999996 
>>> np.nansum(a)/np.sum(~np.isnan(a)) 
3.0

我意識到浮點表示的侷限性。只是好奇爲什麼更笨拙的表達似乎給「更好」的結果。爲什麼scipy.stats.nanmean會給出numpy.nansum的不同結果？

來源

2013-01-23 herrlich10

首先，這裏是scipy.nanmean()讓我們知道我們正在比較：

def nanmean(x, axis=0): 
    x, axis = _chk_asarray(x,axis) 
    x = x.copy() 
    Norig = x.shape[axis] 
    factor = 1.0-np.sum(np.isnan(x),axis)*1.0/Norig 

    x[np.isnan(x)] = 0 
    return np.mean(x,axis)/factor

在數學上，這兩種方法是等效的。在數字上，它們是不同的。

你的方法涉及單個分割，並且它恰巧：

分子（1. + 2. + 4. + 5.）可以準確地被表示爲float;和
分母（4.）是兩個冪。

這意味着劃分的結果是準確的，3.。

stats.nanmean()涉及首先計算平均值[1., 2., 0., 4., 5.]，然後調整它以計入NaNs。碰巧，這意味着（2.4）不能完全表示爲float，所以從這一點來說，計算是不精確的。

我還沒有給它很多想法，但是可能構建一個角色將被顛倒的例子，並且stats.nanmean()會比另一個方法給出更準確的結果。

最讓我驚訝的是，stats.nanmean()不是簡單地做一些事情，如：

In [6]: np.mean(np.ma.MaskedArray(a, np.isnan(a))) 
Out[6]: 3.0

在我看來這是一個優越的方法是什麼目前並。

來源

2013-01-23 09:39:42 NPE

確實，'b = np.r_ [1，2，np.nan，4，8.]'對np.mean更友好，。但我發現很難構建一個反向的例子:) – herrlich10

屏蔽數組很慢（在純Python中實現），所以我猜測提問者提出了什麼（'np.nansum（a）/np.sum（〜np.isnan（a ）''實際上比'np.mean（np.ma.MaskedArray（a，np.isnan（a））''更快'有人應該試試:) –

所以是的，用一個很長的一維數組'dat'， 'np.nansum（dat）/ np.sum（〜np.isnan（dat））'比'np.mean（np.ma.masked_array（dat，np.isnan（dat））''執行速度快10％'。瓶頸的南，，然而，執行10倍的速度。 –

答案是在stats.nanmean代碼：

x, axis = _chk_asarray(x,axis) 
x = x.copy() 
Norig = x.shape[axis] 
factor = 1.0-np.sum(np.isnan(x),axis)*1.0/Norig 
x[np.isnan(x)] = 0 
return np.mean(x,axis)/factor

我相信它有事情做與1.0 - np.sum，總和的減法。

來源

2013-01-23 08:51:16 eumiro

謝謝你指出來源。 – herrlich10

正如@eumiro提到，stats.nanmean計算出的平均值在從circumlocutions你做

從相同的參考代碼的簡單一個襯墊的方式方法不同，

np.sum(np.isnan(x),axis)回報numpy.int32這時候乘* 1.0，得到浮點近似值，而不是當結果是整數時導致結果差異的結果

>>> numpy.int32(1)*1.0/5 
0.20000000000000001 
>>> int(numpy.int32(1))*1.0/5 
0.2 
>>> type(np.sum(np.isnan(x),axis)) 
<type 'numpy.int32'>

來源

2013-01-23 09:03:52 Abhijit

我想你指出了另一個有趣的numpy行爲：'a = np.int_（1）/5.0; np.float_（a） - > 0.20000000000000001; float（a） - > 0.2' – herrlich10

爲什麼scipy.stats.nanmean會給出numpy.nansum的不同結果？

回答

相關問題