Date: Wed, 11 Mar 2009 09:50:39 -0700
Reply-To: stringplayer_2@yahoo.com
Sender: "SAS(r) Discussion" <SAS-L@LISTSERV.UGA.EDU>
From: Dale McLerran <stringplayer_2@YAHOO.COM>
Subject: Re: NLMixed Convergence
In-Reply-To: <893d59bc-1be1-4334-9b10-122a2ec88ec1@m4g2000vbp.googlegroups.com>
Content-Type: text/plain; charset=us-ascii
Because the gamma distribution requires a response variable with
positive values, the zero-inflated gamma distribution can be
fit as a two-stage model:
1) extract the observations with a positive value for the response
and model those observations with the gamma distribution using
the procedure GENMOD
2) construct a modified response variable which is zero when the
original response (NEXTREVENUE) is zero and one when the
original response is positive. Model the probability of the
modified response being zero. Again, you can use the GENMOD
procedure to fit this model.
The models in both 1) and 2) will be exactly correct because the
gamma distribution does not support zero values. Therefore,
there is no mixture of gamma distribution zeros with the zero-
inflation process. In this, it may be incorrect to make statements
about a zero-inflated gamma distribution. Rather, it is probably
better to talk just of a mixture of a gamma response with a zero
response.
If you want to report the "zero-inflated gamma" model likelihood,
you will need to use NLMIXED. You can initialize NLMIXED from the
parameter estimates returned by your two GENMOD stages. Those
parameters should maximize the likelihood - at least if you are
using the same gamma and logistic regression model parameterizations.
Dale
---------------------------------------
Dale McLerran
Fred Hutchinson Cancer Research Center
mailto: dmclerra@NO_SPAMfhcrc.org
Ph: (206) 667-2926
Fax: (206) 667-5977
---------------------------------------
--- On Mon, 3/9/09, Bminer <b_miner@LIVE.COM> wrote:
> From: Bminer <b_miner@LIVE.COM>
> Subject: Re: NLMixed Convergence
> To: SAS-L@LISTSERV.UGA.EDU
> Date: Monday, March 9, 2009, 8:25 AM
> On Mar 9, 10:22 am, Bminer <b_mi...@live.com> wrote:
> > On Mar 9, 9:55 am, Bminer <b_mi...@live.com>
> wrote:
> >
> >
> >
> > > Hi All-
> >
> > > I'm using nlmixed to run a simple zero
> inflated gamma with one
> > > continuous variable and one dummy coded variable.
> >
> > > proc nlmixed data=China;
> > > parms b0_f=0 b1_f=0 b2_f=0
> > > b0_h=0 b1_h=0 b2_h=0
> > > log_theta=0;
> >
> > > eta_f = b0_f + b0_f*hassku +
> b2_f*REVORDONE ;
> > > p_yEQ0 = 1 / (1 + exp(-eta_f));
> >
> > > eta_h = b0_h + b0_h*hassku + b2_h*REVORDONE
> ;
> > > mu = exp(eta_h);
> > > theta = exp(log_theta);
> > > r = mu/theta;
> >
> > > if NEXTREVENUE=0 then
> > > ll = log(p_yEQ0);
> > > else
> > > ll = log(1 - p_yEQ0)
> > > - lgamma(theta) +
> (theta-1)*log(NEXTREVENUE) - theta*log
> > > (r) - NEXTREVENUE/r;
> >
> > > model NEXTREVENUE ~ general(ll);
> > > predict (1 - p_yEQ0)*mu out=expect_zig;
> > > estimate "scale" theta;
> >
> > > run;
> >
> > > I get the following log warnings:
> >
> > > Does anyone have a suggestion on what to try /
> look for to get
> > > convergence?
> >
> > > NOTE: GCONV convergence criterion satisfied.
> > > NOTE: At least one element of the (projected)
> gradient is greater than
> > > 1e-3.
> > > NOTE: Moore-Penrose inverse is used in covariance
> matrix.
> > > WARNING: The final Hessian matrix is not positive
> definite, and
> > > therefore the
> > > estimated covariance matrix is not full
> rank and may be
> > > unreliable.
> > > The variance of some parameter estimates
> is zero or some
> > > parameters
> > > are linearly related to other
> parameters.
> >
> > Typo on the code. s/b:
> >
> > proc nlmixed data=China;
> >
> > > parms b0_f=0 b1_f=0 b2_f=0
> > > b0_h=0 b1_h=0 b2_h=0
> > > log_theta=0;
> >
> > > eta_f = b0_f + b1_f*hassku +
> b2_f*REVORDONE ;
> > > p_yEQ0 = 1 / (1 + exp(-eta_f));
> >
> > > eta_h = b0_h + b1_h*hassku + b2_h*REVORDONE
> ;
> > > mu = exp(eta_h);
> > > theta = exp(log_theta);
> > > r = mu/theta;
> >
> > > if NEXTREVENUE=0 then
> > > ll = log(p_yEQ0);
> > > else
> > > ll = log(1 - p_yEQ0)
> > > - lgamma(theta) +
> (theta-1)*log(NEXTREVENUE) - theta*log
> > > (r) - NEXTREVENUE/r;
> >
> > > model NEXTREVENUE ~ general(ll);
> > > predict (1 - p_yEQ0)*mu out=expect_zig;
> > > estimate "scale" theta;
> >
> > > run;
>
> I changed the tech= option to NEWRAP and got
> convergence......should I
> be happy to have found the right algorithm or skpetical of
> the results
> if other tech options failed?