Romance project: Study 1

Imports

Scores extraction

General love score

# General love score

full_movies$gen_love_score <- as.numeric(str_extract(full_movies$gen_love, "(?<=\\/ SCORE = )\\d+|NA"))

Warning: NAs introduits lors de la conversion automatique

# current_annotation_level <- length(full_movies$gen_love_score)
#                      
# print(current_annotation_level)

Romantic love score

rom_movies$rom_love_score <- as.numeric(str_extract(rom_movies$rom_love, "(?<=\\/ SCORE = )\\d+|NA"))

Warning: NAs introduits lors de la conversion automatique

# Check scores and explanations
# print(rom_movies[1:10, c("rom_love_score", "rom_love")])

Data cleaning: removing outliers

Removing 2020 datapoints:

Validity checks

Comparing average of General Love Score of romance movies vs non romance movies:

# Compare the average of both scores between works categorized under the romance IMDb genre and non-Romance genres. The expectation is that this difference will be significant with a t-test.
# General Love Score
# Filter
gen_love_score_romance <- full_movies$gen_love_score[full_movies$romance == 1]
gen_love_score_non_romance <- full_movies$gen_love_score[full_movies$romance == 0]

# t test
t.test(gen_love_score_romance, gen_love_score_non_romance, mu=0)


    Welch Two Sample t-test

data:  gen_love_score_romance and gen_love_score_non_romance
t = 71.376, df = 2654.5, p-value < 2.2e-16
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 4.082294 4.312931
sample estimates:
mean of x mean of y 
 7.084494  2.886881

# Plotting
mean_gen_love_score_romance <- mean(gen_love_score_romance, na.rm = TRUE)
mean_gen_love_score_non_romance <- mean(gen_love_score_non_romance, na.rm = TRUE)
print(mean_gen_love_score_romance)

[1] 7.084494

print(mean_gen_love_score_non_romance)

[1] 2.886881

scores <- cbind(mean_gen_love_score_romance, mean_gen_love_score_non_romance)
bar_names <- c("Romance", "Non romance")
barplot(scores, beside=TRUE, col=c("pink", "lightblue"), main="Movie Ratings", ylab="Scores", names.arg=bar_names, ylim=c(0,10))

For the General Love Score, their is a significant scores difference between romance and non romance movies (p < 2.2e-16).

Checking annotations

movie_0 <- full_movies[full_movies$gen_love_score == 0, c("gen_love", "gen_love_score")]
movie_2 <- full_movies[full_movies$gen_love_score == 2, c("gen_love", "gen_love_score")]
movie_4 <- full_movies[full_movies$gen_love_score == 4, c("gen_love", "gen_love_score")]
movie_6 <- full_movies[full_movies$gen_love_score == 6, c("gen_love", "gen_love_score")]
movie_8 <- full_movies[full_movies$gen_love_score == 8, c("gen_love", "gen_love_score")]
movie_10 <- full_movies[full_movies$gen_love_score == 10, c("gen_love", "gen_love_score")]

print(movie_8[550,"gen_love"])

# A tibble: 1 × 1
  gen_love                                                                      
  <chr>                                                                         
1 Grease 2 is a musical romantic comedy, and the plot revolves around the love …

Statistical analysis

General Love Score

Models

# A function to include AIC in the upcoming model summaries
summary_AIC <- function(model) {
  summary_output <- summary(model)
  aic <- AIC(model)
  summary_output$aic <- aic
  print(summary_output)
}

Model 1

One linear model with the year as the explanatory variable:

lm1 <- lm(gen_love_score ~ year, data=full_movies)
summary_AIC(lm1)$aic


Call:
lm(formula = gen_love_score ~ year, data = full_movies)

Residuals:
    Min      1Q  Median      3Q     Max 
-3.8634 -1.7753 -0.7019  1.4156  6.8414 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept) -25.781547   3.495211  -7.376 1.73e-13 ***
year          0.014683   0.001747   8.403  < 2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.636 on 12486 degrees of freedom
  (921 observations effacées parce que manquantes)
Multiple R-squared:  0.005623,  Adjusted R-squared:  0.005543 
F-statistic: 70.61 on 1 and 12486 DF,  p-value: < 2.2e-16

[1] 59649.4

Model 2

One linear model with the year as the explanatory variable and the number of movies published that year as the control:

# Number of movies published per year
full_movies_noNA <- full_movies %>%
  filter(gen_love_score != "NA") %>%
  group_by(year) %>%
  mutate(movies_per_year = n()) %>%
  ungroup()

lm2 <- lm(gen_love_score ~ year + movies_per_year, data=full_movies_noNA)
summary_AIC(lm2)$aic


Call:
lm(formula = gen_love_score ~ year + movies_per_year, data = full_movies_noNA)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.1503 -1.7475 -0.7003  1.3998  6.9655 

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)    
(Intercept)     -5.924e+01  8.230e+00  -7.198 6.46e-13 ***
year             3.177e-02  4.187e-03   7.588 3.49e-14 ***
movies_per_year -2.421e-03  5.392e-04  -4.490 7.20e-06 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.634 on 12485 degrees of freedom
Multiple R-squared:  0.007226,  Adjusted R-squared:  0.007067 
F-statistic: 45.44 on 2 and 12485 DF,  p-value: < 2.2e-16

[1] 59631.26

Model 3

One quadratic model with year squared as the explanatory variable:

full_movies$year_sqrt <- full_movies$year^2

lm3 <- lm(gen_love_score ~ year + year_sqrt, data=full_movies)
summary_AIC(lm3)$aic


Call:
lm(formula = gen_love_score ~ year + year_sqrt, data = full_movies)

Residuals:
    Min      1Q  Median      3Q     Max 
-3.8145 -1.7748 -0.7425  1.3505  7.4188 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept) -4.765e+03  5.520e+02  -8.633   <2e-16 ***
year         4.761e+00  5.528e-01   8.613   <2e-16 ***
year_sqrt   -1.188e-03  1.384e-04  -8.586   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.628 on 12485 degrees of freedom
  (921 observations effacées parce que manquantes)
Multiple R-squared:  0.01146,   Adjusted R-squared:  0.0113 
F-statistic: 72.37 on 2 and 12485 DF,  p-value: < 2.2e-16

[1] 59577.88

Model 4

One quadratic model with year squared as the explanatory variable and the number of movies published that year as the control:

full_movies_noNA$year_sqrt <- full_movies_noNA$year^2
lm4 <- lm(gen_love_score ~ year + year_sqrt + movies_per_year, data=full_movies_noNA)
summary_AIC(lm4)$aic


Call:
lm(formula = gen_love_score ~ year + year_sqrt + movies_per_year, 
    data = full_movies_noNA)

Residuals:
    Min      1Q  Median      3Q     Max 
-3.9392 -1.7623 -0.7388  1.3802  7.4428 

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)    
(Intercept)     -4.434e+03  5.699e+02  -7.780 7.83e-15 ***
year             4.421e+00  5.717e-01   7.732 1.14e-14 ***
year_sqrt       -1.101e-03  1.434e-04  -7.677 1.75e-14 ***
movies_per_year -1.299e-03  5.574e-04  -2.331   0.0198 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.628 on 12484 degrees of freedom
Multiple R-squared:  0.01189,   Adjusted R-squared:  0.01165 
F-statistic: 50.07 on 3 and 12484 DF,  p-value: < 2.2e-16

[1] 59574.44

Models comparison

Comparing regression output:

stargazer(lm1, lm2, lm3, lm4, type="html", digits=2, title="Regression models comparison for General Love Score", style="qje", column.labels=c("OLS", "OLS + control", "OLS quadratic", "OLS quadratic + control"), dep.var.labels = "General Love score")

**Regression models comparison for General Love Score**

	General Love score
	OLS	OLS + control	OLS quadratic	OLS quadratic + control
	(1)	(2)	(3)	(4)

year	0.01^***	0.03^***	4.76^***	4.42^***
	(0.002)	(0.004)	(0.55)	(0.57)

movies_per_year		-0.002^***		-0.001^**
		(0.001)		(0.001)

year_sqrt			-0.001^***	-0.001^***
			(0.0001)	(0.0001)

Constant	-25.78^***	-59.24^***	-4,765.06^***	-4,433.68^***
	(3.50)	(8.23)	(551.97)	(569.89)

N	12,488	12,488	12,488	12,488
R²	0.01	0.01	0.01	0.01
Adjusted R²	0.01	0.01	0.01	0.01
Residual Std. Error	2.64 (df = 12486)	2.63 (df = 12485)	2.63 (df = 12485)	2.63 (df = 12484)
F Statistic	70.61^*** (df = 1; 12486)	45.44^*** (df = 2; 12485)	72.37^*** (df = 2; 12485)	50.07^*** (df = 3; 12484)

Notes:	^***Significant at the 1 percent level.
	^**Significant at the 5 percent level.
	^*Significant at the 10 percent level.

Comparing AIC scores:

[1] "Model 1 AIC: 59649.4013988007"

[1] "Model 2 AIC: 59631.2553243927"

[1] "Model 3 AIC: 59577.8770675681"

[1] "Model 4 AIC: 59574.4436038523"

Romantic Love Score

Models

Model 1

One linear model with the year as the explanatory variable:

lm1_rom <- lm(rom_love_score ~ year, data=rom_movies)
summary_AIC(lm1_rom)$aic


Call:
lm(formula = rom_love_score ~ year, data = rom_movies)

Residuals:
    Min      1Q  Median      3Q     Max 
-3.6632 -1.5431 -0.5792  0.9131  6.9131 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept) -20.579694   3.265953  -6.301 3.05e-10 ***
year          0.012007   0.001633   7.354 2.04e-13 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.473 on 12584 degrees of freedom
  (823 observations effacées parce que manquantes)
Multiple R-squared:  0.004279,  Adjusted R-squared:  0.0042 
F-statistic: 54.08 on 1 and 12584 DF,  p-value: 2.045e-13

[1] 58512.03

Model 2

One linear model with the year as the explanatory variable and the number of movies published that year as the control:

# Number of movies published per year
rom_movies_noNA <- rom_movies %>%
  filter(rom_love_score != "NA") %>%
  group_by(year) %>%
  mutate(movies_per_year = n()) %>%
  ungroup()

lm2_rom <- lm(rom_love_score ~ year + movies_per_year, data=rom_movies_noNA)
summary_AIC(lm2_rom)$aic


Call:
lm(formula = rom_love_score ~ year + movies_per_year, data = rom_movies_noNA)

Residuals:
    Min      1Q  Median      3Q     Max 
-3.9511 -1.5247 -0.5598  1.0400  7.0400 

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)    
(Intercept)     -5.373e+01  7.608e+00  -7.062 1.73e-12 ***
year             2.894e-02  3.870e-03   7.476 8.17e-14 ***
movies_per_year -2.382e-03  4.938e-04  -4.823 1.43e-06 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.471 on 12583 degrees of freedom
Multiple R-squared:  0.006116,  Adjusted R-squared:  0.005958 
F-statistic: 38.72 on 2 and 12583 DF,  p-value: < 2.2e-16

[1] 58490.79

Model 3

One quadratic model with year squared as the explanatory variable:

rom_movies$year_sqrt <- rom_movies$year^2

lm3_rom <- lm(rom_love_score ~ year + year_sqrt, data=rom_movies)
summary_AIC(lm3_rom)$aic


Call:
lm(formula = rom_love_score ~ year + year_sqrt, data = rom_movies)

Residuals:
   Min     1Q Median     3Q    Max 
-3.644 -1.582 -0.598  1.356  7.450 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept) -4.426e+03  5.159e+02  -8.579   <2e-16 ***
year         4.424e+00  5.167e-01   8.562   <2e-16 ***
year_sqrt   -1.105e-03  1.294e-04  -8.539   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.466 on 12583 degrees of freedom
  (823 observations effacées parce que manquantes)
Multiple R-squared:  0.01002,   Adjusted R-squared:  0.009859 
F-statistic: 63.65 on 2 and 12583 DF,  p-value: < 2.2e-16

[1] 58441.31

Model 4

One quadratic model with year squared as the explanatory variable and the number of movies published that year as the control:

rom_movies_noNA$year_sqrt <- rom_movies_noNA$year^2
lm4_rom <- lm(rom_love_score ~ year + year_sqrt + movies_per_year, data=rom_movies_noNA)
summary_AIC(lm4_rom)$aic


Call:
lm(formula = rom_love_score ~ year + year_sqrt + movies_per_year, 
    data = rom_movies_noNA)

Residuals:
   Min     1Q Median     3Q    Max 
-3.780 -1.570 -0.593  1.250  7.477 

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)    
(Intercept)     -4.069e+03  5.324e+02  -7.642 2.29e-14 ***
year             4.057e+00  5.341e-01   7.596 3.26e-14 ***
year_sqrt       -1.010e-03  1.340e-04  -7.542 4.93e-14 ***
movies_per_year -1.379e-03  5.104e-04  -2.701  0.00691 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.465 on 12582 degrees of freedom
Multiple R-squared:  0.01059,   Adjusted R-squared:  0.01035 
F-statistic: 44.89 on 3 and 12582 DF,  p-value: < 2.2e-16

[1] 58436.01

Models comparison

Comparing regression outputs:

stargazer(lm1_rom, lm2_rom, lm3_rom, lm4_rom, type="html", digits=2, title="Regression models comparison for Romantic Love Score", style="qje", column.labels=c("OLS", "OLS + control", "OLS quadratic", "OLS quadratic + control"), dep.var.labels = "Romantic Love score")

**Regression models comparison for Romantic Love Score**

	Romantic Love score
	OLS	OLS + control	OLS quadratic	OLS quadratic + control
	(1)	(2)	(3)	(4)

year	0.01^***	0.03^***	4.42^***	4.06^***
	(0.002)	(0.004)	(0.52)	(0.53)

movies_per_year		-0.002^***		-0.001^***
		(0.0005)		(0.001)

year_sqrt			-0.001^***	-0.001^***
			(0.0001)	(0.0001)

Constant	-20.58^***	-53.73^***	-4,425.99^***	-4,069.12^***
	(3.27)	(7.61)	(515.92)	(532.44)

N	12,586	12,586	12,586	12,586
R²	0.004	0.01	0.01	0.01
Adjusted R²	0.004	0.01	0.01	0.01
Residual Std. Error	2.47 (df = 12584)	2.47 (df = 12583)	2.47 (df = 12583)	2.47 (df = 12582)
F Statistic	54.08^*** (df = 1; 12584)	38.72^*** (df = 2; 12583)	63.65^*** (df = 2; 12583)	44.89^*** (df = 3; 12582)

Notes:	^***Significant at the 1 percent level.
	^**Significant at the 5 percent level.
	^*Significant at the 10 percent level.

Comparing AIC scores:

[1] "Model 1 AIC: 58512.0317501282"

[1] "Model 2 AIC: 58490.7858700247"

[1] "Model 3 AIC: 58441.3072774767"

[1] "Model 4 AIC: 58436.0096787003"

Differential love score

rom_scores <- rom_movies$rom_love_score
print(length(rom_scores))
gen_scores <- full_movies$gen_love_score
print(length(gen_scores))
year <- full_movies$year
print(length(year))
scores_df <- data.frame(year, rom_scores, gen_scores)
print(nrow(scores_df))
scores_df$diff_love_score <- scores_df$gen_scores - scores_df$rom_scores

Modelling differential love score

Plots

General Love scores plot

full_movies_noNA_plot <- full_movies_noNA %>%
  group_by(year) %>%
  mutate(mean_gen_love_score = (mean(gen_love_score))) %>%
  ungroup()
  
gen_love_plot <- full_movies_noNA_plot %>%
  ggplot(aes(x=as.numeric(year), y=mean_gen_love_score)) + geom_point(color="#800080", size=2, alpha=0.4) +
  labs(title="Evolution of General Love scores") +
    xlab("Year") +
    ylab("Mean score") +
  theme_minimal() +
  theme(
        plot.title = element_text(face="bold", size=14, hjust = 0.5),  # Adjust title appearance
    axis.title.x = element_text(face="bold", size=12),  # Adjust x-axis label appearance
    axis.title.y = element_text(face="bold", size=12),  # Adjust y-axis label appearance
    axis.text.x = element_text(size=10),  # Adjust x-axis tick labels appearance
    axis.text.y = element_text(size=10),  # Adjust y-axis tick labels appearance
    legend.position = "none"
  )
  
gen_love_plot

Romantic love scores plot

rom_movies_noNA_plot <- rom_movies_noNA %>%
  group_by(year) %>%
  mutate(mean_rom_love_score = (mean(rom_love_score))) %>%
  ungroup()
  
rom_love_plot <- rom_movies_noNA_plot %>%
  ggplot(aes(x=as.numeric(year), y=mean_rom_love_score)) + geom_point(color="#800080", size=2, alpha=0.4) +
  labs(title="Evolution of Romantic Love scores") +
    xlab("Year") +
    ylab("Mean score") +
  theme_minimal() +
  theme(
        plot.title = element_text(face="bold", size=14, hjust = 0.5),  # Adjust title appearance
    axis.title.x = element_text(face="bold", size=12),  # Adjust x-axis label appearance
    axis.title.y = element_text(face="bold", size=12),  # Adjust y-axis label appearance
    axis.text.x = element_text(size=10),  # Adjust x-axis tick labels appearance
    axis.text.y = element_text(size=10),  # Adjust y-axis tick labels appearance
    legend.position = "none"
  )
  
rom_love_plot

Differential love scores plot

# Remove NA values
scores_df_noNA <- scores_df %>%
  filter(gen_scores != "NA") %>%
  filter(rom_scores != "NA")
  
scores_df_noNA_plot <- scores_df_noNA %>%
  group_by(year) %>%
  mutate(mean_gen_scores = (mean(gen_scores))) %>%
  mutate(mean_rom_scores = (mean(rom_scores))) %>%
  mutate(mean_diff_love_score = (mean(diff_love_score))) %>%
  ungroup()
  
scores_df_plot_loess <- scores_df_noNA_plot %>%
  ggplot(aes(x=as.numeric(year), y=mean_diff_love_score)) + geom_point(color="#FF0066", size=2, alpha=0.4) +
  geom_smooth(method="loess", aes(x=as.numeric(year), y=mean_diff_love_score), se=TRUE, color="#B5121E") +
  labs(title="Evolution of Differential Love scores (General Love - Romantic Love) - LOESS method") +
    xlab("Year") +
    ylab("Mean score difference") +
  theme_minimal() +
  theme(
        plot.title = element_text(family="Bookman", size=14, hjust = 0.5),  # Adjust title appearance
    axis.title.x = element_text(family="Bookman", size=12),  # Adjust x-axis label appearance
    axis.title.y = element_text(family="Bookman", size=12),  # Adjust y-axis label appearance
    axis.text.x = element_text(size=10),  # Adjust x-axis tick labels appearance
    axis.text.y = element_text(size=10),  # Adjust y-axis tick labels appearance
    legend.position = "none"
  )
  
scores_df_plot_loess

scores_df_plot_lm <- scores_df_noNA_plot %>%
  ggplot(aes(x=as.numeric(year), y=mean_diff_love_score)) + geom_point(color="#FF0066", size=2, alpha=0.4) +
  geom_smooth(method="lm", aes(x=as.numeric(year), y=mean_diff_love_score), se=TRUE, color="#B5121E") +
  labs(title="Evolution of Differential Love scores (General Love - Romantic Love) _ lm method") +
    xlab("Year") +
    ylab("Mean score difference") +
  theme_minimal() +
  theme(
        plot.title = element_text(family="Bookman", size=14, hjust = 0.5),  # Adjust title appearance
    axis.title.x = element_text(family="Bookman", size=12),  # Adjust x-axis label appearance
    axis.title.y = element_text(family="Bookman", size=12),  # Adjust y-axis label appearance
    axis.text.x = element_text(size=10),  # Adjust x-axis tick labels appearance
    axis.text.y = element_text(size=10),  # Adjust y-axis tick labels appearance
    legend.position = "none"
  )
  
scores_df_plot_lm

Comparing mean general scores and mean romantic scores on the same plot:

# I use scores_df_noNA_plot because it contains the mean scores per year
library(tidyr)

Warning: le package 'tidyr' a été compilé avec la version R 4.3.1

scores_df_long <- scores_df_noNA_plot %>%
  filter(year != 2020) %>%
  pivot_longer(cols=c(mean_rom_scores, mean_gen_scores), names_to="mean_score_type", values_to="mean_score")

colors <- c("mean_rom_scores" = "#83222C", "mean_gen_scores"="#E48896")
all_scores_plot <- scores_df_long %>%
  ggplot(aes(x=year, y=mean_score, color=mean_score_type)) +
  geom_point() +
  geom_smooth(method="loess", aes(group=mean_score_type, color=mean_score_type)) +
                labs(title="Evolution of general and romantic love scores", x="Year", y="Score", color="Type of love score") +
  scale_color_manual(values=colors) +
                theme_minimal()
all_scores_plot