Точность, почему Matlab и Python numpy дают такие разные результаты?

Я знаю об основных типах данных и о том, что типы с плавающей запятой (float, double) не могут точно содержать некоторые числа.

Однако при переносе некоторого кода из Matlab в Python (Numpy) я обнаружил некоторые существенные различия в вычислениях, и я думаю, что он возвращается к точности.

Возьмите следующий код, нормализующий по оси z 500-мерный вектор, только первые два элемента которого имеют ненулевое значение.

Матлаб:

Z = repmat(0,500,1); Z(1)=3;Z(2)=1;
Za = (Z-repmat(mean(Z),500,1)) ./ repmat(std(Z),500,1);
Za(1)
>>> 21.1694

Питон:

from numpy import zeros,mean,std
Z = zeros((500,))
Z[0] = 3
Z[1] = 1
Za = (Z - mean(Z)) / std(Z)
print Za[0]
>>> 21.1905669677

Помимо того, что форматирование показывает немного больше цифр в Python, разница огромная (имхо), более 0,02

И Python, и Matlab используют 64-битный тип данных (афаик). Python использует «numpy.float64» и «двойной» Matlab.

Почему такая огромная разница? Какой из них более правильный?


person Peter Smit    schedule 20.09.2011    source источник
comment
Возможно, должно подойти для Computational Science SE, если его спросят в наши дни.   -  person gerrit    schedule 16.02.2013


Ответы (3)


Возможно, разница возникает из-за вызовов mean и std. Сравните их в первую очередь.

Есть несколько определений для std, некоторые используют квадратный корень из

1 / n * sum((xi - mean(x)) ** 2)

другие используют

1 / (n - 1) * sum((xi - mean(x)) ** 2)

вместо.

С математической точки зрения: эти формулы являются оценкой дисперсии нормально распределенной случайной величины. Распределение имеет два параметра sigma и mu. Если вы точно знаете mu, оптимальная оценка для sigma ** 2 будет

1 / n * sum((xi - mu) ** 2)

Если вам нужно оценить mu по данным, используя mu = mean(xi), оптимальная оценка для sigma**2 будет

1 / (n - 1) * sum((xi- mean(x))**2)
person rocksportrocker    schedule 20.09.2011

Отвечая на ваш вопрос, нет, это не проблема точности. Как @rocksportrocker, существуют два популярных метода оценки стандартного отклонения. В MATLAB std доступны оба варианта, но в качестве стандарта используется другой из то, что вы использовали в Python.

Попробуйте std(Z,1) вместо std(Z):

Za = (Z-repmat(mean(Z),500,1)) ./ repmat(std(Z,2),500,1);Za(1)
sprintf('%1.10f', Za(1))

приводит к

Za(1) = 21.1905669677

в МАТЛАБ. Прочитайте ответ rockpotrocker о том, какой из двух результатов больше подходит для того, что вы хотите сделать ;-).

person Jonas Heidelberg    schedule 20.09.2011
comment
А, я только что увидел, что @rocksportrocker рассказывает вам о математике по этому поводу :-). - person Jonas Heidelberg; 20.09.2011

Согласно документации std в SciPy, параметр с именем ddof:

ddof : int, необязательный
Означает разность степеней свободы. В расчетах используется делитель N - ddof, где N представляет количество элементов. По умолчанию ddof равен нулю.

В numpy ddof по умолчанию равно нулю, а в MATLAB - единице. Итак, я думаю, что это может решить проблему:

std(Z,ddof=1)
person cartoonist    schedule 03.01.2014