Date: Sat, 30 Apr 2005 14:32:43 -0700
Reply-To: Dale McLerran <stringplayer_2@YAHOO.COM>
Sender: "SAS(r) Discussion" <SAS-L@LISTSERV.UGA.EDU>
From: Dale McLerran <stringplayer_2@YAHOO.COM>
Subject: Re: ZINB(k) model
Content-Type: text/plain; charset=us-ascii
--- Alexander Cavallo <acavallo@NAVIGANTCONSULTING.COM> wrote:
> SAS-L (and in particular to Dale)
> I read with interest Dale's posts on SAS-L about estimating ZINB
> using PROC NLMIXED.
> I am interested in estimating a ZINB(k) model.
Alex placed his code for the ZINB(k) model at the end of his post,
but I have moved it up to the top in order to respond to questions
> proc nlmixed;
> parms b0_zer0=0 b1_zero=0 b0_nb=0 b1_nb=0 sigma=1 k=0 alpha=0;
> bounds alpha>0;
> eta_zero = b0_zero + b1_zero * x1 ;
> p_zero = exp(eta_zero) / (1 + exp(eta_zero);
> eta_nb = b0_nb + b1_nb * x1;
> lambda = exp(eta_nb);
> alpha = (lambda**(1-k))/(sigma*sigma);
> c = (1 + sigma * sigma * (lambda**k) )**(-alpha);
> p0 = p_zero + (1-p_zero) * c ;
> pelse = (1-p_zero) * c *
> exp(lgamma(alpha + y) - lgamma(alpha) - lgamma(y)) *
> (lambda / (lambda + alpha))**y;
> if y=0 then loglike = log(p0);
> else loglike = log(pelse);
> model y ~ general(loglike);
> expected = (1-p0) * lambda;
> predict expected out=expect;
> predict eta_zero out=zero;
> predict eta_nb out=nb;
> predict p0 out=prob;
> 1. Does this code look correct?
I've got to be honest. I don't know what the ZINB(k) distribution
is. Have you a reference for the ZINB(k) distribution? How is
ZINB(k) related to ZINB? It is probably easier to discuss the
relationship of NB(k) to NB.
> 2. I want to output datasets with the ilnear predictors for zero
> and NB, plus predicted prob for zero inflation - does my predict
> make sense? I will merge the four datasets by my id var later.
Your PREDICT statements look OK.
> 3. I need compute predicted probs and expected values for additional
> points. Will I have to compute these manually outside NLMIXED, or
> can I
> trick NLMIXED into predicting by including these obs in the dataset
> with Y (dep var) missing?
I don't know. I would certainly try to generate the predicted
values for the additional data points by including observations
in your input data set which have nonmissing predictors and
missing response. I don't see any reason why that should not
work, but I have not tested the idea myself.
> --Alex Cavallo
> Managing Consultant
> Navigant Consulting, Inc.
Fred Hutchinson Cancer Research Center
Ph: (206) 667-2926
Fax: (206) 667-5977
Do You Yahoo!?
Tired of spam? Yahoo! Mail has the best spam protection around