%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% SCRIPT ANMT_Lec2_binomialLikelihoodTable
%
% 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 simple outputted table.
%
% 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')
logLikelihoodTable