## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  message = FALSE,
  warning = FALSE
)

## ----setup--------------------------------------------------------------------
library(mariposa)
library(dplyr)
data(survey_data)

## -----------------------------------------------------------------------------
# Without weights (describes the sample)
unweighted <- survey_data %>%
  summarise(
    mean_age = mean(age, na.rm = TRUE),
    mean_income = mean(income, na.rm = TRUE)
  )

# With weights (represents the population)
age_w <- w_mean(survey_data, age, weights = sampling_weight)
income_w <- w_mean(survey_data, income, weights = sampling_weight)

comparison <- data.frame(
  Variable = c("Age", "Income"),
  Unweighted = c(unweighted$mean_age, unweighted$mean_income),
  Weighted = c(age_w$results$weighted_mean, income_w$results$weighted_mean)
)
print(comparison)

## -----------------------------------------------------------------------------
# Weighted mean
w_mean(survey_data, income, weights = sampling_weight)

# Weighted median
w_median(survey_data, income, weights = sampling_weight)

# Weighted mode
w_modus(survey_data, education, weights = sampling_weight)

## -----------------------------------------------------------------------------
# Weighted standard deviation
w_sd(survey_data, income, weights = sampling_weight)

# Weighted variance
w_var(survey_data, income, weights = sampling_weight)

# Weighted interquartile range
w_iqr(survey_data, income, weights = sampling_weight)

## -----------------------------------------------------------------------------
# Weighted skewness
w_skew(survey_data, income, weights = sampling_weight)

# Weighted kurtosis
w_kurtosis(survey_data, income, weights = sampling_weight)

## -----------------------------------------------------------------------------
# Weighted standard error
w_se(survey_data, income, weights = sampling_weight)

# Weighted quantiles
w_quantile(survey_data, income,
           probs = c(0.25, 0.5, 0.75),
           weights = sampling_weight)

# Weighted range
w_range(survey_data, income, weights = sampling_weight)

## -----------------------------------------------------------------------------
w_mean(survey_data, age, income, life_satisfaction,
       weights = sampling_weight)

## -----------------------------------------------------------------------------
actual_n <- nrow(survey_data)
age_w <- w_mean(survey_data, age, weights = sampling_weight)
effective_n <- age_w$results$effective_n

cat("Actual sample size:", actual_n, "\n")
cat("Effective sample size:", round(effective_n), "\n")
cat("Design effect:", round(actual_n / effective_n, 2), "\n")

## -----------------------------------------------------------------------------
# Descriptive statistics
survey_data %>%
  describe(income, life_satisfaction, weights = sampling_weight)

## -----------------------------------------------------------------------------
# Frequency tables
survey_data %>%
  frequency(education, weights = sampling_weight)

## -----------------------------------------------------------------------------
# t-test
survey_data %>%
  t_test(income, group = gender, weights = sampling_weight)

## -----------------------------------------------------------------------------
# ANOVA
survey_data %>%
  oneway_anova(life_satisfaction, group = education,
               weights = sampling_weight)

## -----------------------------------------------------------------------------
# Chi-square
survey_data %>%
  chi_square(education, employment, weights = sampling_weight)

## -----------------------------------------------------------------------------
# Correlation
survey_data %>%
  pearson_cor(age, income, weights = sampling_weight)

## -----------------------------------------------------------------------------
# Regression
survey_data %>%
  linear_regression(life_satisfaction ~ age + income,
                    weights = sampling_weight)

## -----------------------------------------------------------------------------
survey_data %>%
  group_by(region) %>%
  describe(age, income, life_satisfaction,
           weights = sampling_weight)

## -----------------------------------------------------------------------------
survey_data %>%
  group_by(region) %>%
  t_test(income, group = gender, weights = sampling_weight)

## -----------------------------------------------------------------------------
weight_stats <- survey_data %>%
  summarise(
    min = min(sampling_weight, na.rm = TRUE),
    max = max(sampling_weight, na.rm = TRUE),
    mean = mean(sampling_weight, na.rm = TRUE),
    sd = sd(sampling_weight, na.rm = TRUE)
  )
print(weight_stats)

## -----------------------------------------------------------------------------
missing <- sum(is.na(survey_data$sampling_weight))
total <- nrow(survey_data)
cat("Missing weights:", missing, "/", total,
    "(", round(missing / total * 100, 1), "%)\n")

## -----------------------------------------------------------------------------
trimmed <- survey_data %>%
  mutate(
    weight_trimmed = case_when(
      sampling_weight > 5 ~ 5,
      sampling_weight < 0.2 ~ 0.2,
      TRUE ~ sampling_weight
    )
  )

original <- w_mean(survey_data, income, weights = sampling_weight)
trimmed_result <- w_mean(trimmed, income, weights = weight_trimmed)

cat("Original weighted mean:", round(original$results$weighted_mean, 1), "\n")
cat("Trimmed weighted mean:", round(trimmed_result$results$weighted_mean, 1), "\n")

## -----------------------------------------------------------------------------
# 1. Inspect weight properties
survey_data %>%
  summarise(
    n = n(),
    mean_weight = mean(sampling_weight, na.rm = TRUE),
    sd_weight = sd(sampling_weight, na.rm = TRUE),
    min_weight = min(sampling_weight, na.rm = TRUE),
    max_weight = max(sampling_weight, na.rm = TRUE)
  )

# 2. Compare weighted vs unweighted
vars <- c("age", "income", "life_satisfaction")
for (var in vars) {
  unw <- mean(survey_data[[var]], na.rm = TRUE)
  w_result <- w_mean(survey_data, !!sym(var), weights = sampling_weight)
  w <- w_result$results$weighted_mean
  cat(sprintf("%s: unweighted = %.2f, weighted = %.2f (diff = %.1f%%)\n",
              var, unw, w, (w - unw) / unw * 100))
}

# 3. Weighted group comparison
survey_data %>%
  group_by(region) %>%
  describe(income, weights = sampling_weight)

# 4. Weighted hypothesis test
survey_data %>%
  t_test(income, group = gender, weights = sampling_weight)