10

numpyで使用される数値が非常に少ないため、いくつかの問題があります。関数にfloatを合計すると、float64の精度が失われるという事実に、数値積分に関する私の絶え間ない問​​題をさかのぼるのに数週間かかりました。合計の代わりに積を使用して数学的に同一の計算を実行すると、問題のない値が得られます。

コードサンプルと結果のプロットを次に示します。

from matplotlib.pyplot import *
from numpy import vectorize, arange
import math

def func_product(x):
    return math.exp(-x)/(1+math.exp(x))

def func_sum(x):
    return math.exp(-x)-1/(1+math.exp(x))

#mathematically, both functions are the same

vecfunc_sum = vectorize(func_sum)
vecfunc_product = vectorize(func_product)

x = arange(0.,300.,1.)
y_sum = vecfunc_sum(x)
y_product = vecfunc_product(x)

plot(x,y_sum,    'k.-', label='sum')
plot(x,y_product,'r--',label='product')

yscale('symlog', linthreshy=1E-256)
legend(loc='lower right')
show()

ここに画像の説明を入力してください

ご覧のとおり、非常に低い合計値はゼロの周りに散在しているか、正確にゼロですが、乗算された値は問題ありません...

誰か助けて/説明してもらえますか?どうもありがとう!

4

3 に答える 3

6

浮動小数点の精度は、丸め誤差による加算/減算にかなり敏感です。最終的に1+exp(x)は非常に大きくなるため、exp(x)に1を追加すると、exp(x)と同じものになります。倍精度では、それはどこかにありexp(x) == 1e16ます:

>>> (1e16 + 1) == (1e16)
True
>>> (1e15 + 1) == (1e15)
False

それはおよそ37であることに注意してくださいmath.log(1e16)-これはおおよそあなたのプロットで物事が狂っているところです。

同じ問題が発生する可能性がありますが、規模は異なります。

>>> (1e-16 + 1.) == (1.)
True
>>> (1e-15 + 1.) == (1.)
False

あなたの政権のポイントの大部分について、あなたfunc_productは実際に計算しています:

exp(-x)/exp(x) == exp(-2*x)

これが、グラフの傾きが-2である理由です。

もう一方の極端に言えば、あなたは他のバージョンが(少なくともおおよそ)計算しています:

exp(-x) - 1./exp(x) 

これはおおよそです

exp(-x) - exp(-x)
于 2013-01-11T14:10:17.163 に答える
4

これは壊滅的なキャンセルの例です。

計算がうまくいかない最初のポイントを見てみましょう。x = 36.0

In [42]: np.exp(-x)
Out[42]: 2.3195228302435691e-16

In [43]: - 1/(1+np.exp(x))
Out[43]: -2.3195228302435691e-16

In [44]: np.exp(-x) - 1/(1+np.exp(x))
Out[44]: 0.0

を使用した計算func_productでは、ほぼ等しい数が減算されないため、壊滅的なキャンセルが回避されます。


ちなみに、に変更math.expすると、 (遅い)をnp.exp取り除くことができます:np.vectorize

def func_product(x):
    return np.exp(-x)/(1+np.exp(x))

def func_sum(x):
    return np.exp(-x)-1/(1+np.exp(x))

y_sum = func_sum_sum(x)
y_product = func_product_product(x)
于 2013-01-11T14:14:55.203 に答える
3

問題は、2つの非常に近い値の間の減算を伴うため、数値的に不安定func_sumであるということです。

func_sum(200)たとえば、の計算では、に追加しても効果がないため、同じ値にmath.exp(-200)なります。これは、64ビット浮動小数点の精度の範囲外であるためです。1/(1+math.exp(200))1math.exp(200)

math.exp(200).hex()
0x1.73f60ea79f5b9p+288

(math.exp(200) + 1).hex()
0x1.73f60ea79f5b9p+288

(1/(math.exp(200) + 1)).hex()
0x1.6061812054cfap-289

math.exp(-200).hex()
0x1.6061812054cfap-289

これはなぜfunc_sum(200)ゼロを与えるのかを説明していますが、x軸から外れている点はどうですか?これらは、浮動小数点の不正確さによっても発生します。;math.exp(-x)と等しくないことが時々起こります。1/math.exp(x)理想的には、math.exp(x)はに最も近い浮動小数点値e^xであり、1/math.exp(x)はによって計算される浮動小数点数の逆数に最も近い浮動小数点値であり、math.exp(x)必ずしも。にではありませんe^-x。確かに、math.exp(-100)そして1/(1+math.exp(100))は非常に近く、実際には最後のユニットでのみ異なります。

math.exp(-100).hex()
0x1.a8c1f14e2af5dp-145

(1/math.exp(100)).hex()
0x1.a8c1f14e2af5cp-145

(1/(1+math.exp(100))).hex()
0x1.a8c1f14e2af5cp-145

func_sum(100).hex()
0x1.0000000000000p-197

したがって、実際に計算したのは、との間の差math.exp(-x)です1/math.exp(x)。関数の行をトレースしmath.pow(2, -52) * math.exp(-x)て、の正の値を通過することを確認できますfunc_sum(52は64ビット浮動小数点の仮数のサイズであることを思い出してください)。

于 2013-01-11T14:14:12.197 に答える