Fitting a random effect model for overdispersed grouped binary data

[S_pera]

I need to fit a random effect model for the following SAS dataset. (I have already fitted the standard linear logistic full model.) Now how can I fit the linear logistic model with random effect? 

data exa;

input vty ext y n;

cards;

1 1 10 39

1 1 22 62

1 1 23 81

1 1 26 51

1 1 17 39

1 2 5 6

1 2 53 74

1 2 55 72

1 2 32 51

1 2 46 79

1 2 10 13

2 1 8 16

2 1 10 30

2 1 8 28

2 1 25 45

2 1 0 4

2 2 3 12

2 2 22 41

2 2 15 30

2 2 32 51

2 2 3 7

;

proc logistic;

class vty ext/param=reference ref=first;

model y/n=vty ext vty*ext /scale=none ;

run;

[StatDave]

You probably will want to consider using the binomial cluster model or the beta binomial model, both of which are mentioned in the Overdispersion section of this note. As noted there, these models are discussed and illustrated in the example titled "Modeling Mixing Probabilities" in the PROC FMM documentation.

[S_pera]

Dear StatDave

I am just interested in fitting the random effects model. Thank you StatDave.

rgds

S_pera

[StatDave · Accepted solution]

Then you need to add an observation identifier in the data by including this in the DATA step:

id=_n_;

You can then fit the random effects model in GLIMMIX as follows which also fits the beta binomial, binomial cluster, and GEE models - all of which are ways to deal with the overdispersion. The cumulated results in the PREDS data set show that they all produce similar predicted probabilities.

proc gee data=exa;
   class vty ext id;
   model y/n = vty ext vty*ext / dist=bin;
   repeated subject=id;
   output out=preds pred=pgee;
run;
proc glimmix data=preds;
   class vty ext id;
   model y/n = vty ext vty*ext / dist=bin s;
   random intercept / subject=id;
   output out=preds pred(ilink)=pglim;
run;
proc fmm data=preds;
   class vty ext;
   model y/n =  / dist=binomcluster;
   probmodel vty ext vty*ext;
   output out=preds pred(overall)=pbclus;
run;
proc fmm data=preds;
   class vty ext;
   model y/n = vty ext vty*ext / dist=betabin;
   output out=preds pred(overall)=pbbin;
run;

[S_pera]

Dear StatDave

Thank you very much. It works.

rgds

S_pera