Dlaczego np.dot jest nieprecyzyjny? (tablice n-dim)


15

Załóżmy, że bierzemy np.dotdwie 'float32'tablice 2D:

res = np.dot(a, b)   # see CASE 1
print(list(res[0]))  # list shows more digits
[-0.90448684, -1.1708503, 0.907136, 3.5594249, 1.1374011, -1.3826287]

Liczby. Z wyjątkiem, że mogą zmienić:


PRZYPADEK 1 : plastereka

np.random.seed(1)
a = np.random.randn(9, 6).astype('float32')
b = np.random.randn(6, 6).astype('float32')

for i in range(1, len(a)):
    print(list(np.dot(a[:i], b)[0])) # full shape: (i, 6)
[-0.9044868,  -1.1708502, 0.90713596, 3.5594249, 1.1374012, -1.3826287]
[-0.90448684, -1.1708503, 0.9071359,  3.5594249, 1.1374011, -1.3826288]
[-0.90448684, -1.1708503, 0.9071359,  3.5594249, 1.1374011, -1.3826288]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]

Wyniki różnią się, mimo że wydrukowany wycinek pochodzi z dokładnie tych samych pomnożonych liczb.


PRZYPADEK 2 : spłaszcz a, weź wersję 1D b, a następnie pokrój a:

np.random.seed(1)
a = np.random.randn(9, 6).astype('float32')
b = np.random.randn(1, 6).astype('float32')

for i in range(1, len(a)):
    a_flat = np.expand_dims(a[:i].flatten(), -1) # keep 2D
    print(list(np.dot(a_flat, b)[0])) # full shape: (i*6, 6)
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]

PRZYPADEK 3 : silniejsza kontrola; ustaw wszystkie niezaangażowane wejścia na zero : dodaj a[1:] = 0do kodu PRZYPADKU 1. Wynik: utrzymują się rozbieżności.


PRZYPADEK 4 : sprawdź wskaźniki inne niż [0]; podobnie jak w przypadku [0], wyniki zaczynają stabilizować ustaloną liczbę powiększeń tablicy od momentu ich utworzenia. Wynik

np.random.seed(1)
a = np.random.randn(9, 6).astype('float32')
b = np.random.randn(6, 6).astype('float32')

for j in range(len(a) - 2):
    for i in range(1, len(a)):
        res = np.dot(a[:i], b)
        try:    print(list(res[j]))
        except: pass
    print()

Dlatego w przypadku 2D * 2D wyniki są różne - ale są spójne dla 1D * 1D. Z niektórych moich odczytów wynika, że ​​wynika to z 1D-1D przy użyciu prostego dodawania, podczas gdy 2D-2D używa „bardziej zaawansowanego”, zwiększającego wydajność dodawania, który może być mniej precyzyjny (np. Dodawanie par robi odwrotnie). Niemniej jednak nie jestem w stanie zrozumieć, dlaczego rozbieżności znikają w przypadku, gdy 1 raz azostanie przekroczony przez ustalony „próg”; im większy ai bim później ten próg wydaje się kłamać, ale zawsze istnieje.

Wszyscy powiedzieli: dlaczego jest np.dotnieprecyzyjny (i niespójny) w przypadku tablic ND-ND? Odpowiedni Git


Dodatkowe informacje :

  • Środowisko : Win-10 OS, Python 3.7.4, Spyder 3.3.6 IDE, Anaconda 3.0 2019/10
  • Procesor : i7-7700HQ 2,8 GHz
  • Numpy v1.16.5

Możliwa biblioteka sprawców : Numpy MKL - także biblioteki BLASS; dzięki Bi Rico za odnotowanie


Kod testu warunków skrajnych : jak wspomniano, rozbieżności nasilają się w / w większych tablicach; jeśli powyższe nie jest powtarzalne, poniżej powinno być (jeśli nie, wypróbuj większe przyciemnienia). Moja produkcja

np.random.seed(1)
a = (0.01*np.random.randn(9, 9999)).astype('float32') # first multiply then type-cast
b = (0.01*np.random.randn(9999, 6)).astype('float32') # *0.01 to bound mults to < 1

for i in range(1, len(a)):
    print(list(np.dot(a[:i], b)[0]))

Istotność problemu : pokazane rozbieżności są „małe”, ale już nie tak, gdy działają w sieci neuronowej z miliardami liczb pomnożonymi przez kilka sekund i trylionami przez cały czas działania; podana dokładność modelu różni się o całe 10 procent w danym wątku .

Poniżej znajduje się gif tablic wynikających z karmienia do modelu, co w zasadzie a[0]w / len(a)==1vs len(a)==32.:


INNE PLATFORMY wyniki, zgodnie z podziękowaniami i dzięki testami Paula :

Przypadek 1 odtworzony (częściowo) :

  • Google Colab VM - Intel Xeon 2.3 G-Hz - Jupyter - Python 3.6.8
  • Win-10 Pro Docker Desktop - Intel i7-8700K - jupyter / scipy-notebook - Python 3.7.3
  • Ubuntu 18.04.2 LTS + Docker - AMD FX-8150 - jupyter / scipy-notebook - Python 3.7.3

Uwaga : dają one znacznie mniejszy błąd niż pokazano powyżej; dwa wpisy w pierwszym rzędzie są wyłączone o 1 w najmniej znaczącej cyfrze w stosunku do odpowiadających wpisów w innych wierszach.

Przypadek 1 nie został odtworzony :

  • Ubuntu 18.04.3 LTS - Intel i7-8700K - IPython 5.5.0 - Python 2.7.15+ i 3.6.8 (2 testy)
  • Ubuntu 18.04.3 LTS - Intel i5-3320M - IPython 5.5.0 - Python 2.7.15+
  • Ubuntu 18.04.2 LTS - AMD FX-8150 - IPython 5.5.0 - Python 2.7.15rc1

Uwagi :

  • Te związane Colab notebooków i jupyter środowisk wykazują znacznie mniejszą różnicę (i tylko dla pierwszych dwóch rzędach) niż obserwuje się w moim systemie. Ponadto przypadek 2 nigdy (jeszcze) nie wykazywał niedokładności.
  • W ramach tej bardzo ograniczonej próbki obecne (zadokowane) środowisko Jupyter jest bardziej podatne niż środowisko IPython.
  • np.show_config()zbyt długo, aby opublikować, ale w skrócie: Env IPython są oparte na BLAS / LAPACK; Colab jest oparty na OpenBLAS. W środowisku IPython Linux biblioteki BLAS są instalowane systemowo - w Jupyter i Colab pochodzą z / opt / conda / lib

AKTUALIZACJA : zaakceptowana odpowiedź jest dokładna, ale szeroka i niepełna. Pytanie pozostaje otwarte dla każdego, kto potrafi wyjaśnić zachowanie na poziomie kodu - mianowicie dokładny algorytm zastosowany przez np.doti jak wyjaśnia „spójne niespójności” zaobserwowane w powyższych wynikach (patrz także komentarze). Oto niektóre bezpośrednie implementacje poza moim odszyfrowaniem: sdot.c - arraytypes.c.src


Komentarze nie są przeznaczone do rozszerzonej dyskusji; ta rozmowa została przeniesiona do czatu .
Samuel Liew

Ogólne algorytmy ndarrayszwykle pomijają liczbową utratę precyzji. Ponieważ dla uproszczenia są one reduce-sumwzdłuż każdej osi, kolejność operacji może nie być optymalna ... Pamiętaj, że jeśli masz coś przeciwko błędowi precyzji, równie dobrze możesz użyćfloat64
Vitor SRG

Mogę nie mieć czasu na jutrzejszą recenzję, więc przyznanie nagrody już teraz.
Paul

@Paul W każdym razie zostanie automatycznie przyznana najwyższa głosowana odpowiedź - ale w porządku, dzięki za powiadomienie
OverLordGoldDragon

Odpowiedzi:


7

To wygląda na nieuniknioną niedokładność liczbową. Jak wyjaśniono tutaj , NumPy wykorzystuje wysoce zoptymalizowaną, dokładnie dostrojoną metodę BLAS do mnożenia macierzy . Oznacza to, że prawdopodobnie sekwencja operacji (suma i iloczyn) pomnożyła 2 macierze, zmienia się, gdy zmienia się rozmiar macierzy.

Starając się być jaśniejszymi, wiemy, że matematycznie każdy element wynikowej macierzy można obliczyć jako iloczyn kropkowy dwóch wektorów (ciągi liczb o równej długości). Ale nie w ten sposób NumPy oblicza element macierzy wynikowej. W rzeczywistości istnieją bardziej wydajne, ale złożone algorytmy, takie jak algorytm Strassena , które uzyskują ten sam wynik bez bezpośredniego obliczania iloczynu kropka-kolumna.

Podczas korzystania z takich algorytmów, nawet jeśli element C ij uzyskanej macierzy C = AB jest matematycznie zdefiniowany jako iloczyn iloczynu i-tego rzędu A z j-tą kolumną B , jeżeli pomnożymy macierz A2 o ten sam i-ty wiersz jak A z macierzą B2 mającą tę samą j-tą kolumnę co B , element C2 ij zostanie faktycznie obliczony zgodnie z inną sekwencją operacji (zależy to od całej A2 i B2 macierze), co może prowadzić do różnych błędów numerycznych.

Dlatego, nawet jeśli matematycznie C ij = C2 ij (jak w przypadek 1), różna kolejność operacji następuje przez algorytm w obliczeniach (ze względu na zmianę rozmiaru matrycy) prowadzi do różnych błędów numerycznych. Błąd numeryczny wyjaśnia również nieco inne wyniki w zależności od środowiska oraz fakt, że w niektórych przypadkach w niektórych środowiskach błąd numeryczny może być nieobecny.


2
Dzięki za link, wydaje się, że zawiera istotne informacje - twoja odpowiedź może być jednak bardziej szczegółowa, ponieważ do tej pory jest to parafrazowanie komentarzy pod pytaniem. Na przykład połączone SO pokazuje Ckod bezpośredni i zawiera wyjaśnienia na poziomie algorytmu, więc zmierza we właściwym kierunku.
OverLordGoldDragon,

Nie jest to również „nieuniknione”, jak pokazano na dole pytania - a zakres niedokładności różni się w zależności od środowiska, co pozostaje niewyjaśnione
OverLordGoldDragon

1
@OverLordGoldDragon: (1) Trywialny przykład z dodatkiem: weź numer n , weź liczbę ktaką, że jest poniżej dokładności kostatniej cyfry mantysy. Dla natywnych pływaków Pythona n = 1.0i k = 1e-16działa. Teraz pozwól ks = [k] * 100. Zobacz, że sum([n] + ks) == nchociaż sum(ks + [n]) > nkolejność sumowania ma znaczenie. (2) Nowoczesne procesory mają kilka jednostek do równoległego wykonywania operacji zmiennoprzecinkowych (FP), a kolejność ich a + b + c + dobliczania na CPU nie jest zdefiniowana, nawet jeśli polecenie a + bpojawia się wcześniej c + dw kodzie maszynowym.
9000

1
@OverLordGoldDragon Należy pamiętać, że większość liczb, z którymi program ma się zmierzyć, nie może być reprezentowana przez liczbę zmiennoprzecinkową. Spróbować format(0.01, '.30f'). Jeśli nawet prosta liczba taka jak 0.01nie może być dokładnie reprezentowana przez zmiennoprzecinkowy NumPy, nie trzeba znać głębokich szczegółów algorytmu mnożenia macierzy NumPy, aby zrozumieć punkt mojej odpowiedzi; to znaczy różne macierze początkowe prowadzą do różnych sekwencji operacji , tak że matematycznie równe wyniki mogą się różnić w niewielkim stopniu z powodu błędów numerycznych.
mmj

2
@OverLordGoldDragon re: czarna magia. Jest artykuł, który wymaga przeczytania na kilku CS MOOC. Pamiętam, że nie jest tak świetne, ale myślę, że to jest to: itu.dk/~sestoft/bachelor/IEEE754_article.pdf
Paul
Korzystając z naszej strony potwierdzasz, że przeczytałeś(-aś) i rozumiesz nasze zasady używania plików cookie i zasady ochrony prywatności.
Licensed under cc by-sa 3.0 with attribution required.