Почему здесь происходит клонирование индексации?
Несколько лет назад кто-то опубликовал в Active State Recipes для целей сравнения, три функции python/NumPy; каждый из них принимал те же аргументы и возвращал тот же результат, матрицу расстояний.
Два из них были взяты из опубликованных источников; они оба - или они кажутся мне - идиоматическим кодом. Повторяющиеся вычисления, необходимые для создания матрицы расстояний, управляются синтаксисом элегантного индекса numpy. Здесь один из них:
from numpy.matlib import repmat, repeat
def calcDistanceMatrixFastEuclidean(points):
numPoints = len(points)
distMat = sqrt(sum((repmat(points, numPoints, 1) -
repeat(points, numPoints, axis=0))**2, axis=1))
return distMat.reshape((numPoints,numPoints))
Третий создал матрицу расстояний, используя один цикл (который, очевидно, много циклов, учитывая, что матрица расстояний всего в 1000 двумерных точек имеет один миллион записей). На первый взгляд эта функция выглядела как код, который я использовал для записи, когда я изучал NumPy, и я бы написал код NumPy, сначала написав код Python, а затем переведя его, построчно.
Спустя несколько месяцев после публикации Active State результаты тестов производительности, сравнивающих три, были опубликованы и обсуждены в thread в списке рассылки NumPy.
Функция с циклом фактически значительно превзошла два других:
from numpy import mat, zeros, newaxis
def calcDistanceMatrixFastEuclidean2(nDimPoints):
nDimPoints = array(nDimPoints)
n,m = nDimPoints.shape
delta = zeros((n,n),'d')
for d in xrange(m):
data = nDimPoints[:,d]
delta += (data - data[:,newaxis])**2
return sqrt(delta)
Один участник в теме (Keir Mierle) предложил причину, почему это может быть правдой:
Причина, по которой я подозреваю, что это будет быстрее что он имеет лучшую локальность, полностью заканчивая вычисление на относительно небольшой рабочий набор, прежде чем перейти к следующему. Один вкладыш придется многократно вытягивать потенциально большой массив MxN в процессор.
Под этим плакатом собственный счет, его замечание является лишь подозрением, и, похоже, оно не обсуждалось.
Любые другие мысли о том, как объяснить эти результаты?
В частности, существует ли полезное правило - относительно того, когда цикл и когда индексировать - который можно извлечь из этого примера в качестве руководства при написании numpy-кода?
Для тех, кто не знаком с NumPy или кто не смотрел код, это сравнение не основано на случае с краем - это, конечно, было бы не так интересно, если бы оно было. Вместо этого это сравнение включает функцию, которая выполняет общую задачу при вычислении матрицы (т.е. Создание массива результатов с учетом двух предшествующих); кроме того, каждая функция, в свою очередь, состоит из наиболее распространенных встроенных модулей numpy.
Ответы
Ответ 1
TL; DR Второй код, приведенный выше, представляет собой только цикл по количеству измерений точек (3 раза по циклу for для трехмерных точек), поэтому петли не так много. Реальное ускорение во втором коде выше - это то, что он лучше использует силу Numpy, чтобы избежать создания некоторых дополнительных матриц при поиске различий между точками. Это уменьшает использование памяти и вычислительные усилия.
Дольше Объяснение
Я думаю, что функция calcDistanceMatrixFastEuclidean2 обманывает вас, возможно, своим циклом. Это только цикл по числу измерений точек. Для 1D точек цикл выполняется один раз, для 2D, дважды, и для 3D - трижды. Это действительно не так много циклов.
Проанализируйте код немного, чтобы понять, почему он быстрее, чем другой. calcDistanceMatrixFastEuclidean Я буду называть fast1, а calcDistanceMatrixFastEuclidean2 будет fast2.
fast1 основан на способе Matlab делать вещи, о чем свидетельствует функция repmap. Функция repmap создает массив в этом случае, это только оригинальные данные, повторяющиеся снова и снова. Однако, если вы посмотрите на код функции, это очень неэффективно. Для этого используется множество функций Numpy (3 reshape и 2 repeat s). Функция repeat также используется для создания массива, который содержит исходные данные, причем каждый элемент данных повторяется много раз. Если наши входные данные [1,2,3], то мы вычитаем [1,2,3,1,2,3,1,2,3] из [1,1,1,2,2,2,3,3,3]. Numpy пришлось создавать множество дополнительных матриц между запущенным кодом Numpy C, которого можно было бы избежать.
fast2 использует больше тяжелого подъема Numpy, не создавая столько матриц между вызовами Numpy. fast2 проходит через каждый размер точек, выполняет ли вычитание и сохраняет общее количество квадратов различий между каждым измерением. Только в конце делается квадратный корень. Пока что это может показаться не столь эффективным, как fast1, но fast2 избегает делать вещи repmat, используя индексирование Numpy. Рассмотрим простоту для случая 1D. fast2 создает 1D-массив данных и вычитает его из массива 2D (N x 1) данных. Это создает разностную матрицу между каждой точкой и всеми другими точками без использования repmat и repeat и тем самым обходит создание множества дополнительных массивов. Вот где, по моему мнению, реальная разница в скорости. fast1 создает много лишних промежутков между матрицами (и они создаются дорого вычислительно), чтобы найти различия между точками, а fast2 лучше использует силу Numpy, чтобы избежать этого.
Кстати, вот немного более быстрая версия fast2:
def calcDistanceMatrixFastEuclidean3(nDimPoints):
nDimPoints = array(nDimPoints)
n,m = nDimPoints.shape
data = nDimPoints[:,0]
delta = (data - data[:,newaxis])**2
for d in xrange(1,m):
data = nDimPoints[:,d]
delta += (data - data[:,newaxis])**2
return sqrt(delta)
Разница в том, что мы больше не создаем delta как матрицу нулей.
Ответ 2
dis для удовольствия:
dis.dis(calcDistanceMatrixFastEuclidean)
2 0 LOAD_GLOBAL 0 (len)
3 LOAD_FAST 0 (points)
6 CALL_FUNCTION 1
9 STORE_FAST 1 (numPoints)
3 12 LOAD_GLOBAL 1 (sqrt)
15 LOAD_GLOBAL 2 (sum)
18 LOAD_GLOBAL 3 (repmat)
21 LOAD_FAST 0 (points)
24 LOAD_FAST 1 (numPoints)
27 LOAD_CONST 1 (1)
30 CALL_FUNCTION 3
4 33 LOAD_GLOBAL 4 (repeat)
36 LOAD_FAST 0 (points)
39 LOAD_FAST 1 (numPoints)
42 LOAD_CONST 2 ('axis')
45 LOAD_CONST 3 (0)
48 CALL_FUNCTION 258
51 BINARY_SUBTRACT
52 LOAD_CONST 4 (2)
55 BINARY_POWER
56 LOAD_CONST 2 ('axis')
59 LOAD_CONST 1 (1)
62 CALL_FUNCTION 257
65 CALL_FUNCTION 1
68 STORE_FAST 2 (distMat)
5 71 LOAD_FAST 2 (distMat)
74 LOAD_ATTR 5 (reshape)
77 LOAD_FAST 1 (numPoints)
80 LOAD_FAST 1 (numPoints)
83 BUILD_TUPLE 2
86 CALL_FUNCTION 1
89 RETURN_VALUE
dis.dis(calcDistanceMatrixFastEuclidean2)
2 0 LOAD_GLOBAL 0 (array)
3 LOAD_FAST 0 (nDimPoints)
6 CALL_FUNCTION 1
9 STORE_FAST 0 (nDimPoints)
3 12 LOAD_FAST 0 (nDimPoints)
15 LOAD_ATTR 1 (shape)
18 UNPACK_SEQUENCE 2
21 STORE_FAST 1 (n)
24 STORE_FAST 2 (m)
4 27 LOAD_GLOBAL 2 (zeros)
30 LOAD_FAST 1 (n)
33 LOAD_FAST 1 (n)
36 BUILD_TUPLE 2
39 LOAD_CONST 1 ('d')
42 CALL_FUNCTION 2
45 STORE_FAST 3 (delta)
5 48 SETUP_LOOP 76 (to 127)
51 LOAD_GLOBAL 3 (xrange)
54 LOAD_FAST 2 (m)
57 CALL_FUNCTION 1
60 GET_ITER
>> 61 FOR_ITER 62 (to 126)
64 STORE_FAST 4 (d)
6 67 LOAD_FAST 0 (nDimPoints)
70 LOAD_CONST 0 (None)
73 LOAD_CONST 0 (None)
76 BUILD_SLICE 2
79 LOAD_FAST 4 (d)
82 BUILD_TUPLE 2
85 BINARY_SUBSCR
86 STORE_FAST 5 (data)
7 89 LOAD_FAST 3 (delta)
92 LOAD_FAST 5 (data)
95 LOAD_FAST 5 (data)
98 LOAD_CONST 0 (None)
101 LOAD_CONST 0 (None)
104 BUILD_SLICE 2
107 LOAD_GLOBAL 4 (newaxis)
110 BUILD_TUPLE 2
113 BINARY_SUBSCR
114 BINARY_SUBTRACT
115 LOAD_CONST 2 (2)
118 BINARY_POWER
119 INPLACE_ADD
120 STORE_FAST 3 (delta)
123 JUMP_ABSOLUTE 61
>> 126 POP_BLOCK
8 >> 127 LOAD_GLOBAL 5 (sqrt)
130 LOAD_FAST 3 (delta)
133 CALL_FUNCTION 1
136 RETURN_VALUE
Я не эксперт по dis, но, похоже, вам придется больше смотреть на функции, которые первый зовет, чтобы узнать, почему они занимают некоторое время. Также есть инструмент профилирования производительности с Python, cProfile.
