clear all %This clears all of the variables before running.
%Let's consider the binomial partitioning of transcription factors in an E.
%coli cell. Let's say there are N transcription factors. Upon a cell
%division, what is the distribution of transcription factors in the
%daughter cell? Let's first define some parameters.
NTrials = 2^10; %Number of coin flips.
p = 0.5; %Probability of a "success".
%We'll use another for loop to cycle through the trials.
for i = 1:NTrials
RandomNumber = rand; %Generate the random number.
if RandomNumber < p % A success
Success(i) = 1;
else %A failure
Success(i) = 0;
end %End the coin flip.
end %End the loop.
Success; %This shows a vector of 1's and 0's. There are long streaks as well.
%Let's look at the average number of successes.
sum(Success) / NTrials;
%Given a certain number of trials, how often do we have streaks? How long
%is the average streak? We'll keep a counter of the maximum streak length
%as we parse the success vector. First set a few parameters.
MaxStreak = 0; %Maximum length of the streak.
CurrentCount = 1; %Current number of items in the streak.
%Lets run the the loop
for i=1:NTrials - 1 %Since we iterate from success(i+1), we go to NTrials - 1
%We'll test what the next value in the list is.
if Success(i) == 1
if Success(i+1) == 1
CurrentCount = CurrentCount + 1;
else
MaxStreak = max([MaxStreak, CurrentCount]);
CurrentCount = 1;
end
end
end
%Let's go through this many times.
NSimulations = 1000;
NTrials = 128; %Number of coin flips.
p = 0.5; %Probability of a "success".
for j = 1:NSimulations
for i = 1:NTrials
RandomNumber = rand; %Generate the random number.
if RandomNumber < p % A success
Success(i) = 1;
else %A failure
Success(i) = 0;
end %End the coin flip.
end
%We'll reinitialize the parameters.
MaxStreak(j) = 0;
CurrentCount = 1;
%Lets run the the loop
for i=1:(NTrials - 1) %Since we iterate from success(i+1), we go to NTrials - 1
%We'll test what the next value in the list is.
if Success(i) == 1
if Success(i+1) == 1
CurrentCount = CurrentCount + 1;
else
MaxStreak(j) = max([MaxStreak(j), CurrentCount]);
CurrentCount = 1;
end
end
end
end
[mean(MaxStreak), std(MaxStreak)]
%Let's plot the binomial distribution using MATLABs definition of the
%binomial distribution, binopdf. The first thing we need to choose is the
%range of x values. We can do this by setting a vector
x = 0:10
NTrials = 10
BinProb = binopdf(x, NTrials, p); %Takes arguments of X, N, and P.
%We can plot this to look at the distribution. Let's do this using MATLABs
%"plot" command. We'll plot them as dots with connecting lines. We'll do
%this as "Figure 1"
figure(1)
plot(x, BinProb, '.-')
xlabel('x values');
ylabel('Probability');
%We can also do a bar graph.
figure(2)
bar(x, BinProb)
xlabel('x values');
ylabel('Probability');
%What does this look like with low probability and large trials?
NTrials = 1000;
p = 0.002;
BinProb = binopdf(x, NTrials, p);
figure(3)
bar(x, BinProb)
xlabel('x values')
ylabel('Probability')