Randomized complete block, repeated measures, poisson, count data, nlmixed, error

[nkoerth]

SAS 9.4
I am getting a "Termination due to floating point exception" error with the following program. Any ideas on what I'm doing wrong? (note: could not cut and paste so had to type in program, hopefully no typos)

blk = 4 blocks (random)

trt = 5 treatments (coded as dummy: trt1 trt2 trt3 trt4)

cmp = 20 sites, 1 per trt per block

time = time since trt (1, 3, 5, 8, 12, 16)

abun=count

poisson distribution (and eventually ZIP)

proc sort;

by blk;

proc nlmixed;

parms intercept=0 b1=0 b2=0 b3=0 b4=0 b5=0 b6=0 b7=0 b8=0 b9=0 varblk=0 varcomp=0;

array blkE {4};

array cmpE {20};

eta = intercept + varblk*blkE{blk} + varcmp*cmpE{cmp} + b1*time + (b2 + b6*time)*trt1 + (b3 + b7*time)*trt2 + (b4 + b8*time)*trt3 + (b5 + b9*time)*trt4;

mu=exp(eta);

model abun ~ poisson(mu);

random blkE1 blkE2 blkE3 blkE4 cmpE1 cmpE2 cmpE3 cmpE4 cmpE5 cmpE6 cmpE7 cmpE8 cmpE9 cmpE10 cmpE11 cmpE12 cmpE13 cmpE14 cmpE15

                                                            cmpE16 cmpE17 cmpE18 cmpE19 cmpE20 ~

normal([0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0], [varblk,

0, varblk,

0, 0, varblk,

0, 0, 0, varblk,

0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp,

0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, varcmp])

subject=blk;

run;

Thanks for the help!

[SteveDenham]

I get this error most often after running several times in interactive mode.  The only cure I have found is to reboot SAS.  One thing that might be causing this is mu=exp(eta).  Just plugging in for time=16, if b6-b9 were 1 after some iterations, you would have exp(64), which is starting to get to the point that machine precision can cause problems.

So, search the SAS-L archives for some posts by David Cassell (and a lot of others), where they reparameterize for both Poisson and ZIP models to avoid this overflow problem.

Steve Denham

[nkoerth]

Thanks Steve!

I have tried the program after a fresh start of SAS. No luck.

Now that you mention it, I do remember reading in the archives about specifying the poisson differently to avoid problems, but didn't realize that included my problem.  I will go back and review.  Thanks again.

[nkoerth]

I'm finally back to this project; sorry for the delay. I'm also unable to cut and paste, so bear with me. I'll not retype the whole program from above.

My program now works with dummy data, but I still have 2 questions:

1. Is Laplace approximation (qpoints=1) appropriate?

2. I am not certain how to calculate degrees of freedom for different contrasts. Any advice?

These lines changed:

proc nlmixed df=75 instep=1E-4 qpoints=1;

.

.

.

ll = -mu+abun*log(mu)-lgamma(abun+1);

model abund ~ general(ll);

.

.

*Added after random statement;

contrast 'trt effects' intercept, b2, b3, b4, b5 df=15;

contrast 'CC vs SW' b2, b3 df=1;

contrast 'CC vs CO' b2, intercept df=1;

.

.

.

contrast 'time effect' b1 df=75;

run;

Thanks!

[SteveDenham]

Laplace as a method is OK.

As far as the degrees of freedom, df= denominator degrees of freedom.  I am really never very good at figuring these out, but... I think the contrasts for specific comparisons should have the same denominator degrees of freedom, so I would change the two df=1 to df=15.  This assumes 15 is correct.  My skeleton ANOVA has df=16 (5 treatments - 1)*(4 blocks), but this might be off.  A good way to check would be to run GLIMMIX, realizing that the only thing of interest is the degree of freedom output.

Steve Denham

[nkoerth]

Oops, yes 15 is correct.

I was hoping for more "peer reviews" so also posted to sas-l (hope that wasn't unintentionally rude).

Thanks for you help, Steve!

[Fugue]

Re: your problem with copy/paste . . . not sure if your problem is the same as mine . . .

The workaround for me is to copy from SAS into a text editor (I use Microsoft Word), then copy from Word and paste into the SAS discussion forum.