[Scilab-users] Is there a way to do it with Matrix Operations?

[Heinz Nabielek]

Works like a charm, but I would have never found out myself.
Congratulations to your SciLab skills.
Heinz

PS: In some 20 textbooks and a dozen publications, I could not find such a simple recipe for generating random deviates of a user-defined discrete function............
On 25.09.2017, at 20:43, Rafael Guerra wrote:

> Hi Heinz,
>  
> Check this solution:
>  
> // START OF CODE
> p = [0.1176471 0.2352941 0.0588235 0.3882353 0.2 ];
> P = cumsum(p)';
> np = length(P);
> nj = 7;
> N = 100;
> X = grand(nj,N,'def');
> C = zeros(np,nj);
> P = [0;P];
> for i=1:np
>     C(i,:) = sum((X>=P(i) & X<P(i+1)),2)';
> end
> // END OF CODE
>  
> Regards,
> Rafael
>  
> From: users [mailto:users-bounces at lists.scilab.org] On Behalf Of Heinz Nabielek
> Sent: Monday, September 25, 2017 2:03 PM
> To: Users mailing list for Scilab <users at lists.scilab.org>
> Subject: [Scilab-users] Is there a way to do it with Matrix Operations?
>  
> Dear colleagues:
> 
> in an attempt to code the generation of random deviates for a user-defined probability function p=[0.1176471    0.2352941    0.0588235    0.3882353    0.2 ], I spent only a few minutes to write the Scilab code below and it gives me all the solutions (frequency distribution of random numbers) that I need.
> 
> N=100;X=grand(7,N,'def');
> C=[];for j=1:7;Count(1:5)=0;for k=1:N;i=1;while X(j,k)>P(i);i=i+1;end;Count(i)=Count(i)+1;end;C=[C Count];end;
> and one typical sample run yields this
> C  =
>    15.    9.     6.     12.    8.     12.    10.  
>    20.    26.    38.    20.    23.    26.    24.  
>    6.     7.     4.     7.     5.     10.    4.   
>    38.    39.    32.    37.    48.    30.    39.  
>    21.    19.    20.    24.    16.    22.    23. 
> 
> However, the for and while loops will be terribly inefficient and this is not good for large scale Monte-Carlo simulations.
> 
> Is there a way to do it with Matrix Operations?