1
$\begingroup$

When running a regression with a categorical independent variable, we get results for each level of the variable except for the base, which we can choose.

Now I've always had a hard time on how to interpret these results.

Say we have a study of aneurysm locations. They can be located in, say, 10 different areas.

We want to see if smokers develop aneurysms in other areas than non-smokers. We have our dependent variable (smoker, no/yes) and our independent variable of location with 10 levels.

If we run the regression we might get a significant hit on 3 locations. But this is compared to the base location which let's say is level "5".

So yes, smokers are significantly more likely to get aneurysms in location 1, 2 and 3 compared to location 5. But this doesn't answer my research question of "which areas are smokers more likely to develop aneurysms in?".

What I would like to do is to maybe make a "dummy level" to my categorical variable in which half the patients have that location and half do not and then use that as base to see if ANY of the 10 true levels have a significantly higher risk for smokers. I don't want to compare the levels to each other and I want every level included and not for one of them to be used as a base.

Is there some way to do this or am I using the wrong model to answer my research question?

I assume splitting the categorical variable into 10 dummy variables is probably not so smart.

|cite|improve this question
$\endgroup$
  • $\begingroup$ Do you have observations over time for individuals or one per individual? Also can same individual get aneurysms in several places? Usually for categorial variables the categories are mutually exclusive. (also the title of the question could perhaps be changed, doesn't sound like you are really interested in knowing how to interpret coefficient for categorial variable in regression) $\endgroup$ – Jesper Hybel Oct 29 at 12:43
  • $\begingroup$ One per individual. They can, but it's highly unlikely for my purposes, in the study we can assume they are mutually exclusive. $\endgroup$ – Paze Oct 29 at 13:11
0
$\begingroup$

You're not using the wrong model; you're just not looking at an output that tells you what you want know....

My suspicion is that you are using R, just because this question comes up a lot with newer R users, because so often examples in R rely on the output from the summary function, which provides model coefficients, and usually t tests or other relevant tests for the coefficients.

Instead, what it sounds like what you want, is 1) an anova-like table of effects, and b) (pairwise) comparisons among group means or estimated marginal means. b) will tell you if group 1 is different than group 2, and so on.

In R, a) is shown with e.g. library(car); Anova(model), where 'model' is your model object. b) is shown using the emmeans package. Something like library(emmeans); marginal = emmeans(model, ~ Group); pairs(marginal, adjust="tukey").

|cite|improve this answer
$\endgroup$
  • $\begingroup$ I'm actually using STATA but I get your point. Can you tell me how this is different from dichotomizing each level of the categorical variable and then running a logistic regression on each level by itself? $\endgroup$ – Paze Oct 29 at 21:23
  • $\begingroup$ I'm not sure about if the equivalent is available in STATA, or for which models. There is pwcompare and margins, but I really don't know what these commands do. $\endgroup$ – Sal Mangiafico Oct 29 at 23:18
  • $\begingroup$ The idea is that it is better to fit a single model that includes all your data, and then a post-hoc approach that reflects that model. Using e.g. pairwise chi-square tests or pairwise logistic regressions considers only the data in those pairs and ignores the rest of the data. Pairwise approaches --- at least in simple examples, like a Kruskal-Wallis omnibus test --- can result in incommensurate results. $\endgroup$ – Sal Mangiafico Oct 29 at 23:23
0
$\begingroup$

A way I use to help me interpret is to create a table with the two variables I'm looking at - the dependent and the independent I'm focusing at that moment. In a regression, the other independent variables may change this relation a bit, but usually I can get a good idea regarding what's going on.

Following your example, I could have a table with conditional frequencies by row like: $$ \begin{array}{c|cc} Area \backslash Smoke & Yes\ (1) & No\ (0)\\ \hline 1 & 60\% & 40\% \\ 2 & 50\% & 50\% \\ 3 & 30\% & 70\% \\ \ldots & \ldots & \ldots \\ 10 & 55\% & 45% \\ \end{array} $$ Of course that we could also use conditional frequencies by column.

By looking at the table, we could say that the smokers were more likely to have an aneurysm in areas $1$ and $10$ and that the non-smokers were more likely to have an aneurysm in area $3$.

Edit:

If you have a logistic regression, you could look at the odds ratio ($OR$) or at the coefficients ($B$) to reach a similar conclusion from the reverse path. By the way, $OR = e^B$.

For instance, if the 2nd category is the reference category for the area variable (I admit that I chose the $50\% - 50\%$ category to make my life easier), and non-smoking is the reference category for the smoker variable, we could observe that $OR_1=1.4$ for area $1$, $OR_3=0.44$ for area $3$, and $OR_{10}=1.1$ for area $10$ (I've made up some numbers). This would mean that the smoking individuals, when compared to the non-smokers, were more likely to have aneurysms in areas $1$ and $10$ and less likely to have aneurysms in area $3$. By looking at the $p$-values for each coefficient, you could also say in which cases the $OR$ is significantly different from $1$.

We compare the value of the $OR$ with $1$ because $OR_j=1$ would mean that the odds of finding a smoker in category $j$ would be the same as the odds of finding a smoker in category $2$, i.e., being a smoker would have the same odds of having an aneurysm in area $j$ or in area $2$.

|cite|improve this answer
$\endgroup$
  • $\begingroup$ Right. But if you were to publish this work and your research question results, you would have to then run a logistic regression (or chi2) on each area by itself to be able to say "smokers were more likely to have aneurysms in areas 1 and 10 and non smokers were more likely to have aneurysms in area 3", as you would also be interested in p values and confidence intervals. Then you could add "furthermore, smokers were x times more likely to have aneurysms in area 1 than area 3" by running a regression with the entire variable (all levels). Would that be the correct way to do it? $\endgroup$ – Paze Oct 29 at 13:30
  • $\begingroup$ @Paze , rather than split the analysis into pairwise logistic regressions or chi-square tests, it's better to use one model covering all groups, and use an estimated marginal means procedure to compare the groups pairwise (or each to the control, or whatever is appropriate). Each of SAS, SPSS, and R has something called e.g. EMMEANS or emmeans, at least for some model types. $\endgroup$ – Sal Mangiafico Oct 29 at 13:42
  • $\begingroup$ @Paze: you may also use a multinomial logistic regression with area as the dependent variable and smoker as the independent variable. $\endgroup$ – Ertxiem - reinstate Monica Oct 29 at 13:44
  • $\begingroup$ I deleted my previous comment as it was confusing. I tried running a multinomial log regression swapping the dependent and independent variables but it does the same thing as running a logistic regression. It uses the first as a base. I get the same results at least. $\endgroup$ – Paze Oct 29 at 18:10
0
$\begingroup$

Given the development of aneurysm in one area mutually excludes developments of aneurysm in another it is possible to let the dependent variable be $Y_{ij}=1$ if individual $i$ has developed aneurysm in the area $j$ and $Y_{ij}=0$ otherwise. For each individual you then have a vector of measurements $Y_i =(Y_{i1},...,Y_{iJ})$ as a multivariate output. The vector follows a multinomial distribution the probabilities of which are simply estimated by sample relative frecquencies.

You can fit a multinomial logit model and do testing comparing the relative frecquencies of smoker with non-smokers.

$$H_o : \pi(smoker) = \pi(non-smoker)$$

where $\pi(smoker)$ is a vector of probabilities $$\pi_1,...,\pi_J$$ with $\pi_j$ being the probability that a smoker develops aneurysm in area $j$.

Here is some R code simulating such data and running a test

    library(data.table)
library(mnlogit)

true_difference <- 0
N <- 1000
J <- 10
smoker <- as.numeric(runif(N)<0.5)

# select area effects (parameters to generate probabilities)
AE <- 0.5*rnorm(J)
AE[1] <- AE[1] + true_difference
AE_smoker <- rep(0,J)
AE_smoker[3] <- AE_smoker[3]+true_difference
AE_smoker[4] <- AE_smoker[4]+true_difference
AE_smoker[5] <- AE_smoker[5]+true_difference


AE
AE_smoker


p <- exp(AE)
p <- p/sum(p)
p_smoker <- exp(AE + AE_smoker)
p_smoker <- p_smoker/sum(p_smoker)
p
p_smoker


area <- rep(NA,N)
for (i in 1:N)
    {
        i_probs <- smoker[i]*p_smoker + (1-smoker[i])*p
        area[i] <- sample(1:10,1,prob=i_probs)
    } 


dt <- data.table(id=1:N,area=area,smoker=smoker)
agg_dt <- dt[, .(count = .N), by = .(area,smoker)]
setkey(agg_dt,area,smoker)
agg_dt[,.(difference=diff(count)),by=area]


# Transform data to long format
choice <- rep(0,J*N)
ii <- 1
for (i in 1:N)
    {
        choice[ii:(ii+J-1)][dt$area[i]] <- 1
        ii <- ii + J    
    }

mydata <- data.frame(id=sort(rep(dt$id,J)),
choice=as.logical(choice),
smoker=rep(dt$smoker,each=J),
alt=rep(1:J,N))

# Specify model
# https://arxiv.org/pdf/1404.3177.pdf
fm <- formula(choice ~  1 | smoker | 1)
fit <- mnlogit(fm,data=mydata,choiceVar="alt")


fm.c <- formula(choice ~  1 | 1 | 1)
fit.c <- mnlogit(fm.c,data=mydata,choiceVar="alt")
lrtest(fit, fit.c)
waldtest(fit, fit.c)
# Scoretest do not seem to be working ...
scoretest(fit, fit.c)
|cite|improve this answer
$\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.