# Section 7.1.1: Counting problems with Poisson distribution

Jump to: Source code Text output Plots Library index

    % Boyd & Vandenberghe "Convex Optimization"
    % JoÃ«lle Skaf - 04/24/08
    %
    % The random variable y is nonnegative and integer valued with a Poisson
    % distribution with mean mu > 0. In a simple statistical model, the mean mu
    % is modeled as an affine function of a vector u: mu = a'*u + b.
    % We are given a number of observations which consist of pairs (u_i,y_i),
    % i = 1,..., m, where y_i is the observed value of y for which the value of
    % the explanatory variable is u_i. We find a maximum likelihood estimate of
    % the model parameters a and b from these data by solving the problem
    %           maximize    sum_{i=1}^m (y_i*log(a'*u_i + b) - (a'*u_i + b))
    % where the variables are a and b.

    % Input data
    rand('state',0);
    n = 10;
    m = 100;
    atrue = rand(n,1);
    btrue = rand;

    u = rand(n,m);
    mu = atrue'*u + btrue;

    % Generate random variables y from a Poisson distribution
    % (The distribution is actually truncated at 10*max(mu) for simplicity)
    L  = exp(-mu);
    ns = ceil(max(10*mu));
    y  = sum(cumprod(rand(ns,m))>=L(ones(ns,1),:));

    % Maximum likelihood estimate of model parameters
    cvx_begin
        variables a(n) b(1)
        maximize sum(y.*log(a'*u+b) - (a'*u+b))
    cvx_end

    Successive approximation method to be employed.
       For improved efficiency, SDPT3 is solving the dual problem.
       SDPT3 will be called several times to refine the solution.
       Original size: 276 variables, 103 equality constraints
       92 exponentials add 736 variables, 460 equality constraints
    -----------------------------------------------------------------
     Cones  |             Errors              |
    Mov/Act | Centering  Exp cone   Poly cone | Status
    --------+---------------------------------+---------
     92/ 92 | 6.764e-01  3.304e-02  0.000e+00 | Solved
     90/ 92 | 5.602e-02  2.250e-04  0.000e+00 | Solved
     55/ 88 | 1.386e-03  1.261e-07  0.000e+00 | Solved
      0/ 34 | 1.389e-04  1.143e-09  0.000e+00 | Solved
    -----------------------------------------------------------------
    Status: Solved
    Optimal value (cvx_optval): +102.57