Hi,
I am trying to fit a ordinal GEE model using ordgee {geepack}. In order to check
the validity of the function, I specified the correlation structure as
independence (i.e. constr = "independence") and compared the result
with that using polr {MASS}.
Because a GEE model with an independent working correlation structure is
equivalent to an ordinary GLM model, we would expect the same results from both
functions (there might be small differences due to dissimilarity in
approximation methods). However, I got a very different results from the two
approaches.
Here is my R codes and ouputs for the two approaches:
1.)
library(geepack)summary(ordgee(UACR~HR+SP+bmi+diabetes+age+HDL+ht.med+male+SCHV+SCRV+smoke+statin+AP,id=x18,mean.link="logit",corstr="independence"))
Call:
ordgee(formula = UACR ~ HR + SP + bmi + diabetes + age + HDL + ht.med + male +
SCHV + SCRV + smoke + statin + AP, id = x18, mean.link = "logit",
corstr = "independence")
Mean Model:
Mean Link: logit
Variance to Mean Relation: binomial
Coefficients:
estimate san.se wald p
Inter:0 -1.5843906608 2.408716706 4.326663e-01 0.510683514
Inter:1 -2.6094525748 2.401539294 1.180645e+00 0.277224675
Inter:2 -3.6006918946 2.410172287 2.231905e+00 0.135187061
HR 0.0491831800 0.018913783 6.762020e+00 0.009311831
SP -0.0287228998 0.011465183 6.276168e+00 0.012237244
bmi 0.0363585224 0.029345236 1.535102e+00 0.215348334
diabetes 0.2557417872 0.464391207 3.032742e-01 0.581837057
age 0.0188629344 0.022226239 7.202555e-01 0.396060127
HDL -0.0210835675 0.014822770 2.023156e+00 0.154916842
ht.med 0.0089636059 0.383824776 5.453804e-04 0.981368392
male 1.2806017705 0.473921757 7.301547e+00 0.006889526
SCHV -0.0008670472 0.005862042 2.187699e-02 0.882414803
SCRV 0.5667065765 0.879731084 4.149699e-01 0.519457682
smoke -0.6064323366 0.352556881 2.958737e+00 0.085414655
statin -0.0008867941 0.408462933 4.713465e-06 0.998267753
AP 0.1159595043 0.036258436 1.022809e+01 0.001383178
Scale is fixed.
Correlation Model:
Correlation Structure: independence
Returned Error Value: 0
Number of clusters: 285 Maximum cluster size: 6
2.)
library(MASS)
summary(polr(UACR~HR+SP+bmi+diabetes+age+HDL+ht.med+male+SCHV+SCRV+smoke+statin+AP,method="logistic"))
Re-fitting to get Hessian
Call:
polr(formula = UACR ~ HR + SP + bmi + diabetes + age + HDL + ht.med + male +
SCHV + SCRV + smoke + statin + AP, method = "logistic")
Coefficients:
Value Std. Error t value
HR 0.0004674288 0.009690505 0.048235757
SP -0.0106449514 0.005305501 -2.006398816
bmi 0.0298158386 0.015501658 1.923396794
diabetes 0.6563205176 0.251591908 2.608671009
age 0.0034051085 0.010467630 0.325298895
HDL -0.0153783914 0.006227563 -2.469407699
ht.med -0.1146770998 0.194130303 -0.590722303
male 0.1798842955 0.227476463 0.790782014
SCHV 0.0022770371 0.002675591 0.851040793
SCRV 0.0010954808 0.440261143 0.002488252
smoke 0.0207174113 0.165287870 0.125341390
statin 0.4320833418 0.206676105 2.090630370
AP 0.0419551861 0.015799215 2.655523401
Intercepts:
Value Std. Error t value
0|1 -0.9883 1.1258 -0.8778
1|2 0.1712 1.1246 0.1522
2|3 1.3306 1.1258 1.1819
Residual Deviance: 1376.931
AIC: 1408.931
Is there anything wrong with my code? I guess my real question should be, could
I safely use the ordgee program to fit a ordinal GEE model with exchangeable
correlation structure?
Thanks very much for your time and professional opinion,
High Seng
[[alternative HTML version deleted]]