I am attemtping to write a function that produces m random simulations of 5 Bernoulli Trials. I create a histogram showing the distribution of the number of successes across the m simulations.
I then need to also plot a line showing the theoretical / normalised distribution around the theoretical mean number of successes.
Here is my function as of now:
function x = generate_binomial_bernoulli(n,p,m)
% generate Bi(n, p) outcomes m times
emperical = zeros(1,m); % allocate array for m simulations
for i = 1:m % iterate over m simulations
successes = 0; % count the number of successful trials per simualtion (0-5)
for j = 1:n % iterate through the n trials
u = rand; % generate random nuumber from 0-1
if (u <= p) % if random number is <= p
successes++; % count it as a success
endif
end
emperical(i) = successes; % store the number of successful trials in this simulation
end
close all; % close any existing graphs
x_values = [0:n]; % array of x-axis values
hist(emperical, x_values, "facecolor", "r"); % plot empirical data
xlim([-0.5 (n + 0.5)]); % set x-axis to allow for histogram bar widths
hold on; % hold current graph
mean = n * p; % theoretical mean
norm = normpdf(x_values, mean, 1); % normalised y values
plot(x_values, norm, "color", "b"); % plot theoretical distribution
legend('Emprical', 'Theoretical');
end
When the function is called as
generate_binomial_bernoulli(5, 0.2, 100)
I expect to see a red histogram showing the results of 100 simulations (emperical results), and then a blue line plot normalised around the mean of 1 success (theoretical result).
Here is the graph produced:
The emperical results are displayed correctly, but the theoretical plot is only extending to the height of a very low value on the y-axis.
Where is my function going wrong?