Почему np.dot неточен? (n-тусклые массивы)

Допустим, мы берем np.dot из двух 'float32' 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]

Числа. Кроме того, они могут измениться:

СЛУЧАЙ 1: срезa

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]

Результаты отличаются, даже если напечатанный фрагмент получен из одинаковых умноженных чисел.

СЛУЧАЙ 2: сплющитьaвозьмем 1D версию b, затем нарезать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]

СЛУЧАЙ 3: усиление контроля; установить для всех не задействованных входов ноль: добавитьa[1:] = 0к коду CASE 1. Результат: расхождения сохраняются.

СЛУЧАЙ 4: проверьте индексы, отличные от[0]; как для[0], результаты начинают стабилизировать фиксированное количество увеличений массива с момента их создания. Выход

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

Следовательно, для случая 2D * 2D результаты различаются, но согласуются для случая 1D * 1D. Судя по некоторым из моих чтений, это, похоже, происходит от 1D-1D, использующего простое добавление, тогда как 2D-2D использует "более причудливое" добавление, повышающее производительность, которое может быть менее точным (например, попарное сложение дает противоположное). Тем не менее, я не могу понять, почему расхождения исчезают, если 1 разaпроходит мимо установленного "порога"; большийa а также b, чем позже этот порог кажется ложным, но он существует всегда.

Все сказали: почему np.dotнеточные (и непоследовательные) для массивов ND-ND? Соответствующий Git

Дополнительная информация:

  • Среда: ОС Win-10, Python 3.7.4, IDE Spyder 3.3.6, Anaconda 3.0 2019/10
  • Процессор: i7-7700HQ 2,8 ГГц
  • Numpy v1.16.5

Возможная библиотека виновников: Numpy MKL - также библиотеки BLASS; спасибо Bi Rico за то, что отметили

Код стресс-теста: как уже отмечалось, расхождения усугубляются по частоте с большими массивами; если выше не воспроизводится, то должно быть ниже (если нет, попробуйте более крупные тусклости). Мой результат

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

Серьезность проблемы: показанные расхождения "небольшие", но уже не так при работе в нейронной сети с миллиардами чисел, умноженными за несколько секунд, и триллионами за все время выполнения; заявленная точность модели отличается на целые десятки процентов в этой ветке.

Ниже приведена гифка с массивами, полученными в результате подачи в модель того, что в основном a[0], ж / len(a)==1 vs. len(a)==32:

Результаты ДРУГИХ ПЛАТФОРМ, согласно и благодаря тестированию Paul:

Случай 1 воспроизведен (частично):

  • Google Colab VM - Intel Xeon 2.3 G-Hz - Jupyter - Python 3.6.8
  • Настольный ПК Docker Win-10 Pro - 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

Примечание: они дают гораздо меньшую ошибку, чем показано выше; две записи в первой строке отличаются на 1 в младшем разряде от соответствующих записей в других строках.

Случай 1 не воспроизведен:

  • Ubuntu 18.04.3 LTS - Intel i7-8700K - IPython 5.5.0 - Python 2.7.15+ и 3.6.8 (2 теста)
  • 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

Примечания:

  • В связанных Colab ноутбуков и jupyter сред показывают гораздо меньшее расхождение (и только для первых двух строк), чем наблюдается в моей системе. Кроме того, случай 2 никогда (пока) не показал неточности.
  • В этом очень ограниченном примере текущая (Dockerized) среда Jupyter более восприимчива, чем среда IPython.
  • np.show_config()слишком долго для публикации, но в итоге: env IPython основаны на BLAS/LAPACK; Colab основан на OpenBLAS. В среде IPython Linux библиотеки BLAS устанавливаются системой - в Jupyter и Colab они берутся из /opt/conda/lib.

ОБНОВЛЕНИЕ: принятый ответ точный, но общий и неполный. Вопрос остается открытым для всех, кто может объяснить поведение на уровне кода, а именно точный алгоритм, используемыйnp.dot, и как это объясняет "постоянные несоответствия", наблюдаемые в приведенных выше результатах (также см. комментарии). Вот несколько прямых реализаций, выходящих за рамки моего расшифровки: sdot.c - arraytypes.c.src

Источник

1 ответ

Решение

Это похоже на неизбежную неточность числовых значений. Как объясняется здесь, NumPy использует оптимизированный, тщательно настроенный метод BLAS для умножения матриц. Это означает, что, вероятно, последовательность операций (сумма и произведения), выполняемых для умножения двух матриц, изменяется при изменении размера матрицы.

Мы знаем, что математически каждый элемент результирующей матрицы может быть вычислен как скалярное произведение двух векторов (последовательности чисел равной длины). Но это не то, как NumPy вычисляет элемент результирующей матрицы. На самом деле существуют более эффективные, но сложные алгоритмы, такие как алгоритм Штрассена, которые получают тот же результат без прямого вычисления скалярного произведения строка-столбец.

При использовании таких алгоритмов, даже если элемент C ij результирующей матрицы C = AB математически определен как скалярное произведение i-й строки A на j-й столбец B, если вы умножите матрицу A2, имеющую та же i-я строка, что и A с матрицей B2, имеющей тот же j-й столбец, что и B, элемент C2 ij будет фактически вычисляться после другой последовательности операций (которая зависит от всего A2 и B2 матрицы), что может приводить к различным численным ошибкам.

Вот почему, даже если математически C ij = C2 ij (как в вашем СЛУЧАЕ 1), другая последовательность операций, которой следует алгоритм в вычислениях (из-за изменения размера матрицы), приводит к различным численным ошибкам. Числовая ошибка объясняет также несколько разные результаты в зависимости от среды и тот факт, что в некоторых случаях для некоторых сред числовая ошибка может отсутствовать.

Источник
Другие вопросы по тегам