# Setup
if (!require("dplyr")) install.packages("dplyr")
## Loading required package: dplyr
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
if (!require("ggplot2")) install.packages("ggplot2")
## Loading required package: ggplot2
if (!require("gt")) install.packages("gt")
## Loading required package: gt
if (!require("gtExtras")) install.packages("gtExtras")
## Loading required package: gtExtras
library(dplyr)
library(ggplot2)
library(gt)
library(gtExtras)

options(scipen = 999)
# Load Data
mydata <- read.csv("SandwichAd.csv") #Edit: TYPE YOUR DATA FILE'S NAME

mydata <- mydata %>%
  mutate(
    DV = Calories, #Edit: TYPE YOUR DEPENDENT VARIABLE'S NAME
    IV = Actor  #Edit: TYPE YOUR INDEPENDENT VARIABLE'S NAME
  )
# Histograms of DV per IV group with group means
# Calculate group means
group_means <- mydata %>%
  group_by(IV) %>%
  summarise(mean_DV = mean(DV, na.rm = TRUE), .groups = "drop")

Graphic <- ggplot(mydata, aes(x = DV)) +
  geom_histogram(binwidth = diff(range(mydata$DV, na.rm = TRUE)) / 30,
                 color = "black", fill = "#1f78b4", alpha = 0.7) +
  # Add group-specific mean lines
  geom_vline(data = group_means, aes(xintercept = mean_DV),
             color = "red", linetype = "dashed", linewidth = 1) +
  facet_grid(IV ~ .) +  # one histogram per group, stacked vertically
  labs(
    title = "DV Distributions by IV Group",
    x = "Dependent Variable",
    y = "Count"
  ) +
  theme_minimal()

# Show graphic
Graphic

# Descriptive Statistics
mydata %>%
  group_by(IV) %>%
  summarise(
    count = n(),
    mean  = mean(DV, na.rm = TRUE),
    sd    = sd(DV, na.rm = TRUE),
    min   = min(DV, na.rm = TRUE),
    max   = max(DV, na.rm = TRUE)
  )
## # A tibble: 2 × 6
##   IV      count  mean    sd   min   max
##   <chr>   <int> <dbl> <dbl> <int> <int>
## 1 Buff       50  584.  128.   316   827
## 2 Dad bod    50  656.  141.   358   946
# Normality Check (Shapiro-Wilk)
mydata %>%
  group_by(IV) %>%
  summarise(
    W_statistic = shapiro.test(DV)$statistic,
    p_value     = shapiro.test(DV)$p.value
  )
## # A tibble: 2 × 3
##   IV      W_statistic p_value
##   <chr>         <dbl>   <dbl>
## 1 Buff          0.976   0.385
## 2 Dad bod       0.985   0.790
# Inferential Tests
# Run Welch's t-test (default for two groups)
t_res <- t.test(DV ~ IV, data = mydata, var.equal = FALSE)

# Run Wilcoxon rank-sum test if group sizes < 40 and 
# distributions are non-normal
wilcox.test(mydata$DV ~ mydata$IV)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  mydata$DV by mydata$IV
## W = 882, p-value = 0.01129
## alternative hypothesis: true location shift is not equal to 0
# Create a tidy summary of results
t_summary <- tibble(
  Group1 = levels(as.factor(mydata$IV))[1],
  Group2 = levels(as.factor(mydata$IV))[2],
  Mean1  = t_res$estimate[1],
  Mean2  = t_res$estimate[2],
  t      = t_res$statistic,
  df     = t_res$parameter,
  p      = t_res$p.value,
  CI_low = t_res$conf.int[1],
  CI_high= t_res$conf.int[2]
)
# Present Results as a gt Table (APA style)
Table <- t_summary %>%
  gt() %>%
  # Round group means and CI to 2 decimals
  fmt_number(columns = c(Mean1, Mean2, CI_low, CI_high), decimals = 2) %>%
  # Round test statistics, df, and p-value to 3 decimals
  fmt_number(columns = c(t, df, p), decimals = 3) %>%
  tab_header(
    title = "Independent Samples t-Test Results",
    subtitle = "Welch's t-test (unequal variances assumed)"
  ) %>%
  cols_label(
    Group1 = "Group 1",
    Group2 = "Group 2",
    Mean1  = "Mean (Group 1)", 
    Mean2  = "Mean (Group 2)",
    t      = "t Statistic",
    df     = "Degrees of Freedom",
    p      = "p-value",
    CI_low = "95% CI (Lower)",
    CI_high= "95% CI (Upper)"
  )
# Visual output
# Show the graphic
Graphic

# Show the table
Table
Independent Samples t-Test Results
Welch's t-test (unequal variances assumed)
Group 1 Group 2 Mean (Group 1) Mean (Group 2) t Statistic Degrees of Freedom p-value 95% CI (Lower) 95% CI (Upper)
Buff Dad bod 584.26 656.28 −2.671 97.069 0.009 −125.53 −18.51