2013-01-11 14 views
9

Numpy ile kullanılan düşük sayılar nedeniyle bazı problemlerim var. Sayısal entegrasyon ile sürekli problemlerimi geri takip etmek bana birkaç hafta sürdü, float64 hassasiyeti bir fonksiyona kayarken, hassasiyeti kaybeder. Matematiksel olarak özdeş hesaplamanın bir ürün yerine toplam yerine gerçekleştirilmesi, tamam olan değerlere yol açar. Gördüğünüz gibiSayılar eklerken python/numpy'de float hassas dökümü

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() 

enter image description here

sıfıra etrafa dağılmış oldukça düşüktür ya da özetlenebilir değerleri ise tam sıfır gibidir: Burada

bir kod örneği ve sonuçların bir komplodur çoğaltılan değerler gayet iyi ...

Lütfen birisi yardım edebilir mi? Çok teşekkürler!

+2

ise' np var .float128'. Aynı sorunları var mı? – eumiro

+0

Yep, float96 için aynı sorun ... –

+0

Kesin problemden kaçınmak için neden günlük alanında bu toplamları yapmıyorsunuz? http://lingpipe-blog.com/2012/02/16/howprevent-overflow-underflow-logistic-regression/ – jeff7

cevap

5

Kayan nokta hassasiyeti, dönüş hatası nedeniyle ekleme/çıkarma işlemine oldukça hassastır. Sonunda, 1+exp(x) o kadar büyük oluyor ki, 1'e exp ekliyor (x), exp (x) ile aynı şeyi veriyor. exp(x) == 1e16 etrafında bir yerlerde çift hassasiyet içinde:

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

Not math.log(1e16) o yaklaşık 37 - şeyler sizin arsa üzerinde çılgın nereye kabaca hangisi.

Aynı problem olabilir, ancak farklı ölçeklerde: En rejimde puan büyük çoğunluğu için

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

, senin func_product aslında hesaplıyor:

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

Hangi nedenle Grafiğinizin -2 eğimli güzel bir eğimi var.

Diğer uçta onu alarak, diğer versiyonu (en azından yaklaşık olarak) hesaplanması duyuyorsunuz:

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

yaklaşık sorun senin func_sumnumerically unstable çünkü olmasıdır

exp(-x) - exp(-x) 
2

olan iki çok yakın değer arasında bir çıkarma içerir.

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 

: math.exp(200) için 1 ekleyerek etkisi yoktur çünkü 64 bit kayan nokta hassas dışında olduğundan func_sum(200) hesaplanmasında

, örneğin, math.exp(-200) ve 1/(1+math.exp(200)), aynı değere sahip Bu, neden func_sum(200)'un sıfır verdiğini açıklıyor, ancak x ekseninden ayrılan noktalar nedir? Bunlar aynı zamanda kayan noktaların yanlış olmasından kaynaklanır; bazen math.exp(-x)'un 1/math.exp(x)'a eşit olmadığı; İdeal olarak, math.exp(x), e^x'a en yakın kayan nokta değeridir ve 1/math.exp(x), math.exp(x) tarafından hesaplanan kayan nokta sayısının karşılıklı olarak en yakın kayan nokta değeridir, e^-x için zorunlu değildir.Nitekim math.exp(-100) ve 1/(1+math.exp(100)) sadece son ünitesinde farklılık çok yakın ve aslında: Varsa

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 

Peki gerçekte hesapladık math.exp(-x) ve 1/math.exp(x) arasında, farktır. math.pow(2, -52) * math.exp(-x) işlevinin satırını, func_sum pozitif değerlerinden geçtiğini görmek için (52 değerinin, 64 bit kayan noktadaki büyüklüğün boyutu olduğunu hatırlayın) izleyebilirsiniz.

4

Bu, catastrophic cancellation örneğidir. hesaplama ters gider ilk noktada

atalım, ne zaman 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 

o felaket iptal önler, böylece neredeyse eşit sayıda çıkarma yapmaz func_product kullanarak hesaplama. Eğer math.expnp.exp geçtiğinizdeyse arada


, sen (yavaş olan) np.vectorize kurtulabilirsiniz: `np.float64` yeterli değildir

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)