Estimate poverty and mean in CHN using new files

CHN_update.R

#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# libraries   ---------
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~


library(data.table)
library(collapse)
library(ggplot2)

options(pipload.verbose = FALSE)

#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Load data   ---------
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
yrs   <- c(2016, 2017, 2018)
yrs   <- c(1984:2018)
yrs   <- c(1984:2020)


povlines <- c(2.15, 3.65, 6.85)

# ld   <- pipload::pip_load_cache("CHN", yrs, type = "list")
ld   <- pipload::pip_load_cache("CHN", 
                                type = "list", 
                                version = "20230626_2017_01_02_TEST")

#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
## population --------

pop <- pipload::pip_load_aux("pop")
pop <- pop[country_code == "CHN"]

#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
## Group data means in LCU --------

gdm <- pipload::pip_load_aux("gdm")
means <- 
  gdm[country_code == "CHN" 
  ][, 
    mean := survey_mean_lcu* (12/365)
  ][, 
    .(surveyid_year, pop_data_level, mean )]

means <- split(means, by = "surveyid_year")

means <- 
  purrr::map(.x = means, 
             .f = ~{
               y <-       .x[, mean]
               names(y) <- .x[, pop_data_level]
               y
             })

years     <- names(means)
names(ld) <- years

#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Functions   ---------
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~


gd_povstats <- function(df, mn) {
  
  levels <- df[, unique(as.character(reporting_level))]
  
  y <- purrr::map_df(
    .x = levels, 
    .f = ~{
      
      dfx <- df[reporting_level == .x]
      mnx <- mn[.x]
      
      cpi  <- dfx[, unique(cpi)]
      ppp  <- dfx[, unique(ppp)]
      
      mnx <- wbpip::deflate_welfare_mean(
        welfare_mean = mnx, 
        ppp = ppp, 
        cpi = cpi)
      
      st <- 
        purrr::map(
          .x = povlines, 
          .f = ~{
            rs <- wbpip:::gd_compute_poverty_stats(
              welfare = dfx$welfare, 
              population = dfx$weight, 
              povline = .x,
              requested_mean = mnx)
            
            rs <- as.data.table(rs)
          }
        ) |> 
        rbindlist(use.names = TRUE) |> 
        ftransform(data_level = .x, 
                   mean       = mnx) 
      
      return(st)
    }
  )
  return(y)
}


#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Execution   ---------
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~


dd <- purrr::map_df(.x = years, 
                    .f = ~{
                      mn <- means[[.x]]
                      df <- ld[[.x]]
                      z <- gd_povstats(df, mn)
                      z <- copy(z)
                      z[, year := as.numeric(.x)]
                      return(z)
                    })

df <- joyn::merge(dd, 
                  pop, 
                  by         = c("year", 
                                 "data_level = pop_data_level"), 
                  match_type = "m:1", 
                  keep       = "left", 
                  reportvar  = FALSE)

dfn <- 
  df |> 
  fgroup_by(year, poverty_line) |> 
  get_vars(c("headcount", 
             "poverty_gap", 
             "poverty_severity", 
             "watts", 
             "pop")) |> 
  fmean(w = pop, keep.w = FALSE) |> 
  ftransform(data_level = "national")


vars <- names(dfn)
dff  <- rbindlist(list(df[, ..vars], dfn), use.names = TRUE) |>
  setorder(year, poverty_line, data_level) # It invisibly returns the data


filename <- fs::path(tdirp, "CHN_pov_mean", ext = "dta")
haven::write_dta(dff,filename)


# chart
ggplot(dff[year >= 2000], 
       aes(x = year, 
           y = headcount, 
           color = as.factor(poverty_line))
) +
  geom_line() +
  facet_wrap(~data_level) +
  theme_minimal() +
  theme(legend.title = element_blank(), 
        legend.position = "bottom")