## ----setup, message = FALSE, warning = FALSE----------------------------------
library(rankingQ)
data(identity)

library(dplyr)
library(ggplot2)
library(tidyr)

## -----------------------------------------------------------------------------
round(prop.table(table(identity$anc_correct_identity)) * 100, digits = 1)

## -----------------------------------------------------------------------------
identity$random_identity <- case_when(
  identity$anc_correct_identity == 1 ~ 0,
  TRUE ~ 1
)

unbiased_correct_prop(
  sum(identity$random_identity == 0) / sum(!is.na(identity$random_identity)),
  J = 4
)

## -----------------------------------------------------------------------------
## party, religion, gender, and race hold the marginal rank of each item

# Perform bias correction
out_direct <- imprr_direct(
  data = identity,
  ## Not strictly necessary: when `main_q` lists the ranking columns,
  ## J is inferred from their number
  J = 4,
  main_q = c("party", "religion", "gender", "race"),
  anc_correct = "anc_correct_identity",
  # setting to 10 only for our vignette
  n_bootstrap = 10
)

## -----------------------------------------------------------------------------
# Bias correction for the entire population with the uniform assumption

out_direct_uniform <- imprr_direct(
  data = identity,
  J = 4,
  main_q = c("party", "religion", "gender", "race"),
  anc_correct = "anc_correct_identity",
  population = "all",
  assumption = "uniform",
  n_bootstrap = 10
)

## -----------------------------------------------------------------------------
# Bias correction for the entire population with contaminated sampling

out_direct_contaminated <- imprr_direct(
  data = identity,
  J = 4,
  main_q = c("party", "religion", "gender", "race"),
  anc_correct = "anc_correct_identity",
  population = "all",
  assumption = "contaminated",
  n_bootstrap = 10
)

## -----------------------------------------------------------------------------
# Estimated proportion of random responses with a 95% CI
out_direct$est_p_random

## -----------------------------------------------------------------------------
# View the results based on the quantity of interest
out_direct$results %>%
  filter(qoi == "average rank")

# View the results based on the item
out_direct$results %>%
  filter(item == "party")

## -----------------------------------------------------------------------------
# Plot the result
out_direct$results %>%
  mutate(
    item = factor(
      item,
      levels = c("party", "religion", "gender", "race")
    )
  ) %>%
  plot_avg_ranking()

## -----------------------------------------------------------------------------
# Perform bias correction
out_weights <- imprr_weights(
  data = identity,
  J = 4,
  main_q = c("party", "religion", "gender", "race"),
  anc_correct = "anc_correct_identity"
)

## -----------------------------------------------------------------------------
# Perform bias correction with the uniform preference assumption
out_weights_uniform <- imprr_weights(
  data = identity,
  J = 4,
  main_q = c("party", "religion", "gender", "race"),
  anc_correct = "anc_correct_identity",
  population = "all",
  assumption = "uniform"
)

## -----------------------------------------------------------------------------
# Perform bias correction with the uniform preference assumption
out_weights_contaminated <- imprr_weights(
  data = identity,
  J = 4,
  main_q = c("party", "religion", "gender", "race"),
  anc_correct = "anc_correct_identity",
  population = "all",
  assumption = "contaminated"
)

## -----------------------------------------------------------------------------
# View the estimated weights
out_weights$rankings %>%
  select(ranking, weights)

## -----------------------------------------------------------------------------
# View the estimated corrected PMF
out_weights$rankings %>%
  select(ranking, prop_obs, prop_bc)

## -----------------------------------------------------------------------------
identity_w <- out_weights$results
head(identity_w)
# save(identity_w, file = "data/identity_w.rda")