%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% SCRIPT ANMT_Lec2_binomialLikelihoodTable_Fancy
%
% Simulates Binomial random samples from Binomial(.3), and determines the
% log-likelihood of each given a vector of guesses 0.2:0.1:0.9.
%
% Log likelihoods are reported in a fancy table using the program
% displaytable.m available in Mathworks and on the class website.
%
% INPUT: Random samples from Binomial(.3). Created within program.
%
% OUTPUT: likelihoodTable(Grows dinamically): Table with likelihood values
%
%                                                 Cesar Zamudio, 8.31.2013
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Clean house
clear; clc; format bank;

%Simulate required random variables in three vectors.
%Remember that the whole sample is not needed when dealing with the
%binomial log-likelihood: only the SUM of the successes is needed. This is
%what is reported by the binornd function and is labeled "x" in Assig. 2.
bin50=binornd(50,.3);
bin100=binornd(100,.3);
bin1000=binornd(1000,.3);

%Log-likelihood table
logLikelihoodTable=zeros(3,5);

%%%%%%%%%%%%%%%%%%%%%%%%%
%       ALGORITHM       %
%%%%%%%%%%%%%%%%%%%%%%%%%
%Obtain log-likelihoods for first sample (N=50)
N=50;
%The cycle's "i" index will be used to represent the parameter p.

%The only problem with using "i" as p is that we cannot use i
%directly to access our likelihood table matrix cells (because there's
%no such thing as a location in, say, (1, 1.5). So, let's create an
%index that will grow along with i, to know exactly where we are.
index=1;

%Finally, understand that "index" is started on the OUTSIDE of the cycle.
%This is so it will not be reset WITHIN the cycle, thus causing the index
%to never advance. On the other hand, p is started INSIDE of the
%cycle, because we want it to change as i changes.

for i=0.2:0.1:0.9
    
    p=i;
    
    %For each value of p, we must calculate the log-likelihood. The
    %log-likelihood has three parts, so we will calculate them separately,
    %and then sum them to form the final log-likelihood.
    
    %IMPORTANT: In general, the constant is irrelevant for the
    %log-likelihood calculation. We will calculate the first part of the
    %log-likelihood, but because of its size, it's best to omit it.
    
    %PART 1. Obtain the log of constant k as defined in Assignment 2.
    %This is a straightforward calculation.
    LOGLIKPART1=factorial(N)/  (  factorial(bin50)*factorial(N-bin50)  )  ;
    
    %PART 2. x*log(p). Also straightforward.             
    LOGLIKPART2=bin50*log(p);
    
    %PART 3. (n-x)*log(1-p). Also a simple calculation.
    LOGLIKPART3=(N-bin50)*log(1-p);
    
    %Finally, form the log-likelihood.
    LOGLIKELIHOOD=LOGLIKPART2+LOGLIKPART3;             
    
    %Store the log-likelihood in the corresponding table cell.
    logLikelihoodTable(1,index)=LOGLIKELIHOOD;
    
    %MAKE SURE TO INCREASE YOUR INDEX, so you update the position where you
    %are.
    index=index+1;
end

%Obtain log-likelihoods for second sample (N=100)
N=100;
index=1;

for i=0.2:0.1:0.9
    
    p=i;
    
    %Just omit Part 1
           
    LOGLIKPART2=bin100*log(p);
        
    LOGLIKPART3=(N-bin100)*log(1-p);
        
    LOGLIKELIHOOD=LOGLIKPART2+LOGLIKPART3;             
    
    %REMEMBER TO STORE VALUES IN THE NEXT ROW OF YOUR LIKELIHOOD TABLE
    logLikelihoodTable(2,index)=LOGLIKELIHOOD;
    
    %MAKE SURE TO INCREASE YOUR INDEX, so you update the position where you
    %are.
    index=index+1;
end

%Obtain log-likelihoods for second sample (N=1000)
N=1000;
index=1;

for i=0.2:0.1:0.9
    
    p=i;
    
    %Just omit Part 1
           
    LOGLIKPART2=bin1000*log(p);
        
    LOGLIKPART3=(N-bin1000)*log(1-p);
        
    LOGLIKELIHOOD=LOGLIKPART2+LOGLIKPART3;             
    
    %REMEMBER TO STORE VALUES IN THE NEXT ROW OF YOUR LIKELIHOOD TABLE
    logLikelihoodTable(3,index)=LOGLIKELIHOOD;
    
    %MAKE SURE TO INCREASE YOUR INDEX, so you update the position where you
    %are.
    index=index+1;
end

disp('Binomial LLs from p=0.2 to p=0.9. True p=.3 (starred)') %Display a comment
disp(' ') %This command creates one white space.

%Display summary of results using "displaytable" function from Mathworks.
displaytable(logLikelihoodTable,{'N' 'p'},{{'N=50' 'N=100' 'N=1,000'}, ...
              {'p=0.2' 'p=0.3*' 'p=0.4' 'p=0.5' 'p=0.6' 'p=0.7' 'p=0.8' 'p=0.9'}},10,2)




