Rationale

The Cultivation Theory states that with heavy or prolonged viewership to television or other forms of electronic media that there is a perceived influence on how people view the world around them.

With this in mind, the influence that people receive from their viewing habits from television could impact how they view different facets of their life. This could range from how they view their immediate household or their neighborhood, and range all the way up to the city or country that they live in.

Hypothesis

The average amount of weekly hours a participant watches television is linear and positively related to how much of the U.S. population the participant estimates works full-time in law enforcement/ criminal justice, medicine, or emergency response services.

Variables & Method

400 volunteer study participants from a random sample of U.S. adults had their television habits monitored for six months with a device tracking only their watching habits and not others in the household. This devices worked on televisions and other devices used for television content consumption. Following data being collected for six months, participants were given a questionnaire. This questionnaire gathered the data on how the participants estimate the U.S. population that works full-time in law enforcement/ criminal justice, medicine, or emergency response services.

The dependent variable of the analysis was the quantitative measure of how participants estimated the U.S. population’s careers in law enforcement/ criminal justice, medicine, or emergency response services. The independent variable was the quantitative measure of how many hours a week participants watched television.

A bivariate linear regression was conducted to examine whether the association between weekly television viewing and the participants’ estimated percentage of people working in these fields was statistically significant. 

Results & Discussion

The scatterplot and regression table below summarize the relationship and association between the independent and dependent variable.

## `geom_smooth()` using formula = 'y ~ x'

## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.

## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.

Regression Analysis Results
Coefficient Estimates
Term Estimate Std. Error t p-value
(Intercept) 14.9525 1.4656 10.2026 0.0000
IV 0.3686 0.0275 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 6.4347
Leverage estimates for 10 largest outliers
Row # Leverage
359 0.0275
346 0.0236
108 0.0212
198 0.0168
236 0.0168
191 0.0158
388 0.0158
39 0.0130
333 0.0130
392 0.0130

These results support the hypothesis, as the regression showed a weak, yet positive correlation. It shows a statistically significant relationship between hours viewing television each week and the estimated percentage of the U.S. population in law enforcement/ criminal justice, medicine, or emergency response services.  

Code

##################################################
# 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$video
mydata$IV <- mydata$pct


##################################################
# 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