Statistics needed for whole population and groups (by time and run) important to the study of inequalities (education, households type, gender, age group, income quintile)

• mean: ghqcase, ghqlike, employment, hours worked, income, poverty
• median: ghqlike, hours worked, income
• count: ghq cases, employed, hours worked categories, poverty
• iqr: ghqlike, hours worked, income
• 25th percentile: ghqlike, hours worked, income
• 75th percentile: ghqlike, hours worked, income
• 10th percentile: income
• 90th percentile: income
• 90/10 ratio: income

We have 6 outcomes, some of which must be constructed from the LABSim output:

1. mean of continuous GHQ12 score (variable dhm)
2. prevalence (mean of a dummy) of GHQ12 caseness (dummy=1 if dhm <= 24)
3. prevalence (mean of a dummy) of employment (dummy=1 if les_c4 == "EmployedOrSelfEmployed")
4. mean hours worked (convert laboursupplyweekly to numerical)
5. mean household income (variable equivaliseddisposableincomeyearl)
6. prevalence (mean of a dummy) of poverty (dummy =1 if atriskofpoverty == 1 || atriskofpoverty == null)

We construct 13 groups, and also include results for the whole population:

1. male (dgn == "Male")
2. female (dgn == "Female")
3. age 25-44 (dag >= 25 && dag < 45)
4. age 45-64 (dag >= 45 && dag < 65)
5. household with children (one of n_children_1-17 != 0)
6. household without children (all of n_children_1-17 == 0 or missing)
7. employed (les_c4 == "EmployedOrSelfEmployed")
8. unemployed (les_c4 == "NotEmployed")
9. in-work poverty (grp_emp == 1 && out_poverty == 1)
10. out-of-work poverty (grp_emp == 0 && out_poverty == 1)
11. Low education (deh_c3 == "Low")
12. medium education (deh_c3 == "Medium")
13. high education (deh_c3 == "High")

We want:

1. aggregate mean outcomes for the baseline and reform, in each year, and year-run pairs.
2. Effect estimates (difference in aggregate outcomes between the baseline and reform), again for each year, and each year-run pair.
3. Ranking of outcomes and effect estimates, again for each year, and each year-run pair.

Table structure:

• Write a document outlining general rules for contributions.

This should make any code reviews much easier.

As a quick initial run of creating graphs, here's the running of the code currently in R/outputting graphs.R:

library(readr)
library(tidyverse)
library(SPHSUgraphs)

out_data <-
  read_csv("C:/Programming/covid19_effect_estimates/data/new_data.csv",
           show_col_types = FALSE)


# tidying dataset ---------------------------------------------------------


compare_results <- out_data |>
  filter(grp_all == TRUE, !is.na(run)) |>
  select(-contains("eff"), -starts_with("grp")) |>
  pivot_longer(
    -c(scenario, run, time),
    names_to = c("metric", "outcome", "policy"),
    values_to = "val",
    names_pattern = "(.*)_(.*)_(baseline|reform)"
  ) |>
  pivot_wider(
    c(scenario, run, time, outcome, policy),
    names_from = metric,
    values_from = val
  )


# faceted graph -----------------------------------------------------------

compare_results |>
  ggplot(aes(time, out, colour = policy, fill = policy)) +
  geom_vline(xintercept = 2019, colour = "red") +
  stat_summary(
    fun.data = mean_se,
    geom = "ribbon",
    alpha = 0.5,
    colour = NA
  ) +
  stat_summary(fun.data = mean_se, geom = "line") +
  stat_summary(fun.data = mean_se, geom = "point") +
  facet_wrap(~ outcome, scales = "free_y") +
  scale_fill_manual(
    "Policy",
    aesthetics = c("fill", "colour"),
    labels = c("Baseline", "Covid policy"),
    values = sphsu_cols("University Blue", "Thistle", names = FALSE)
  ) +
  theme(legend.position = "bottom")

As a small initial point - these outputs currently have a very small range (accidently put just as standard error of means across 50 runs in file). Should intervals combine the sd's of the means of each run, rather than taking the variance between the mean outputs?

Suggested edits so far:

• Output separate graphs
• Update ribbons to 95% CIs (rather than default se)
• Correct variable names and titles
• Standardise colour and display

^{Created on 2022-07-15 by the reprex package (v2.0.1)}

This code needs some data validation.

@dkopasker I'd like you to describe here what you expect from every variable in the raw data files, that includes their range, possible NA or NaN, etc. In addition, we need to clearly state how the code processes such values. Common options include dropping such entries, asking the aggregate functions to ignore them, or replacing with some imputed values (mean of some sort, median).
This approach should make the data analysis much more reproducible.

We should also consider LABsim output as potentially corrupted as the code itself is not tested properly. Constant changes in the code do not help here either. That means this script must notify every user in the case any input value is out of expected range.

Hi @vkhodygo,

The R code used to aggregate the results from 1,000 runs of the simulation has an unusual amount of ties up to the eighth decimal place. For example, two observations for out_ghqcase_baseline in grp_age25 in 2020 have a value of 0.54180604. This happens multiple times across various outcomes, groups, and years. Could you please review the code to ensure there is not an error?

• Reduce the size
• Introduce more general functions, no hardcoded names
• Add docstrings
• Add tests