WARNING: The final Hessian matrix is full rank but has at least one negative eigenvalue.

[Yeganeh]

Hi  everybody,

I am trying to run zip model in proc nlmixed (SAS 9.4):

proc nlmixed data=migraine  maxiter=10000;                                                                                        

parms b0 =1 b1 = 0 b2= 0 c0 =1 c1 = 0 c2= 0 z0= 0 z1= 0 z2= 0 z3= 0 rho=-0.2 s2u =1 s2v=1;
                                                                                                                      

/*linear predictor for mixture probability of zero part*/                                                                          

linp_pi = b0 + b1*time+ b2*group*time+ u;                                                                                                   

pi = exp(linp_pi)/(1+exp(linp_pi));                       

/*linear predictor for mean of counts (poisson part)*/                                                                                              

linp_mu = c0 + c1*time+ c2*group*time+ v;                                                                                                  

mu = exp(linp_mu);   

logsigu2 = z0 + z1*group;                                                                                                                    

logsigv2 = z2 + z3*group;

s2u=exp(logsigu2);                                                                                                                           

suv=rho*sqrt(exp(logsigu2)*exp(logsigv2));                                                                                                   

s2v=exp(logsigv2);       

if attack=0 then  ll = log((pi) + (1-pi)*exp(-mu));                                                                                                                           

else ll = log((1-pi)) + attack*log(mu) - lgamma(attack+1)- mu;                                                                               

model attack ~ general(ll);                                                                                                                                                                

random u v ~ normal([0,0],[s2u,suv,s2v]) subject=id;                                                                                                                        run;

I am attaching longitudinal data. While running my model for "cov hess, tech=newrap, method= gauss, maxiter, qmax noad, seed=12345", Unfortunately, I always face this warning and have some large SE or no SE .

WARNING1: The final Hessian matrix is full rank but has at least one negative eigenvalue. Second-order
optimality condition violated.

WARNING2: 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.

I changed initial values, method, qmax...again and again. Sometimes I have Error: Optimization cannot be completed. I don't know what to do.

I have a principal question:

How do I find appropriate initial values?

Thank you for your time.

Parameter Estimates
Parameter	Estimate	Standard Error	DF	t Value	Pr > |t|	95% Confidence Limits	Gradient
b0	-2.6089	0.4103	64	-6.36	<.0001	-3.4285	-1.7893	6.81751
b1	-1.5772	14.2856	64	-0.11	0.9124	-30.1160	26.9617	0.071280
b2	-2.5369	12.5514	64	-0.20	0.8405	-27.6111	22.5373	0.078330
c0	0.7799	0.1335	64	5.84	<.0001	0.5133	1.0466	-25.8045
c1	-0.4147	0.09284	64	-4.47	<.0001	-0.6002	-0.2292	19.9297
c2	0.09872	0.05501	64	1.79	0.0775	-0.01118	0.2086	-37.4677
z0	0.2445	2.2688	64	0.11	0.9145	-4.2881	4.7770	4.09923
z1	0.4000	1.5060	64	0.27	0.7914	-2.6087	3.4087	6.45820
z2	-0.9555	1.6735	64	-0.57	0.5700	-4.2988	2.3877	-1.27785
z3	-1.3047	1.0613	64	-1.23	0.2234	-3.4249	0.8155	-3.62090
rho	0.07197	2.0419	64	0.04	0.9720	-4.0072	4.1512	0.76599
s2u	2.8419	.	64	.	.	.	.	1.73494
s2v	0.02830	0.06871	64	0.41	0.6818	-0.1090	0.1656	-73.7197

[sbxkoenk]

Hello @Yeganeh ,

You ask: "How do I find appropriate initial values?"

There is no golden rules in terms of specifying starting values for the parameters for PARMS statement in NLMIXED. You can always start with using the default starting values (parameters not listed in the PARMS statement are assigned an initial value of 1).

Using a grid search (possible in PARMS statement) or specifying a guess based upon subject matter knowledge or previous studies might help.

But before going further in that direction, I would like to ask why you turn to PROC NLMIXED for your zip model.

There are several procedures that can deal with zip models (zip = Zero-inflated Poisson regression).

Have you tried :

• GENMOD procedure
• (HP)COUNTREG Procedure
• (HP)FMM Procedure
• GLIMMIX procedure
  ??

Maybe your model has a little twist that above procedures cannot accomplish, but then you can maybe use above procedures to generate appropriate starting values for at least some of your parameters.

Kind regards,

Koen

[Yeganeh]

Dear Koen,

Thank you for sharing your valuable advice with me. 

When I used proc genmod, parameters of zero model had a high SE.

proc genmod data = migraine;
model attack= Time Group*Time/ dist=zip;
zeromodel Time  Group*Time/link = logit;
run;

Also, I ran proc glimmix and proc countreg but, the initial values didn't work correctly.(unable to estimate SE or high SEs)

proc glimmix data = migraine noclprint method=laplace;
class id;
model attack= Time Group*Time / solution dist=poisson;
random intercept / subject = id;
run;

proc countreg data=migraine;
model attack= Time Time *Group / dist=zip;
zeromodel attack~Time Time *Group/ link=logistic;
run;

I appreciate your guidance.

[StatDave]

Several comments here... First, I suspect what you want is for the model to allow for the groups to have separate intercepts and separate slopes. If so, then the model would be GROUP GROUP*TIME. If you use that model, there seems to be no real evidence of overdispersion which would be expected if the observed zeros were excessive. This can be seen by fitting the simple Poisson model ignoring the repeated measures.

proc genmod;
class group id;
model attack=group group*time / d=p;
effectplot;
run;

Note that the Pearson/DF and deviance/DF values are less than 1 suggesting no overdispersion. This might be why you have trouble fitting the zero-inflated model - in fact, trying to fit that model (without random effects) fails in GENMOD, and fails even if only an intercept is included in the model for extra zeros.

proc genmod;
class group;
model attack=group group*time / dist=zip;
zeromodel group group*time;
run;

So, probably the best model is the simple Poisson model above. You can account for the repeated measures be adding the REPEATED statement in the above code:

repeated subject=id / type=exch;

[Yeganeh]

Hi,

I'm so grateful for your suggestions. I check it.

[StatDave · Accepted solution]

I don't understand the concern. As I mentioned, the ZIP model seems inappropriate because the simple Poisson model shows no evidence of overdispersion and because of the fitting problems when trying to fit the ZIP model. The Poisson GEE model seems to provide a reasonable fit and the parameter estimates and test results seem consistent with the plot of the fitted model. So, the only code needed to fit a reasonable model is this: 

proc genmod;
class group id;
model attack=group group*time / d=p;
repeated subject=id / type=exch;
effectplot;
run;

[Yeganeh]

Ok, I get it now. That's clear.

Thank you for your immediate reply.