Model: decision, baseline
The reference profiled individual has the following characteristics:
- Romanian
- Female
- No German proficiency
- Doctor
- 40 years old
Therefore, the values of the Intercept correspond to this profile.
# Use a reference dataset different from the original,
# in case we want to apply filters or perform transformations.
d <- E %>%
mutate(Observation = 1:n())
# Define the outcome variable and calculate its size.
Y <- d$Decision
nO <- length(Y)
# The profile factors are cleaned, establishing reference categories and
# selecting the relevant variables that will be passed onto the final matrix.
# Store their names to be used later in the figures.
X <- d %>%
mutate(`Nationality Netherlands` = ifelse(Nationality == "Netherlands", 1, 0)) %>%
mutate(`Gender Male` = ifelse(Gender == "Male", 1, 0)) %>%
mutate(`Profession Nurse` = ifelse(Profession == "Nurse", 1, 0)) %>%
mutate(`Language Full` = ifelse(Language == "Full", 1, 0)) %>%
mutate(`Age 25` = ifelse(Age == "25", 1, 0)) %>%
mutate(`Age 55` = ifelse(Age == "55", 1, 0)) %>%
mutate(`First shown` = ifelse(ProfileOrder == 1, 1, 0)) %>%
select(`Nationality Netherlands`,
`Gender Male`,
`Profession Nurse`,
`Language Full`,
`Age 25`,
`Age 55`,
`First shown`) %>%
as.matrix()
X <- cbind("(Intercept)" = 1, X)
nF <- dim(X)[[2]]
feature.label <- dimnames(X)[[2]]
feature.label.order <- c("Nationality Netherlands",
"Language Full",
"Profession Nurse",
"Age 25",
"Age 55",
"Gender Male")
# Define a matrix for later use in multivariate normal priors.
b0 <- rep(0, nF)
B0 <- diag(nF)
diag(B0) <- 2.5^-2
# Identify outcome, exercise, population and treatment with numbers,
# count their size and keep labels for later.
id.outcome <- as.numeric(d$Outcome)
outcome.label <- levels(d$Outcome)
nOutcome <- length(outcome.label)
id.exercise <- d$Exercise
exercise.label <- sort(unique(d$Exercise))
nE <- length(exercise.label)
id.population <- as.numeric(d$Population)
population.label <- levels(d$Population)
nP <- length(population.label)
id.treatment <- as.numeric(d$Treatment)
treatment.label <- levels(d$Treatment)
nT <- length(treatment.label)
first.shown.profile <- ifelse(d$ProfileOrder == 1, 1, 0)
# Store all the variables needed by JAGS in a single list
D <- list(
nO = nO,
id_outcome = id.outcome, nOutcome = nOutcome,
id_exercise = id.exercise , nE = nE,
id_population = id.population, nP = nP,
id_treatment = id.treatment, nT = nT,
first_shown_profile = first.shown.profile,
X = unname(X), nF = nF, b0 = b0, B0 = B0,
Y = Y)Model description:
\[ \begin{aligned} Y_{o} & \sim & \mathcal{B}(\pi_{o}) && \text{Main data component} \\ \pi_{o} & = & \text{logit} (\theta_{O,p,t,f} F_{o,f}) && \text{Linear relationship} \\ % \theta_{O,p,t,f} & \sim &\mathcal{N}(\Theta_{O,f}, \sigma_{\theta_{O,f}}) && \text{Priors for explanatory variables} \\ % \Theta_{O,f} & \sim &\mathcal{N}(\Omega_{f}, \sigma_{\Theta_{f}}) && \text{Prior for the effects shared by outcome} \\ % \sigma_{\theta_{O,f}} & \sim & \mathcal{IG}(1, 1) && \text{Prior for SD, pooling of population and treatment} \\ % \sigma_{\Theta_{f}} & \sim & \mathcal{IG}(1, 1) && \text{Prior for SD, pooling of outcome} \\ % \Omega_{f} & \sim &\mathcal{N}(0, 1) && \text{Prior for the effects of every feature} \\ % \end{aligned} \]
Where:
- \(Y\): Outcome variable capturing whether a profile has been prioritized (1) or not (0).
- \(o\): Observation
- \(o\): Outcomes (Welfare/ Rights)
- \(p\): Population (General population / Public administration)
- \(t\): Treatment (Shown / Not shown)
- \(F\): Matrix with the observations of Features (the discrimination sources, plus intercept and first shown profile), for each experimental data point.
- \(f\): Feature
- \(\theta_{O,p,t,f}\): Main parameters of interest capturing the discrimination effects by outcome, population, treatment and profiles’ feature.
- \(\Theta_{O,f}\): Hyper-parameters capturing the shared effect of outcomes and features over population and treatment.
- \(\Omega_{f}\): Hyper-parameters capturing the shared effect of features over outcomes, population and treatment.
- \(\sigma_{\theta_{O,f}}\): Between population/treatment and within outcome/feature standard deviations.
- \(\sigma_{\Theta_{f}}\): Between outcome and within feature standard deviations.
The \(\sigma_{\theta_{O,f}}\) can also be seen as a pooling factor. The higher, the more freely are the effects of population and treatment to vary from the overall \(\Theta_{O,f}\) shared by outcome and feature. If it is restricted to be very close to 0 then it assumes that there should be a lot of variation from the two populations and two treatments to make them different. Otherwise, they “borrow strength” from the overall means.
So far we have used a very relaxed pooling (almost none). A model with a restriction very close to zero is the “Cauchy” model.
This model estimates approximately the following number of parameters:
# theta
(nOutcome * nP * nT * nF) +
# Theta
(nOutcome * nF) +
# Omega
(nF) +
# sigma_theta
(nOutcome * nF) +
# sigma_Theta
(nF)-> [1] 112M <- "Decision, Baseline"
M.lab <- "decision-baseline"
# Define the JAGS model (Bugs language) in a character vector
# and write it in a file which name is taken from the labels that we just
# defined the model will be called.
m <- "
model {
for (o in 1:nO) {
Y[o] ~ dbern(p[o])
logit(p[o]) <- inprod(theta[id_outcome[o],id_population[o],id_treatment[o],1:nF], X[o,1:nF])
}
# Priors for effects
for (f in 1:nF) {
for (ocm in 1:nOutcome) {
for (p in 1:nP) {
for (t in 1:nT) {
#theta[ocm,p,t,f] ~ dnorm(Theta[ocm,f], sigma_theta[ocm,f]^-2)
theta[ocm,p,t,f] ~ dnorm(Theta[ocm,f], tau_theta[ocm,f])
}
}
#sigma_theta[ocm,f] ~ dt(0, 0.5^-2, 1)T(0,)
tau_theta[ocm,f] ~ dgamma(1, 1)
sigma_theta[ocm,f] <- 1 / sqrt(tau_theta[ocm,f])
#Theta[ocm,f] ~ dnorm(Omega[f], sigma_Theta[f]^-2)
Theta[ocm,f] ~ dnorm(Omega[f], tau_Theta[f])
}
#sigma_Theta[f] ~ dt(0, 0.5^-2, 1)T(0,)
tau_Theta[f] ~ dgamma(1, 1)
sigma_Theta[f] <- 1 / sqrt(tau_Theta[f])
}
Omega ~ dmnorm(b0, B0)
}
"
write(m, file = paste("models/model-", M.lab, ".bug", sep = ""))
# Define which parameters we want to get from the model.
par <- NULL
par <- c(par, "theta")
par <- c(par, "Theta", "sigma_theta")
par <- c(par, "Omega", "sigma_Theta")
#if (run.pcp) par <- c(par, "p")
par.pcp <- c("p")
# Set the sampling characteristics
adapt <- 5e2
burnin <- 1e3#1e4
run <- 2e3
n.pcp <- 50
chains <- 3
method <- "parallel"
thin <- 5# Run the sampler and save its results
# This chunk only runs when required (eval = FALSE), so that we do not have to
# run the model every time
t0 <- proc.time()
#rj <- run.jags(model = paste("models/model-", M.lab, ".bug", sep = ""),
# data = dump.format(D, checkvalid = FALSE),
# inits = list(list(.RNG.name = "base::Super-Duper", .RNG.seed = 1),
# list(.RNG.name = "base::Super-Duper", .RNG.seed = 2),
# list(.RNG.name = "base::Wichmann-Hill", .RNG.seed = 3)),
# modules = "glm",
# n.chains = chains,
# adapt = adapt,
# burnin = burnin, sample = run,
# thin = thin,
# monitor = par, method = method, summarise = FALSE)
my.jags <- function(seed, model, data, inits = NULL, n.adapt, n.burnin, n.samples, thin) {
require(rjags)
load.module('glm')
load.module('lecuyer')
load.module('dic')
jm <- jags.model(model, n.chains = 1, n.adapt = n.adapt,
inits = list(list(.RNG.name = "base::Super-Duper", .RNG.seed = seed)),
data = data)
update(jm, n.iter = burnin)
s <- coda.samples(jm, variable.names = par, n.iter = run * thin, thin = thin)
s.pcp <- coda.samples(jm, variable.names = par.pcp, n.iter = n.pcp, thin = 1)
return(list(s, s.pcp))
}
s.both <- mclapply(
X = seq_len(chains),
FUN = my.jags,
model = paste("models/model-", M.lab, ".bug", sep = ""),
data = D,
n.adapt = adapt, n.burnin = burnin, n.samples = run, thin = thin,
mc.preschedule = FALSE,
mc.cores = chains)
# Not very elegant, but works now. Must be improved
s <- as.mcmc.list(list(s.both[[1]][[1]][[1]],
s.both[[2]][[1]][[1]],
s.both[[3]][[1]][[1]]))
s.p <- as.mcmc.list(list(s.both[[1]][[2]][[1]],
s.both[[2]][[2]][[1]],
s.both[[3]][[2]][[1]]))
save(s, file = paste("samples-", M.lab, ".RData", sep = ""))
save(s.p, file = paste("samples-p-", M.lab, ".RData", sep = ""))
proc.time() - t0load(file = paste("samples-", M.lab, ".RData", sep = ""))
load(file = paste("samples-p-", M.lab, ".RData", sep = ""))ggmcmc(ggs(s, family = "^theta|^Theta|^Omega"),
file = paste0("ggmcmc-", M.lab, ".pdf"),
param_page = 8,
plot = c("traceplot", "geweke", "Rhat"))Part-worth contributions
Figures 1 and 4 show the discrimination on all the characteristics. The parameters are the result of a logistic transformation, and therefore must be interpreted in terms of the odds (by exponentiating them) or predicted probabilities. They represent the relative importances of each feature (except for age, that the relative importances must be roughly divided by 2).
The caterpillar plot shows the distribution of the values for the parameters, where the dot represents the median expected effect, and the band covers 95 percent of the credible interval (also called the highest posterior density, HPD). In Bayesian statistics, however, there is not such thing as “statistically significant” for a parameter, as a clear cut indication, but results are presented in terms of probabilities of interest (such as the probability that the parameter is positive). See a Bayesian interpretation below, in the “Differential discrimination” output.
Comparing general population against public administration allows us to see the effect of having a smaller sample in the later group, by the sizes of the uncertainty bands. At this stage, I would avoid making any sort of interpretation of the public administration results (except, maybe, the strong deviation of the parameter regarding language effects in rights when treatment is not shown, that is remarkably far away from the general population).
For both outcomes, language is the most important feature to discriminate. In terms of magnitude, though, it is more important in rights than in welfare, but specially when the treatment is shown. For welfare, there are also signs of profession being relevant for discrimination, which is not at all suggested in rights.
# Use ggmcmc's features to get all the parameters in a tidy format
# Save some of them for later, to make model comparisons.
L.theta <- plab("theta", list(Outcome = outcome.label,
Population = population.label,
Treatment = treatment.label,
Feature = feature.label))
S.theta <- ggs(s, family = "^theta\\[", par_labels = L.theta) %>%
mutate(Model = M)
save(S.theta, file = paste("samples-", M.lab, "-theta", ".RData", sep = ""))
ci.theta <- ci(S.theta)Report baseline probabilities of picking the reference individual.
tb <- ci.theta %>%
filter(Feature == "(Intercept)") %>%
arrange(desc(abs(median))) %>%
mutate(`Odds Ratio` = exp(median)) %>%
mutate(`Expected odds ratio` = or(`Odds Ratio`)) %>%
mutate(`Expected probability` = inv.logit(median)) %>%
select(Outcome, Treatment, Population, `Expected odds ratio`, `Expected probability`)
tc <- "Expected odds ratio and probability of prioritization for the reference profile."
if (knitr::is_latex_output()) {
kable(tb, format = "latex", caption = tc, longtable = TRUE, booktabs = TRUE) %>%
kable_styling(font_size = 10)
} else {
kable(tb, format = "html", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10, position = "center", bootstrap_options = "striped", full_width = T)
}| Outcome | Treatment | Population | Expected odds ratio | Expected probability |
|---|---|---|---|---|
| Rights | Not shown | Public administration | ▽ 51% | 0.3288 |
| Welfare | Not shown | Public administration | ▽ 49% | 0.3364 |
| Rights | Shown | Public administration | ▽ 48% | 0.3441 |
| Welfare | Shown | Public administration | ▽ 46% | 0.3515 |
| Rights | Not shown | General population | ▽ 46% | 0.3519 |
| Rights | Shown | General population | ▽ 44% | 0.3589 |
| Welfare | Not shown | General population | ▽ 42% | 0.3670 |
| Welfare | Shown | General population | ▽ 41% | 0.3698 |
tb <- ci.theta %>%
filter(Feature == "(Intercept)") %>%
arrange(rev(Population), Outcome, Treatment) %>%
mutate(`Odds Ratio` = exp(median)) %>%
mutate(`Expected odds ratio` = or(`Odds Ratio`)) %>%
mutate(`Expected probability` = inv.logit(median)) %>%
select(Population, Outcome, Treatment, `Expected odds ratio`, `Expected probability`)
tc <- "Expected odds ratio and probability of prioritization for the reference profile."
if (knitr::is_latex_output()) {
kable(tb, format = "latex", caption = tc, longtable = TRUE, booktabs = TRUE) %>%
kable_styling(font_size = 10) %>%
pack_rows("Public administration", 5, 8) %>%
pack_rows("General population", 5, 8)
} else {
kable(tb, format = "html", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10, position = "center", bootstrap_options = "striped", full_width = T) %>%
pack_rows("Public administration", 1, 4) %>%
pack_rows("General population", 5, 8)
}| Population | Outcome | Treatment | Expected odds ratio | Expected probability |
|---|---|---|---|---|
| Public administration | ||||
| Public administration | Rights | Not shown | ▽ 51% | 0.3288 |
| Public administration | Rights | Shown | ▽ 48% | 0.3441 |
| Public administration | Welfare | Not shown | ▽ 49% | 0.3364 |
| Public administration | Welfare | Shown | ▽ 46% | 0.3515 |
| General population | ||||
| General population | Rights | Not shown | ▽ 46% | 0.3519 |
| General population | Rights | Shown | ▽ 44% | 0.3589 |
| General population | Welfare | Not shown | ▽ 42% | 0.3670 |
| General population | Welfare | Shown | ▽ 41% | 0.3698 |
ggplot(ci.theta, aes(ymin = low, ymax = high,
y = median, x = Feature,
color = Outcome)) +
coord_flip() +
geom_point(position = position_dodge(width = 0.4)) +
geom_linerange(position = position_dodge(width = 0.4)) +
geom_linerange(aes(ymin = Low, ymax = High), size = 1, position = position_dodge(width = 0.4)) +
geom_hline(aes(yintercept = 0), lty = 3) +
xlab("Parameter") + ylab("HPD") +
facet_grid(Treatment ~ Population) +
scale_color_discrete_qualitative(palette = "Dark 2")
Figure 1: Prioritization of features, by population and treatment.
ggplot(ci.theta, aes(ymin = low, ymax = high,
y = median, x = Feature,
color = Treatment)) +
coord_flip() +
geom_point(position = position_dodge(width = 0.4)) +
geom_linerange(position = position_dodge(width = 0.4)) +
geom_linerange(aes(ymin = Low, ymax = High), size = 1, position = position_dodge(width = 0.4)) +
geom_hline(aes(yintercept = 0), lty = 3) +
xlab("Parameter") + ylab("HPD") +
facet_grid(Outcome ~ Population) +
scale_color_discrete_qualitative(palette = "Harmonic")
Figure 2: Prioritization of features, by population and outcome.
ci.theta %>%
filter(!(Feature %in% c("(Intercept)", "First shown"))) %>%
droplevels() %>%
mutate(Feature = factor(Feature, rev(feature.label.order))) %>%
mutate(Outcome = factor(Outcome, rev(levels(Outcome)))) %>%
mutate(low = exp(low), high = exp(high),
median = exp(median),
Low = exp(Low), High = exp(High)) %>%
ggplot(aes(ymin = low, ymax = high,
y = median, x = Feature,
color = Treatment)) +
coord_flip() +
geom_point(position = position_dodge(width = 0.4)) +
geom_linerange(position = position_dodge(width = 0.4)) +
geom_linerange(aes(ymin = Low, ymax = High), size = 1, position = position_dodge(width = 0.4)) +
geom_hline(aes(yintercept = 1), lty = 3) +
xlab("Parameter") +
ylab("Odds Ratio (HPD)") +
# ylab("HPD") +
facet_grid(Outcome ~ Population) +
scale_color_grey(start = 0.6, end = 0.2) +
theme(strip.text.x = element_text(face = "bold"),
strip.text.y = element_text(face = "bold")) +
scale_y_continuous(trans = log_trans(), breaks = c(0.7, 1, 1.5, 2)) +
theme(panel.spacing = unit(2, "lines"))
Figure 3: Prioritization of features, by population and outcome.
ggplot(ci.theta, aes(ymin = low, ymax = high,
y = median, x = Feature,
color = Population)) +
coord_flip() +
geom_point(position = position_dodge(width = 0.4)) +
geom_linerange(position = position_dodge(width = 0.4)) +
geom_linerange(aes(ymin = Low, ymax = High), size = 1, position = position_dodge(width = 0.4)) +
geom_hline(aes(yintercept = 0), lty = 3) +
xlab("Parameter") + ylab("HPD") +
facet_grid(Treatment ~ Outcome) +
scale_color_discrete_qualitative(palette = "Dynamic")
Figure 4: Prioritization of features, by treatment and outcome.
ci.theta %>%
filter(!(Feature %in% c("(Intercept)", "First shown"))) %>%
droplevels() %>%
mutate(Feature = factor(Feature, rev(feature.label.order))) %>%
mutate(Outcome = factor(Outcome, rev(levels(Outcome)))) %>%
mutate(low = exp(low), high = exp(high),
median = exp(median),
Low = exp(Low), High = exp(High)) %>%
ggplot(aes(ymin = low, ymax = high,
y = median, x = Feature,
color = Population)) +
coord_flip() +
geom_point(position = position_dodge(width = 0.4)) +
geom_linerange(position = position_dodge(width = 0.4)) +
geom_linerange(aes(ymin = Low, ymax = High), size = 1, position = position_dodge(width = 0.4)) +
geom_hline(aes(yintercept = 1), lty = 3) +
xlab("Parameter") +
ylab("Odds Ratio (HPD)") +
# ylab("HPD") +
facet_grid(Treatment ~ Outcome) +
scale_color_grey(start = 0.8, end = 0.2) +
scale_y_continuous(trans = log_trans(), breaks = c(0.7, 1, 1.5, 2)) +
theme(strip.text.x = element_text(face = "bold"),
strip.text.y = element_text(face = "bold")) +
theme(panel.spacing = unit(2, "lines"))
Figure 5: Prioritization of features, by treatment and outcome.
L.Theta <- plab("Theta", list(Outcome = outcome.label,
Feature = feature.label))
S.Theta <- ggs(s, family = "^Theta\\[", par_labels = L.Theta)
ggplot(ci(S.Theta), aes(ymin = low, ymax = high,
y = median, x = Feature,
color = Outcome)) +
coord_flip() +
geom_point(position = position_dodge(width = 0.3)) +
geom_linerange(position = position_dodge(width = 0.3)) +
geom_linerange(aes(ymin = Low, ymax = High), size = 1, position = position_dodge(width = 0.3)) +
geom_hline(aes(yintercept = 0), lty = 3) +
xlab("Parameter") + ylab("HPD") +
scale_color_discrete_qualitative(palette = "Dark 2")Figure 6: Prioritization of features (grand mean).
L.Omega <- plab("Omega", list(Feature = feature.label))
S.Omega <- ggs(s, family = "^Omega\\[", par_labels = L.Omega)
ggs_caterpillar(S.Omega) +
geom_vline(xintercept = 0, lty = 3)Figure 7: Prioritization of features (grand-grand mean).
List median expected effects.
tb <- S.Omega %>%
filter(Feature != "(Intercept)") %>%
ci() %>%
arrange(desc(abs(median))) %>%
mutate(`Odds Ratio` = exp(median)) %>%
mutate(`Expected effect` = or(`Odds Ratio`)) %>%
select(Feature, `Odds Ratio`, `Expected effect`)
tc <- "Odds ratio of expected grand-effects, sorted by magnitude."
if (knitr::is_latex_output()) {
kable(tb, format = "latex", caption = tc, longtable = TRUE, booktabs = TRUE) %>%
kable_styling(font_size = 10)
} else {
kable(tb, format = "html", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10, position = "center", bootstrap_options = "striped", full_width = T)
}| Feature | Odds Ratio | Expected effect |
|---|---|---|
| Language Full | 2.1741 | △ 120% |
| Nationality Netherlands | 1.4403 | △ 44% |
| Gender Male | 0.7355 | ▽ 26% |
| First shown | 1.3551 | △ 36% |
| Age 55 | 0.8653 | ▽ 13% |
| Profession Nurse | 1.1448 | △ 14% |
| Age 25 | 0.9244 | ▽ 7.6% |
On average, overall, individuals with full German are 2.17 times more likely to be prioritized. Dutch citizens are 44 percent more likely to be prioritized. Nurses are 14 percent more likely. Compared to individuals aged 40, younger and older ones are 7.6 and 14 percent less likely. Men are 26 percent less likely to be prioritized.
Attribute importance
Figures 8 and 9 show the relative attribute importance. Language accounts for approximately half the importance that respondents give to the prioritization (more in rights than in welfare). In public administration population the relative importances are more distributed across attributes. But at the same time, the relative importances between profession and nationality change more radically: they base their decision more on nationality when rights are involved, and profession when welfare is involved.
importance <- function(x) {
if (length(x) == 1) {
return(abs(x))
} else {
return(max(x) - min(x))
}
}
S.importance <- S.theta %>%
filter(!Feature %in% c("First shown", "(Intercept)")) %>%
ci() %>%
select(Outcome, Population, Treatment, Feature, value = median) %>%
separate(Feature, " ", into = c("Attribute", "Level")) %>%
group_by(Outcome, Population, Treatment, Attribute) %>%
summarize(Importance = importance(value))
S.relative.importance <- S.importance %>%
group_by(Outcome, Population, Treatment) %>%
mutate(`Relative importance` = Importance / sum(Importance)) %>%
mutate(Attribute = factor(Attribute, rev(c("Nationality",
"Language",
"Profession",
"Age",
"Gender"))))
save(S.relative.importance, file = paste("relative-importance-", M.lab, "-theta", ".RData", sep = ""))ggplot(S.relative.importance,
aes(x = Outcome, y = `Relative importance`, fill = Attribute)) +
geom_bar(position = "fill", stat = "identity") +
facet_grid(Population ~ Treatment) +
scale_fill_brewer(type = "qual", palette = "Set1")Figure 8: Relative importances, by outcome, population and treatment. Focus on comparing outcomes.
ggplot(S.relative.importance,
aes(x = Treatment, y = `Relative importance`, fill = Attribute)) +
geom_bar(position = "fill", stat = "identity") +
facet_grid(Population ~ Outcome) +
scale_fill_brewer(type = "qual", palette = "Set1")Figure 9: Relative importances, by outcome, population and treatment. Focus on comparing treatment.
S.relative.importance %>%
ungroup() %>%
mutate(Outcome = factor(Outcome, rev(levels(Outcome)))) %>%
ggplot(aes(x = Treatment, y = `Relative importance`, fill = Attribute)) +
geom_bar(position = "fill", stat = "identity") +
facet_grid(Population ~ Outcome) +
theme(strip.text.x = element_text(face = "bold"),
strip.text.y = element_text(face = "bold")) +
scale_fill_grey(start = 0.8, end = 0.2)Figure 10: Relative importances, by outcome, population and treatment. Focus on comparing treatment.
ggplot(S.relative.importance,
aes(x = Population, y = `Relative importance`, fill = Attribute)) +
geom_bar(position = "fill", stat = "identity") +
facet_grid(Outcome ~ Treatment) +
scale_fill_brewer(type = "qual", palette = "Set1")Figure 11: Relative importances, by population, treatment and outcome. Focus on comparing population.
S.relative.importance %>%
ungroup() %>%
mutate(Outcome = factor(Outcome, rev(levels(Outcome)))) %>%
mutate(Population = str_replace(Population, " ", "\n")) %>%
ggplot(aes(x = Population, y = `Relative importance`, fill = Attribute)) +
geom_bar(position = "fill", stat = "identity") +
facet_grid(Outcome ~ Treatment) +
theme(strip.text.x = element_text(face = "bold"),
strip.text.y = element_text(face = "bold")) +
scale_fill_grey(start = 0.8, end = 0.2)Figure 12: Relative importances, by population, treatment and outcome. Focus on comparing population.
Differential discrimination
Figure 13 shows the differential discrimination that welfare outcomes have over rights. Histograms show the distribution of the values of the differential discrimination in each simulation. Therefore, a histogram that clearly overlaps 0.5 indicates that there is no trace of differential discrimination. When the histogram barely overlaps zero it indicates evidence of a differential discrimination, either welfare effects being greater than rights (right to the zero) or rights effects being greater than welfare (left to the zero).
For the differential effect on discrimination based on Nationality, it can be read as follows: evidence that there is a differential effect of higher discrimination in welfare than in rights when the treatment is shown is 100 percent (also when the treatment is not shown). More specifically, the evidence is > 0.99983 (one in 6,000).
# With every iteration and every chaing a value of interest has been
# generated, a sample from the target distribution of interest is retrieved.
# Now simply pivot the cleaned version so that two columns with the different
# values of Outcome are generated, and can be compared to assess the
# differential discrimination.
dd.outcome <- S.theta %>%
select(Iteration, Chain, Outcome, Population, Treatment, Feature, value) %>%
pivot_wider(names_from = Outcome, values_from = value) %>%
mutate(`Differential discrimination on welfare` = Rights - Welfare)
# Prepare the table
tb <- dd.outcome %>%
group_by(Population, Treatment, Feature) %>%
summarize(`Prob Welfare > Rights` = length(which(`Differential discrimination on welfare` > 0)) / n()) %>%
arrange(desc(`Prob Welfare > Rights`))
# Improve the readability of the numbers that will appear in the figures
fig.labels <- tb %>%
mutate(P = paste0(100 * signif(`Prob Welfare > Rights`, 2), "%")) %>%
mutate(P = str_replace(P, "^0\\.", "\\."))
ggplot(filter(dd.outcome, Population == "General population"),
aes(x = `Differential discrimination on welfare`)) +
geom_histogram(binwidth = 0.01) +
facet_grid(Treatment ~ Feature) +
geom_text(data = filter(fig.labels, Population == "General population"),
aes(x = -Inf, y = Inf, label = P), size = rel(0.8),
hjust = -0.5, vjust = +1.5) +
geom_vline(xintercept = 0, lty = 3) +
xlim(c(-0.75, 0.75)) +
ggtitle("Differential discrimination on welfare (General population)")
Figure 13: Differential discrimination on welfare over rights. Only general population is shown. Values indicate the probability of having found a differential discrimination favouring welfare.
tb <- dd.outcome %>%
group_by(Population, Treatment, Feature) %>%
summarize(`Prob Welfare > Rights` = length(which(`Differential discrimination on welfare` > 0)) / n()) %>%
arrange(desc(`Prob Welfare > Rights`))
tc <- "Probability of finding a differential discrimination on welfare over rights."
if (knitr::is_latex_output()) {
kable(tb, format = "latex", caption = tc, longtable = TRUE, booktabs = TRUE) %>%
kable_styling(font_size = 8)
} else {
kable(tb, format = "html", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 8, position = "center", bootstrap_options = "striped", full_width = T)
}| Population | Treatment | Feature | Prob Welfare > Rights |
|---|---|---|---|
| General population | Not shown | Nationality Netherlands | 1.0000 |
| General population | Not shown | Language Full | 1.0000 |
| General population | Shown | Language Full | 1.0000 |
| Public administration | Not shown | Language Full | 1.0000 |
| Public administration | Shown | Nationality Netherlands | 0.9995 |
| General population | Shown | Nationality Netherlands | 0.9987 |
| Public administration | Not shown | Nationality Netherlands | 0.9982 |
| Public administration | Shown | Language Full | 0.9942 |
| General population | Not shown | First shown | 0.9828 |
| Public administration | Not shown | Age 25 | 0.9482 |
| Public administration | Shown | Age 25 | 0.9108 |
| General population | Not shown | Age 25 | 0.9097 |
| General population | Shown | Gender Male | 0.8113 |
| General population | Not shown | Gender Male | 0.8012 |
| General population | Shown | Age 25 | 0.7692 |
| Public administration | Shown | First shown | 0.6980 |
| Public administration | Shown | Gender Male | 0.6838 |
| Public administration | Not shown | Gender Male | 0.6595 |
| Public administration | Not shown | First shown | 0.4970 |
| Public administration | Shown | (Intercept) | 0.3920 |
| Public administration | Not shown | (Intercept) | 0.3902 |
| Public administration | Shown | Age 55 | 0.2735 |
| General population | Shown | First shown | 0.2265 |
| General population | Shown | (Intercept) | 0.2183 |
| General population | Not shown | (Intercept) | 0.1432 |
| Public administration | Not shown | Age 55 | 0.0068 |
| General population | Shown | Age 55 | 0.0002 |
| General population | Not shown | Profession Nurse | 0.0000 |
| General population | Not shown | Age 55 | 0.0000 |
| General population | Shown | Profession Nurse | 0.0000 |
| Public administration | Not shown | Profession Nurse | 0.0000 |
| Public administration | Shown | Profession Nurse | 0.0000 |
Figure 14 shows the differential discrimination that treated individuals (shown) have over not treated (not shown).
dd.treatment <- S.theta %>%
select(Iteration, Chain, Outcome, Population, Treatment, Feature, value) %>%
pivot_wider(names_from = Treatment, values_from = value) %>%
mutate(`Differential discrimination on treatment shown` = `Not shown` - Shown)
tb <- dd.treatment %>%
group_by(Population, Outcome, Feature) %>%
summarize(`Prob Shown > Not shown` = length(which(`Differential discrimination on treatment shown` > 0)) / n()) %>%
arrange(desc(`Prob Shown > Not shown`))
fig.labels <- tb %>%
mutate(P = paste0(100 * signif(`Prob Shown > Not shown`, 2), "%")) %>%
mutate(P = str_replace(P, "^0\\.", "\\."))
ggplot(filter(dd.treatment, Population == "General population"),
aes(x = `Differential discrimination on treatment shown`)) +
geom_histogram(binwidth = 0.01) +
facet_grid(Outcome ~ Feature) +
geom_text(data = filter(fig.labels, Population == "General population"),
aes(x = -Inf, y = Inf, label = P), size = rel(0.8),
hjust = -0.5, vjust = +1.5) +
geom_vline(xintercept = 0, lty = 3) +
xlim(c(-0.75, 0.75)) +
ggtitle("Differential discrimination on treatment shown")
Figure 14: Differential discrimination on treated individuals over non treated. Only general population is shown. Values indicate the probability of having found a differential discrimination favouring treated individuals.
tb <- dd.treatment %>%
group_by(Population, Outcome, Feature) %>%
summarize(`Prob Shown > Not shown` = length(which(`Differential discrimination on treatment shown` > 0)) / n()) %>%
arrange(desc(`Prob Shown > Not shown`))
tc <- "Probability of finding a differential discrimination on treatment shown over non shown."
if (knitr::is_latex_output()) {
kable(tb, format = "latex", caption = tc, longtable = TRUE, booktabs = TRUE) %>%
kable_styling(font_size = 8)
} else {
kable(tb, format = "html", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 8, position = "center", bootstrap_options = "striped", full_width = T)
}| Population | Outcome | Feature | Prob Shown > Not shown |
|---|---|---|---|
| General population | Welfare | Language Full | 0.9995 |
| Public administration | Welfare | Age 55 | 0.9640 |
| General population | Rights | Nationality Netherlands | 0.9253 |
| General population | Rights | First shown | 0.8888 |
| Public administration | Rights | Language Full | 0.7623 |
| Public administration | Welfare | First shown | 0.7493 |
| General population | Welfare | Age 55 | 0.6408 |
| Public administration | Rights | Profession Nurse | 0.6375 |
| Public administration | Welfare | Nationality Netherlands | 0.6228 |
| Public administration | Welfare | Gender Male | 0.6215 |
| Public administration | Rights | Gender Male | 0.5838 |
| General population | Rights | Age 25 | 0.5788 |
| Public administration | Rights | Age 55 | 0.5408 |
| Public administration | Rights | First shown | 0.5362 |
| Public administration | Rights | Nationality Netherlands | 0.5060 |
| General population | Rights | Age 55 | 0.5030 |
| Public administration | Rights | Age 25 | 0.4575 |
| General population | Welfare | Gender Male | 0.4557 |
| General population | Rights | Gender Male | 0.4418 |
| Public administration | Welfare | Language Full | 0.4300 |
| General population | Welfare | (Intercept) | 0.4128 |
| Public administration | Welfare | Age 25 | 0.3642 |
| General population | Welfare | Age 25 | 0.3405 |
| General population | Welfare | Profession Nurse | 0.3265 |
| General population | Rights | (Intercept) | 0.3133 |
| Public administration | Welfare | (Intercept) | 0.2925 |
| Public administration | Rights | (Intercept) | 0.2878 |
| General population | Rights | Language Full | 0.2827 |
| General population | Welfare | Nationality Netherlands | 0.2688 |
| Public administration | Welfare | Profession Nurse | 0.1748 |
| General population | Rights | Profession Nurse | 0.1167 |
| General population | Welfare | First shown | 0.0430 |
Figure 15 shows the differential discrimination that public administration have over general population.
dd.population <- S.theta %>%
select(Iteration, Chain, Outcome, Population, Treatment, Feature, value) %>%
pivot_wider(names_from = Population, values_from = value) %>%
mutate(`Differential discrimination on public administration` = `Public administration` - `General population`)
tb <- dd.population %>%
group_by(Outcome, Treatment, Feature) %>%
summarize(`Prob PA > General` = length(which(`Differential discrimination on public administration` > 0)) / n()) %>%
arrange(desc(`Prob PA > General`))
fig.labels <- tb %>%
mutate(P = paste0(100 * signif(`Prob PA > General`, 2), "%")) %>%
mutate(P = str_replace(P, "^0\\.", "\\."))
ggplot(filter(dd.population, Treatment == "Shown"),
aes(x = `Differential discrimination on public administration`)) +
geom_histogram(binwidth = 0.01) +
facet_grid(Outcome ~ Feature) +
geom_text(data = filter(fig.labels, Treatment == "Shown"),
aes(x = -Inf, y = Inf, label = P), size = rel(0.8),
hjust = -0.5, vjust = +1.5) +
geom_vline(xintercept = 0, lty = 3) +
xlim(c(-0.75, 0.75)) +
ggtitle("Differential discrimination on public administration")
Figure 15: Differential discrimination on public administration over general population. Only treatment is shown. Values indicate the probability of having found a differential discrimination favouring public administration.
tb <- dd.population %>%
group_by(Outcome, Treatment, Feature) %>%
summarize(`Prob PA > General` = length(which(`Differential discrimination on public administration` > 0)) / n()) %>%
arrange(desc(`Prob PA > General`))
tc <- "Probability of finding a differential discrimination on public administration over general population."
if (knitr::is_latex_output()) {
kable(tb, format = "latex", caption = tc, longtable = TRUE, booktabs = TRUE) %>%
kable_styling(font_size = 8)
} else {
kable(tb, format = "html", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 8, position = "center", bootstrap_options = "striped", full_width = T)
}| Outcome | Treatment | Feature | Prob PA > General |
|---|---|---|---|
| Welfare | Shown | Profession Nurse | 0.9920 |
| Welfare | Shown | Language Full | 0.9867 |
| Welfare | Not shown | First shown | 0.9655 |
| Welfare | Not shown | Profession Nurse | 0.9295 |
| Rights | Shown | Age 25 | 0.9263 |
| Welfare | Not shown | Age 55 | 0.9140 |
| Rights | Not shown | Profession Nurse | 0.9110 |
| Rights | Shown | Nationality Netherlands | 0.9040 |
| Rights | Not shown | Age 25 | 0.8860 |
| Rights | Not shown | Age 55 | 0.8645 |
| Rights | Shown | First shown | 0.8513 |
| Rights | Shown | Age 55 | 0.8262 |
| Rights | Not shown | First shown | 0.6592 |
| Rights | Not shown | Nationality Netherlands | 0.6588 |
| Welfare | Shown | Age 25 | 0.5808 |
| Welfare | Not shown | Nationality Netherlands | 0.5700 |
| Rights | Shown | Profession Nurse | 0.5510 |
| Welfare | Not shown | Age 25 | 0.5157 |
| Welfare | Not shown | Language Full | 0.4967 |
| Welfare | Shown | First shown | 0.4640 |
| Rights | Not shown | Language Full | 0.3567 |
| Welfare | Shown | Nationality Netherlands | 0.2820 |
| Welfare | Not shown | Gender Male | 0.2618 |
| Welfare | Shown | Age 55 | 0.2558 |
| Rights | Shown | (Intercept) | 0.2548 |
| Rights | Not shown | Gender Male | 0.2493 |
| Welfare | Shown | (Intercept) | 0.1993 |
| Rights | Shown | Gender Male | 0.1545 |
| Rights | Not shown | (Intercept) | 0.1428 |
| Welfare | Shown | Gender Male | 0.1347 |
| Welfare | Not shown | (Intercept) | 0.0790 |
| Rights | Shown | Language Full | 0.0530 |
Model fit
Model fit is assessed using the percent of cases correctly predicted (PCP). On average, we would find that simply by chance we can predict 0.5 of the prioritizations. The PCP is the percentage of cases correctly predicted (those predicted as being prioritized indeed being prioritized plus those predicted not being prioritized indeed not being prioritized, divided by all cases).
Figure 16 shows the posterior distribution of PCP, and the following table its median expected value.
threshold <- d %>%
group_by(Outcome, Population, Treatment) %>%
summarize(Threshold = length(which(Decision == 1)) / n())
L.p <- plab("p", list(Observation = 1:nO)) %>%
mutate(Observation = as.numeric(as.character(Observation)))
S.full <- ggs(s.p, family = "^p\\[", par_labels = L.p, sort = FALSE) %>%
select(Iteration, Chain, Observation, value) %>%
left_join(select(d, Observation, Population, Treatment, Outcome, Decision)) %>%
left_join(threshold) %>%
mutate(Correct = if_else( (value < Threshold & Decision == 0) |
(value > Threshold & Decision == 1),
TRUE, FALSE))
S <- S.full %>%
group_by(Outcome, Population, Treatment, Iteration, Chain) %>%
summarize(PCP = length(which(Correct)) / n())
ggplot(S, aes(x = PCP, color = Population)) +
geom_density() +
facet_grid(Treatment ~ Outcome) +
expand_limits(x = c(0, 1)) +
geom_vline(xintercept = 0, lty = 3) +
scale_color_discrete_qualitative(palette = "Dynamic")Figure 16: Model fit: Percent correctly predicted.
t.pcp <- S %>%
group_by(Population, Treatment, Outcome) %>%
summarize(`Median PCP` = mean(PCP))
t.pcp.m <- mutate(t.pcp, Model = M)
save(t.pcp.m, file = paste0("table-pcp-", M.lab, ".RData"))
tc <- "Posterior median percent correctly predicted, by outcome, treatment and
population."
if (knitr::is_latex_output()) {
kable(t.pcp, format = "latex", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10)
} else {
kable(t.pcp, format = "html", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10, position = "center", bootstrap_options = "striped", full_width = F)
}| Population | Treatment | Outcome | Median PCP |
|---|---|---|---|
| General population | Not shown | Rights | 0.6258 |
| General population | Not shown | Welfare | 0.5938 |
| General population | Shown | Rights | 0.6248 |
| General population | Shown | Welfare | 0.5849 |
| Public administration | Not shown | Rights | 0.6264 |
| Public administration | Not shown | Welfare | 0.5989 |
| Public administration | Shown | Rights | 0.6155 |
| Public administration | Shown | Welfare | 0.6010 |
List of the 10 worstly predicted cases, on average.
wpo <- S.full %>%
group_by(Observation) %>%
summarize(Distance = mean(Decision - value)) %>%
arrange(desc(abs(Distance))) %>%
slice(1:10) %>%
select(Observation) %>%
unlist(use.names = FALSE)
rm(S.full)
invisible(gc()) # Save memory
E[wpo,] %>%
select(id, Population, Treatment, Outcome,
Nationality, Gender, Language, Profession, Age,
Decision) %>%
arrange(id) %>%
data.frame()-> id Population Treatment Outcome Nationality Gender Language
-> 1 10 General population Not shown Rights Romania Male No
-> 2 477 General population Not shown Rights Romania Male No
-> 3 572 General population Not shown Rights Romania Male No
-> 4 582 General population Not shown Rights Romania Male No
-> 5 658 General population Not shown Rights Romania Male No
-> 6 850 General population Not shown Rights Romania Male No
-> 7 862 General population Not shown Rights Romania Male No
-> 8 1144 General population Not shown Rights Romania Male No
-> 9 1182 General population Not shown Rights Romania Male No
-> 10 1248 General population Not shown Rights Romania Male No
-> Profession Age Decision
-> 1 Nurse 55 1
-> 2 Nurse 55 1
-> 3 Nurse 55 1
-> 4 Nurse 55 1
-> 5 Nurse 55 1
-> 6 Nurse 55 1
-> 7 Nurse 55 1
-> 8 Nurse 55 1
-> 9 Nurse 55 1
-> 10 Nurse 55 1Model fit contrasting respondent quality
load("poor_quality_individual.RData")
threshold <- d %>%
group_by(Outcome, Population, Treatment) %>%
summarize(Threshold = length(which(Decision == 1)) / n())
L.p <- plab("p", list(Observation = 1:nO)) %>%
mutate(Observation = as.numeric(as.character(Observation)))
S.full <- ggs(s.p, family = "^p\\[", par_labels = L.p, sort = FALSE) %>%
select(Iteration, Chain, Parameter, Observation, value) %>%
left_join(select(d, id, Observation, Population, Treatment, Outcome, Decision)) %>%
left_join(threshold) %>%
mutate(Correct = if_else( (value < Threshold & Decision == 0) |
(value > Threshold & Decision == 1),
TRUE, FALSE))
S <- S.full %>%
mutate(Respondent = ifelse(id %in% id.poor, "Low quality", "Regular")) %>%
group_by(Outcome, Population, Treatment, Respondent, Iteration, Chain) %>%
summarize(PCP = length(which(Correct)) / n())
ggplot(S, aes(x = PCP, y = Treatment,
color = Respondent)) +
geom_boxplot(position = position_dodge(width = 0.3)) +
facet_grid(Population ~ Outcome) +
scale_color_discrete_qualitative(palette = "Dynamic")Figure 17: Model fit: Percent correctly predicted.
t.pcp <- S %>%
group_by(Population, Treatment, Outcome, Respondent) %>%
summarize(`Median PCP` = mean(PCP))
t.pcp.m.respondent <- mutate(t.pcp, Model = M)
save(t.pcp.m.respondent, file = paste0("table-pcp-respondent-", M.lab, ".RData"))
tc <- "Posterior median percent correctly predicted, by outcome, treatment, population and respondent quality."
if (knitr::is_latex_output()) {
kable(t.pcp, format = "latex", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10)
} else {
kable(t.pcp, format = "html", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10, position = "center", bootstrap_options = "striped", full_width = F)
}| Population | Treatment | Outcome | Respondent | Median PCP |
|---|---|---|---|---|
| General population | Not shown | Rights | Low quality | 0.5886 |
| General population | Not shown | Rights | Regular | 0.6284 |
| General population | Not shown | Welfare | Low quality | 0.5515 |
| General population | Not shown | Welfare | Regular | 0.5969 |
| General population | Shown | Rights | Low quality | 0.5683 |
| General population | Shown | Rights | Regular | 0.6292 |
| General population | Shown | Welfare | Low quality | 0.6159 |
| General population | Shown | Welfare | Regular | 0.5826 |
| Public administration | Not shown | Rights | Low quality | 0.5938 |
| Public administration | Not shown | Rights | Regular | 0.6282 |
| Public administration | Not shown | Welfare | Low quality | 0.6623 |
| Public administration | Not shown | Welfare | Regular | 0.5954 |
| Public administration | Shown | Rights | Low quality | 0.5893 |
| Public administration | Shown | Rights | Regular | 0.6179 |
| Public administration | Shown | Welfare | Low quality | 0.5956 |
| Public administration | Shown | Welfare | Regular | 0.6015 |
Hypothesis 1: Dutch are treated more favourably
tb <- S.theta %>%
filter(Feature == "Nationality Netherlands") %>%
group_by(Outcome, Population, Treatment) %>%
summarize(`Prob (H1)` = length(which(value > 0)) / n()) %>%
arrange(desc(`Prob (H1)`))
tc <- "Evidence for Hypothesis 1."
if (knitr::is_latex_output()) {
kable(tb, format = "latex", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10)
} else {
kable(tb, format = "html", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10, position = "center", bootstrap_options = "striped", full_width = F)
}| Outcome | Population | Treatment | Prob (H1) |
|---|---|---|---|
| Rights | General population | Not shown | 1 |
| Rights | General population | Shown | 1 |
| Rights | Public administration | Not shown | 1 |
| Rights | Public administration | Shown | 1 |
| Welfare | General population | Not shown | 1 |
| Welfare | General population | Shown | 1 |
| Welfare | Public administration | Not shown | 1 |
| Welfare | Public administration | Shown | 1 |
Hypothesis 2: Fluent german speakers are treated more favourably
tb <- S.theta %>%
filter(Feature == "Language Full") %>%
group_by(Outcome, Population, Treatment) %>%
summarize(`Prob (H2)` = length(which(value > 0)) / n()) %>%
arrange(desc(`Prob (H2)`))
tc <- "Evidence for Hypothesis 2."
if (knitr::is_latex_output()) {
kable(tb, format = "latex", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10)
} else {
kable(tb, format = "html", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10, position = "center", bootstrap_options = "striped", full_width = F)
}| Outcome | Population | Treatment | Prob (H2) |
|---|---|---|---|
| Rights | General population | Not shown | 1 |
| Rights | General population | Shown | 1 |
| Rights | Public administration | Not shown | 1 |
| Rights | Public administration | Shown | 1 |
| Welfare | General population | Not shown | 1 |
| Welfare | General population | Shown | 1 |
| Welfare | Public administration | Not shown | 1 |
| Welfare | Public administration | Shown | 1 |
Hypothesis 3: Bureaucrats discriminate less
This hypothesis only considers absolute discrimination (absolute values), not the direction).
tb <- S.theta %>%
pivot_wider(c(Iteration, Chain, Outcome, Treatment, Feature),
names_from = Population, values_from = value) %>%
group_by(Outcome, Treatment, Feature) %>%
summarize(`Prob (H3)` = length(which(abs(`Public administration`) < abs(`General population`))) / n()) %>%
arrange(desc(`Prob (H3)`))
tc <- "Evidence for Hypothesis 3."
if (knitr::is_latex_output()) {
kable(tb, format = "latex", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10)
} else {
kable(tb, format = "html", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10, position = "center", bootstrap_options = "striped", full_width = F)
}| Outcome | Treatment | Feature | Prob (H3) |
|---|---|---|---|
| Rights | Shown | Language Full | 0.9470 |
| Rights | Not shown | Profession Nurse | 0.9015 |
| Rights | Not shown | Age 55 | 0.8645 |
| Rights | Shown | Age 55 | 0.8262 |
| Rights | Shown | Age 25 | 0.7183 |
| Welfare | Shown | Nationality Netherlands | 0.7180 |
| Rights | Not shown | Age 25 | 0.6742 |
| Rights | Not shown | Language Full | 0.6433 |
| Welfare | Shown | Age 25 | 0.5793 |
| Welfare | Shown | First shown | 0.5360 |
| Rights | Shown | Profession Nurse | 0.5277 |
| Welfare | Not shown | Age 25 | 0.5157 |
| Welfare | Not shown | Language Full | 0.5033 |
| Welfare | Not shown | Age 55 | 0.4592 |
| Welfare | Not shown | Nationality Netherlands | 0.4300 |
| Rights | Not shown | Nationality Netherlands | 0.3412 |
| Rights | Not shown | First shown | 0.3408 |
| Welfare | Not shown | Gender Male | 0.2618 |
| Rights | Shown | (Intercept) | 0.2548 |
| Rights | Not shown | Gender Male | 0.2493 |
| Welfare | Shown | Age 55 | 0.2487 |
| Welfare | Shown | (Intercept) | 0.1993 |
| Rights | Shown | Gender Male | 0.1545 |
| Rights | Shown | First shown | 0.1487 |
| Rights | Not shown | (Intercept) | 0.1428 |
| Welfare | Shown | Gender Male | 0.1347 |
| Rights | Shown | Nationality Netherlands | 0.0960 |
| Welfare | Not shown | (Intercept) | 0.0790 |
| Welfare | Not shown | Profession Nurse | 0.0705 |
| Welfare | Not shown | First shown | 0.0345 |
| Welfare | Shown | Language Full | 0.0133 |
| Welfare | Shown | Profession Nurse | 0.0080 |
Hypothesis 4: respondents treat more favourably applicants that are similar to them
Not possible to assess at the aggregate level.
Hypothesis 5: Managerial accountability reduces discrimination
In other words, the existence of the treatment (shown) reduces discrimination.
This hypothesis only considers absolute discrimination (absolute values), not the direction).
tb <- S.theta %>%
pivot_wider(c(Iteration, Chain, Outcome, Population, Feature),
names_from = Treatment, values_from = value) %>%
group_by(Outcome, Population, Feature) %>%
summarize(`Prob (H5)` = length(which(
(`Not shown` < 0 & `Shown` > `Not shown`) |
(`Not shown` > 0 & `Shown` < `Not shown`))) / n()) %>%
arrange(desc(`Prob (H5)`))
tc <- "Evidence for Hypothesis 5."
if (knitr::is_latex_output()) {
kable(tb, format = "latex", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10)
} else {
kable(tb, format = "html", caption = tc, booktabs = TRUE) %>%
kable_styling(font_size = 10, position = "center", bootstrap_options = "striped", full_width = F)
}| Outcome | Population | Feature | Prob (H5) |
|---|---|---|---|
| Welfare | General population | Language Full | 0.9995 |
| Rights | General population | Nationality Netherlands | 0.9253 |
| Rights | General population | First shown | 0.8888 |
| Rights | General population | Profession Nurse | 0.8833 |
| Welfare | Public administration | Age 55 | 0.7962 |
| Rights | Public administration | Language Full | 0.7623 |
| Welfare | Public administration | First shown | 0.7493 |
| Rights | Public administration | Age 25 | 0.7280 |
| Rights | Public administration | (Intercept) | 0.7122 |
| Welfare | Public administration | (Intercept) | 0.7075 |
| Rights | General population | (Intercept) | 0.6867 |
| Welfare | General population | Age 25 | 0.6595 |
| Welfare | Public administration | Age 25 | 0.6550 |
| Welfare | Public administration | Nationality Netherlands | 0.6228 |
| Rights | Public administration | Profession Nurse | 0.6207 |
| Welfare | General population | (Intercept) | 0.5872 |
| Rights | General population | Gender Male | 0.5582 |
| Welfare | General population | Gender Male | 0.5443 |
| Rights | Public administration | First shown | 0.5362 |
| Rights | Public administration | Nationality Netherlands | 0.5060 |
| Rights | General population | Age 55 | 0.4970 |
| Rights | Public administration | Age 55 | 0.4652 |
| Rights | General population | Age 25 | 0.4370 |
| Welfare | Public administration | Language Full | 0.4300 |
| Rights | Public administration | Gender Male | 0.4162 |
| Welfare | General population | Age 55 | 0.4122 |
| Welfare | Public administration | Gender Male | 0.3785 |
| Welfare | General population | Profession Nurse | 0.3265 |
| Rights | General population | Language Full | 0.2827 |
| Welfare | General population | Nationality Netherlands | 0.2688 |
| Welfare | Public administration | Profession Nurse | 0.1748 |
| Welfare | General population | First shown | 0.0430 |