Ordinal regression model vs. logistic regression model
Elias Eythorsson
November 26th, 2023
1 Ordinal regression model vs. logistic regression model
This script assumes that you have already run mgus_bm_prediction_20231125. This analysis has been completed in response to the revision request by the Annals of Internal Medicine, and will test the hypothesis that an ordinal regression is more precise than logistic regression, as the underlying statistical model for the iStopMM prediction model. First we fit the ordinal regression model, that ‘borrows’ information from more common outcome categories (0-4% with n=455 and 5-9% with n = 430) to increase statistical precision on less common outcome categories (10-14% with n = 91 and >15% with n = 67) that are of clinical interest.
fit_pred <- lrm(outcome ~ rcs(igg, 4) + rcs(iga, 4) + rcs(igm, 4) + rcs(m_protein, 4) +
rcs(flc_ratio, 4) + isotype, data = df_pred, x = TRUE, y = TRUE)As seen below the likelihood ratio chi-square test is 401.11 on 18 degrees of freedom, with an r squared of 0.358 and a concordance index of 0.766. Importantly, the coefficient estimate of the intercept for y >= 10-14% is 1.5889 with a standard error of 0.8686
fit_pred## Logistic Regression Model
##
## lrm(formula = outcome ~ rcs(igg, 4) + rcs(iga, 4) + rcs(igm,
## 4) + rcs(m_protein, 4) + rcs(flc_ratio, 4) + isotype, data = df_pred,
## x = TRUE, y = TRUE)
##
##
## Frequencies of Responses
##
## 0-4% 5-9% 10-14% >15%
## 455 430 91 67
##
## Model Likelihood Discrimination Rank Discrim.
## Ratio Test Indexes Indexes
## Obs 1043 LR chi2 401.11 R2 0.358 C 0.766
## max |deriv| 4e-08 d.f. 18 R2(18,1043)0.307 Dxy 0.533
## Pr(> chi2) <0.0001 R2(18,882.4)0.352 gamma 0.533
## Brier 0.085 tau-a 0.335
##
## Coef S.E. Wald Z Pr(>|Z|)
## y>=5-9% 4.2561 0.8788 4.84 <0.0001
## y>=10-14% 1.5889 0.8686 1.83 0.0674
## y>=>15% 0.1966 0.8704 0.23 0.8213
## igg -0.3544 0.0937 -3.78 0.0002
## igg' 1.5930 0.4062 3.92 <0.0001
## igg'' -4.4209 1.1912 -3.71 0.0002
## iga -0.1330 0.2827 -0.47 0.6380
## iga' 2.7083 1.7328 1.56 0.1181
## iga'' -6.0916 3.9008 -1.56 0.1184
## igm 1.0890 0.7754 1.40 0.1602
## igm' -16.7238 8.5061 -1.97 0.0493
## igm'' 34.2726 17.1540 2.00 0.0457
## m_protein -0.6602 0.1832 -3.60 0.0003
## m_protein' 46.0506 11.3215 4.07 <0.0001
## m_protein'' -57.9203 14.3889 -4.03 <0.0001
## flc_ratio -2.1258 0.3555 -5.98 <0.0001
## flc_ratio' 38.5632 6.5415 5.90 <0.0001
## flc_ratio'' -90.7608 15.5093 -5.85 <0.0001
## isotype=IgA -0.2766 0.2610 -1.06 0.2892
## isotype=IgG -0.5350 0.2038 -2.63 0.0087
## isotype=Light-chain -0.5640 0.3257 -1.73 0.0833
## It should be noted that the concordance index for the ordinal model shown above does not represent the same statistic as for the logistic model, shown below. It is instead the ‘global’ C-statistic for the ordinal model. A comparable C-statistic to the logistic model, , i.e. the proportion of all pairs of participants, of which the participant who had ≥10-14% BMPC, had a higher predicted probability of ≥10-14% BMPC than the participant who had <10-14% BMPC, was calculated in mgus_bm_prediction_20231125 and is contained in the object ‘pred_C2’. It is 0.8529 and is essentially indistinguishable from the C-statistic from the binary logistic model below. Having said this, the R-squared from the ordinal model is directly comparable.
pred_C2## [1] 0.8546772 0.8185687 0.8842890Next we fit the same model as a dichotomous logistic regression model
fit_pred_2 <- lrm(
formula = outcome_cont >= 2 ~ rcs(igg, 4) + rcs(iga, 4) + rcs(igm, 4) + rcs(m_protein, 4) +
rcs(flc_ratio, 4) + isotype, data = df_pred, x = TRUE, y = TRUE
)As seen below the likelihood ratio chi-square test is 292.62 on 18 degrees of freedom, with an R-squared of 0.427 and a concordance index of 0.870. The R-squared is considerably higher than the ordinal model, but as stated above the C-statistic is indistinguishable.
fit_pred_2## Logistic Regression Model
##
## lrm(formula = outcome_cont >= 2 ~ rcs(igg, 4) + rcs(iga, 4) +
## rcs(igm, 4) + rcs(m_protein, 4) + rcs(flc_ratio, 4) + isotype,
## data = df_pred, x = TRUE, y = TRUE)
##
## Model Likelihood Discrimination Rank Discrim.
## Ratio Test Indexes Indexes
## Obs 1043 LR chi2 292.62 R2 0.427 C 0.870
## FALSE 885 d.f. 18 R2(18,1043)0.231 Dxy 0.739
## TRUE 158 Pr(> chi2) <0.0001 R2(18,402.2)0.495 gamma 0.739
## max |deriv| 2e-05 Brier 0.082 tau-a 0.190
##
## Coef S.E. Wald Z Pr(>|Z|)
## Intercept 2.4808 1.3398 1.85 0.0641
## igg -0.3375 0.1494 -2.26 0.0239
## igg' 1.1877 0.6772 1.75 0.0794
## igg'' -2.9731 1.9908 -1.49 0.1353
## iga -0.3882 0.4146 -0.94 0.3491
## iga' 1.8762 2.7364 0.69 0.4929
## iga'' -3.2413 6.2620 -0.52 0.6047
## igm 0.9693 1.2309 0.79 0.4310
## igm' -16.2027 14.9225 -1.09 0.2776
## igm'' 32.7951 30.6593 1.07 0.2848
## m_protein -0.7035 0.3831 -1.84 0.0663
## m_protein' 59.3497 22.3347 2.66 0.0079
## m_protein'' -75.3620 28.3033 -2.66 0.0078
## flc_ratio -2.6218 0.5007 -5.24 <0.0001
## flc_ratio' 47.7156 9.4546 5.05 <0.0001
## flc_ratio'' -112.3930 22.4743 -5.00 <0.0001
## isotype=IgA -0.2954 0.4720 -0.63 0.5314
## isotype=IgG -0.8259 0.3610 -2.29 0.0222
## isotype=Light-chain -0.2312 0.5461 -0.42 0.6720
## The argument for preferring an ordinal model over the logistic model, is the added statistical precision obtained by ‘borrowing’ information from the more frequent outcome categories. This can be quantified by comparing how much tighter the 95% confidence intervals are for predictions made by the ordinal model are compared to the logistic model, (the variance of the log-odds). The ratio of the variance of the log-odds between the two represents the relative efficiency of the two models for a given sample size, and can be re-stated as the approximate increase in sample size that would be needed so that the logistic model is as precise as the ordinal model for the current sample size.
First we must show that the proportional odds assumption is met. Firstly we examine the unadjusted means of the predictor variables over different levels of the outcome. Even in this unadjusted comparison, proportional odds are shown to be highly likely.
plot.xmean.ordinaly(
outcome ~ igg + iga + igm + m_protein + flc_ratio,
data = df_pred
)
Now we do a more extensive check. We fit the exact same model as the full model. We then produce three dichotomous models and compare the results of these. This was shown in supplementary figures 1 and 2. For brevity of this document, we do not include all the code instead choosing to show the figures alone.
supplementary_figure_1
Predictions of the proportional odds model (the prediction model) are compared with binary submodels of each of the specific outcomes: 5-9% bone marrow plasma cells (BMPC) or higher, 10-14% BMPC or higher, and 15% BMPC or higher. The top row illustrates the correlation between predictions made by the prediction model (X-axis) and the binary submodels (Y-axis). The mean correlation is depicted as a blue line. A theoretical perfect correlation is represented by an oblique dotted line. The middle row depicts the risk difference between the prediction model and the binary submodels as a function of predicted probabilities of the prediction model. The mean risk difference is shown as a blue line. The minimum and maximum mean differences are shown as red dotted lines. The bottom row contains histograms illustrating the distributions of the differences in predicted risks between the prediction model and the binary submodels. This figure provides strong evidence that the proportional odds assumptions hold.
supplementary_figure_2
2 Precision of the ordinal compared to logistic model
First we calculate the predicted log-odds of ≥10-14% bone marrow plasma cells for each person in the sample using the ordinal model.
l_ordinal_logodds <- list()
for (i in 1:nrow(df_pred)) {
ordinal_logodds <- Predict(
fit_pred, kint = 2,
isotype = as.character(df_pred[i, "isotype"]),
m_protein = as.numeric(df_pred[i, "m_protein"]),
flc_ratio = as.numeric(df_pred[i, "flc_ratio"]),
igg = as.numeric(df_pred[i, "igg"]),
iga = as.numeric(df_pred[i, "iga"]),
igm = as.numeric(df_pred[i, "igm"])
)
ordinal_logodds$num <- i
l_ordinal_logodds[[i]] <- ordinal_logodds
}
ordinal_logodds <- do.call(rbind, l_ordinal_logodds)Then we calculate the width of the confidence interval of the predicted log-odds for each person in the sample.
ordinal_logodds$width <- ordinal_logodds$upper - ordinal_logodds$lowerNow we do all of the above steps for the logistic model, calculating the predicted log-odds for each participant.
l_logistic_logodds <- list()
for (i in 1:nrow(df_pred)) {
logistic_logodds <- Predict(
fit_pred_2,
isotype = as.character(df_pred[i, "isotype"]),
m_protein = as.numeric(df_pred[i, "m_protein"]),
flc_ratio = as.numeric(df_pred[i, "flc_ratio"]),
igg = as.numeric(df_pred[i, "igg"]),
iga = as.numeric(df_pred[i, "iga"]),
igm = as.numeric(df_pred[i, "igm"])
)
logistic_logodds$num <- i
l_logistic_logodds[[i]] <- logistic_logodds
}
logistic_logodds <- do.call(rbind, l_logistic_logodds)And calculating the width of the confidence interval for the predicted log-odds
logistic_logodds$width <- logistic_logodds$upper - logistic_logodds$lowerNow we can calculate the ratio of the widths of the confidence intervals for each persons predicted log-odds of ≥10-14% bone marrow plasma cells between the ordinal and logistic model. The widths form a distribution which is shown below. The widths produced by the ordinal model are always tighter, by a median of 40% (range: 18-80%)
quantile(ordinal_logodds$width/logistic_logodds$width)## 0% 25% 50% 75% 100%
## 0.1977347 0.5574244 0.6031247 0.6727922 0.8195177These substantive gains in statistical efficiency can be restated as how much larger a sample size would need to be using the logistic model, to equal the precision of the ordinal model using our sample size of 1043 persons. We show that, to match the current efficiency of the ordinal model, we would require at a minimum of 1231 persons, most likely 1456, and up to 1879.
quantile((1+ (1-ordinal_logodds$width/logistic_logodds$width))*1043)## 0% 25% 50% 75% 100%
## 1231.243 1384.278 1456.941 1504.606 1879.763Again restated, choosing the logistic over the ordinal model would be equivalent to discarding between 188-836 persons worth of information.
We have previously shown that there is no discordance between predictions made by the two models and that the proportional odds assumption holds, in supplemental figures 1 and 2. # Summary We have shown that the ordinal model is considerably more precise than the logistic model, and we have shown that all the assumptions of the ordinal model hold. Although not as commonly used, the ordinal model is fairly intuitive to explain to non-statistician readers, especially how it is used in the current manuscript. After careful consideration the author’s feel that the gains in efficiency more than outweigh the slight added complexity and argue that the ordinal model should be used.