1
$\begingroup$

I am trying to tie the odds ratio from a 2x2 cross classification table to the intercepts of a logistic regression on those 2 variables. I have a cross classification table that produces 2 odds ratios and the results of a logistic regression of PLACE3 ~ VIOL should produce intecepts should match the odds ratio of the contingency table. i.e. Odds ratio = exp(intercepts) BUT the POLR package is not producing the correct intercepts.

Here is the data. In the logistic regression PLACE3 is the outcome and VIOl is the independent variable. You can see the PLACE3 vs. VIOL contingency table below and the logistic regression of PLACE3 ~ VIOL. The odds ratios in the contingency table 1.79 and 3.1 are correct but the polr function seems off. Any thoughts on why exp(summary(m)$zeta) does not produce 1.79 and 3.1?

For reference this is from Lemeshow's Applied Logisitic Regression book page 274.

library(data.table)
aps <- fread('http://www.umass.edu/statdata/statdata/data/aps.dat')
colnames(aps) = c("ID","PLACE","PLACE3","AGE","RACE","GENDER","NEURO","EMOT","DANGER","ELOPE","LOS","BEHAV","CUSTD",
                    "VIOL")
head(aps)

Here is a cross classification table of PLACE3 vs. VIOl variables

table(aps$PLACE3,aps$VIOL) 
      0   1
  0  80 179
  1  26 104
  2  15 104

using PLACE3 = 0 as the reference the 2 odds ratios from the contingency table are

(104*80)/(179*26)  #1.79
(104*80)/(179*15)  #3.10

These odds ratios should be the same as exponentiating the slope coefficients from a logistic model PLACE3 ~ VIOL which is below

aps$constant = rep(1,dim(aps)[1])
m <- polr(as.factor(PLACE3) ~ constant + as.factor(VIOL), data = aps, Hess=TRUE,model=TRUE,method = c("logistic"))
summary(m)

> summary(m)
Call:
polr(formula = as.factor(PLACE3) ~ constant + as.factor(VIOL), 
    data = aps, Hess = TRUE, model = TRUE, method = c("logistic"))

Coefficients:
                  Value Std. Error t value
as.factor(VIOL)1 0.8454     0.2112   4.003

Intercepts:
    Value  Std. Error t value
0|1 0.6869 0.1884     3.6464 
1|2 1.8608 0.2032     9.1557 

Residual Deviance: 1031.75 
AIC: 1037.75 

But you can see the exponentiation of the zeta vector is not 1.79 and 3.10

exp(summary(m)$zeta)

> exp(summary(m)$zeta)
     0|1      1|2 
1.987495 6.429049 
$\endgroup$
2
  • $\begingroup$ Why do you think those odds ratios are appropriate for a proportional odds model? Also you are ignoring the coefficient for VIOL. Also you do not need to add a constant and as you can see polr ignored it. $\endgroup$
    – mdewey
    Jul 13, 2016 at 16:18
  • $\begingroup$ I am not sure what you mean when you say appropriate. For a single dichotomous independent variable the odds ratio calculated manually should match that of the logistic regression. I am not doing model fitting or anything here. Just pure mechanics. I see that the polr does not include the constant even though one is provided. Can you advise how to make polr recognize the constant? i have constant + factor(VIOL). $\endgroup$ Jul 13, 2016 at 16:29

1 Answer 1

1
$\begingroup$

The relevant odds ratios for a proportional model are (80 * (104 + 104)) / ((26 + 15) * 179) and ((80 + 26) * 104) / ((179 + 104) * 15) so whatever model you are trying to fit here it is not the proportional odds model.

The documentation for polr is perhaps a little sparse so you might want to look at something like the article by Bender and Benner here entitled "Calculating ordinal regression models in SAS and S-Plus". Despite the title it is relevant for R and also more generally.

But as I mention whatever you are trying to do the proportional odds model does not seem to be what you want.

Are you by any chance wanting to model the other way round? glm(VIOL ~ PLACE3, family = binomial, data = aps)

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.