In FPRgroup/FPEMcountry: fpemlocal

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
knitr::opts_chunk$set(echo = TRUE)
library(FPEMlocal)
devtools::load_all()

Introduction

In this vignette we will fit FPET to multiple countries and aggregate the samples to obtain results for aggregate levels.

Table of Contents

1. Fit models
2. Aggregate
3. Calculate results

Fit models

Let's task ourselves with obtaining results for each sub-region within Africa.

divisions %>%
  dplyr::filter(name_region == "Africa") %>%
  dplyr::select(sub_region_numeric_code, name_sub_region) %>%
  unique()

First, we need fits for all countries in Africa. Start by pulling a vector of country codes.

divs <- divisions %>%
  dplyr::filter(name_region == "Africa") %>%
  dplyr::filter(division_numeric_code %in% contraceptive_use$division_numeric_code) %>% # Filtering on the countries available in `contraceptive_use` because we are using this package data
  dplyr::pull(division_numeric_code)

Next, map the country codes to the function fit_fp_c. Note: Supply country codes as a list to be mapped.

first_year <- 1970
last_year <- 2030
union <- "Y"
fits <- purrr::pmap(list(division_numeric_code = divs %>% as.list(),
                         is_in_union = union,
                         first_year = first_year,
                         last_year = last_year
                         ),
                    fit_fp_c)

Aggregate

Within out list of fits we have more complex lists. We only need the samples from these fits. This function plucks the posterior samples and combines each array of samples along index 1 making one large array. The first index of this large array corresponds to the country of the samples. Use dimnames() on the resulting list to see how to index a country. This is the identical format of a posterior samples array from FPEMglobal.

posterior_samples <- pluck_abind_fp_c(fits, divs)

To aggregate these samples we need to construct two tibbles. First, create a tibble with the aggregation level and corresponding division codes. The aggregation level column must be named division_level and contain division codes. In this example we will obtain aggregate results for all sub-regions in Africa. To aggregate samples to each sub-region, use the sub_region_numeric_code column of the data set division. The use may construct a custom division level column.

## Get divisions data
division_level_data <- divisions %>%
   dplyr::rename(division_level = sub_region_numeric_code) %>%
   dplyr::rename(division_level_name = name_sub_region) %>%
   dplyr::select(division_numeric_code, division_level, division_level_name) %>%
   dplyr::filter(division_numeric_code %in% divs)
division_level_data

Lastly, we need the population counts for these countries filtered by the time frame and union status.

## Get population data
population_data <- population_counts %>%
  dplyr::filter(division_numeric_code %in% divs) %>%
  dplyr::filter(is_in_union == union) %>%
  dplyr::filter(mid_year >= first_year) %>%
  dplyr::filter(mid_year <= last_year)
population_data

Now supply these two tibbles of data and the large array of posterior samples to the function aggregate_fp. This function returns a named list of aggregated samples.

devtools::load_all()
posterior_samples_aggregated <- aggregate_fp_multilevel(division_level_data,
                                         population_data,
                                         posterior_samples =  posterior_samples)
names(posterior_samples_aggregated)

Calculate results

We can map the list of samples to calc_fp to obtain results for each sub-region

results_list <- purrr::pmap(list(posterior_samples = posterior_samples_aggregated,
                            country_population_counts = list(population_data), #need to wrap any complex inputs in another list()
                            first_year = first_year),
                       calc_fp)
results_list$`Southern Africa`$contraceptive_use_modern

FPRgroup/FPEMcountry documentation built on April 24, 2023, 4:32 p.m.