[R] GEE with own link function

[Johanna BRANDT]

Hello,

I want to fit a GEE with a user-defined link function.

For the user-defined link-function I still read 
http://finzi.psych.upenn.edu/R/Rhelp01/archive/6555.html and 
http://finzi.psych.upenn.edu/R/Rhelp02a/archive/25727.html.

Only for testing purposes I added a new link function 
(corlogit) in make.link (as well as in binomial) with 
exactly the same code as logit before using my code.

I tried it with glm() and it works without any problem (I 
get the same results using 'binomial(link="logit")' or 
'binomial(link="corlogit")').


But with library(geepack) (Version: 0.2-10) I get different 
results for logit and test.logit as link functions. First 
there was an error, so that I had to modify 'geese.fit' (I 
only added test.logit where I found logit).

With library(gee) (Version: 4.13-10) R exits without error 
when fitting a GEE with binomial(link="corlogit"). Because 
of the error 'unkown link' I modified 'gee' the same way as 
for geepack().


I never changed the original application files, but I wrote 
"own programmes" (using the original code and just adding 
corlogit in the list of the link functions). Are there any 
other functions where I had to add or modify something? What 
else can I do?

Thank you for your help!
Johanna Brandt


(I'm using R 2.0.1 under Windows 2000)


## Example for geese() from the R-Help #####################


I took the example from the help:
 > data(ohio)
 > summary(geese(resp ~ age + smoke + age:smoke, id=id, 
data=ohio,
+              family=binomial(link="logit"), corstr="exch", 
scale.fix=TRUE))

Call:
geese(formula = resp ~ age + smoke + age:smoke, id = id, 
data = ohio,
     family = binomial(link = "logit"), scale.fix = TRUE, 
corstr = "exch")

Mean Model:
  Mean Link:                 logit
  Variance to Mean Relation: binomial

  Coefficients:
                estimate     san.se        wald          p
(Intercept) -1.90049529 0.11908698 254.6859841 0.00000000
age         -0.14123592 0.05820089   5.8888576 0.01523698
smoke        0.31382583 0.18575838   2.8541747 0.09113700
age:smoke    0.07083184 0.08852946   0.6401495 0.42365667

Scale is fixed.

Correlation Model:
  Correlation Structure:     exch
  Correlation Link:          identity

  Estimated Correlation Parameters:
       estimate     san.se     wald p
alpha 0.354531 0.03582698 97.92378 0

Returned Error Value:    0
Number of clusters:   537   Maximum cluster size: 4

 >
 > ## Korrigiert
 > summary(fit.korr <- geese(resp ~ age + smoke + age:smoke, 
id=id, data=ohio,
+              family=binomial(link="corlogit"), 
corstr="exch", scale.fix=TRUE))

Call:
geese(formula = resp ~ age + smoke + age:smoke, id = id, 
data = ohio,
     family = binomial(link = "corlogit"), scale.fix = TRUE, 
corstr = "exch")

Mean Model:
  Mean Link:                 corlogit
  Variance to Mean Relation: binomial

  Coefficients:
                estimate     san.se        wald          p
(Intercept) -1.12581067 0.06344341 314.8891093 0.00000000
age         -0.07680433 0.03128947   6.0252497 0.01410264
smoke        0.17083868 0.10162807   2.8258236 0.09275930
age:smoke    0.03672858 0.04872412   0.5682249 0.45096515

Scale is fixed.

Correlation Model:
  Correlation Structure:     exch
  Correlation Link:          identity

  Estimated Correlation Parameters:
        estimate     san.se    wald p
alpha 0.3545883 0.03583136 97.9315 0

Returned Error Value:    0
Number of clusters:   537   Maximum cluster size: 4



## Example for gee() from the R-Help #######################

if(require(MASS)) {
data(OME)
## not fully appropriate link for these data.
(fm.korr <- gee(cbind(Correct, Trials-Correct) ~ Loud + Age 
+ OME, id = ID,
            data = OME, family = binomial(link="corlogit"), 
corstr = "exchangeable"))
summary(fm.korr)
}