1.
Conduct the same analysis comparing action movies versus romantic movies using the median rating instead of the mean rating. What was different and what was the same?
Section 9.6 Case study: are action or romance movies rated higher?
Let’s apply our knowledge of hypothesis testing to answer the question: "Are action or romance movies rated higher on IMDb?" IMDb is a database on the internet providing information on movie and television show casts, plot summaries, trivia, and ratings. We will investigate if, on average, action or romance movies get higher ratings on IMDb.
Subsection 9.6.1 IMDb ratings data
The
movies dataset in the ggplot2movies package contains information on 58,788 movies that have been rated by users of IMDb.com.
movies
# A tibble: 58,788 × 24 title year length budget rating votes r1 r2 r3 <chr> <int> <int> <int> <dbl> <int> <dbl> <dbl> <dbl> 1 $ 1971 121 NA 6.4 348 4.5 4.5 4.5 2 $1000 a Touchdown 1939 71 NA 6.0 20 0.0 14.5 4.5 3 $21 a Day Once a Month 1941 7 NA 8.2 5 0.0 0.0 0.0 4 $40,000 1996 70 NA 8.2 6 14.5 0.0 0.0 5 $50,000 Climax Show, The 1975 71 NA 3.4 17 24.5 4.5 0.0 6 $pent 2000 91 NA 4.3 45 4.5 4.5 4.5 7 $windle 2002 93 NA 5.3 200 4.5 0.0 4.5 8 '15' 2002 25 NA 6.7 24 4.5 4.5 4.5 9 '38 1987 97 NA 6.6 18 4.5 4.5 4.5 10 '49-'17 1917 61 NA 6.0 51 4.5 0.0 4.5 # ℹ 58,778 more rows # ℹ 15 more variables: r4 <dbl>, r5 <dbl>, r6 <dbl>, r7 <dbl>, r8 <dbl>, # r9 <dbl>, r10 <dbl>, mpaa <chr>, Action <int>, Animation <int>, # Comedy <int>, Drama <int>, Documentary <int>, Romance <int>, Short <int>
We will focus on a random sample of 68 movies that are classified as either "action" or "romance" movies but not both. We disregard movies that are classified as both so that we can assign all 68 movies into either category. Furthermore, since the original
movies dataset was a little messy, we provide a pre-wrangled version of our data in the movies_sample data frame included in the moderndive package. If you are curious, you can look at the necessary data-wrangling code to do this on GitHub.
movies_sample
# A tibble: 68 × 4 title year rating genre <chr> <int> <dbl> <chr> 1 Underworld 1985 3.1 Action 2 Love Affair 1932 6.3 Romance 3 Junglee 1961 6.8 Romance 4 Eversmile, New Jersey 1989 5 Romance 5 Search and Destroy 1979 4 Action 6 Secreto de Romelia, El 1988 4.9 Romance 7 Amants du Pont-Neuf, Les 1991 7.4 Romance 8 Illicit Dreams 1995 3.5 Action 9 Kabhi Kabhie 1976 7.7 Romance 10 Electric Horseman, The 1979 5.8 Romance # ℹ 58 more rows
The variables include the
title and year the movie was filmed. Furthermore, we have a numerical variable rating, which is the IMDb rating out of 10 stars, and a binary categorical variable genre indicating if the movie was an Action or Romance movie. We are interested in whether Action or Romance movies got a higher rating on average.
Let’s perform an exploratory data analysis of this data. Recall from Chapter 2 that a boxplot is a visualization we can use to show the relationship between a numerical and a categorical variable. Another option you saw in Chapter 2 would be to use a faceted histogram. However, in the interest of brevity, let’s only present the boxplot in Figure 9.6.1.
ggplot(data = movies_sample, aes(x = genre, y = rating)) +
geom_boxplot() +
labs(y = "IMDb rating")
Eyeballing Figure 9.6.1, romance movies have a higher median rating. Do we have reason to believe, however, that there is a significant difference between the mean
rating for action movies compared to romance movies? It is hard to say just based on this plot. The boxplot does show that the median sample rating is higher for romance movies. However, there is a large amount of overlap between the boxes. Recall that the median is not necessarily the same as the mean either, depending on whether the distribution is skewed.
Let’s calculate some summary statistics split by the binary categorical variable
genre: the number of movies, the mean rating, and the standard deviation split by genre. We will do this using dplyr data wrangling verbs. Notice in particular how we count the number of each type of movie using the n() summary function.
movies_sample |>
group_by(genre) |>
summarize(n = n(), mean_rating = mean(rating), std_dev = sd(rating))
# A tibble: 2 × 4 genre n mean_rating std_dev <chr> <int> <dbl> <dbl> 1 Action 32 5.28 1.36 2 Romance 36 6.32 1.61
Observe that we have 36 movies with an average rating of 6.322 stars and 32 movies with an average rating of 5.275 stars. The difference in these average ratings is thus 6.322 - 5.275 = 1.047. So there appears to be an edge of 1.047 stars in favor of romance movies. The question is, however, are these results indicative of a true difference for all romance and action movies? Or could we attribute this difference to chance sampling variation?
Subsection 9.6.2 Sampling scenario
Let’s now revisit this study in terms of terminology and notation related to sampling we studied in Subsection 7.2.1. The study population is all movies in the IMDb database that are either action or romance (but not both). The sample from this population is the 68 movies included in the
movies_sample dataset.
Since this sample was randomly taken from the population
movies, it is representative of all romance and action movies on IMDb. Thus, any analysis and results based on movies_sample can generalize to the entire population. What are the relevant population parameter and point estimates? We introduce the fourth sampling scenario:
| Scenario | Population parameter | Notation | Point estimate | Notation |
|---|---|---|---|---|
| 1 | Population proportion | \(p\) | Sample proportion | \(\widehat{p}\) |
| 2 | Population mean | \(\mu\) | Sample mean | \(\overline{x}\) or \(\widehat{\mu}\) |
| 3 | Difference in population proportions | \(p_1 - p_2\) | Difference in sample proportions | \(\widehat{p}_1 - \widehat{p}_2\) |
| 3 | Difference in population proportions | \(\mu_1 - \mu_2\) | Difference in sample proportions | \(\widehat{x}_1 - \widehat{x}_2\) or \(\widehat{\mu}_1 - \widehat{\mu}_2\) |
So, whereas the sampling bowl exercise in Section 7.1 concerned proportions, the almonds exercise in Subsection 8.2.1 concerned means, the case study on whether yawning is contagious in Section 8.4 and the music genre activity in Section 9.2 concerned differences in proportions, we are now concerned with differences in means.
In other words, the population parameter of interest is the difference in population mean ratings \(\mu_a - \mu_r\text{,}\) where \(\mu_a\) is the mean rating of all action movies on IMDb and similarly \(\mu_r\) is the mean rating of all romance movies. Additionally the point estimate/sample statistic of interest is the difference in sample means \(\overline{x}_a - \overline{x}_r\text{,}\) where \(\overline{x}_a\) is the mean rating of the \(n_a = 32\) movies in our sample and \(\overline{x}_r\) is the mean rating of the \(n_r = 36\) in our sample. Based on our earlier exploratory data analysis, our estimate \(\overline{x}_a - \overline{x}_r\) is \(5.275 - 6.322 = -1.047\text{.}\)
So there appears to be a slight difference of \(-1.047\) in favor of romance movies. The question is, however, could this difference of \(-1.047\) be merely due to chance and sampling variation? Or are these results indicative of a true difference in mean ratings for all romance and action movies on IMDb? To answer this question, we will use hypothesis testing.
Subsection 9.6.3 Conducting the hypothesis test
We will be testing:
\begin{align*}
H_0 \amp: \mu_a - \mu_r = 0\\
\text{vs } H_A \amp: \mu_a - \mu_r \neq 0
\end{align*}
In other words, the null hypothesis \(H_0\) suggests that both romance and action movies have the same mean rating. This is the "hypothesized universe" we will assume is true. On the other hand, the alternative hypothesis \(H_A\) suggests that there is a difference. Unlike the one-sided alternative we used in the popularity exercise \(H_A: p_m - p_d \gt 0\text{,}\) we are now considering a two-sided alternative of \(H_A: \mu_a - \mu_r \neq 0\text{.}\)
Furthermore, we will pre-specify a low significance level of \(\alpha = 0.001\text{.}\) By setting this value low, all things being equal, there is a lower chance that the \(p\)-value will be less than \(\alpha\text{.}\) Thus, there is a lower chance that we will reject the null hypothesis \(H_0\) in favor of the alternative hypothesis \(H_A\text{.}\) In other words, we will reject the hypothesis that there is no difference in mean ratings for all action and romance movies, only if we have quite strong evidence. This is known as a "conservative" hypothesis testing procedure.
Subsubsection 9.6.3.1 1. specify variables
Let’s now perform all the steps of the
infer workflow. We first specify() the variables of interest in the movies_sample data frame using the formula rating ~ genre. This tells infer that the numerical variable rating is the outcome variable, while the binary variable genre is the explanatory variable. Note that unlike previously when we were interested in proportions, since we are now interested in the mean of a numerical variable, we do not need to set the success argument.
movies_sample |>
specify(formula = rating ~ genre)
Response: rating (numeric)
Explanatory: genre (factor)
# A tibble: 68 × 2
rating genre
<dbl> <fct>
1 3.1 Action
2 6.3 Romance
3 6.8 Romance
4 5 Romance
5 4 Action
6 4.9 Romance
7 7.4 Romance
8 3.5 Action
9 7.7 Romance
10 5.8 Romance
# ℹ 58 more rows
Observe at this point that the data in
movies_sample has not changed. The only change so far is the newly defined Response: rating (numeric) and Explanatory: genre (factor) meta-data.
Subsubsection 9.6.3.2 2. hypothesize the null
We set the null hypothesis \(H_0: \mu_a - \mu_r = 0\) by using the
hypothesize() function. Since we have two samples, action and romance movies, we set null to be "independence" as we described in Section 9.4.
movies_sample |>
specify(formula = rating ~ genre) |>
hypothesize(null = "independence")
Response: rating (numeric)
Explanatory: genre (factor)
Null Hypothesis: independence
# A tibble: 68 × 2
rating genre
<dbl> <fct>
1 3.1 Action
2 6.3 Romance
3 6.8 Romance
4 5 Romance
5 4 Action
6 4.9 Romance
7 7.4 Romance
8 3.5 Action
9 7.7 Romance
10 5.8 Romance
# ℹ 58 more rows
Subsubsection 9.6.3.3 3. generate replicates
After we have set the null hypothesis, we generate "shuffled" replicates assuming the null hypothesis is true by repeating the shuffling/permutation exercise you performed in Section 9.2. We will repeat this resampling without replacement of
type = "permute" a total of 1000 times. Feel free to run the code to check out what the generate() step produces.
movies_sample |>
specify(formula = rating ~ genre) |>
hypothesize(null = "independence") |>
generate(reps = 1000, type = "permute") |>
View()
Subsubsection 9.6.3.4 4. calculate summary statistics
Now that we have 1000 replicated "shuffles" assuming the null hypothesis \(H_0\) that both
Action and Romance movies on average have the same ratings on IMDb, let’s calculate() the appropriate summary statistic for these 1000 replicated shuffles. Since the unknown population parameter of interest is the difference in population means \(\mu_{a} - \mu_{r}\text{,}\) the test statistic of interest here is the difference in sample means \(\overline{x}_{a} -
\overline{x}_{r}\text{.}\)
For each of our 1000 shuffles, we can calculate this test statistic by setting
stat = "diff in means". Furthermore, since we are interested in \(\overline{x}_{a} - \overline{x}_{r}\text{,}\) we set order = c("Action", "Romance"). Let’s save the results in a data frame called null_distribution_movies:
null_distribution_movies <- movies_sample |>
specify(formula = rating ~ genre) |>
hypothesize(null = "independence") |>
generate(reps = 1000, type = "permute") |>
calculate(stat = "diff in means", order = c("Action", "Romance"))
null_distribution_movies
# A tibble: 1,000 × 2
replicate stat
<int> <dbl>
1 1 0.511
2 2 0.346
3 3 -0.327
4 4 -0.209
5 5 -0.433
6 6 -0.103
7 7 0.387
8 8 0.169
9 9 0.257
10 10 0.334
# ℹ 990 more rows
Observe that we have 1000 values of
stat, each representing one instance of \(\overline{x}_{a} - \overline{x}_{r}\text{.}\) The 1000 values form the null distribution, which is the technical term for the sampling distribution of the difference in sample means \(\overline{x}_{a} - \overline{x}_{r}\) assuming \(H_0\) is true. What happened in real life? What was the observed difference in popularity rates? What was the observed test statistic \(\overline{x}_{a} - \overline{x}_{r}\text{?}\) Recall from our earlier data wrangling, this observed difference in means was \(5.28 - 6.32 = -1.04\text{.}\) We can also compute the observed test statistic \(\overline{x}_{a} - \overline{x}_{r}\) using the code that constructed the null distribution but with the hypothesize() and generate() steps removed. We save this in obs_diff_means:
obs_diff_means <- movies_sample |>
specify(formula = rating ~ genre) |>
calculate(stat = "diff in means", order = c("Action", "Romance"))
obs_diff_means
# A tibble: 1 × 1 stat <dbl> 1 -1.04
Subsubsection 9.6.3.5 5. visualize the p-value
Lastly, in order to compute the \(p\)-value, we have to assess how "extreme" the observed difference in means of \(-1.04\) is. We do this by comparing \(-1.04\) to our null distribution, which was constructed in a hypothesized universe of no true difference in movie ratings. We visualize both the null distribution and the \(p\)-value in Figure 9.6.3. Unlike our example in Subsection 9.4.1 involving music popularity, since we have a two-sided \(H_A: \mu_a - \mu_r \neq 0\text{,}\) we have to allow for both possibilities for more extreme, so we set
direction = "both".
visualize(null_distribution_movies, bins = 10) +
shade_p_value(obs_stat = obs_diff_means, direction = "both")
Let’s go over the elements of this plot. First, the histogram is the null distribution. Second, the solid line is the observed test statistic, or the difference in sample means we observed in real life of \(5.28 - 6.32 = -1.04\text{.}\) Third, the two shaded areas of the histogram form the \(p\)-value, or the probability of obtaining a test statistic just as or more extreme than the observed test statistic assuming the null hypothesis \(H_0\) is true.
What proportion of the null distribution is shaded? In other words, what is the numerical value of the \(p\)-value? We use the
get_p_value() function to compute this value:
null_distribution_movies |>
get_p_value(obs_stat = obs_diff_means, direction = "both")
# A tibble: 1 × 1
p_value
<dbl>
1 0.004
This \(p\)-value of 0.004 is very small. In other words, there is a very small chance that we would observe a difference of \(5.28 - 6.32 = -1.04\) in a hypothesized universe where there was truly no difference in ratings.
But this \(p\)-value is larger than our (even smaller) pre-specified \(\alpha\) significance level of 0.001. Thus, we are inclined to fail to reject the null hypothesis \(H_0: \mu_a - \mu_r = 0\text{.}\) In non-statistical language, the conclusion is: we do not have the evidence needed in this sample of data to suggest that we should reject the hypothesis that there is no difference in mean IMDb ratings between romance and action movies. We, thus, cannot say that a difference exists in romance and action movie ratings, on average, for all IMDb movies.
Exercises 9.6.3.6 Exercises
2.
What conclusions can you make from viewing the faceted histogram looking at
rating versus genre that you could not see when looking at the boxplot?
3.
Describe in a paragraph how we used Allen Downey’s diagram to conclude if a statistical difference existed between mean movie ratings for action and romance movies.
4.
Why are we relatively confident that the distributions of the sample ratings will be good approximations of the population distributions of ratings for the two genres?
5.
Using the definition of \(p\)-value, write in words what the \(p\)-value represents for the hypothesis test comparing the mean rating of romance to action movies.
6.
What is the value of the \(p\)-value for the two-sided hypothesis test comparing the mean rating of romance to action movies?
7.
Test your data-wrangling knowledge and EDA skills:
-
Use
dplyrandtidyrto create the necessary data frame focused on only action and romance movies (but not both) from themoviesdata frame in theggplot2moviespackage. -
Make a boxplot and a faceted histogram of this population data comparing ratings of action and romance movies from IMDb.
-
Discuss how these plots compare to the similar plots produced for the
movies_sampledata.
