class: center, middle, inverse, title-slide # Ordinal Regression ### Yue Jiang ### STA 210 / Duke University / Spring 2023 --- ### Licorice and post-op sore throat <img src="img/licorice.jpg" width="80%" style="display: block; margin: auto;" /> --- ### Licorice and post-op sore throat <img src="ordinal_files/figure-html/unnamed-chunk-2-1.png" style="display: block; margin: auto;" /> --- ### Licorice and post-op sore throat <img src="ordinal_files/figure-html/unnamed-chunk-3-1.png" style="display: block; margin: auto;" /> --- ### Licorice and post-op sore throat <img src="ordinal_files/figure-html/unnamed-chunk-4-1.png" style="display: block; margin: auto;" /> --- ### Ordinal data In this case, our outcome is *ordered* (and categorical). Although we have "numbers" on a pain scale from 0 to 10, these numbers don't have a linear relationship - a 4 isn't necessarily twice as painful as a 4; a 3 isn't three times as painful as 1. However, these data *are* ordered. We do know that a 4 is more painful than a 3, which is more painful than a 2, etc. .question[ What are some potential pitfalls of using an ordinary least squares regression? How about using a logistic regression approach (how would we even do this)? ] --- ### A cumulative link model We might consider an outcome `\(Y\)` that looks at all outcomes at once. For `\(j\)` total ordered categories, we might model the **cumulative probability** for observation `\(i\)`: `\begin{align*} \gamma_{ij} &= P(Y_i \le j)\\ &= P(Y_i = 1) + P(Y_i = 2) + \cdots + P(Y_i = j). \end{align*}` Note that `\(\gamma_{ij}\)` is limited to values from 0 to 1, as it is a probability. We might consider a model `\begin{align*} logit(\gamma_{ij}) = \beta_{0;j} + \beta_1x_{i1} + \cdots + \beta_px_{ip} \end{align*}` .question[ Take a look at this model - what do the `\(\beta_{0;j}\)` terms mean/represent? ] --- ### A cumulative link model `\begin{align*} \gamma_{ij} &= P(Y_i \le j)\\ &= P(Y_i = 1) + P(Y_i = 2) + \cdots + P(Y_i = j). \end{align*}` The `\(\beta_{0;j}\)` terms are constants representing the "baseline" value for each category (on the logit scale). `\begin{align*} logit(\gamma_{ij}) = \beta_{0;j} + \beta_1x_{i1} + \cdots + \beta_px_{ip} \end{align*}` .question[ What is the interpretation of `\(\beta_1\)` through `\(\beta_p\)`? What are we implicitly assuming? ] --- ### A cumulative link model In this case, we have the same covariate relationships `\(\beta_1, \cdots, \beta_p\)` across **all** of the categories. This means that `\(\beta_k\)` represents the conditional change in (transformed) cumulative probabilities given a 1 unit difference in `\(X_{ik}\)`. --- ### Ordered logistic regression The ordered logistic regression model is a cumulative link model that assumes a logit transformation of the cumulative probabilities: `\begin{align*} logit(\gamma_{ij}) &= \beta_{0;j} + \beta_1x_{i1} + \cdots + \beta_px_{ip}\\ log\left(\frac{\gamma_{ij}}{1 - \gamma_{ij}}\right) &= \beta_{0;j} + \beta_1x_{i1} + \cdots + \beta_px_{ip}\\ log\left(\frac{P(Y_i \le j)}{P(Y_i > j)}\right) &= \beta_{0;j} + \beta_1x_{i1} + \cdots + \beta_px_{ip} \end{align*}` Exponentiating, we have `\begin{align*} \frac{P(Y_i \le j)}{P(Y_i > j)} &= \beta_{0;j} + \beta_1x_{i1} + \cdots + \beta_px_{ip} \end{align*}` .question[ What is the outcome here? How might we interpret `\(\exp(\beta_{0;j})\)`? How might we interpret the `\(\beta\)` terms here? ] --- ### The proportional odds assumption Remember that we have only one `\(\beta_k\)` term for each predictor across *all* categories. This implies that changes in `\(X_k\)` have the same conditional relationship with odds of being in category 1 vs. 2, 6 vs. 7, or any `\(j-1\)` vs. `\(j\)`. .question[ When might this be a reasonable assumption? When might this assumption be violated? How might we modify the model in the case that this assumption does not hold? How might we gut-check this assumption using the data? ] --- ### Fitting an ordinal logistic model One quirk about the function used in R is that it calculates `\begin{align*} logit(\gamma_{ij}) &= \beta_{0;j} - x_{i1}\eta_1 - \cdots - x_{ip}\eta_p \end{align*}` with `\(\eta_k = -\beta_k\)`. --- ### Fitting an ordinal logistic model ```r library(MASS) m1 <- polr(factor(postOp4hour_throatPain) ~ treat + preOp_calcBMI, data = dat) summary(m1) ``` ``` ## Call: ## polr(formula = factor(postOp4hour_throatPain) ~ treat + preOp_calcBMI, ## data = dat) ## ## Coefficients: ## Value Std. Error t value ## treat -1.11716 0.28973 -3.8558 ## preOp_calcBMI 0.01853 0.03272 0.5662 ## ## Intercepts: ## Value Std. Error t value ## 0|1 0.7007 0.8624 0.8125 ## 1|2 1.5582 0.8695 1.7921 ## 2|3 2.5771 0.8901 2.8952 ## 3|4 3.6860 0.9485 3.8862 ## 4|5 4.1043 0.9904 4.1439 ## 5|6 4.8097 1.1081 4.3406 ## 6|7 5.5089 1.3139 4.1928 ## ## Residual Deviance: 476.6836 ## AIC: 494.6836 ## (2 observations deleted due to missingness) ``` --- ### Fitting an ordinal logistic model ```r exp(coef(m1)) ``` ``` ## treat preOp_calcBMI ## 0.3272071 1.0187007 ``` .question[ How might we interpret the `\(\exp(\eta)\)` term corresponding to treatment? ] --- ### Fitting an ordinal logistic model ```r exp(coef(m1)) ``` ``` ## treat preOp_calcBMI ## 0.3272071 1.0187007 ``` Patients who receive licorice treatment have approximately 1/3 the odds of having the next higher pain category (e.g., 6 vs. 5, or 2 vs. 1, etc.) compared to patients receiving placebo, while controlling for BMI. For every `\(kg/m^2\)` increase in BMI, the odds of being in the next higher pain category is multiplied by approx. 1.02, while controlling for treatment. --- ### Fitting an ordinal logistic model ```r exp(confint(m1)) ``` ``` ## 2.5 % 97.5 % ## treat 0.1830655 0.5719223 ## preOp_calcBMI 0.9554557 1.0866260 ``` --- ### Fitting an ordinal logistic model ```r round(head(m1$fitted.values), 3) ``` ``` ## 0 1 2 3 4 5 6 7 ## 1 0.770 0.118 0.069 0.029 0.005 0.005 0.002 0.002 ## 2 0.799 0.105 0.059 0.025 0.004 0.004 0.002 0.002 ## 3 0.789 0.109 0.062 0.026 0.005 0.004 0.002 0.002 ## 4 0.784 0.111 0.064 0.027 0.005 0.005 0.002 0.002 ## 5 0.778 0.114 0.066 0.028 0.005 0.005 0.002 0.002 ## 6 0.761 0.121 0.072 0.030 0.005 0.005 0.003 0.003 ```