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

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.

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

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.