Deux principes à optimiser

Version imprimable

23/11/2012, 17h41
TheReveller

Deux principes à optimiser
Bonjour,

Je cherche à rendre plus rapide deux trucs.

Premièrement, j'aimerais pouvoir générer aussi rapidement que possible un nombre aléatoire X entier de 1 à M excluant N.

Présentement, je le fais ainsi :
Code:

1 2 choix = [1:N-1,N+1:M]; X = choix(ceil(rand()*(M-1)));
Sauriez-vous faire plus rapide ? Cela prend 4e-6 secondes sur ma machine.

Deuxièmement, un autre cas un peu plus important. J'ai des N valeurs et je leur attribue une probabilité d'être choisie à l'aide de la relation N_i / sum(N_i). J'aimerais choisir à nouveau N valeurs, mais cette fois selon les probabilités. Donc, une valeur pourrait ne pas être choisie et une autre pourrait être choisie plusieurs fois, selon sa probabilité, en autant que la longueur finale des valeurs choisies est N.

Exemple :

Code:

[1, 10, 3, 1, 15]

La somme est 30, les probabilités sont :

Code:

[1/30, 1/3, 1/10, 1/30, 1/2]

Le vecteur résultant pourrait être :

Code:

[15, 10, 15, 10, 15]

Merci.

PS : Je posterai mon algorithme pour ce cas prochainement. Dans mon algorithme, c'est mon utilisation de repmat qui pose problème.

Pour le second problème, je te propose la fonction toto :

Code:

1
2
3
4
5
6
7
8
9
10
11
12
function r = toto(x)
    arrP = [1/30, 1/3, 1/10, 1/30, 1/2];
    arrVal = [1, 10, 3, 1, 15];
    if x<=arrP(1)
        r = arrVal(1);
    elseif x<=sum(arrP(1:2))
        r = arrVal(2);
    elseif x<=sum(arrP(1:3))
        r = arrVal(3);    
        % continuer jusqu'a 5
    end
end

toto(rand) devrait te renvoyer ce que tu veux.

23/11/2012, 18h58
Invité

Bonsoir,

Pour le premier, randi serait je pense utile.
Pour le deuxième, si tu as la Communications System Toolbox, tu peux utiliser la fonction randsrc.

Merci FLB pour la proposition. Malheureusement, j'avais oublié de mentionner que mon cas réel implique un vecteur de longueur variable qui pourrait même avoir une longueur de 100 éléments, voire plus.

Merci Winjerome pour l'information, je ne connaissais pas randsrc. Malheureusement, cette fonction n'est pas optimisée pour mon cas. Elle est fonctionnelle, mais elle est plus lourde qu'elle ne le paraît.

Effectivement, comme pour l'utilisation de randi, il faut parfois évaluer les performances de certains outils et fonctions que MATLAB nous propose lorsqu'on veut optimiser en vitesse. MATLAB est un outil extrêmement puissant pour l'utilisation de matrices, mais on doit éviter autant que possible les boucles for. Aussi, certaines fonctions MATLAB possèdent tellement d'options et de vérifications qu'elles peuvent devenir lourdes. Toutes ces options les obligent parfois à passer par des boucles for.

D'ailleurs, je viens déjà d'analyser l'utilisation de la fonction repmat et je préfère utiliser un produit matriciel.

Voici quelques petits exemples de ce que j'avance :
Code:

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 clc; close all; clear all; pause(1); tic; for i=1:1e6 randi([1,100]); end toc; clear all; pause(1); tic; for i=1:1e6 ceil(rand()*100); end toc; clear all; vecteur = 1:100; pause(1); tic; for i=1:1e6 repmat(vecteur,100,1); end toc; clear all; vecteur = 1:100; pause(1); tic; for i=1:1e6 ones(100,1)*vecteur; end toc;
Code:

1 2 3 4 Elapsed time is 1.489855 seconds. Elapsed time is 0.072176 seconds. Elapsed time is 28.582096 seconds. Elapsed time is 14.823746 seconds.
Bref, on s'aperçoit dans cet exemple que générer un nombre entier aléatoirement de 1 à 100 prends 20 fois moins de temps. Aussi, copier un vecteur ligne de 100 éléments en une matrice de 100 lignes prend 2 fois moins de temps.

Voici ce que j'ai comme algorithme jusqu'à maintenant :
Code:

1 2 3 4 5 6 7 8 9 10 11 clear all; close all; clc; vecteur = [1, 10, 3, 1, 15]'; longueur = numel(vecteur); probs = vecteur/sum(vecteur); pickprobs = ones(longueur,1) * (cumsum(probs) - probs)'; randpick = rand(longueur,1) * ones(1,longueur); vecteur(sum(randpick >= pickprobs,2))
PS : Les transposées sont uniquement pour mieux visualiser mes besoins personnels dans l'algorithme réel entier.

Merci de votre aide,

Éric

24/11/2012, 01h51
TheReveller

Pour le premier, je viens de trouver une meilleure façon qui me permet de diminuer le temps par un facteur 4 ou plus. Au lieu d'une moyenne de 4e-6 seconde par itération, j'ai une moyenne de 1e-6 seconde par itération.

Voici l'astuce toute bête :

Code:

X = mod(N+ceil(rand()*(M-1))-1,M)+1;

force brute

Salut,
je trouve que la force brute est plus rapide :
Code:

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 clc % ta méthode clear M = 100 ; N = 25 ; tic for i=1:1e6 X = mod(N+ceil(rand()*(M-1))-1,M)+1 ; end toc % Elapsed time is 0.170726 seconds. % méthode brutale essai-erreur tic clear M = 100 ; N = 25 ; Mm = M-1 ; for i=1:1e6 X=N ; while X==N, X = ceil(rand()*Mm) ; end end toc % Elapsed time is 0.097532 seconds.
Note : effacer la mémoire entre les essais est nécessaire pour les rendre comparables.

Si la vitesse est critique, pourquoi ne pas pré-calculer tes nombres aléatoires ? C'est à peu près aussi rapide et, surtout, ce temps de calcul peut être externalisé (mis hors de ta section critique).
Code:

1 2 3 4 5 6 7 8 9 10 % méthode par lot clear M = 100 ; N = 25 ; tic for i=1:1e4 X = ceil(rand(1,100).*M) ; X = X(X ~= N) ; end toc % Elapsed time is 0.131254 seconds.

04/12/2012, 21h09
TheReveller

Merci pour la proposition.

J'aurais probablement dû expliquer mon cas réel afin d'illustrer tous les petits détails qui pourraient accrocher dans une manipulation matricielle.

Je dois en fait modifier les valeurs de plusieurs vecteurs de longueur N. Dans mon cas réel, N ne depasse pas 101, mais on pourrait aller jusqu'à 501 selon la limite que j'ai jugée pour mon code.

Disons que tous les M vecteurs de longueur N contiennent des valeurs entières aléatoires. La borne inférieure de chacune des valeurs des vecteurs est 1 tandis que les bornes supérieures pour chacune des valeurs sont données par les éléments du vecteur C contenant donc N bornes.

Je veux modifier aléatoirement chacune des valeurs de chacun des vecteurs en respectant lleurs bornes supérieures et en ne reprenant pas la même valeur.

J'espère être clair, sinon je m'expliquerai mieux dès que j'aurai du temps et que je serai sur un ordinateur.

Merci,

Éric

--------------------------------------------------------------------------

PS : Dans mon astuce toute bête, N peut être un vecteur de données et M un vecteur de bornes supérieures et le vecteur X retourné sera valable.

X = mod(N+ceil(rand()*(M-1))-1,M)+1;
04/12/2012, 22h45
ol9245
en effet je ne comprend pas tout :aie:
- ou est C dans ta formule ?
- si tu dois tirer M vecteurs de longueur fixe 101, pourquoi je ne vois pas quelque-part un rand (M, 101) ?
disons quelque chose qui s'inspire de ça :
Code:

1 2 3 M = [10 100 1000] ; C = ceil(bsxfun(@times, rand(101, length(M)), M)) ;
04/12/2012, 23h27
TheReveller

Bon, j'ai maintenant un ordinateur, alors je m'explique et en même temps j'apporte quelques remarques et détails de plus.

Premièrement, le problème de la "force brute" de ta proposition est qu'elle est probabiliste. Je parlais dans mon exemple d'un nombre borné par [1,100] généré aléatoirement différent de 25. Il n'y a donc que 1% de chances de retomber sur ce même nombre. Cependant, si mon exemple était plutôt avec un nombre borné par [1,3] et différent de 2, alors il y aurait 33% de chances de retomber sur ce même nombre.

Bon, oublions tous les noms de variables mentionnés précédemment.

J'ai des vecteurs de longueur numel(N), soient N, NMAX et X.

Par exemple : N = [2,3,1,2,4] et NMAX = [2,5,1,3,8]Je veux obtenir un vecteur X de valeurs entières aléatoires tel que ses valeurs ne dépassent pas les valeurs NMAX correspondantes (element-wise) et soient différentes des valeurs correspondantes de N.

Une des valeurs maximales est 1, c'est le seul cas où la valeur aléatoire sera reprise puisqu'il n'y a qu'une seule possibilité : 1.

Bref X <= NMAX et X(NMAX ~= 1) ~= N.
Nous aurions par exemple, une réponse possible pour ce cas :

Code:

X = [1,4,1,3,6]

Maintenant, il faut faire en sorte que pour un cas où numel(N) est très grand (par exemple 101) et qu'on doit répéter cet exercice un grand nombre de fois (par exemple 1e6 fois), le temps soit minimal.

Merci,

Éric

sans spécialement chercher la performance, c'est comme ça que je ferais.
après, si tu veux vraiment aller vite, il te reste le C dans un mex-file.
Code:

1 2 3 4 5 6 L = 101 ; N = [2,3,1,2,4] ; NMAX = [2,5,1,3,8] ; X = ceil(bsxfun(@times, rand(L, length(M)), NMAX-1)) ; ajouter_un = (X>=repmat(N,L,1)) | (X==0) ; X(ajouter_un) = X(ajouter_un) + 1 ;

Une meilleur solution concernant le second problème :

Code:

1
2
3
4
 
arrP = [1/30, 1/3, 1/10, 1/30, 1/2];
arrVal = [1, 10, 3, 1, 15];
arrVal(min(find(rand<=cumsum(arrP)))) %renvoi une des valeurs de arrVal avec une probabilité arrP

05/12/2012, 16h53
TheReveller

Bonjour,

Je n'ai pas pu tester cette proposition, mais est-ce que le vecteur renvoyé est de meme longueur ?

J'ai un vecteur d'une certaine longueur contenant des valeurs auxquelles je calcule les probabilités correspondantes et je veux obtenir un vecteur de même longueur contenant des valeurs choisies selon les probabilités.

Merci,

Éric

Code:

1 2 3 4 5 6 arrP = [1/30, 1/3, 1/10, 1/30, 1/2]; arrVal = [1, 10, 3, 1, 15]; res = zeros(size(arrVal)); for k = 1:length(arrVal) res(k) = arrVal(min(find(rand<=cumsum(arrP)))); end
Maintenant oui ;)
Il est peut être possible de vectoriser tout ça, mais le code va commencer à devenir illisible...

Bonjour,

Merci pour vos suggestions. Malheureusement, je ne peux les tester actuellement puisque mon MATLAB roule un code pendant 2x 5h.

J'ai fait ce bout de code pour comparer. J'espère ne pas avoir fait d'erreur dans mon template puisque je l'ai fait dans le Bloc-Notes.

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
close all;
clc;
 
% Début - Initialisation (rien à modifier)
clear all;
L = 101; % Longueur désirée
NMAX = ceil(rand(1,L)*8); % NMAX aléatoires de 1 à 8
overallmax = max(NMAX);
N = mod(ceil(rand(1,L)*overallmax)-1,NMAX)+1; % Vecteur N valide selon NMAX
pause(1);
% Fin - Initialisation (rien à modifier)
 
% Début - Comparaison (ma solution, rien à modifier)
tic;
for i=1:1e5
	X = mod(N+ceil(rand(1,L).*(NMAX-1))-1,NMAX)+1;
	if any(X < 1) || any(N > NMAX) || any(X(NMAX ~= 1) == N(NMAX ~= 1))
		disp('Solution invalide');
	end
end
toc;
% Fin - Comparaison (ma solution, rien à modifier)
 
% Début - Initialisation (rien à modifier)
clear all;
L = 101; % Longueur désirée
NMAX = ceil(rand(1,L)*8); % NMAX aléatoires de 1 à 8
overallmax = max(NMAX);
N = mod(ceil(rand([1,L])*overallmax)-1,NMAX)+1; % Vecteur N valide selon NMAX
pause(1);
% Fin - Initialisation (rien à modifier)
 
% Début - Comparaison (votre solution)
tic;
for i=1:1e5
 
	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
	% Insérez ici votre suggestion %
	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
	if any(X < 1) || any(N > NMAX) || any(X(NMAX ~= 1) == N(NMAX ~= 1))
		disp('Solution invalide');
	end
end
toc;
% Fin - Comparaison (votre solution)

Merci de votre aide,

Éric

PS : S'il y a un problème avec ce code ou si ma solution donne des solutions invalides, je vous reviens là-dessus plus tard parce que je sais que c'est normalement fonctionnel.

13/12/2012, 19h33
ol9245

quand je dois faire toutner matlab pendant des heures, je compile le code et je le fais tourner à côté. Comme ça j'ai toujours Matlab dispo pour travailler sur la prochaine étape.