Dans Matlab, quand est-il optimal d'utiliser bsxfun?

Ma question: J'ai remarqué que beaucoup de bonnes réponses aux questions Matlab sur SO utilisent fréquemment la fonction bsxfun. Pourquoi?

Motivation: Dans la documentation Matlab pour bsxfun, l'exemple suivant est fourni:

A = magic(5);
A = bsxfun(@minus, A, mean(A))

Bien sûr, nous pourrions faire la même opération en utilisant:

A = A - (ones(size(A, 1), 1) * mean(A));

Et en fait, un simple test de vitesse démontre que la deuxième méthode est environ 20% plus rapide. Alors pourquoi utiliser la première méthode? Je suppose qu'il y a des circonstances où l'utilisation bsxfunsera beaucoup plus rapide que l'approche "manuelle". Je serais vraiment intéressé de voir un exemple d'une telle situation et une explication pour expliquer pourquoi c'est plus rapide.

Aussi, un dernier élément à cette question, toujours de la documentation Matlab pour bsxfun: "C = bsxfun (fun, A, B) applique l'opération binaire élément par élément spécifiée par la fonction handle fun aux tableaux A et B, avec singleton extension activée. ". Que signifie l'expression «avec l'expansion singleton activée»?

J'utilise trois raisons bsxfun( documentation , lien de blog )

1. bsxfunest plus rapide que repmat(voir ci-dessous)
2. bsxfun nécessite moins de frappe
3. Utiliser bsxfun, comme utiliser accumarray, me fait me sentir bien dans ma compréhension de Matlab.

bsxfunrépliquera les tableaux d'entrée selon leurs "dimensions singleton", c'est-à-dire les dimensions le long desquelles la taille du tableau est 1, de sorte qu'ils correspondent à la taille de la dimension correspondante de l'autre tableau. C'est ce que l'on appelle "expasion singleton". En passant, les dimensions singleton sont celles qui seront supprimées si vous appelez squeeze.

Il est possible que pour de très petits problèmes, l' repmatapproche soit plus rapide - mais à cette taille de baie, les deux opérations sont si rapides qu'elles ne feront probablement aucune différence en termes de performances globales. Il y a deux raisons importantes bsxfunest plus rapide: (1) le calcul se produit dans le code compilé, ce qui signifie que la réplication réelle du tableau ne se produit jamais, et (2) bsxfunest l'une des fonctions Matlab multithread.

J'ai effectué une comparaison de vitesse entre repmatet bsxfunavec R2012b sur mon ordinateur portable décemment rapide.

Pour moi, bsxfunc'est environ 3 fois plus rapide que repmat. La différence devient plus prononcée si les tableaux s'agrandissent

Le saut dans l'exécution de repmatse produit autour d'une taille de tableau de 1 Mo, ce qui pourrait avoir quelque chose à voir avec la taille de mon cache de processeur - bsxfunne devient pas aussi mauvais qu'un saut, car il n'a besoin que d'allouer le tableau de sortie.

Vous trouverez ci-dessous le code que j'ai utilisé pour le chronométrage:

n = 300;
k=1; %# k=100 for the second graph
a = ones(10,1);
rr = zeros(n,1);
bb=zeros(n,1);
ntt=100;
tt=zeros(ntt,1);
for i=1:n;
   r = rand(1,i*k);
   for it=1:ntt;
      tic,
      x=bsxfun(@plus,a,r);
      tt(it)=toc;
   end;
   bb(i)=median(tt);
   for it=1:ntt;
      tic,
      y=repmat(a,1,i*k)+repmat(r,10,1);
      tt(it)=toc;
   end;
   rr(i)=median(tt);
end

Dans mon cas, je l'utilise bsxfuncar cela m'évite de penser aux problèmes de colonne ou de ligne.

Afin d'écrire votre exemple:

A = A - (ones(size(A, 1), 1) * mean(A));

Je dois résoudre plusieurs problèmes:

1) size(A,1)ousize(A,2)

2) ones(sizes(A,1),1)ouones(1,sizes(A,1))

3) ones(size(A, 1), 1) * mean(A)oumean(A)*ones(size(A, 1), 1)

4) mean(A)oumean(A,2)

Quand j'utilise bsxfun, je dois juste résoudre le dernier:

a) mean(A)oumean(A,2)

Vous pourriez penser que c'est paresseux ou quelque chose du genre, mais quand j'utilise bsxfun, j'ai moins de bogues et je programme plus rapidement .

De plus, il est plus court, ce qui améliore la vitesse de frappe et la lisibilité .

Question très intéressante! Je suis récemment tombé sur exactement une telle situation en répondant à cette question. Considérez le code suivant qui calcule les indices d'une fenêtre glissante de taille 3 via un vecteur a:

a = rand(1e7,1);

tic;
idx = bsxfun(@plus, [0:2]', 1:numel(a)-2);
toc

% equivalent code from im2col function in MATLAB
tic;
idx0 = repmat([0:2]', 1, numel(a)-2);
idx1 = repmat(1:numel(a)-2, 3, 1);
idx2 = idx0+idx1;
toc;

isequal(idx, idx2)

Elapsed time is 0.297987 seconds.
Elapsed time is 0.501047 seconds.

ans =

 1

Dans ce cas, bsxfunc'est presque deux fois plus rapide! Il est utile et rapide car il évite l'allocation explicite de mémoire pour les matrices idx0et idx1, les enregistrer dans la mémoire, puis les relire juste pour les ajouter. Étant donné que la bande passante mémoire est un atout précieux et souvent le goulot d'étranglement sur les architectures d'aujourd'hui, vous voulez l'utiliser à bon escient et réduire les besoins en mémoire de votre code pour améliorer les performances.

bsxfunvous permet de faire exactement cela: créer une matrice basée sur l'application d'un opérateur arbitraire à toutes les paires d'éléments de deux vecteurs, au lieu d'opérer explicitement sur deux matrices obtenues en répliquant les vecteurs. C'est une expansion singleton . Vous pouvez également le considérer comme le produit extérieur de BLAS:

v1=[0:2]';
v2 = 1:numel(a)-2;
tic;
vout = v1*v2;
toc
Elapsed time is 0.309763 seconds.

Vous multipliez deux vecteurs pour obtenir une matrice. Juste que le produit externe n'effectue que des multiplications et bsxfunpeut appliquer des opérateurs arbitraires. En remarque, il est très intéressant de voir que bsxfunc'est aussi rapide que le produit extérieur BLAS. Et BLAS est généralement considéré comme offrant la performance.

Edit Grâce au commentaire de Dan, voici un excellent article de Loren discutant exactement de cela.

À partir de R2016b, Matlab prend en charge l' expansion implicite pour une grande variété d'opérateurs, de sorte que dans la plupart des cas, il n'est plus nécessaire d'utiliser bsxfun:

Auparavant, cette fonctionnalité était disponible via la bsxfunfonction. Il est maintenant recommandé de remplacer la plupart des utilisations de bsxfunpar des appels directs aux fonctions et aux opérateurs qui prennent en charge l' expansion implicite . Par rapport à l'utilisation bsxfun, l' expansion implicite offre une vitesse plus rapide , une meilleure utilisation de la mémoire et une meilleure lisibilité du code .

Il y a une discussion détaillée de l' expansion implicite et de ses performances sur le blog de Loren. Pour citer Steve Eddins de MathWorks:

Dans R2016b, l' expansion implicite fonctionne aussi vite ou plus rapidement que bsxfundans la plupart des cas. Les meilleurs gains de performances pour l' expansion implicite sont avec de petites tailles de matrice et de tableau. Pour de grandes tailles de matrice, l'expansion implicite a tendance à être à peu près la même vitesse que bsxfun.

Les choses ne sont pas toujours cohérentes avec les 3 méthodes courantes repmat:, extension par indexation, et bsxfun. Cela devient plus intéressant lorsque vous augmentez encore la taille du vecteur. Voir l'intrigue:

bsxfundevient en fait légèrement plus lent que les deux autres à un moment donné, mais ce qui m'a surpris, c'est que si vous augmentez encore plus la taille du vecteur (éléments de sortie> 13E6), bsxfun redevient soudainement plus rapide d'environ 3x. Leurs vitesses semblent sauter par étapes et l'ordre n'est pas toujours cohérent. Je suppose que cela pourrait également dépendre de la taille du processeur / de la mémoire, mais en général, je pense que je m'en tiendrai bsxfunautant que possible.