Mixture Poisson parameter estimates

mixpoissfit.m

clf('reset');

% Definition of pdf
pdf_mixpoiss = @(x, p, lambda1, lambda2) p * poisspdf(x, lambda1) ...
                                          + (1 - p) * poisspdf(x, lambda2);

% Make a demo data
rng(218, 'twister');
p = 0.4;
n =  10000;
lambda = [0.86 11.2];
n1 = n * p;
n2 = n * (1 - p);

x1 = poissrnd(lambda(1), n1, 1);
x2 = poissrnd(lambda(2), n2, 1);
x = cat(1, x1, x2);

% Setting initial value for parameters
pStart = .5;
lambdaStart = quantile(x,[.25 .75]);
start = [pStart lambdaStart];

% Setting boundaries for parmaters
lb = [0 0 0];
ub = [1 Inf Inf];

% Estimatation of parameters
%statset('mlecustom')
options = statset('MaxIter',4000, 'MaxFunEvals',4000);
paramEsts = mle(x, 'pdf', pdf_mixpoiss, 'start', start, ...
  'lower', lb, 'upper', ub, 'options', options)

% Plotting a demo data and estimation
figure(1);

bins = 0:1:max(x)+1;
h = bar(bins, histc(x ,bins) / length(x), 'histc');
xlim([0, max(x)]+1)
set(h,'FaceColor',[.9 .9 .9]);

pdfgrid = pdf_mixpoiss(bins, paramEsts(1), paramEsts(2), paramEsts(3));
hold on; plot(bins, pdfgrid,'-'); hold off;
xlim([0, max(x)]+1)
xlabel('x'); ylabel('Probability Density');

paramEsts =
0.3982 0.8372 11.1374

Requirements: Matlab, Statistics Toolbox
Link: Fitting Custom Univariate Distributions