ERROR in using NLPNRA and PROC IML for MLE of custom function

[roystriker77r]

Dear community, 

I am trying to follow Rick Wicklin's example of using NLPNRA with PROC IML for maximum likelihood estimation. Unfortunately, I am encountering an error with NLPNRA.

I am trying to use this script to find the MLE of a two parameter, two variable custome probability mass function:

/*Define Density Function*/

data multivariate_data;
call streaminit(5);
%let theta1 = 5;
%Let theta2 = 3;
	array Xi[2];
    do i = 1 to 1000;
        do j=1 to 2;
        Xi[j]= rand('Bernoulli',.5);
        end;
  	N=exp(&theta1*(Xi[1]+Xi[2])+&theta2*Xi[1]*Xi[2]); /*numerator of density function*/
  	D=(1+2*exp(&theta1)+exp(2*&theta1+&theta2)); /*denominator of density function*/
  	P=N/D;
        output;
        
    end;
    drop i j;
run;

proc print data= multivariate_data;
var P;
RUN;


proc iml;
use multivariate_data;
read all var {"Xi1","Xi2"} into x;
close multivariate_data;
/* write the log-likelihood function the distribution of interest */
start LogLik(param) global (x);
   theta1 = param[1];
   theta2 = param[2];
   n = nrow(x);
   return ( theta1*(sum(Xi1)+sum(Xi2)) + theta2*sum(Xi1)*sum(Xi2)-n*log(1+2*exp(theta1)+exp(2*theta1+theta2)) );
finish LogLik;



p = {2 4};/* initial guess for solution */
con = {., .}; /* Upper bounds for parameters (if any) */
optn = {1,  /* find max of function       */
        3}; /* print a lot of output      */
call nlpnra(rc,xres,"LogLik",p,optn,con);

However, after running this code I receive an error code following the nlpnra call that says error in argument BLC. After reviewing the documentation I thought that
BLC was an option argument so I am confused as to why I am receiving this error. Anyt help would be much appreciated.
Thanks!

[Rick_SAS · Accepted solution]

Yes, the BLC parameter is optional, and you opted to specify it by using the CON argument. But the CON argument is not the correct format. If you don't want a constrained optimization, then omit the CON argument. Otherwise, specify it correctly.

I don't know what distribution you are trying to fit, but I want to point out that your input data set only has three unique pieces of information:

• When (x1,x2)=(0,0), the value P(x1,x2)=2.259E-6
• When (x1,x2)=(1,1), the value P(x1,x2)=0.999
• When (x1,x2)=(1,0) or (0,1), the value P(x1,x2)=0.00033

I also do not think you can perform an optimization in this function because it is unbounded. For example, look at the cross section for theta1=3 as you let theta2 vary:

/* visualize the LL when theta1=3 */
theta1 = 3;
theta2 = T( do(0,8,0.1) );
LL = j(nrow(theta2), 1);
do i = 1 to nrow(theta2);
   LL[i] = LogLik( theta1 || theta2[i] );
end;
title "LL(3, theta2)";
call series(theta2, LL) grid={x y};

As you can see, the objective function does not have a maximum. I suggest you check the computations that you used to construct the LogLik function.

[roystriker77r]

Dear community,  I am trying to follow Rick Wicklin's example of using NLPNRA with PROC IML for maximum likelihood estimation. Unfortunately, I am encountering an error with NLPNRA. I am trying to use this script to find the MLE of a two parameter, two variable custome probability mass function: 

/*Define Density Function*/

data multivariate_data;

call streaminit(5);

%let theta1 = 5;

%Let theta3 = 3;

array Xi[2];

    do i = 1 to 1000;

        do j=1 to 2;

        Xi[j]= rand('Bernoulli',.5);

        end;

  N=exp(&theta1*(Xi[1]+Xi[2])+&theta3*Xi[1]*Xi[2]);

  D=(1+2*exp(&theta1)+exp(2*&theta1+&theta3));

  P=N/D;

        output;

    end;

    drop i j;

run;

proc print data= multivariate_data;

var P;

RUN;

proc iml;

use multivariate_data;

read all var {"Xi1","Xi2"} into x;

close multivariate_data;

/* write the log-likelihood function the distribution of interest */

start LogLik(param) global (x);

   theta1 = param[1];

   theta3 = param[2];

   n = nrow(x);

   return ( theta1*(sum(Xi1)+sum(Xi2)) + theta3*sum(Xi1)*sum(Xi2)-n*log(1+2*exp(theta1)+exp(2*theta1+theta3)) );

finish LogLik;

p = {2 4};/* initial guess for solution */

con = {., .}; /* Upper bounds for parameters (if any) */

optn = {1,  /* find max of function       */

        3}; /* print a lot of output      */

call nlpnra(rc,xres,"LogLik",p,optn,con);

However, after running this code I receive an error code following the nlpnra call that says error in argument BLC. After reviewing the documentation I thought that BLC was an option argument so I am confused as to why I am receiving this error. Anyt help would be much appreciated. Thanks!