[R] gee() and geeglm() errors

[Lilia Verchinina]

Dear R users,

I have a question regarding the gee() (from package gee) and geeglm() (from
package geepack). Here is the head() from my data:

R> head(compl.samples_long)
    UserID Intervention PHQ_base compl_bin Sex EFE Neuro depr0 subject time
PHQ
1.1  10004            1        6         2   2  20    23     2       1    1
  5
2.1  10009            1        0         2   2   6    13     2       2    1
 11
3.1  10013            1        2         2   2  30    31     1       3    1
  1
4.1  10015            2        4         1   2  48    31     1       4    1
  3
5.1  10016            1        1         2   2  10    32     2       5    1
 12
7.1  10020            1        5         2   1  47    23     1       7    1
  5
    stress hours errors id
1.1      1    70      1  1
2.1      0    50      1  2
3.1      1    60      2  3
4.1      1    55      1  4
5.1      0    70      1  5
7.1      0    85      2  7



When I try to use geeglm(), I get error:

R>
geeInt=(geeglm(PHQ~as.factor(compl_bin)+Neuro+PHQ_base+as.factor(depr0)+EFE+as.factor(Sex)+hours+stress+as.factor(errors),family
= poisson,data=compl.samples_long ,id=subject, corst="exchangeable"))


 Error in geese.fit(xx, yy, id, offset, soffset, w, waves = waves, zsca,  :

  nrow(zsca) and length(y) not match

The text of the error message does not seem helpfull, since whatever typo
you do, you get the same message.

However, I do not get any complains if I use gee() function from gee
package, but I get a strange working correlation matrix:

R>
geeInt=(gee(PHQ~as.factor(compl_bin)+Neuro+PHQ_base+as.factor(depr0)+EFE+as.factor(Sex)+hours+stress+as.factor(errors),family
= poisson,data=compl.samples_long ,id=subject, corst="exchangeable"))

R> summary(geeInt)

Coefficients:
                              Estimate  Naive S.E.    Naive z Robust S.E.
(Intercept)               0.7538052712 0.181793256  4.1464974 0.194869157
as.factor(Intervention)2 -0.1305672837 0.070750465 -1.8454619 0.069618104
Neuro                     0.0225507686 0.004412580  5.1105631 0.004767788
PHQ_base                  0.0344067686 0.017511477  1.9648125 0.015846200
as.factor(depr0)2        -0.0456275666 0.071731191 -0.6360910 0.072515391
EFE                       0.0004936796 0.003356067  0.1471006 0.003139989
as.factor(Sex)2           0.0571359598 0.071662335  0.7972941 0.073569880
hours                     0.0096543634 0.001757295  5.4938782 0.001958358
stress                   -0.0763068749 0.060170447 -1.2681786 0.054164383
as.factor(errors)2       -0.3436635241 0.080389082 -4.2750025 0.079987789
                           Robust z
(Intercept)               3.8682636
as.factor(Intervention)2 -1.8754789
Neuro                     4.7298174
PHQ_base                  2.1712946
as.factor(depr0)2        -0.6292122
EFE                       0.1572234
as.factor(Sex)2           0.7766216
hours                     4.9298254
stress                   -1.4088017
as.factor(errors)2       -4.2964499


Working Correlation
     [,1] [,2] [,3] [,4]
[1,]    1    0    0    0
[2,]    0    0    0    0
[3,]    0    0    0    0

I also tried to use SAS proc genmode and did not get any complains but the
p-values are different from that of based on gee() output. If you have any
ideas why geeglm does not work and why gee() gives such strange Working
correlation, it will be very helpfull!!

Thank you!