Glimmix with repeated measurement: how to apply the random statement and a time variable?

[SteveDenham]

OK, I will try not to throw SAS/STAT 13.1, or even 12.1 syntax and methods at you since you are at SAS9.2, and both of my machines are 9.3 and higher.

I wouldn't go into the spatial coordinate part yet, but there is definitely a problem with factor1 and factor2 being completely confounded with datenum.  Could you compare the results from the following two analyses:

proc glimmix data=hen;

  class datenum treatment;

  model count/total = treatment datenum/ link=logit dist=binomial solution ddfm=kr(firstorder);

  random datenum/residual subject=treatment type=ar(1);

  Lsmeans treatment/pdiff;

  run ;

proc glimmix data=hen; /*For this to work, the data must be sorted by factor1 and factor2 within treatment*/

  class treatment;

  model count/total = treatment factor1|factor2 / link=logit dist=binomial solution ddfm=kr(firstorder);

  random factor1 factor2/ subject=treatment ;

  Lsmeans treatment/pdiff;

  run ;

Now, I don't know that this will accomplish much at this point.  My hope is that, due to the confounding, we can end up replacing factor1 and factor2 with datenum, and get a reasonable result.  Or vice versa.

Steve Denham

[mgerritsen]

Hi Steve,

So you first model works. I have the subject-observation per subject structure i expect. Furthermore both datenum and treatment ar specified as having a significant effect. Though the DF values are a bit strange, like 69,8 or 17.23, etc. I have also attempted this model using method laplace and removing bothd the ddfm and residual options. The output is nearly equal only the DF values seem better, 125.

For the other model i sorted the data using:

  proc sort data=hen;

  by treatment factor1 factor2;

  run;

I hope this was what you intended. The model itself provided the expected structure furthermore all variables had an output but were non-significant. However, the DF for all tests was 1.

Marielle

[SteveDenham]

Drop the ddfm=kr(firstorder) option.  The correlation between time points/factors is eating up the degrees of freedom.  However, that correlation needs to be accounted for in calculating the standard errors, so add EMPIRICAL as an option in the PROC GLIMMIX statement.

Steve Denham

[mgerritsen]

Steve,

I believe you did it! My output describes a model using g-side cov. parameters only  with 6 subjects- 26 observations per subject. Furthermore, all variables have an output with treatment, factor2 and the factor1*factor2 being significant.

So if i understand the application of emperical correctly (in a laymans description): because factor1, factor2 and datenum were confounded both factors could replace datenum in the random statement to indicate the repetitions. By sorting both factors in treatment  these repetitions were physically grouped in the dataset already, creating the correct covaraite structure for applying the emperical option. Therefore there was no need to indicate factor1 and factor2 in the class statement or to specify type=ar(1).  The emperical statement corrects the standard errors for the bias caused by the loss of the (i'm not sure how to call this: covariance structure) specified by type=ar(1) or  the ddfm options?

What I was wondering in  Stroup's Generalized Linear Mixed Models they indicate that earlier emperical functions could cause a bias themselves if the dataset is to small and they also discuss an correction presented nearly ten years ago. My dataset it not that large and i wonder if this correction is default in SAS or if i should account for it somehow?

Mariëlle

[SteveDenham]

I am glad it is working and is interpretable.

One thing about the sandwich/empirical estimators is that they can be biased, but there are several options.  In particular, the empirical=mbn option goes a long way to eliminate that bias, so you might try that.

Steve Denham

[mgerritsen]

Hi Steve,

I wanted to thank you again for your help. I'm referencing your input in my final report and i wanted to know if you would like to be specified in a certain manner? Currently i have used your name (username you use on the website) and a link to this thread.

please let me know,

Marielle

[SteveDenham]

That is fine with me.  Let me know where it gets published (if it is in a journal).

Steve Denham