Can I avoid double forcycle by bsxfun?

I am wondering, if it is possible to use double bsxfun, or something similar. I have this piece of code:

N = 5;
prob = [0.1 0.2 0 0.1; 0 0 0.05 0.1; 0.2 0.2 0 0.1];
r = rand(size(prob,1),N);

P = zeros(N, size(prob,1));
csm = cumsum(normP([zeros(size(prob,1),1),prob]),2);

for i = 1:N
 P(i,:) = sum(bsxfun(@ge,r(:,i),csm),2);
end

Prob matrix contains rows with elements between 0 and 1, each row is probability distribution (after normalisation done by normP). First row of prob matrix is used to generate first element of vector P (values1,2 or 4), second row for the second element (values 3 or 4), and so on.

e.g.: P =

       2     4     2
       4     4     1
       2     4     2
       2     4     4
       2     3     1

I have already vectorised generating elements for one vector P, but I need to generate several (N) vectors. There must be way to avoid the for cycle.

In the attachment, there is the normP function. I will be glad for help, thank you.

Michal Roubalik

P.S. Code of normP is here:

function nP = normP(P)
%
% Probability matrix normalization to hold sum of rows be equal to one
% i.e. sum(nP,2) = ones(N,1)
%
% No dependencies
srowP = sum(P,2);
good  = srowP>0;
bad   = ~good;
nP    = zeros(size(P));

% good case
if any(good)
    nP(good,:) = bsxfun(@rdivide, P(good,:), srowP(good));
end
% bad case
if any(bad)
    nP(bad,:) = nan(size(P(bad,:)));
end

Solution

Here's one way to complete vectorization with permute -

P = squeeze(sum(bsxfun(@ge,permute(r,[1 3 2]),csm),2)).'

Explanation

1) Push dim-2 of r to dim-3 position bringing in singleton-dim at dim=2. So, when paired with csm for bsxfun(@ge operation, we would have a 3D array as output, whose :

dim-1 : r, csm 's dim-1
dim-2 : csm's dim-2
dim-3 : r's dim-2

2) The original operation had iterative output with dim-2 still representing csm's dim-2 and sum-reduction along it. So, in our 3D array output, we need to sum-reduce along dim-2 as well.

3) Final steps involve squeeze-ing and transpose to correspond to the way we are iteratively saving output P.