The binary model will give probabilities that depend on the frequency distribution of individual outcome categories within the coarsened category combinations. And it depends on the single threshold being ‘magical’, i.e., one every expert would
The partial prop. Odds model would be better than either, but may require knowing ahead of time which predictors are likely to not act in prop. Odds.
This study tested the effect of Fluoxetine on functional recovery after stroke, which was measured at 6 months using the modified Rankin Scale (mRS; 0-5). The authors used ordinal (PO) logistic regression for the analysis. Their model was the following:
\ mRS_{6m} \text ~ \ Group \ + \ Time \ + \ Predicted \ Outcome \ + \ Motor \ Deficit \ + \ Aphasia
where,
Time = time since stroke (2-8 days vs 9-15 days)
Predicted Outcome = probability of mRS of 0-2 at 6 months (≤0.15 vs >0.15), obtained from a prediction model
Motor Deficit = presence of motor deficit at randomization (yes vs no)
Aphasia = presence of aphasia at randomization (yes vs no)
I see the dichotomization of Time and Predicted Outcome as a problem. The very use of a prediction/ probability as a covariate, in place of the baseline mRS for instance, is problematic in my view. It seems that they didn’t have the baseline mRS available to use in the model. Also, they didn’t assess whether the PO assumption was met/ appropriate.
Any thoughts/ opinions?
To do a study of mRS without collecting mRS at baseline is a key design flaw. And the PO assumption is the least of their worries. The flaws you pointed out, especially dichotomization of a continuous risk measure, are huge. A spline in logit risk should have been used. Dichtomizing time is pretty bad also.
I saw that you mentioned on a different post that when we have an ordinal covariate (mRS) it is better to model it as nominal instead of as ordinal. Does that apply here too?
I didn’t mean to say that. I was saying that if there’s only say 4 levels it’s not worth the trouble of specially modeling it as ordinal.
I’m fitting this model but the bootstrap is not working:
m <- lrm(y ~ age + sex + lead73, x = T, y = T, data = lead)
b <- bootcov(m, B = 500, coef.reps = T, group = lead$y)
Warning: at least one resample excluded highest Y values, invalidating bootstrap covariance matrix estimate
The model has 39 intercepts.
@f2harrell, any suggestion?
Also, would the code below be correct for getting the contrast on the mean scale?
M <- Mean(m)
c <- contrast(b,
list(lead73 = quantile(lead$lead73, 0.75)),
list(lead73 = quantile(lead$lead73, 0.25)))
print(c, fun = M)
Thank you
fun= goes with contrast not print. And contrast will use the bootcov output to get approximate CLs on differences in means. The problem with some bootstrap samples not including some Y levels is a very tricky one. Sometimes it’s OK to give bootcov a group=Y type of argument to balance the bootstrap but performance of this has not been studied. At other times a bit of rounding of Y will suffice. This is another reason to go Bayesian.