Cultivation theory is the idea that television consumption and exposure can shape people’s perceptions of the real world. Any content a person may see on TV over a long period of time can affect how they view certain aspects of reality even if those perceptions are not true. Assumptions and societal norms about crime, race, gender, and even professions may be affected depending on how much TV a person watches per week.
There is a positive correlation between how many hours per week viewers watch TV and their estimated percentage of the U.S. population who is employed full-time in either law enforcement/criminal justice, medicine, or emergency response services.
A random sample of 400 volunteer U.S. adult participants were monitored for their weekly hours spent watching television over a six month span. They agreed to connect a monitoring device to their household televisions, which was configured to only count the hours the participant spent watching. At the end of the six month period, the participants completed a questionnaire. In it, they were asked about the estimated percentages of the U.S population who were employed full time in one of these three worker categories: law enforcement/criminal justice, medicine, and emergency services. Afterwards, the percentages were summed.
The dependent variable in the analysis is a continuous measure of the participants’ estimated percentage. The independent variable is a continuous measure of the number of weekly hours the participants had watched television.
A bivariate regression was used to test the correlation between the two variables and whether or not that correlation is statistically significant.
The graph below shows the correlation between the two variables with the regression line. The regression table shows a more detail regarding the results of the scatter plot and that the relationship is statistically significant. The model fit table shows the correlation’s strength and the overall residual.
## `geom_smooth()` using formula = 'y ~ x'
| Regression Analysis Results | ||||
| Coefficient Estimates | ||||
| Term | Estimate | Std. Error | t | p-value |
|---|---|---|---|---|
| (Intercept) | 23.2076 | 2.2026 | 10.5363 | 0.0000 |
| IV | 0.8440 | 0.0630 | 13.4056 | 0.0000 |
| Model Fit Statistics | |||||
| Overall Regression Performance | |||||
| R-squared | Adj. R-squared | F-statistic | df (model) | df (residual) | Residual Std. Error |
|---|---|---|---|---|---|
| 0.3111 | 0.3093 | 179.7107 | 1.0000 | 398.0000 | 9.7373 |
The results provided here support the hypothesis. There is a positive correlation between the two variables; however, the correlation is not strong. The R-square according to the model fit table is 31.1%, which indicates a weak positive correlation.
Below is the code used to conduct the analysis.
##################################################
# 1. Install and load required packages
##################################################
if (!require("tidyverse")) install.packages("tidyverse")
if (!require("gt")) install.packages("gt")
if (!require("gtExtras")) install.packages("gtExtras")
library(tidyverse)
library(gt)
library(gtExtras)
##################################################
# 2. Read in the dataset
##################################################
# Replace "YOURFILENAME.csv" with the actual filename
mydata <- read.csv("Cultivation.csv")
# ################################################
# # (Optional) 2b. Remove specific cases by row number
# ################################################
# # Example: remove rows 10 and 25
# rows_to_remove <- c(10, 25) # Edit and uncomment this line
# mydata <- mydata[-rows_to_remove, ] # Uncomment this line
##################################################
# 3. Define dependent variable (DV) and independent variable (IV)
##################################################
# Replace YOURDVNAME and YOURIVNAME with actual column names
mydata$DV <- mydata$pct
mydata$IV <- mydata$video
##################################################
# 4. Explore distributions of DV and IV
##################################################
# Make a histogram for DV
DVGraph <- ggplot(mydata, aes(x = DV)) +
geom_histogram(color = "black", fill = "#1f78b4")
# Make a histogram for IV
IVGraph <- ggplot(mydata, aes(x = IV)) +
geom_histogram(color = "black", fill = "#1f78b4")
##################################################
# 5. Fit and summarize initial regression model
##################################################
# Suppress scientific notation
options(scipen = 999)
# Fit model
myreg <- lm(DV ~ IV, data = mydata)
# Model summary
summary(myreg)
##################################################
# 6. Visualize regression and check for bivariate outliers
##################################################
# Create scatterplot with regression line as a ggplot object
RegressionPlot <- ggplot(mydata, aes(x = IV, y = DV)) +
geom_point(color = "#1f78b4") +
geom_smooth(method = "lm", se = FALSE, color = "red") +
labs(
title = "Scatterplot of DV vs IV with Regression Line",
x = "Independent Variable (IV)",
y = "Dependent Variable (DV)"
) +
theme_minimal()
##################################################
# 7. Check for potential outliers (high leverage points)
##################################################
# Calculate leverage values
hat_vals <- hatvalues(myreg)
# Rule of thumb: leverage > 2 * (number of predictors + 1) / n may be influential
threshold <- 2 * (length(coef(myreg)) / nrow(mydata))
# Create table showing 10 largest leverage values
outliers <- data.frame(
Obs = 1:nrow(mydata),
Leverage = hatvalues(myreg)
) %>%
arrange(desc(Leverage)) %>%
slice_head(n = 10)
# Format as a gt table
outliers_table <- outliers %>%
gt() %>%
tab_header(
title = "Leverage estimates for 10 largest outliers"
) %>%
cols_label(
Obs = "Row #",
Leverage = "Leverage"
) %>%
fmt_number(
columns = Leverage,
decimals = 4
)
##################################################
# 8. Create nicely formatted regression results tables
##################################################
# --- Coefficient-level results ---
reg_results <- as.data.frame(coef(summary(myreg))) %>%
tibble::rownames_to_column("Term") %>%
rename(
Estimate = Estimate,
`Std. Error` = `Std. Error`,
t = `t value`,
`p-value` = `Pr(>|t|)`
)
reg_table <- reg_results %>%
gt() %>%
tab_header(
title = "Regression Analysis Results",
subtitle = "Coefficient Estimates"
) %>%
fmt_number(
columns = c(Estimate, `Std. Error`, t, `p-value`),
decimals = 4
)
# --- Model fit statistics ---
reg_summary <- summary(myreg)
fit_stats <- tibble::tibble(
`R-squared` = reg_summary$r.squared,
`Adj. R-squared` = reg_summary$adj.r.squared,
`F-statistic` = reg_summary$fstatistic[1],
`df (model)` = reg_summary$fstatistic[2],
`df (residual)` = reg_summary$fstatistic[3],
`Residual Std. Error` = reg_summary$sigma
)
fit_table <- fit_stats %>%
gt() %>%
tab_header(
title = "Model Fit Statistics",
subtitle = "Overall Regression Performance"
) %>%
fmt_number(
columns = everything(),
decimals = 4
)
##################################################
# 9. Final print of key graphics and tables
##################################################
DVGraph
IVGraph
RegressionPlot
outliers_table
reg_table
fit_table