MATLAB: Probability of binary sequence with Monte Carlo

binary sequencemonte carlo

Hello everyone !
Context:
I'm trying to compute the probability that a binary sequence "stop" under a certain constraint: if the sequence has a majority of "1" we stop.
ex: Let's say n is the length of the sequence; At n=2 we have : {(00),(01),(10),(11)} possibilities and so at n=3 {(000),(001),(010),(011),(100),(101)} possibilities (note: we didn't took (1,1) car we already stop at n=2) and so on… the probability to stop at n=2 is P(2)=1/4 and at n=3, P(3)=2/6
Question:
My goal is to know what is the probability that at "step n" we stop (ie, how many sequence at n-step "stop" over how many possibilities, taking in account tthat some possibility already stop before)?
Tentative:
The probleme is the more n will grow up the bigger possibilities we have and im not sure the simulation can support for n tend to big number.
So I heard about "Monte carlo simulation"; my idea was to reapeat a certain time the experience. Like I generate a random binary sequence of length-n and i look the majority (ie if sum(sequence)>=n/2 , it's important that this appear in code) and try to generate statistic. I don't know how much "monte carlo simulation" could help me in my problem.
I did this:
rng('shuffle')
nb_mort=0;  %number of sequence who stop
n=0;       %length of sequence
maj=0;     %majority
M=1000;  %nb of run
for n_run=1:M;
    
    while (n <= 100) | (maj <= n/2 )
        
        
        count=zeros(n);
        n=n+1
        
        seq = randi([0 1],1,n)
        
        maj=sum(seq)
                       
    end
    
    clear n;    %I want he start again the loop but renitialize his length
    n=0;
    
       nb_mort=nb_mort+1;
                       
    %count(n)=count(n)+1;
    
    %proba_n(n)=count(n)/M;
    
    %nb_mort_mean(n_run)=nb_mort
    
end
Don't take in account, the last line in comment Proba, because the proba is wrong, i didn't found…
I don't know why also at the end of my "while" he doesnt renialize the n??
Im not super familiar with Monte Carlo and im beginner on Matlab
Thank you for your help

Best Answer

First point is that probability of stopping at n=3 is not 2/6, it is (3/4)*(2/6) because you have to include the probablility that you didn't stop at n=2. Regardless ...
Seems like you could simply generate a long string of m digits and then use cumsum on that compared to 1:m to see where you stopped for that string. I.e., first index where the cumsum element is greater than the (1:m)/2 values is the stopping point. If you didn't stop in m digits put all of those results into a single "greater than m" category. Wrap that in a large Monte-Carlo loop and then divide your counts by the number of Monte-Carlo runs to get the individual probabilities.
And it looks like you have a rule to never stop at n=1. True?
E.g., initialize some values
m = 100; % some value for maximum stopping point we will look for
count = zeros(1,m+1);
m2 = (1:m)/2; % the comparison vector
And then each iteration of the Monte-Carlo loop would be something like
seq = randi([0 1],1,m);
cseq = cumsum(seq);
cseq(1) = 0; % force that we won't stop at n=1
x = find(cseq > m2, 1); % find the first stopping point
if isempty(x)
    count(m+1) = count(m+1) + 1;  % we didn't stop, so lump this into special category
else
    count(x) = count(x) + 1; % increment the 1st place we stopped
end
Related Question