Эффективно вычисляйте множество внутренних продуктов в Matlab

Я работаю над проектом, для которого мне нужно вычислить множество внутренних продуктов в больших измерениях. Я знаю, что мы всегда должны пытаться векторизовать операции в матлабе, но я не уверен, как это сделать...

Допустим, у нас есть две матрицы, A и B размера N x d, где N — это количество внутренних продуктов для вычисления в d измерениях.

Это легко реализовать с помощью циклов for, но я подозреваю, что существуют более эффективные способы. Реализация цикла for может выглядеть так:

innerprods=zeros(N,1);
for i=1:N
    innerprods(i)=A(i,:)*B(i,:)';
end

У кого-нибудь есть идеи, как это векторизовать? Я предполагаю, что bsxfun должен вступить в игру в какой-то момент, но я не могу понять, как его использовать для этого...

Заранее спасибо!


matlab vectorization bsxfun
person Jeff Dun    schedule 06.06.2013    source источник

Ответы (2)


arrow_upward
2
arrow_downward

Ваше подозрение относительно bsxfun полностью оправдано! Вы можете эффективно вычислять внутренние продукты, используя следующий однострочный код с bsxfun и старый добрый оператор двоеточия:

innerprods=sum(reshape(bsxfun(@times,A(:),B(:)),N,d),2);

Для более подробного объяснения того, что здесь происходит:

    bsxfun(@times,A(:),B(:)) --> element-wise product, returns Nd x 1 vector
    reshape( ^^^ , N, d)      --> reshape into N x d matrix
    sum( ^^^ , 2)             --> sum over columns to get N x 1 vector

РЕДАКТИРОВАТЬ: я сделал несколько замеров времени, чтобы дать некоторое представление о разнице в производительности. Я использовал matlab R2010b (не спрашивайте .. .). На рисунке показаны некоторые эмпирические результаты с использованием цикла в вашем вопросе и однострочника, который я предлагаю для N=500 и различного количества измерений. Использование циклов от 2 до 8 раз медленнее, чем подход с использованием bsxfun. Разница в скорости может быть больше для более новых версий matlab за счет лучшего распараллеливания.

введите здесь описание изображения

person Marc Claesen    schedule 06.06.2013

arrow_upward
6
arrow_downward

Что уж говорить о простом:

sum(A.*B,2)

Простой тест (R2013a WinVista 32bit Core duo):

N = 8e3;
A = rand(N);
B = rand(N);

tic
out = zeros(N,1);
for ii = 1:N
    out(ii) = A(ii,:)*B(ii,:)';
end
toc

tic
out2 = sum(A.*B,2);
toc

all(out-out2 < 1e5*eps) % note difference in precision

раз

loop     5.6 sec
multsum  0.8 sec

Дополнительные тесты на R2013a Win7 64 Xeon E5

Avg time loop:           2.00906 seconds
Avg time multSum:        0.18114 seconds
Avg time bsxfun:         0.18203 seconds
Avg time reshapeMultSum: 0.18088 seconds

Основные пункты выноса:

  • зацикливание в этом случае чрезвычайно неэффективно (ожидается);
  • bsxfun() полностью избыточен, хотя накладные расходы несущественны (ожидается);
  • Преобразование в столбец, а затем обратно в матрицу, как ожидается, также принесет пользу из-за внутреннего «предпочтения» механизма MATLAB для операций по столбцам (т. е. сначала по строкам);
  • для ясности синтаксиса я по-прежнему рекомендую multSum: sum(A.*B,2).

Набор тестов (среднее время на 100 испытаний, фиксированный размер матрицы 8e3, равенство результатов 1e5*eps):

N = 8e3;
A = rand(N);
B = rand(N);

tic
for r = 1:100
    out = zeros(N,1);
    for ii = 1:N
        out(ii) = A(ii,:)*B(ii,:)';
    end
end
sprintf('Avg time loop: %.5f seconds', toc/100)

tic
for r = 1:100; out2 = sum(A.*B,2);                                  end
sprintf('Avg time multSum: %.5f seconds', toc/100)

tic
for r = 1:100; out3 = sum(reshape(bsxfun(@times,A(:),B(:)),N,N),2); end
sprintf('Avg time bsxfun: %.5f seconds', toc/100)

tic
for r = 1:100; out4 = sum(reshape(A(:).*B(:),N,N),2);               end
sprintf('Avg time reshapeMultSum: %.5f seconds', toc/100)
person Oleg    schedule 06.06.2013
comment
Спасибо за Ваш ответ! Я отметил ответ Марка, так как он лучше работает при больших значениях N и d. - person Jeff Dun; 07.06.2013
comment
@JeffDun Я немного озадачен, потому что bsxfun полностью излишен. Я постараюсь сделать больше тестов и посмотреть, откроется ли что-то интересное. - person Oleg; 07.06.2013
comment
Привет, @OlegKomarov, приятно видеть, что ты делаешь дополнительные тесты. Я согласен, что ваш подход кажется эквивалентным моему. На самом деле, ваш код мне понравился больше, чем мой, когда я его увидел. Однако по какой-то причине мои тайминги расходятся. Для больших матриц подход bsxfun оказывается до 40% быстрее. Я понятия не имею, почему, хотя. - person Marc Claesen; 07.06.2013
comment
Можете ли вы опубликовать свою архитектуру, версию Matlab и тайминги для вышеуказанного теста для сравнения? Кроме того, что вы называете большими матрицами? (Обратите внимание, что 8e3 на 8e3 имеет 64 миллиона элементов и ~ 488 МБ ОЗУ) - person Oleg; 07.06.2013