2

I have a dichotomous outcome (gallstones/no gallstones) and an ordinal predictor variable consisting of four classes (body mass index <25(ref.), 25-30, 30-35, 35-45). I want to perform a test for trend using R.

When running the glm function only the bmi 35-45 is significantly associated to the outcome (gallstones):

glm(formula = gallstones ~ bmi, family = binomial, data = m)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-0.9282  -0.4058  -0.3624  -0.3624   2.3477  

Coefficients:
               Estimate Std. Error z value Pr(>|z|)    
(Intercept)     -2.6902     0.1476 -18.223  < 2e-16 ***

bmi(25,30]   0.2347     0.2131   1.102   0.2706    
bmi(30,35]   0.6108     0.3291   1.856   0.0635 .  
bmi(35,45]   2.0712     0.4915   4.214 2.51e-05 ***

To examine if there is a dose-response relationship with a linear trend with rising bmi I perform the analysis of variance with the chi square test:

anova(inc_test_3, test="Chisq")
Analysis of Deviance Table

Model: binomial, link: logit

Response: gallstones

Terms added sequentially (first to last)


          Df Deviance Resid. Df Resid. Dev Pr(>Chi)   
NULL                      1489     810.05            
bmi       3   15.892      1486     794.16 0.001193 **

It seems that there is significant difference between the obtained estimates from the ordinal bmi predictor.

Is this the correct way to do a test of trend on an ordinal variable?

Scortchi - Reinstate Monica
  • 27,560
  • 8
  • 81
  • 248
Daniel
  • 21
  • 1
  • 3
  • 1
    In your analysis of deviance you haven't tested for trend at all, but simply compared the null model with one having a four-level categorical predictor. [Orthogonal polynomial coding](http://www.ats.ucla.edu/stat/r/library/contrast_coding.htm#ORTHOGONAL) is useful for examining trends, & in fact is the default when you tell R the predictor is ordinal (see `?ordered`). – Scortchi - Reinstate Monica Apr 13 '15 at 11:59
  • Thank you - I was not aware of the ordered function in R. However, what if we concider the variable to be a factor? I am not sure if it is all right to assume that the jumps in this BMI variable are completely ordinal. I have created the BMI variable with the default cut function in R. – Daniel Apr 13 '15 at 13:31
  • (1) Not sure what you mean by "consider the variable to be a factor". Orthogonal polynomial coding doesn't produce a substantively different model from e.g. reference level coding. (2) It doesn't seem at all likely that the chance of having a gallstone jumps at BMIs of 25, 30, 35, & 45, or that it's constant between those BMIs. See [here](http://stats.stackexchange.com/questions/68834/). If you have BMI then use it as a predictor rather than cutting it, allowing for curvilinearity with polynomials or (probably better) splines. – Scortchi - Reinstate Monica Apr 13 '15 at 13:41
  • Sorry for jumping to conclusions to soon:) What i mean is the scenario where the continuous variable has no association on the outcome but I assume it to have an association based on previous analyses. Therefore, I would like to cut the variable into a multi-class categorical variable. Would that be totally wrong? – Daniel Apr 13 '15 at 14:38
  • Sorry - I don't understand the "therefore". Are you saying that when you treated the predictor as having a *linear* relationship to the response (on the logit scale) you found no (or only a weak) association, & you therefore binned the predictor to allow for a curvilinear relationship? If so, that's not *totally* wrong, but not a very good way of modelling it - read the post I linked to. – Scortchi - Reinstate Monica Apr 13 '15 at 14:58
  • @Scortchi The orthogonal polynomial coding worked perfect! You want to post it as an answer so i can rate it? – Daniel Apr 18 '15 at 08:42

0 Answers0