--- title: "Bootstrap weights" output: rmarkdown::html_vignette vignette: > %\VignetteIndexEntry{Bootstrap weights} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r, include = FALSE} knitr::opts_chunk$set( collapse = TRUE, comment = "#>", eval = nzchar(Sys.getenv("COMPILE_VIG_CANPUMF")) ) ``` ```{r setup, warning=FALSE} library(dplyr) library(canpumf) ``` PUMF data are a synthetic sample, altered from the original survey responses for privacy reasons. As with any survey sample, point estimates carry sampling uncertainty, and it is good practice to quantify it. Some PUMF ship with bootstrap weights for exactly this purpose; when they do not, **canpumf** can generate them with `add_bootstrap_weights()`. This vignette explains how `add_bootstrap_weights()` works: the method it uses, where the weights live, how stratification works, and — importantly — what happens when you call it again on a dataset that already has weights. ## The method: a resampling bootstrap For each of `n_replicates` replicates, a sample of `n` rows is drawn **with replacement** from the `n` rows of the (stratum of the) survey. The bootstrap weight for row $i$ in replicate $b$ is $$ w^{*}_{i,b} = w_i \cdot c_{i,b}, $$ where $w_i$ is the original survey weight and $c_{i,b}$ is the number of times row $i$ was selected in replicate $b$. Because each row is selected once on average, $\mathbb{E}[c_{i,b}] = 1$ and the replicate weights are unbiased for the original weights. The spread of an estimate across the replicates is what gives you its sampling variance (see [Estimating uncertainty](#estimating-uncertainty)). The default is `n_replicates = 500`. Pass `seed` for reproducibility. ## Two modes `add_bootstrap_weights()` works on two kinds of input and behaves differently for each: * **A DuckDB-backed lazy table** returned by `get_pumf()` (the typical case). The replicate weights are written into the DuckDB file as a persistent, separate table and exposed through a VIEW that joins the survey table to the weights. The returned `tbl` points at that view. Because the weights are persisted, **subsequent calls reuse them** instead of recomputing. * **An in-memory `data.frame` / `tibble`** (e.g. after `collect()`). The weights are generated in memory and the augmented data frame is returned. Nothing is persisted. ```{r} cis <- get_pumf("CIS", "2019") cis_bsw <- cis |> add_bootstrap_weights(weight_col = "FWEIGHT", n_replicates = 200, seed = 42) # 200 replicate columns CPBSW1 … CPBSW200 are now available grep("^CPBSW", colnames(cis_bsw), value = TRUE) |> head() ``` ## Where the weights are stored (DuckDB path) The replicate weights are **not** copied into the main survey table. Instead: * They are written to a separate DuckDB table, auto-named `paste0("pumf_bsw_", tolower(weight_col))` (e.g. `pumf_bsw_pweight`). The auto-naming means generating weights for two different weight columns (e.g. a household weight and a person weight) produces two independent BSW tables that never collide. * A companion VIEW (e.g. `eng_bsw_pweight`) joins the survey table to the BSW table on a row identifier, so the returned `tbl` exposes every survey column plus the replicate columns. `bsw_info()` summarises what is stored: ```{r} bsw_info(cis_bsw) ``` This separation keeps the main table untouched and lets you add, inspect, or remove replicate weights independently of the survey data. ## Identifying rows To join replicate weights back to the survey rows a stable row identifier is needed. With `id_col = NULL` (the default): * the survey registry's `bsw_join_key` is used when one is defined (e.g. `PEFAMID` for the SFS), requiring no change to the table; otherwise * a `pumf_row_id` column (DuckDB `rowid`) is added to the survey table once. You can always pass `id_col` explicitly to use a natural key. ## Stratified bootstrap weights Many survey designs are stratified, and the resampling should respect that: each replicate should resample **within** each stratum, preserving the stratum sample sizes. Pass `strata_cols`: ```{r, eval = FALSE} cis_bsw <- cis |> add_bootstrap_weights(weight_col = "FWEIGHT", strata_cols = "PROV", n_replicates = 200, seed = 42) ``` Resolution order for the strata: 1. `strata_cols` if you pass it, 2. otherwise the registry `bsw_strata` for the survey, 3. otherwise, for the LFS, the default `c("SURVYEAR", "SURVMNTH")` so each month is resampled separately, 4. otherwise unstratified. Pass `strata_cols = character(0)` to force unstratified weights even when a default exists (e.g. for the LFS). ## Estimating uncertainty {#estimating-uncertainty} To get a bootstrap standard error for an estimate, compute the estimate once with the original weights and once with each replicate weight; the spread of the replicate estimates is the sampling variability. The example below estimates the total population represented by the survey (the sum of the weights), which depends only on the weight columns, but the same pattern applies to any weighted statistic — replace `sum(.x)` with your estimator (a weighted mean, share, total, …) evaluated with each weight column. ```{r} est <- cis_bsw |> summarise(across(c(FWEIGHT, matches("^CPBSW[0-9]+$")), ~ sum(.x))) |> collect() point_estimate <- est$FWEIGHT replicate_estimates <- est |> select(matches("^CPBSW[0-9]+$")) |> unlist() # Bootstrap variance: mean squared deviation of the replicate estimates from # the full-sample estimate; the standard error is its square root. # Confidence intervals can be obtained by taking the appropriate quantiles. confidence_interval <- quantile(replicate_estimates,c(0.025,0.975)) std_error <- sqrt(mean((replicate_estimates - point_estimate)^2)) c(estimate = point_estimate, se = std_error, conf=confidence_interval) ``` The same pattern works after `group_by()` / `summarise()`: carry the `CPBSW*` columns through the summary, then collapse them into a standard error per group. ## Incremental re-runs Calling `add_bootstrap_weights()` again on a survey that already has weights does only the work needed to satisfy the request. The decision is driven by two questions: are there enough replicate columns, and do all rows have weights? ### Reuse — nothing to do If the stored weights already cover every row and there are at least `n_replicates` of them, they are reused without recomputation (requesting fewer than are stored simply returns a subset). This is why repeatedly opening a weighted survey is cheap. ### More replicates requested If you ask for more replicates than are stored (and no rows were added), the additional replicate columns are generated and appended; the existing columns are kept unchanged. When stratified, the new columns are resampled within strata, consistent with the existing ones. ```{r, eval = FALSE} # First 200 replicates … cis_bsw <- add_bootstrap_weights(cis, "FWEIGHT", n_replicates = 200, seed = 42) # … later, extend to 500: only CPBSW201 … CPBSW500 are generated. cis_bsw <- add_bootstrap_weights(cis, "FWEIGHT", n_replicates = 500, seed = 42) ``` ### Rows added to the survey table This is the subtle case. A bootstrap replicate resamples the **full population** (or full stratum), so adding rows invalidates the replicate weights of the affected resampling universe — you cannot simply generate weights for the new rows in isolation. `add_bootstrap_weights()` detects that not all rows have weights and **deletes and regenerates** the affected weights: * **Unstratified:** the whole population changed, so every row's replicate weights are regenerated. * **Stratified:** only the strata that gained rows are regenerated in full; strata that are already complete keep their existing weights. This makes re-weighting an updated multi-year file (for example, appending a new LFS month or year) inexpensive — only the new strata are recomputed. Extending the number of replicates and adding rows can happen in the same call; each is handled on the rows it applies to. > Extending a weighted table after appending rows relies on the new rows having > a valid identifier. Use a natural `id_col` (or the registry `bsw_join_key`) > when you expect to grow the table, so the new rows can be matched. ### Forcing a full regeneration Pass `overwrite = TRUE` to discard the stored weights and regenerate from scratch regardless of what is already there. ```{r, eval = FALSE} cis_bsw <- add_bootstrap_weights(cis, "PWEIGHT", n_replicates = 200, seed = 42, overwrite = TRUE) ``` ## Multiple weight columns Because each BSW table is named after its weight column, hierarchical surveys with several weights are handled by simply calling the function once per weight: By default all replicate columns share the `CPBSW` prefix, so when storing weights for more than one weight column give each a distinct `prefix` to tell the two sets of replicates apart (the column names are survey-specific): ```{r, eval = FALSE} tbl <- get_pumf("", "") |> add_bootstrap_weights(weight_col = "", prefix = "HHBSW") |> add_bootstrap_weights(weight_col = "", prefix = "PPBSW") ``` ## Filtered input tables Bootstrap weights always cover the complete physical survey table, so variance estimates use the full sample. If the input `tbl` carries `filter()` operations, they are captured and re-applied to the returned view, so the visible rows still match your subset. Other operations (`select()`, `mutate()`, …) are not replayed — apply them to the returned tbl instead. ## Connection lifecycle Generating weights on a DuckDB-backed table requires exclusive write access, so `add_bootstrap_weights()` shuts down the connection held by the input `tbl`. **Use the returned tbl** afterwards; the input `tbl` (and any other lazy tables backed by the same DuckDB file) are invalid after the call. ## Inspecting and removing weights ```{r, eval = FALSE} bsw_info(cis_bsw) # list stored BSW tables and replicate counts cis_clean <- remove_bootstrap_weights(cis_bsw, weight_col = "PWEIGHT") ``` `remove_bootstrap_weights()` drops the BSW table and its companion view and returns a clean tbl on the original survey table.