% By David (TA), prelim verson (i.e. could be better commented 
% and cleand up a bit, but does the job well).

function BinomialEx
close all
n_rolls = 30;  % change these 2 values to experiment with binom a bit
p_scss  = 3/8; % in the class problem set n_tosses = 5 and p_sccs = 1/6

% write a little helper f-n 'find_probs' (see below)
prob_values = find_probs(n_rolls, p_scss)
n_succs = [0:n_rolls];

figure;
plot(n_succs, prob_values, 'bo', 'MarkerFaceColor', [.49 1 .63])
title('A Binomial Probability Mass Function (pmf)')

cum_prob_vals = cumsum(prob_values)

%pause(5)
figure;
stairs(n_succs, cum_prob_vals) 
title('A Binomial Cumulative Distribution Function (cdf)')

% Helper function 'find_probs'
% returns a vector of size n+1, k-th element of which (where k runs 0, 1, ..thru n}) 
% is the probability of k successes out of n trials 
function p_values = find_probs(n, p)
p_values = zeros(1,n+1); 
  for i = 1:n+1
    p_values(i) = nchoosek(n,i-1)*p^(i-1)*(1-p)^(n-(i-1));
  end
end

% what follows may or may not be useful (comment out and try out if desired).
%% Toss a coin a bunch of times
% n = 10
% x = rand(1, n) < .5
% 
% %% Toss a coin a bunch of times; count the number of heads
% % 10 tosses of a fair coin. 1 - Heads, 0 - Tails
% 
% n = 10
% t = rand(1, n) < .5
% x = sum(t)
%
% %% Toss a coin a bunch of times; do it again, and again ..., and again...
% clear
% p = .7
% n = 100
% trials = 10000
% x = zeros(1,trials);
% for i = 1:trials
%     t = rand(1,n) < p;
%     x(i) = sum(t);
% end
% hist(x, 0:n)
% 
% %% Plot empirical cdf
% n = 10
% x = sum(ceil(6*rand(2,n)))                    % Roll a pair of dice and look at the sum of each roll
% stairs( [min(x) sort(x)], [0:1/length(x):1] ) % Plot the c.d.f of x
% 
% 
 end
