# Verify relationship between a categorical and numerical variable

I have got some data about movie directors and number of cuts per minutes in their movies. I whould like to verify if there is any relation between this two variables. I have tried looking at the mean values of the data grouped by director, and it seems to me there could be something interesting. However, the number of instances is quite low for some directors. Any suggestions on how to verify the relationship? Note that the variable of the cuts per minute is not distributed normally.

# EDIT

Here is the data:

data = structure(list(Director = c("Director A", "Director A", "Director A",
"Director A", "Director B", "Director A", "Director B", "Director C",
"Director C", "Director C", "Director C", "Director A", "Director A",
"Director A", "Director A", "Director B", "Director C", "Director C",
"Director D", "Director D", "Director B", "Director B", "Director C",
"Director C", "Director C", "Director C", "Director C", "Director B",
"Director B", "Director B", "Director C", "Director C", "Director C",
"Director B", "Director B", "Director B", "Director B", "Director B",
"Director B", "Director B", "Director E", "Director E", "Director F",
"Director F", "Director G", "Director G", "Director C", "Director C",
"Director E", "Director E", "Director E", "Director E", "Director E",
"Director C", "Director C", "Director C", "Director E", "Director C"
), CutsPerMinute = c(10.5140679608446, 15.9323927435056, 12.5931626536961,
11.9933878240714, 10.7534044174676, 10.5722374320877, 10.150267917226,
10.14185419837, 10.2274315628792, 8.597567605187, 11.2610061781615,
11.5530433481127, 8.28255075816311, 8.71178186135948, 8.1589465201578,
7.68692552413314, 6.84104588989501, 9.02599848594941, 8.46017478960381,
5.30101556092088, 7.38173747636346, 7.88623562496654, 6.99593391518614,
6.79296728065352, 6.99907913122217, 7.54379856090484, 8.08708813022303,
7.83580886474942, 9.17641479573382, 6.62232278314871, 9.33009967506452,
5.61068565096716, 8.38412344816077, 8.42626951824825, 7.66648997484084,
7.10413146060059, 4.17657904626811, 7.53333742519325, 7.67502218583697,
7.98868602365391, 6.09848933756457, 6.73728682033344, 8.2223199922461,
7.26367081188495, 5.42732930399057, 7.10515586961886, 7.47580967835688,
6.93978155365154, 9.40494846731573, 7.95957905150531, 8.61859633144952,
7.56455661421347, 7.64150388044662, 9.6414259939393, 8.10302490573573,
8.164055371306, 7.77019990100153, 10.2018773770501)), class = "data.frame", row.names = c(NA,
-58L))

• You’ve tagged this with anova (with which I agree), so it sounds like you have some idea of how to proceed. What do you have in mind for ANOVA?
– Dave
Jun 26, 2022 at 20:55
• @Dave Hello! Thanks for your comment. I do not have a good knowledge of anova, but reading online it seemed like a solution. I tagged it in order to know if it could be a nice approach for my particular problem. I tried applying it and get a very high feature importance for the director variable, but I think this has violated the assumption of normality of the data. Jun 26, 2022 at 21:11
• You do know the movies and the directors? I'd add information about the genre in this analysis (unless you already know that all the movies have the same genre). Jun 27, 2022 at 0:03
• Sadly I don't know that, sorry! Jun 27, 2022 at 0:05

As Dave mentioned above, ANOVA (Analysis of Variance) might be a starting point for what you're looking for. This technique will help you identify differences between the variances by looking at the variation both between and within groups.

To help you get started, I generated some data with three directors (A, B and C) with random cut counts, with some variation among the directors:

> set.seed(123)
> n = 100
> cuts = rpois(n, 2)
> dirs = factor(sample(paste("Director", LETTERS[1:3]), n, replace = TRUE))

## Introduce variation by adding random counts to directors A and B with different expected values.
> flag_dir_a = (dirs == "Director A")
> cuts[flag_dir_a] = cuts[flag_dir_a] + rpois(sum(flag_dir_a), 5)
> flag_dir_b = (dirs == "Director B")
> cuts[flag_dir_b] = cuts[flag_dir_b] + rpois(sum(flag_dir_b), 1)

## Create a dataframe using the simulated data above.
> dat = data.frame(director = dirs, cuts = cuts)

director cuts
1 Director A    6
2 Director B    3
3 Director A    8
4 Director B    4
5 Director A    7
6 Director C    0


If we look at the mean cuts by director, we will observe some variation (by construction):

> with(dat, tapply(cuts, director, mean))

Director A Director B Director C
7.538462   2.846154   2.171429


The standard deviations also vary:

> with(dat, tapply(cuts, director, sd))

Director A Director B Director C
2.453255   1.581779   1.635735


We can check this visually as well:

with(dat, boxplot(cuts ~ director))


You could then run an ANOVA analysis. The output below tells you that a significant effect has been detected (p-value associated with the F-statistic is less than 0.05) - i.e. there is a significant evidence to reject the null hypothesis that all directors' average cut counts are equal (of course, we already knew this).

> aov(cuts ~ director, dat)

Df Sum Sq Mean Sq F value Pr(>F)
director     2  491.6  245.80   70.85 <2e-16 ***
Residuals   97  336.5    3.47
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1


Finally, since ANOVA is just a type of linear model, you could fit a linear regression model to assess the model fit.

> mod_lm = lm(cuts ~ director, data = dat)
> summary(mod_lm)

Call:
lm(formula = cuts ~ director, data = dat)

Residuals:
Min      1Q  Median      3Q     Max
-3.5385 -1.1714 -0.1714  1.1538  5.4615

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept)          7.5385     0.3653   20.64   <2e-16 ***
directorDirector B  -4.6923     0.4716   -9.95   <2e-16 ***
directorDirector C  -5.3670     0.4822  -11.13   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.863 on 97 degrees of freedom
Multiple R-squared:  0.5936,    Adjusted R-squared:  0.5853
F-statistic: 70.85 on 2 and 97 DF,  p-value: < 2.2e-16



Hope this helps.

Edited:

Using your data, we can see:

and the ANOVA output is:

            Df Sum Sq Mean Sq F value   Pr(>F)
director     6  79.87  13.312   4.977 0.000444 ***
Residuals   51 136.40   2.674
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1



Now, whilst the cut/minute isn't distributed normally, I think ANOVA is fairly robust against non-normality as long as there are no gross outliers in the data.

If you're concerned about non-normality then you could consider an alternative, non-parametric test such as Kruskal-Wallis. The null hypothesis (against which we have significant evidence below) here is that cut/minute are identically distributed among the directors.

> kruskal.test(cuts ~ director, data = dat)

data:  cuts by director
Kruskal-Wallis chi-squared = 16.995, df = 6, p-value = 0.0093

• Thank you a lot for your help! Just one question, could the distribution of the data violet some assumptions or be an issue in general? I'll add the real data to the post. Jun 26, 2022 at 21:30