Choosing the right modeling procedure

Question

I am working with a dataset (electronic medical record extract) and trying to find out whether there is a difference in a continuous integer outcome (RxTotal, mean=7.4, range 2-34, skew=1.6, kurtosis=2.3) for a categorical predictor (Race4) after controlling for another continuous integer predictor (TotalRxLen - a complex variable to describe total days of pill supply made available over time, with overlapping/early prescriptions truncated).

So my model is as follows: RxTotal = Race4 + TotalRxLen

The purpose of this is to know whether or not non-white patients had to have more Rx visits (RxTotal) in order to reach the same total pill supply over time (TotalRxLen), as a way of demonstrating medical provider racial bias. I know that this model may sounds a little backwards because RxTotal is the outcome and TotalRxLen is the predictor. I think it is okay because the analysis is seeking to identify adjusted mathematical relationships and not to explore or say anything about causal order. I am only interested in Total Pill Supply (TotalRxLen) as a covariate, not at all as an outcome.

I already know that doctors provided significantly shorter prescriptions to Black and Hispanic patients, factoring in refills, compared to white non-Hispanic patients, and yet Black and Hispanic patients did not have a shorter Total Pill Supply, nor did they have a smaller RxTotal. By adjusting for Total Pill Supply I am trying to "sop up" some of the variation. I am wondering what modeling procedure to use. I ran this Shapiro Wilk test:

proc univariate data=totaldays2 normal;
where rxtotal >1;
var RxTotal;
run;

Shapiro Wilk p <.001 so I reject the normality assumption. I believe that this means I cannot use this PROC GLM?:

proc glm data=totaldays2 ;
where rxtotal >1;
class  race4;
model RxTotal = TotalRxLen race4  / solution ;
run; 

I have been researching this online and am having a hard time identifying what procedure and model type to use. Can anyone provide any guidance on what procedure and model type to use, and/or what further assumption testing I need to perform? Or whether I in fact perhaps actually can use the PROC GLM procedure I included?

Thank you!

Answer 1

With ordered outcomes like your integer values, ordinal regression is a useful choice. It makes no assumptions about the distribution of outcome values (or of residuals around model estimates). It can be thought of as a generalization of some non-parametric methods. There are a couple of different approaches with different underlying assumptions, but in practice results can be fairly robust to violations of the assumptions. Frank Harrell has provided a list of resources on the topic.

Modeling of predictors is as for other types of multiple regression. You probably want to model your Total Pill Supply covariate flexibly, for example with a regression spline or another type of generalized additive model. That allows you to let the data tell you the form of the association with outcome and keeps you from falling into the trap of imposing an unrealistic form yourself.

An additional thought

The above is intended as a recommendation for how to handle an ordinal outcome. I can't say whether this particular model can accomplish what you want. See the answer from @Björn (+1) with respect to that.

Answer 2

These sound like very complex causal questions, where you seem to be ignoring many important aspects (e.g. what diseases do the people have, is the diagnosis probability the same for all of them even if they have the same disease, how old are they etc.) that may be alternative explanations of what you see. Running some regression model adjusting for one or two things most likely won't give you a useful answer to the causal question you really want to address. It is also not necessarily clear that excluding those with 0 visits would be the right thing to do. It's important to think this through (tools like DAGs can be a good way for trying to be rigorous about this).

One would not normally expect that treating a count variable like number of Rx visits as being a continuous variable with (approximately) normally distributed error terms is the right thing to do. Things like negative binomial regression or an ordinal approach might be more appropriate. However, note that in general one would not look at the normality of the variable under analysis, but of the residuals of the regression model. So, the test you conducted is not relevant to the question. It is of course in general debatable whether such a pre-test is a good idea, especially in a case like this where we know that the residuals just cannot be exactly normal (so that it's clear that the only reason why a test for the normality of residuals would not be significant would be if the sample size were too small).