topic Re: Log Likelihood for BETA distribution in Proc NLMIXED in Statistical Procedures

Log Likelihood for BETA distribution in Proc NLMIXED

LadyIDO — Tue, 03 Jan 2023 02:30:15 GMT

Hello,

I'm trying to estimate the beta distribution parameters (alpha, beta) using Proc NLMIXED. Because beta is not one of the distributions with built-in log likelihood, I used the general statement as follows:

proc nlmixed data=test;
  parms a=5 b=5;
  *s=beta(a,b);
  s=gamma(a*b)/(gamma(a)*gamma(b));
  ll=(1-a)*xbeta+(1-b)*(1-xbeta)-log(s);
  model xbeta~general(ll);
run;

I tried using the beta function directly and with gamma function. Both gave me the 'QUANEW optimization cannot be completed' error.

The test data is just 1000 beta variates with alpha=beta=5. I even set the initial values to 5, but still got this error.

Anne

Re: Log Likelihood for BETA distribution in Proc NLMIXED

StatDave — Tue, 03 Jan 2023 03:15:13 GMT

If you want to fit a model to a response variable that is distributed as beta, then you can use the FMM or GLIMMIX procedure. Specify DIST=BETA in the MODEL statement. If your goal is to estimate the parameters of the beta distribution for one population, you can use the HISTOGRAM statement in PROC UNIVARIATE as mentioned in this note.

Re: Log Likelihood for BETA distribution in Proc NLMIXED

Rick_SAS — Tue, 03 Jan 2023 14:18:31 GMT

First, the beta function is

s=gamma(a+b)/(gamma(a)*gamma(b));

not

s=gamma(a*b)/(gamma(a)*gamma(b));

Second, the beta distribution when a=b is somewhat degenerate, in the sense that the mean is exactly 1/2. This might be relevant for your nonconvergence issues. I'm not sure.

But if you want to use PROC NLMIXED, the easiest way is to specify the LOGPDF of the beta distribution directly:

data test;
call streaminit(54321);
do i = 1 to 1000;
   xbeta = rand("Beta", 5, 5);
   output;
end;
run;

proc nlmixed data=test;
  parms a=4 b=5.5;
  ll= logpdf("Beta", xbeta, a, b);;
  model xbeta~general(ll);
run;

Re: Log Likelihood for BETA distribution in Proc NLMIXED

LadyIDO — Tue, 03 Jan 2023 15:32:46 GMT

Thank you so much Rick! Using logpdf did it! I didn't know about logpdf function. This would really help with future work needing log likelihood.

Thanks again!

Anne

Re: Log Likelihood for BETA distribution in Proc NLMIXED

LadyIDO — Tue, 03 Jan 2023 15:39:58 GMT

I was trying to understand how to use nlmixed to estimate parameters of 
less common distribution like the following for beta-binomial. I didn't know about
FMM before. I'll check it out. Thanks.

proc nlmixed data=&data fconv=1E-14 df=&rdf alpha=0.05 ecov cov;
  parms mu=0.5 gamma=0.5;
  bounds mu>=0, mu<=1, gamma>0, gamma<1;
  ll=0;
  theta=gamma/(1-gamma);
  do i=1 to &ntrials;
    if i<=&nsucc
      then ll=ll+log(mu+(i-1)*theta);
        else ll=ll+log((1-mu)+(i-&nsucc-1)*theta);
    ll=ll-log(1+(i-1)*theta);
  end;
  alpha=mu/theta;
  beta=(1-mu)/theta;
  model &nsucc~general(ll);
run;