Numpy：如何避免舍入錯誤

我認爲我有一個Numpy和小數字的問題。Numpy：如何避免舍入錯誤

你能幫我找到以下問題的解決方案：

import numpy as np 

def gaussian(xx, mu=0, sigma=1): 
    return 1./(np.sqrt(2*np.pi)*sigma)*np.exp(-(mu-xx)**2/(2*sigma**2)) 

factors = (1., 0.1, 0.01, 0.001, 0.0001) 
for factor in factors: 
    xx = np.linspace(5000, 5200, 1000) 
    yy = 1.-(factor*xx*(1.+gaussian(xx, 5100))) 
    step = xx[1] - xx[0] 

    print np.sum((1.-yy/(1.-factor*xx))*step)

此代碼應計算爲-1所有不同factors的。但輸出是：

1.0 -1.00019611689 
0.1 -1.00196463662 
0.01 -1.0200000008 
0.001 -1.24390245353 
0.0001 1.04081641153

所以問題是，該係數越小，越我惹上麻煩，因爲我認爲有numpy的/ Python的精度做。

即使對於小因素，該如何評估等式？

非常感謝您的幫助。

來源

2014-02-07 dotcs

我懷疑問題的至少一部分是在數學中（或者更可能是在將數學翻譯成Python + NumPy代碼時）。你爲什麼認爲這應該評估爲-1？這是基於什麼？ –

將'yy'部分重寫爲'1-ax（1 + g（x））'，然後將和中的部分改寫爲1-yy /（1-ax）= 1 - （（1-ax因爲g（x）是一個高斯函數，它的峯值在'5100'處，因此g（x）= g（x）和1的stddev，我們在'200'（從'5000'到'5200'）的範圍內進行整合，'g（x）'將是1，因此結果將是'-1'。認爲數學不是問題，而是Numpy/Python處理小數字的方式 – dotcs

不可以：你的代數有一個錯誤，你的表達式簡化爲'axg（x）/（1-ax）'，而不是' - g（x）'。 –

這不是一個NumPy問題。您的期望結果應該是-1不正確。

您正在有效計算一個函數f(x)的定積分，大約在5100周圍。您正在集成的功能簡化爲factor * x * gaussian(x)/(1 - factor * x)。我們可以很容易地爲您使用的因子計算積分值的回差估計值：數量factor * x/(1 - factor * x)在感興趣的範圍內變化相當緩慢，大約爲5095至5105;對於這個範圍以外的任何東西，高斯將會使貢獻忽略不計。因此，第一次近似我們可以將此數量視爲常數。然後我們剩下的那個常數乘以gaussian(x)的積分，這將足夠接近1。因此預期產量應該在factor * 5100/(1 - factor * 5100)的地區。（雖然這不會工作這麼好當factor接近1/5100。）

因此，例如在factor是0.0001的情況下，factor * 5100/(1 - factor * 5100)大約有1.0408163265306123價值。這與你所看到的答案非常接近，使NumPy做出了或多或少的正確的事情。

來源

2014-02-07 23:40:19

是的，你是完全正確。我想簡化我遇到的問題，並犯了上述錯誤。我對此感到非常抱歉，雖然我檢查了我的簡化，但是我搞砸了。謝謝你的努力。對我感到羞恥，並讚揚你的改正。 – dotcs

Numpy：如何避免舍入錯誤

回答

相關問題