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

## -----------------------------------------------------------------------------
library(xmap)
library(dplyr)
library(stringr)
library(ggplot2)

## -----------------------------------------------------------------------------
split_isiccomb <- function(threefour_df) {
  interim <- list()
  interim$isiccomb.rows <-
    threefour_df |>
    filter(str_detect(isiccomb, "[:alpha:]"))
  interim$isiccomb.avg <-
    interim$isiccomb.rows |>
    group_by(country, year, isiccomb) |>
    mutate(value = tidyr::replace_na(value, 0)) |>
    summarise(
      avg.value = mean(value),
      n_isic = n_distinct(isic),
      n_rows = n(),
      .groups = "drop"
    )
  left_join(
    threefour_df,
    interim$isiccomb.avg,
    by = c("country", "year", "isiccomb")
  ) |>
    rename(value.nosplit = value) |>
    mutate(
      value = coalesce(avg.value, value.nosplit),
      split.isiccomb = !is.na(avg.value)
    ) |>
    select(country, year, isic, isiccomb, value, value.nosplit, split.isiccomb)
}

split_links <- indstat$masked_sample |>
  split_isiccomb() |>
  mutate(weights = value / 1000) |>
  tidyr::drop_na(weights)

crossmaps <- split_links |>
  nest_by(country, year, .key = "links") |>
  mutate(
    valid = validate_as_xmap(links, isiccomb, isic, weights),
    xmap = list(as_xmap_tbl(links, isiccomb, isic, weights))
  ) |>
  ungroup()

## -----------------------------------------------------------------------------
weight_type <- crossmaps |>
  select(country, year, links) |>
  tidyr::unnest(links) |>
  dplyr::left_join(
    indstat$country_lookup[c("code", "iso3c", "name")],
    by = join_by(country == code)
  ) |>
  rename(country_iso3c = "iso3c") |>
  mutate(weight_type = ifelse(weights == 1, "unit", "frac"))

weight_type |>
  ggplot(aes(x = year, fill = weight_type)) +
  facet_wrap(~country_iso3c) +
  geom_bar(position = "stack") +
  scale_fill_manual(values = c("frac" = "#ef8a62", "unit" = "#67a9cf")) +
  theme_minimal() +
  theme(legend.position = "bottom") +
  labs(y = "`isic` codes", x = "Year", fill = "Weight type")

## -----------------------------------------------------------------------------
group_summary <- weight_type |>
  summarise(
    any_isiccomb = any(weight_type == "frac"),
    pct_isiccomb = sum(weight_type == "frac") / n(),
    .by = c(country_iso3c, year)
  ) |>
  mutate(pct_year_isiccomb = sum(any_isiccomb) / n(), .by = c(country_iso3c))

group_summary |>
  ggplot(aes(
    x = as.factor(year),
    fill = pct_isiccomb,
    y = forcats::fct_reorder(country_iso3c, pct_year_isiccomb)
  )) +
  geom_tile(color = "black") +
  ## abbreviate year labels
  scale_x_discrete(labels = ~ str_sub(.x, -2)) +
  scale_fill_stepsn(
    colours = c("white", RColorBrewer::brewer.pal(6, "Oranges")),
    breaks = c(0.0001, seq(0.2, 1, 0.2)),
    limits = c(0, 1),
    labels = scales::label_number(accuracy = 0.01)
  ) +
  # theme_minimal() +
  # move legend to bottom
  theme(legend.position = "bottom") +
  labs(
    x = "Year",
    y = "Country",
    fill = "Fraction of `isic` derived from `isiccomb` split",
    title = "Country, Year crossmaps for splitting reported values across combinations of ISIC Rev. 3 codes",
    subtitle = "Sorted by no. of years with at least one `isic` value derived from an `isiccomb` code"
  )

## -----------------------------------------------------------------------------
group_summary |>
  ggplot(aes(
    x = as.factor(year),
    y = forcats::fct_reorder(country_iso3c, pct_year_isiccomb),
    fill = any_isiccomb
  )) +
  geom_tile(color = "black") +
  scale_x_discrete(labels = ~ str_sub(.x, -2)) +
  scale_fill_manual(values = c("TRUE" = "#ef8a62", "FALSE" = "#67a9cf")) +
  theme_minimal() +
  theme(legend.position = "bottom") +
  labs(
    x = "Year",
    y = "Country",
    fill = "Crossmap involves at least 1 `isiccomb` split",
    title = "Country, Year crossmaps for splitting reported values across combinations of ISIC Rev. 3 codes",
    subtitle = "Sorted by no. of years with at least one `isic` value derived from an `isiccomb` code"
  )

## -----------------------------------------------------------------------------
(isic3_xmap <- split_links |>
  distinct(isic) |>
  mutate(isic = as.character(isic), isic3 = str_sub(isic, 1, 3), weight = 1) |>
  as_xmap_tbl(from = isic, to = isic3, weight_by = weight))

## -----------------------------------------------------------------------------
(isic3_summary <- isic3_xmap |>
  summarise(
    `.from` = glue::glue_collapse(.from, "+"),
    `n.from` = n(),
    .by = c(.to)
  ) |>
  arrange(desc(n.from)))

## -----------------------------------------------------------------------------
isic3_agg <- isic3_summary |>
  mutate(agg = n.from > 1)

isic3_agg |>
  tidyr::uncount(n.from, .id = "unit", .remove = FALSE) |>
  mutate(agg = forcats::fct_relevel(factor(agg), "TRUE", "FALSE")) |>
  ggplot(aes(
    y = forcats::fct_reorder(.to$isic3, n.from),
    x = unit,
    fill = n.from,
  )) +
  facet_wrap(vars(agg), scales = "free_y") +
  geom_tile(width = 0.8, height = 0.8, color = "white") +
  scale_x_continuous(breaks = scales::breaks_width(1)) +
  scale_fill_viridis_c(breaks = scales::breaks_width(1)) + # + coord_equal()
  labs(
    x = "isic code (ordinal position within isic3 group)",
    y = "isic3 code",
    fill = "no. of isic codes aggregated",
    title = "Number of isic codes aggregated into each isic3 code",
    subtitle = "Faceted by whether the isic3 code aggregates more than one isic code"
  )

## -----------------------------------------------------------------------------
isiccomb_to_isic3 <- split_links |>
  mutate(isic = as.character(isic), isiccomb = as.character(isiccomb)) |>
  nest_by(country, year, .key = "links") |>
  mutate(
    xmap1 = list(as_xmap_tbl(links, isiccomb, isic, weights)),
    composed = list(tidyr::unpack(
      compose_xmap(xmap1, isic3_xmap),
      everything()
    ))
  ) |>
  reframe(composed)

isiccomb_to_isic3 |>
  filter(country == "276", year == 1991, isiccomb == "151A")

## -----------------------------------------------------------------------------
comb_values <- indstat$masked_sample |>
  mutate(isiccomb = as.character(isiccomb)) |>
  group_by(country, year, isiccomb) |>
  summarise(value = value[!is.na(value)][1], .groups = "drop") |>
  filter(!is.na(value))

direct_isic3 <- comb_values |>
  nest_by(country, year, .key = "values") |>
  mutate(
    xmap_grp = list(
      isiccomb_to_isic3 |>
        filter(.data$country == .env$country, .data$year == .env$year) |>
        as_xmap_tbl(isiccomb, isic3, weight_by)
    ),
    applied = list(apply_xmap(
      values,
      xmap_grp,
      values_from = value,
      keys_from = isiccomb
    ))
  ) |>
  reframe(applied)

direct_isic3

## -----------------------------------------------------------------------------
two_step_isic3 <- split_links |>
  mutate(isic3 = str_sub(as.character(isic), 1, 3)) |>
  group_by(country, year, isic3) |>
  summarise(value_two_step = sum(value, na.rm = TRUE), .groups = "drop")

direct_isic3 |>
  full_join(two_step_isic3, by = c("country", "year", "isic3")) |>
  mutate(diff = value - value_two_step) |>
  summarise(max_abs_diff = max(abs(diff), na.rm = TRUE))