Setup

# Εγκατάσταση
# install.packages(c("tidyverse", "pwr", "broom", "scales", "janitor"))

library(tidyverse)
library(pwr)
library(broom)
library(scales)
library(janitor)

set.seed(42)
ads <- read_csv("marketing_AB.csv") |>
  clean_names() |>
  mutate(
    group     = factor(test_group, levels = c("psa", "ad")),
    converted = as.integer(converted)
  )

glimpse(ads)
## Rows: 588,101
## Columns: 8
## $ x1            <dbl> 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16…
## $ user_id       <dbl> 1069124, 1119715, 1144181, 1435133, 1015700, 1137664, 11…
## $ test_group    <chr> "ad", "ad", "ad", "ad", "ad", "ad", "ad", "ad", "ad", "a…
## $ converted     <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0,…
## $ total_ads     <dbl> 130, 93, 21, 355, 276, 734, 264, 17, 21, 142, 209, 47, 6…
## $ most_ads_day  <chr> "Monday", "Tuesday", "Tuesday", "Tuesday", "Friday", "Sa…
## $ most_ads_hour <dbl> 20, 22, 18, 10, 14, 10, 13, 18, 19, 14, 11, 13, 20, 13, …
## $ group         <fct> ad, ad, ad, ad, ad, ad, ad, ad, ad, ad, ad, ad, ad, ad, …

Μέρος Α — Simulated A/B Test

TODO 1: Δημιουργία simulated dataset

n_control   <- 8000
n_treatment <- 8000
p_control   <- 0.08
p_treatment <- 0.10

experiment <- tibble(
  user_id   = 1:(n_control + n_treatment),
  group     = rep(c("control", "treatment"), times = c(n_control, n_treatment)),
  converted = c(
    rbinom(n_control,   1, p_control),
    rbinom(n_treatment, 1, p_treatment)
  )
)

head(experiment)
## # A tibble: 6 × 3
##   user_id group   converted
##     <int> <chr>       <int>
## 1       1 control         0
## 2       2 control         1
## 3       3 control         0
## 4       4 control         0
## 5       5 control         0
## 6       6 control         0

TODO 2: Summary statistics ανά ομάδα

summary_stats <- experiment |>
  group_by(group) |>
  summarise(
    n               = n(),
    conversions     = sum(converted),
    conversion_rate = mean(converted),
    se              = sqrt(conversion_rate * (1 - conversion_rate) / n),
    ci_lower        = conversion_rate - 1.96 * se,
    ci_upper        = conversion_rate + 1.96 * se,
    .groups = "drop"
  )

summary_stats
## # A tibble: 2 × 7
##   group         n conversions conversion_rate      se ci_lower ci_upper
##   <chr>     <int>       <int>           <dbl>   <dbl>    <dbl>    <dbl>
## 1 control    8000         677          0.0846 0.00311   0.0785   0.0907
## 2 treatment  8000         764          0.0955 0.00329   0.0891   0.102

TODO 3: Οπτικοποίηση conversion rate

ggplot(summary_stats, aes(x = group, y = conversion_rate, fill = group)) +
  geom_col(width = 0.5) +
  geom_errorbar(aes(ymin = ci_lower, ymax = ci_upper), width = 0.2) +
  scale_y_continuous(labels = percent) +
  labs(
    title = "Conversion Rate ανά ομάδα (Simulated Data)",
    x = "Ομάδα", y = "Conversion Rate"
  ) +
  theme_minimal() +
  theme(legend.position = "none")

TODO 4: Prop.test

# Πίνακας conversions & n ανά ομάδα
tab <- summary_stats |> arrange(desc(group))  # treatment πρώτο

test_result <- prop.test(
  x = tab$conversions,
  n = tab$n,
  correct = FALSE
)

test_result
## 
##  2-sample test for equality of proportions without continuity correction
## 
## data:  tab$conversions out of tab$n
## X-squared = 5.7725, df = 1, p-value = 0.01628
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  0.002005125 0.019744875
## sample estimates:
##   prop 1   prop 2 
## 0.095500 0.084625

Ερμηνεία: Το p-value είναι 0.0163. Αν p < 0.05, απορρίπτουμε την H₀ ότι τα conversion rates είναι ίσα.

TODO 5: Χειρωνακτική επαλήθευση

# Τιμές από τα simulated data
p1 <- summary_stats$conversion_rate[summary_stats$group == "treatment"]
p2 <- summary_stats$conversion_rate[summary_stats$group == "control"]
n1 <- summary_stats$n[summary_stats$group == "treatment"]
n2 <- summary_stats$n[summary_stats$group == "control"]

# (α) Pooled estimate
p_pool <- (p1 * n1 + p2 * n2) / (n1 + n2)

# (β) Pooled SE
se_pool <- sqrt(p_pool * (1 - p_pool) * (1/n1 + 1/n2))

# (γ) Διαφορά
delta <- p1 - p2

# (δ) 95% CI για τη διαφορά
ci_lower_manual <- delta - 1.96 * se_pool
ci_upper_manual <- delta + 1.96 * se_pool

cat("Pooled p̂:", round(p_pool, 4), "\n")
## Pooled p̂: 0.0901
cat("Pooled SE:", round(se_pool, 4), "\n")
## Pooled SE: 0.0045
cat("Διαφορά δ:", round(delta, 4), "\n")
## Διαφορά δ: 0.0109
cat("95% CI (χειρωνακτικό): [", round(ci_lower_manual, 4), ",", round(ci_upper_manual, 4), "]\n")
## 95% CI (χειρωνακτικό): [ 0.002 , 0.0197 ]
cat("95% CI (prop.test):    [", round(test_result$conf.int[1], 4), ",", round(test_result$conf.int[2], 4), "]\n")
## 95% CI (prop.test):    [ 0.002 , 0.0197 ]

Τα δύο CI είναι πολύ κοντά. Μικρές αποκλίσεις οφείλονται στο ότι το prop.test() χρησιμοποιεί ελαφρώς διαφορετικό τρόπο υπολογισμού SE (Wald vs Wilson), αλλά η ουσία παραμένει η ίδια.

TODO 6: Power Analysis

# Cohen's h
h <- ES.h(p1 = 0.10, p2 = 0.08)
cat("Cohen's h:", round(h, 4), "\n")
## Cohen's h: 0.07
# Ελάχιστο n ανά ομάδα για power = 80%
power_result <- pwr.2p.test(h = h, sig.level = 0.05, power = 0.80)
power_result
## 
##      Difference of proportion power calculation for binomial distribution (arcsine transformation) 
## 
##               h = 0.069988
##               n = 3204.715
##       sig.level = 0.05
##           power = 0.8
##     alternative = two.sided
## 
## NOTE: same sample sizes
cat("\nΑπαιτούμενο n ανά ομάδα:", ceiling(power_result$n), "\n")
## 
## Απαιτούμενο n ανά ομάδα: 3205
cat("Τρέξαμε n ανά ομάδα:", n_control, "\n")
## Τρέξαμε n ανά ομάδα: 8000

Σχολιασμός: Αν τρέξαμε περισσότερα άτομα από το ελάχιστο, σημαίνει ότι είχαμε υπερ-αρκετή ισχύ (overpowered). Αν λιγότερα, θα κινδυνεύαμε να μη βρούμε τη διαφορά ακόμα κι αν υπάρχει.

Μέρος Β — Πραγματικά Δεδομένα (Kaggle Marketing AB)

TODO 7: Invariants check — Τυχαιοποίηση

# (α) Αναλογία ad/psa
ads |>
  count(group) |>
  mutate(pct = n / sum(n))
## # A tibble: 2 × 3
##   group      n    pct
##   <fct>  <int>  <dbl>
## 1 psa    23524 0.0400
## 2 ad    564577 0.960
# (β) Κατανομή ανά ημέρα εβδομάδας και ομάδα
ads |>
  count(group, most_ads_day) |>
  ggplot(aes(x = most_ads_day, y = n, fill = group)) +
  geom_col(position = "dodge") +
  labs(
    title = "Κατανομή χρηστών ανά ημέρα εβδομάδας",
    x = "Ημέρα", y = "Πλήθος"
  ) +
  theme_minimal()

Σχολιασμός: Ελέγχουμε αν η κατανομή χρηστών ανά ημέρα είναι παρόμοια μεταξύ των δύο ομάδων, ώστε να επιβεβαιώσουμε ότι η τυχαιοποίηση δούλεψε σωστά.

TODO 8: Conversion rate & prop.test στα πραγματικά δεδομένα

# Summary stats
real_stats <- ads |>
  group_by(group) |>
  summarise(
    n               = n(),
    conversions     = sum(converted),
    conversion_rate = mean(converted),
    se              = sqrt(conversion_rate * (1 - conversion_rate) / n),
    ci_lower        = conversion_rate - 1.96 * se,
    ci_upper        = conversion_rate + 1.96 * se,
    .groups = "drop"
  )

real_stats
## # A tibble: 2 × 7
##   group      n conversions conversion_rate       se ci_lower ci_upper
##   <fct>  <int>       <int>           <dbl>    <dbl>    <dbl>    <dbl>
## 1 psa    23524         420          0.0179 0.000863   0.0162   0.0195
## 2 ad    564577       14423          0.0255 0.000210   0.0251   0.0260
# Prop test (ad vs psa)
real_test <- prop.test(
  x = real_stats$conversions,
  n = real_stats$n,
  correct = FALSE
)

tidy(real_test)
## # A tibble: 1 × 9
##   estimate1 estimate2 statistic  p.value parameter conf.low conf.high method    
##       <dbl>     <dbl>     <dbl>    <dbl>     <dbl>    <dbl>     <dbl> <chr>     
## 1    0.0179    0.0255      54.3 1.71e-13         1 -0.00943  -0.00595 2-sample …
## # ℹ 1 more variable: alternative <chr>

Ερμηνεία: Το p-value είναι 1.705281e-13. Αν p < 0.05, υπάρχει στατιστικά σημαντική διαφορά μεταξύ ad και psa.

TODO 9: Segmentation — Conversion rate ανά ημέρα

day_stats <- ads |>
  group_by(group, most_ads_day) |>
  summarise(
    n               = n(),
    conversion_rate = mean(converted),
    se              = sqrt(conversion_rate * (1 - conversion_rate) / n),
    ci_lower        = conversion_rate - 1.96 * se,
    ci_upper        = conversion_rate + 1.96 * se,
    .groups = "drop"
  )

# Σωστή σειρά ημερών
day_order <- c("Monday", "Tuesday", "Wednesday", "Thursday", "Friday", "Saturday", "Sunday")
day_stats$most_ads_day <- factor(day_stats$most_ads_day, levels = day_order)

ggplot(day_stats, aes(x = most_ads_day, y = conversion_rate, color = group, group = group)) +
  geom_line(linewidth = 1) +
  geom_point(size = 2) +
  geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper, fill = group), alpha = 0.15, color = NA) +
  scale_y_continuous(labels = percent) +
  labs(
    title = "Conversion Rate ανά ημέρα εβδομάδας",
    x = "Ημέρα", y = "Conversion Rate"
  ) +
  theme_minimal()

Σχολιασμός: Παρατηρούμε αν υπάρχει κάποια ημέρα με αισθητά μεγαλύτερη διαφορά μεταξύ ad και psa. Αυτό θα μπορούσε να αξιοποιηθεί στρατηγικά (π.χ. εστίαση διαφημίσεων τη συγκεκριμένη ημέρα).

TODO 10 (BONUS): Επιχειρηματική απόφαση

# Κατώφλι επιχειρηματικής ουσίας
d_min <- 0.005

# Conversion rates
p_ad  <- real_stats$conversion_rate[real_stats$group == "ad"]
p_psa <- real_stats$conversion_rate[real_stats$group == "psa"]

# Lifts
absolute_lift <- p_ad - p_psa
relative_lift <- absolute_lift / p_psa

cat("Conversion rate (ad):", round(p_ad, 4), "\n")
## Conversion rate (ad): 0.0255
cat("Conversion rate (psa):", round(p_psa, 4), "\n")
## Conversion rate (psa): 0.0179
cat("Absolute lift:", round(absolute_lift, 4), "\n")
## Absolute lift: 0.0077
cat("Relative lift:", percent(relative_lift, accuracy = 0.1), "\n")
## Relative lift: 43.1%
cat("Κατώφλι δ_min:", d_min, "\n")
## Κατώφλι δ_min: 0.005
# CI της διαφοράς
ci_diff <- real_test$conf.int
cat("95% CI της διαφοράς: [", round(ci_diff[1], 5), ",", round(ci_diff[2], 5), "]\n")
## 95% CI της διαφοράς: [ -0.00943 , -0.00595 ]

Ανάλυση 6 περιπτώσεων CI:

  • A: Ολόκληρο το CI > δ_min → Ξεκάθαρα σημαντικό, εφαρμόζουμε
  • B: CI περιλαμβάνει δ_min, αλλά > 0 → Στατιστικά σημαντικό, αλλά αβέβαιο αν πρακτικά σημαντικό
  • C: CI > 0, αλλά < δ_min → Στατιστικά σημαντικό, αλλά πρακτικά ασήμαντο
  • D: CI περιλαμβάνει 0 → Δεν απορρίπτουμε H₀
  • E: CI < 0, αλλά > -δ_min → Αρνητικό αποτέλεσμα, αλλά μικρό
  • F: Ολόκληρο το CI < -δ_min → Η διαφήμιση βλάπτει
# Ποια περίπτωση ισχύει;
if (ci_diff[1] > d_min) {
  cat("Περίπτωση A: Ξεκάθαρα θετικό & πρακτικά σημαντικό. Σύσταση: ΕΦΑΡΜΟΓΗ της καμπάνιας.\n")
} else if (ci_diff[1] > 0 & ci_diff[2] > d_min) {
  cat("Περίπτωση B: Στατιστικά σημαντικό, πρακτική σημασία αβέβαιη. Σύσταση: Πιθανή εφαρμογή, αλλά χρειάζεται περαιτέρω αξιολόγηση.\n")
} else if (ci_diff[1] > 0 & ci_diff[2] < d_min) {
  cat("Περίπτωση C: Στατιστικά σημαντικό αλλά πρακτικά ασήμαντο. Σύσταση: ΜΗ εφαρμογή — το κόστος πιθανώς ξεπερνά το όφελος.\n")
} else if (ci_diff[1] < 0 & ci_diff[2] > 0) {
  cat("Περίπτωση D: Δεν απορρίπτουμε H₀. Σύσταση: ΜΗ εφαρμογή — δεν υπάρχει αρκετή απόδειξη.\n")
} else if (ci_diff[1] > -d_min & ci_diff[2] < 0) {
  cat("Περίπτωση E: Αρνητικό αλλά μικρό αποτέλεσμα. Σύσταση: ΜΗ εφαρμογή.\n")
} else {
  cat("Περίπτωση F: Η διαφήμιση βλάπτει σημαντικά. Σύσταση: ΑΠΟΦΥΓΗ.\n")
}
## Περίπτωση F: Η διαφήμιση βλάπτει σημαντικά. Σύσταση: ΑΠΟΦΥΓΗ.

Απαντήσεις στις ερωτήσεις

1. Ποιο είναι το p-value; Απορρίπτουμε την H₀;

Στα simulated data, το p-value του prop.test() είναι 0.0163. Αν είναι μικρότερο του 0.05, απορρίπτουμε την H₀ ότι τα conversion rates είναι ίσα, δηλαδή η νέα διαφήμιση έχει στατιστικά σημαντικό αποτέλεσμα.

2. Συμπίπτει το χειρωνακτικό CI με αυτό του prop.test();

Τα δύο CI είναι πολύ κοντά αλλά όχι ταυτόσημα. Αυτό συμβαίνει γιατί το prop.test() χρησιμοποιεί τη μέθοδο Wilson score interval, ενώ ο χειρωνακτικός υπολογισμός χρησιμοποιεί Wald interval (normal approximation). Η διαφορά είναι αμελητέα σε μεγάλα δείγματα.

3. Πόσα άτομα χρειαζόντουσαν; Τι συνεπάγεται;

Σύμφωνα με την power analysis, χρειαζόμασταν 3205 άτομα ανά ομάδα για power 80%. Τρέξαμε 8000 ανά ομάδα. Αν 8000 > απαιτούμενο n, το πείραμα ήταν overpowered (καλό — αυξάνει την πιθανότητα ανίχνευσης του πραγματικού effect). Αν 8000 < απαιτούμενο n, θα ήταν underpowered (κίνδυνος να χάσουμε πραγματικό effect).