catchment

The goal of catchment is to estimate health facility catchment populations using Bayesian gravity models. The primary purpose of these population estimates is to calculate catchment-level incidence rates using routine reporting data from health information systems.

The model is fit with TMB and provides:

• a learned distance-decay parameter (exponential or power family), estimated from the data rather than fixed a priori;
• uncertainty on catchment populations, propagated through the random effects via TMB::sdreport();
• optional negative-binomial likelihood for overdispersed counts;
• validation tooling — posterior-predictive checks (pp_check()) and leave-one-facility-out cross-validation (loo_facility_cv());
• plot methods for dominant catchments and per-facility access-probability surfaces.

See the getting-started and prior-predictive vignettes for full walkthroughs.

Installation

Fitting the catchment model requires the R-INLA and Template Model Builder (TMB) packages. INLA is not on CRAN, and the installation process may depend on the type of computer you are using. I recommend following the installation instructions for INLA first. Be sure to close out of any active R session before installation!

The catchment package is also not on CRAN (yet!). You can install the development version of catchment from GitHub with:

# install.packages("devtools")
devtools::install_github("PATH-Global-Health/catchment")

Example

# Load libraries
library(terra)
library(sf)
library(catchment)

# 1. Pre-processing --------------------------------------
## Load example data
data("example_shp")   # an sf polygon
data("example_locs")  # a data frame of facility locations
pop <- example_pop()  # a terra SpatRaster

## Get friction surface with the traveltime package
## (https://github.com/idem-lab/traveltime)
fric <- traveltime::get_friction_surface(surface = "walk2020", extent = pop) |>
  terra::mask(terra::vect(example_shp))
fric <- terra::resample(fric, pop, method = "average")

## Create output folder
f <- tempfile()
fs::dir_create(fs::path(f, "tt"))

# 2. Travel time surfaces --------------------------------
## Create individual travel time rasters
create_travel_surface(friction_surface = fric, extent_file = pop,
  points = example_locs, id_col = "label", x_col = "x", y_col = "y",
  output_dir = fs::path(f, "tt"), individual_surfaces = TRUE)

## Organize travel time matrix
tmat <- travel_mat_from_folder(dir = fs::path(f, "tt"), reference = pop)

# 3. Fit catchment model ---------------------------------
## Organize input data. Passing the raw travel-time matrix lets the model
## learn the distance-decay parameter (sparsity controls live here).
catch_dat <- prepare_data(prob_mat_init = tmat, pop_raster = pop,
  location_data = example_locs, minimum_time = 10, force_threshold = 300,
  n_fac_limit = 10, mesh.args = list(cutoff = 0.1, max.edge = c(0.1, 4)))

## Fit catchment model (Poisson + learned exponential decay by default)
mod <- catchment_model(catch_dat)
mod                          # family, decay, convergence, pdHess
mod$decay_param              # estimated decay (tau in minutes)

# 4. Post-processing -------------------------------------
## Estimated catchment populations, with uncertainty
catchment_populations(mod)
catchment_populations(mod, uncertainty = TRUE)   # SEs + 95% CIs

# 5. Validation & visualization --------------------------
pp_check(mod)                # posterior-predictive coverage + dispersion
loo_facility_cv(mod)         # leave-one-facility-out CV
plot(mod)                    # dominant-catchment map
plot_prob_surface(mod, id_label = example_locs$label[1])