Translating PROC GLIMMIX Models#586
Conversation
|
Thank you for this nice work @Franosei ! This is so very useful for our users. |
|
|
||
| This scenario translates a SAS `PROC GLIMMIX` binary random-intercept and random-slope model into R. The outcome is asthma control, coded as `AVAL`, with repeated observations clustered within `USUBJID`. | ||
|
|
||
| The model estimates the effect of treatment arm, scaled time, scaled age, sex, race, baseline severity, and scaled adherence on the odds of asthma control. The SAS model uses `PROC GLIMMIX` with independent subject-level random intercept and random slope terms. The R implementation uses `lme4::glmer()` with the matching `||` random-effects structure. |
There was a problem hiding this comment.
could be: "with independent (uncorrelated) subject-level random intercept and random slope terms."
|
|
||
| This scenario translates a SAS `PROC GLIMMIX` binary random-intercept model for timepoint-specific responder rates into R. The outcome is `AVAL`, where `1` indicates responder status. Repeated visits are clustered within `USUBJID`. | ||
|
|
||
| The SAS model uses `PROC GLIMMIX` with a binary logit model, treatment-by-visit interaction, and subject-level random intercept. The R implementation uses `glmmTMB` with the same conditional GLMM structure. `geepack` is used as a marginal sensitivity analysis. |
There was a problem hiding this comment.
Why did you want the geepack model here? The SAS and glmmTMB models seem pretty similar.
| emm_r <- emmeans::emmeans( | ||
| fit_r, | ||
| ~ ARMCD, | ||
| type = "response" | ||
| ) |
There was a problem hiding this comment.
For this one, to mimic SAS (and get the rate), you could do:
emmeans::emmeans(fit_r, ~ ARMCD, type = "response", offset = 0)
| pairs( | ||
| emm_r, | ||
| reverse = TRUE, | ||
| adjust = "none", | ||
| type = "response" | ||
| ) |
There was a problem hiding this comment.
This one SAS (I think) reports the ratio as the difference on the log scale, so:
emm_r <- emmeans::emmeans(
fit_r, ~ ARMCD)
pairs(
emm_r,
reverse = TRUE,
adjust = "none")
| scenario_09$AVAL_ORD <- ordered( | ||
| scenario_09$AVAL, | ||
| levels = c(1, 2, 3, 4, 5) | ||
| ) |
There was a problem hiding this comment.
I think this is the source of the sign discrepancy with SAS. Presumably, it orders the categories from last to first (like it does with CLASS variables as predictors). So maybe reverse the order of the levels:
levels = c(5,4,3,2,1)
| Odds_Ratio = exp(fixed_r_sas_scale[, "Estimate"]), | ||
| Lower_95_CI = exp( | ||
| fixed_r_sas_scale[, "Estimate"] - | ||
| 1.96 * fixed_r_sas_scale[, "Std. Error"] |
There was a problem hiding this comment.
A bit picky. Instead of 1.96 you could be more precise:
qnorm(p = 1-0.05/2)
| | 4 | 4 | 0.135 | 0.247 | 0.112 | | ||
| | 4 | 5 | 0.079 | 0.212 | 0.132 | | ||
|
|
||
| The likelihood, AIC, random-intercept variance, and fixed-effect estimates were aligned after reversing the R coefficient signs to match the SAS cumulative-logit direction. The active arm shifted probability away from lower `AVAL` categories and towards higher `AVAL` categories. This is consistent with the SAS cumulative-logit odds ratios and the R category probability estimates. |
There was a problem hiding this comment.
See comment above about reversing the response levels before fitting the model.
|
|
||
| - R uses GEE as a marginal sensitivity analysis for a conditional GLIMMIX model. | ||
|
|
||
| - Ordinal cumulative-logit models use different sign conventions. |
There was a problem hiding this comment.
See other comments. This is probably just SAS's default. Presumably, if you don't specify an ordering in R it would be first to last, whereas SAS uses last to first.
3 participants
Summary
The article compares SAS and R implementations across common GLIMMIX use cases, including binary GLMMs, repeated-measures models, structured covariance models, count models, ordinal outcomes, and post-fit inference.
Main changes
R vs SAS: Translating PROC GLIMMIX Models.PROC GLIMMIXand R packages.lme4,glmmTMB,nlme,mmrm,geepack,emmeans, andordinal.The article focuses on whether the same statistical model, estimand, covariance structure, and post-fit interpretation can be reproduced or closely approximated in R.
Some scenarios are direct translations, while others are described as approximate or sensitivity-based where the R implementation does not exactly match the SAS
PROC GLIMMIXspecification.