# How to check if there is a linear relationship between a categorical feature and a continuous feature?

The target feature in my data is numeric, continuous. I have several predictive features, some of them are categorical. Now in order to train a linear regression model, I've used dummy variables and created new binary features in the data. So each categorical feature created new binary features of 1 and 0.

These are the columns in my data, as you can see many of them start with occupation since it's a categorical data that had several values:

Index(['age', 'educational-num', 'hours-per-week', 'occupation_Armed-Forces',
'occupation_Craft-repair', 'occupation_Exec-managerial',
'occupation_Farming-fishing', 'occupation_Handlers-cleaners',
'occupation_Machine-op-inspct', 'occupation_Other-service',
'occupation_Priv-house-serv', 'occupation_Prof-specialty',
'occupation_Protective-serv', 'occupation_Sales',
'occupation_Tech-support', 'occupation_Transport-moving',
'gender_Male'],
dtype='object')


Now, linear regression have several assumptions that should stand. The first and the most important assumption in linear regression is features should have linear and additive relationship with the target variable.

There are plenty of explenations for this and how to check this in python, like this for example.

But this doesn't explain how can this be done when some of my predictive features are categorical. Checking this for continuous features is pretty easy, using scatter plots or with pearson corelation. But for a categorical this is not the case and I'm wondering how can I prove that this linearity assumptions stands.

thanks!

• Why would you do anything different? In other words, what do you not find satisfying about plotting the 0/1 variable on the horizontal axis and the continuum of values on the vertical axis?
– Dave
Commented Dec 24, 2022 at 0:23
• Basically you use the same p-value of this categorical predictor variable to check whether it's significant based on your accepted significance level. Commented Dec 24, 2022 at 0:45
• @mohottnad Care to post that as an answer? There are problems with screening based on p-values, sure, but comments are not the place for answers.
– Dave
Commented Dec 24, 2022 at 0:53
• @Dave How would a plot like that display linear relationship? imagine 0 and 1 on the x-axis with many dots stacked in a vertical manner, that won't do me any good.
– CORy
Commented Dec 24, 2022 at 9:42
• @mohottnad I don't understand your solution. What p-value are you talking about?
– CORy
Commented Dec 24, 2022 at 9:43

This is what various distribution plots are for. Standard plots include histograms and boxplots. I like kernel density estimation (KDE) plots, too. Briefly, KDE can be thought of as a continuous histogram.

Let's simulate some data and look at the distributions.

library(ggplot2)
set.seed(2023)
N <- 100
x <- sample(c("Cat", "Dog", "Elephant"), N, replace = T)
y <- rnorm(N)
d <- data.frame(
Species = x,
Weight = y
)
ggplot(d, aes(y = Weight, fill = Species)) +
geom_boxplot(alpha = 1.0)


The plots show the groups not to be so different, consistent with the simulation setup that creates the weight without regard for the species. If you increase the sample size to 1000, the differences becomes even less.

If you see something like this, then you have some evidence that the distributions are not so different for the thre species.

Even better, however, is to visualize the entire distribution, since boxplots ignore a lot of information (such as multi-modality). The common way to do this is with a histogam, though I would prefer to look at the empirical CDFs or KDE plots.

library(ggplot2)
set.seed(2023)
N <- 1000
x <- sample(c("Cat", "Dog", "Elephant"), N, replace = T)
y <- rnorm(N)
d <- data.frame(
Species = x,
Weight = y
)
ggplot(d, aes(x = Weight, fill = Species)) +
geom_density(alpha = 0.3)
d0 <- data.frame(
Weight = y[x == "Cat"],
Quantile = ecdf(y[x == "Cat"])(y[x == "Cat"]),
Species = "Cat"
)
d1 <- data.frame(
Weight = y[x == "Dog"],
Quantile = ecdf(y[x == "Dog"])(y[x == "Dog"]),
Species = "Dog"
)
d2 <- data.frame(
Weight = y[x == "Elephant"],
Quantile = ecdf(y[x == "Elephant"])(y[x == "Elephant"]),
Species = "Elephant"
)
d <- rbind(d0, d1, d2)
ggplot(d, aes(x = Weight, y = Quantile, col = Species)) +
geom_line(size = 1.5)


The KDE for the three species are all about on top of each other.

The empirical CDFs for the three species are all about on top of each other.

Next, let's look at linear differences. To me, this would mean that the distributions are the same except for a linear shift, such as $$N(0, 1)$$, $$N(1, 1)$$, and $$N(4, 1)$$.

library(ggplot2)
set.seed(2023)
N <- 1000
x <- c(
rep("Cat", N),
rep("Dog", N),
rep("Elephant", N)
)
y <- c(
rnorm(N, 0, 1),
rnorm(N, 1, 1),
rnorm(N, 4, 1)
)
d <- data.frame(
Species = x,
Weight = y
)
ggplot(d, aes(y = Weight, fill = Species)) +
geom_boxplot()

ggplot(d, aes(x = Weight, fill = Species)) +
geom_density(alpha = 0.3)
d0 <- data.frame(
Weight = y[x == "Cat"],
Quantile = ecdf(y[x == "Cat"])(y[x == "Cat"]),
Species = "Cat"
)
d1 <- data.frame(
Weight = y[x == "Dog"],
Quantile = ecdf(y[x == "Dog"])(y[x == "Dog"]),
Species = "Dog"
)
d2 <- data.frame(
Weight = y[x == "Elephant"],
Quantile = ecdf(y[x == "Elephant"])(y[x == "Elephant"]),
Species = "Elephant"
)
d <- rbind(d0, d1, d2)
ggplot(d, aes(x = Weight, y = Quantile, col = Species)) +
geom_line(size = 1.5)


All three of these visualizations suggest that the difference in the distribution for each of the three animals is just in shifting up or down.

Finally, let's consider nonlinear differences. This one is a bit iffy, because the conditional mean is either shifted up, shifted down, or not changed, so the relationship between a categorical feature and the conditional mean is necessarily linear. However, the entire distribution does not have to shift the same way, and this can suggest alternative modeling, such as generalized linear models.

library(ggplot2)
set.seed(2023)
N <- 1000
x <- c(
rep("Cat", N),
rep("Dog", N),
rep("Elephant", N)
)
y <- c(
rexp(N, 5) - 1/5,
rt(N, 5),
rchisq(N, 1) - 1
)
d <- data.frame(
Species = x,
Weight = y
)
ggplot(d, aes(y = Weight, fill = Species)) +
geom_boxplot()

ggplot(d, aes(x = Weight, fill = Species)) +
geom_density(alpha = 0.3)
d0 <- data.frame(
Weight = y[x == "Cat"],
Quantile = ecdf(y[x == "Cat"])(y[x == "Cat"]),
Species = "Cat"
)
d1 <- data.frame(
Weight = y[x == "Dog"],
Quantile = ecdf(y[x == "Dog"])(y[x == "Dog"]),
Species = "Dog"
)
d2 <- data.frame(
Weight = y[x == "Elephant"],
Quantile = ecdf(y[x == "Elephant"])(y[x == "Elephant"]),
Species = "Elephant"
)
d <- rbind(d0, d1, d2)
ggplot(d, aes(x = Weight, y = Quantile, col = Species)) +
geom_line(size = 1.5)


These three plots show that there is much more to the difference between the three distributions than just sliding up and down the real line. They have different skewnesses. They have different variances. (They actually have equal means.)

Depending on what you are modeling, you might be quite interested in these differences.

If you want to quantify the strength of a relationship between a categorical variable and a continuous outcome, an analogous statistic to Pearson correlation between a continuous feature and a categorical outcome would be to use regression and take the square root of the $$R^2$$. If you do this for a continuous feature, you get the (magnitude of) the Pearson correlation, so this seems like a reasonable generalization. I will demonstrate below.

set.seed(2023)
x <- c(
rep("Cat", N),
rep("Dog", N),
rep("Elephant", N)
)
y <- c(
rnorm(N, 1, 1),
rnorm(N, 2, 1),
rnorm(N, 3, 1)
)
L <- lm(y ~ x)
sqrt(summary(L)$r.squared)  I get a fairly strong "correlation" of $$0.6278021$$. Since the feature is categorical, there is not really a notion of direction, so the sign does not matter, though I would go with the positive square root out of convenience. If you need to convey the strength of this relationship, you can use the boxplots, KDEs, and CDFs above (maybe histograms, too), but another option, which I confess I have not used (but I do like the idea), is to simulate some contnuous data with that calculated "correlation". For instance, the "correlation" between the outcome $$y$$ and categorical feature is $$0.6278021$$. Simulate some data with that correlation and graph that bivariate data, such as below. library(MASS) set.seed(2023) X <- MASS::mvrnorm(N, c(0, 0), matrix(c( 1, sqrt(summary(L)$r.squared),
sqrt(summary(L)\$r.squared), 1
), 2, 2))
d <- data.frame(
x = X[, 1],
y = X[, 2]
)
ggplot(d, aes(x = x, y = y)) +
geom_point()


To me, this demonstrates the fairly strong relationship between the categorical feature and continuous outcome, and it does it in way that should be familiar and comfortable to stakeholders.