Margins of error for combined estimates
Source:vignettes/theory-moe-propagation.Rmd
theory-moe-propagation.RmdEach American Community Survey (ACS) estimate for a census tract is
published with a margin of error (MOE), the half-width of its 90 percent
confidence interval. An estimate for a drive-time area combines the
estimates of the tracts that overlap the area, and catchmentACS computes
its margin of error from the margins of error of those tract estimates.
The margins of error are computed by
cacs_intersect_weight(), cacs_propagate_moe(),
and cacs_derive_rates(), whose help pages describe the
arguments and the columns of the results.
What the margins of error assume
Margins of error are not added like estimates. For independent
estimates the variances add, so the sum of two independent estimates
that each have a margin of error M has
a margin of error of \sqrt{2}\,M. The
formulas for sums, proportions, and ratios follow the approximations in
chapter 8 of the Census Bureau’s handbook for ACS data users (U.S. Census Bureau 2020). The package applies
the formula for a sum to the weighted tract estimates, treating the
weights as fixed constants. The margins of error for medians and
per-person values (rows whose estimand_family is
"median_proxy" or
"area_weighted_scalar_proxy") come from an average of tract
values weighted by area, an approximation made by the package rather
than one given in the handbook.
The margins of error are combined as if the tract estimates were independent. The handbook notes that its approximation formulas leave out the covariance between estimates, so they can overstate or understate the margin of error depending on the direction of the correlation. If the tract estimates in an area are positively correlated, the true margin of error of a count, a median, or a per-person value is larger than the reported one. For those three kinds of estimate the reported margins of error should be read as understatements of unknown size. For a rate, errors that move the numerator and the denominator in the same direction partly offset each other in the ratio, and the default ratio formula gives a wider margin than the proportion formula (both are given below). The package does not compute the net effect for a given rate.
According to the handbook, the variance replicate tables that the Census Bureau publishes for some 5-year detailed tables, including estimates for tracts, account for the covariance that the approximation formulas leave out. Folch et al. (2023) compare three ways of computing the margin of error of combined ACS estimates and find that only the one based on these replicate tables accounts for the covariance between the estimates, and that it performs best. The package does not use these tables.
The weights are treated as fixed numbers (the results record this
with weight_uncertainty_propagated = FALSE). The margins of
error therefore do not include error from area weighting itself. One
source they leave out is the assumption that whatever a variable counts
is spread evenly over each tract’s area. The other, for medians and
per-person values, is the weighting of tracts by area. They also take
each drive-time area as given, so they do not include uncertainty in the
area itself, which depends on the routing service and the road data used
to compute it.
Standard errors and the 90 percent level
Because the published margins of error are at the 90 percent
confidence level, the handbook converts a margin of error M to a standard error (SE) with \mathrm{SE} = M / 1.645 (U.S. Census Bureau 2020, chap. 7). A margin of
error at another confidence level is the standard error multiplied by
z, the normal quantile for the level,
qnorm(1 - (1 - level) / 2):
M = z\,\mathrm{SE}, \qquad \mathrm{SE} = \frac{M}{z}.
At the 90 percent level the package uses z
= 1.645, the value in the handbook, rather than
qnorm(0.95), which is 1.644854; the two differ by about
0.00015.
The package combines margins of error through their standard errors. For two independent estimates with margins of error M_1 and M_2 at the 90 percent level, the margin of error of their sum is
M_{\text{sum}} = z \sqrt{\left(\frac{M_1}{1.645}\right)^2 + \left(\frac{M_2}{1.645}\right)^2},
which at the 90 percent level, where z =
1.645, is \sqrt{M_1^2 + M_2^2},
the handbook’s formula for the margin of error of a sum. All the
formulas in this article have this form: each tract margin of error is
divided by 1.645 to give a standard error, the standard errors are
combined, and the combined standard error is multiplied by z. The level is set by the level
argument of cacs_propagate_moe(), 0.9 by default, and
cacs_derive_rates() uses the same level for the rates.
cacs_se_to_moe() and cacs_moe_to_se()
convert between standard errors and margins of error at a given
level.
Counts, medians, and per-person values
The help page of cacs_propagate_moe() names the
package’s four margin-of-error formulas Families A, B, C1, and C2. C1
and C2 also appear in the package’s warnings. The column
moe_formula_effective of the results names the formula used
for each row: "weighted_sum" (A),
"weighted_mean" (B), "proportion_subset" (C1),
or "general_ratio_conservative" (C2).
In the notation of
vignette("methodology", package = "catchmentACS"), a count
for the drive-time area of site s is
\widehat{Y}_s = \sum_j w^{cov}_{sj}\,
Y_j. Here Y_j is the estimate
for tract j and w^{cov}_{sj} its coverage weight, the share
of the tract’s area inside the drive-time area; the sum runs over the
tracts kept for the site, drive time, and variable. With the weights
fixed, the term for tract j has the
standard error w^{cov}_{sj} M_j /
1.645, where M_j is the
published margin of error of Y_j. The
margin of error of the count (Family A) is therefore
M(\widehat{Y}_s) = z \sqrt{\sum_j \left(\frac{w^{cov}_{sj}\, M_j}{1.645}\right)^2},
which at the 90 percent level is \sqrt{\sum_j (w^{cov}_{sj} M_j)^2}, the handbook’s formula for a sum applied to the weighted tract estimates. Each tract enters with its margin of error multiplied by its coverage weight. A tract that lies wholly inside the area enters with its full margin of error, and a tract with a tenth of its area inside with a tenth of it.
A median or a per-person value is estimated by \widehat{\bar X}_s = \sum_j w^{mean}_{sj}\, X_j, the average of the tract estimates X_j weighted by their area shares w^{mean}_{sj}, which sum to one. Its margin of error (Family B) is the same formula with the area shares in place of the coverage weights:
M(\widehat{\bar X}_s) = z \sqrt{\sum_j \left(\frac{w^{mean}_{sj}\, M_j}{1.645}\right)^2}.
This is the margin of error of the weighted average under the same
assumptions, independent tract estimates and fixed weights. The average
of tract medians or per-person values only stands in for the median or
per-person value of the drive-time area. Its margin of error describes
the sampling error of the average itself, and the difference between
that average and the drive-time area’s own value is not part of it.
Unlike sampling error, that difference does not shrink as the ACS
margins of error do (see
vignette("theory-spatial-aggregation", package = "catchmentACS")).
Since the area shares sum to one, \sum_j
w^{mean}_{sj} M_j is a weighted average of the tract margins of
error. As long as no kept tract has a negative margin of error, the
margin of error of the average is never larger than this weighted
average at the 90 percent level. It is smaller whenever two or more of
the kept tracts have a margin of error above zero.
Rates: the proportion and ratio formulas
A rate \widehat{R}_s = \widehat{A}_s / \widehat{B}_s divides two coverage-weighted counts, the numerator \widehat{A}_s and the denominator \widehat{B}_s, whose margins of error M(\widehat{A}_s) and M(\widehat{B}_s) come from the formula for counts above. For a ratio of two estimates the handbook gives two approximations (U.S. Census Bureau 2020, chap. 8). The proportion formula (C1) is for a proportion, a ratio whose numerator is part of its denominator, such as the people below the poverty level among the people for whom poverty status is determined:
M(\widehat{R}_s) = \frac{\sqrt{M(\widehat{A}_s)^2 - \widehat{R}_s^2\, M(\widehat{B}_s)^2}}{\widehat{B}_s}.
The ratio formula (C2) is for a ratio whose numerator is not part of its denominator:
M(\widehat{R}_s) = \frac{\sqrt{M(\widehat{A}_s)^2 + \widehat{R}_s^2\, M(\widehat{B}_s)^2}}{\widehat{B}_s}.
In the first-order (delta method) approximation, the variance of a ratio is
\mathrm{SE}(\widehat{R}_s)^2 \approx \frac{\mathrm{SE}(\widehat{A}_s)^2 - 2\widehat{R}_s\,\mathrm{Cov}(\widehat{A}_s, \widehat{B}_s) + \widehat{R}_s^2\,\mathrm{SE}(\widehat{B}_s)^2}{\widehat{B}_s^2}.
With the covariance of the numerator and the denominator set to zero, z times the square root of this expression is the ratio formula; with the covariance set to \widehat{R}_s\,\mathrm{SE}(\widehat{B}_s)^2, it is the proportion formula.
The two formulas differ by 2\widehat{R}_s^2
M(\widehat{B}_s)^2 under the square root, so the ratio formula
never gives the narrower margin. How much wider it is depends on \widehat{R}_s M(\widehat{B}_s) relative to
M(\widehat{A}_s): with r = \widehat{R}_s^2 M(\widehat{B}_s)^2 /
M(\widehat{A}_s)^2 below 1, the ratio formula gives \sqrt{(1 + r) / (1 - r)} times the margin of
error of the proportion formula. When the numerator is not zero, r is also the square of the ratio of the
denominator’s relative margin of error, M(\widehat{B}_s) / \widehat{B}_s, to the
numerator’s, M(\widehat{A}_s) /
\widehat{A}_s, because \widehat{R}_s =
\widehat{A}_s / \widehat{B}_s. The factor therefore varies from
rate to rate and from area to area, and it grows without bound as r approaches 1. When r is greater than 1, the value under the
square root of the proportion formula is negative; the next section
describes what the package does then. The second calculation by hand in
vignette("theory-derived-rates", package = "catchmentACS")
compares the two formulas for the labor force participation rate of one
drive-time area.
The numerator of each of the five built-in rates is part of its
denominator, so for all five the handbook’s formula is the proportion
formula. By default, cacs_derive_rates() and
cacs_run() nevertheless use the ratio formula for all five
rates.
vignette("theory-derived-rates", package = "catchmentACS")
lists the formula that each rate gets for each value of
formula_dispatch and gives the reason recorded for
unemp_rate.
A negative value under the square root
When the proportion formula is chosen for a rate and the value under its square root is negative, the package uses the ratio formula for that row instead, as the Census Bureau’s handbook advises (U.S. Census Bureau 2020, chap. 8). The value is negative when r > 1, that is, when the relative margin of error of the numerator is smaller than that of the denominator. The value under the square root of the ratio formula is a sum of two squares and cannot be negative, so the ratio formula never fails for this reason. It gives no margin of error only when the denominator is zero or one of its four inputs is missing.
The row records the substitution: moe_formula_requested
keeps "proportion_subset", the formula chosen, while
moe_formula_effective is
"general_ratio_conservative", with
moe_fallback = TRUE and
moe_fallback_reason = "negative_variance". A rate computed
with the formula chosen for it has moe_fallback = FALSE and
moe_fallback_reason = "n/a". These columns are described in
the help page of cacs_run(). When some rows have a
substituted formula or cannot be computed,
cacs_derive_rates() ends with one warning that counts them
by reason.
Zero denominators and missing values
A rate whose denominator is zero, or closer to zero than
sqrt(.Machine$double.eps) (about 1.5e-08), is
NA, as is its margin of error, with
moe_fallback = TRUE and
moe_fallback_reason = "zero_denominator". The check catches
only a denominator of about zero; a small positive denominator gives a
rate and a margin of error as usual.
A missing tract estimate or margin of error is not replaced by zero.
If a tract kept for an area has a missing estimate, the estimate and the
margin of error for the area are NA; if only the tract’s
margin of error is missing, only the margin of error is NA.
A rate is NA when its numerator or denominator, or the
margin of error of either, is missing, and the help page of
cacs_derive_rates() describes the columns that mark such
rates.
The Census Bureau’s data API puts negative annotation codes, such as
-555555555, in place of some estimates and margins of
error. cacs_intersect_weight() sets the six codes
(-222222222, -333333333, -555555555, -666666666, -888888888, and
-999999999) to NA, with a warning that counts them (a
margin-of-error code next to a missing estimate is not counted), and
they are then handled as the missing values described above. A code in
place of a tract margin of error, for example, makes the margin of error
of the variable NA for the areas that include the tract,
and a rate that uses the variable is NA. Any other negative
margin of error is not read as missing: the formulas above square it
like any other value, with no warning. cacs_acs_prefetch()
returns every negative margin of error as NA.
A count: its margin of error from the tract margins
The examples below use data bundled with the package. In these data
the ACS estimates and margins of error are random numbers, drawn
separately for each tract and variable, so the examples show how the
formulas work. One row in twenty instead has a missing estimate and the
code -555555555 in place of a margin of error. Because the
estimate is missing as well, those rows make a combined estimate
NA, and cacs_intersect_weight() gives no
warning about their codes. The draws ignore how the counts of a table
nest, so in some tracts the labor force is larger than the population 16
years and over. How often the proportion formula fails in these data,
for instance, says nothing about real ACS data.
We use the area whose weights are recomputed in
vignette("theory-spatial-aggregation", package = "catchmentACS"),
the 10-minute drive-time area of site AL_SITE_17.
cacs_intersect_weight() combines the tract estimates, and
with keep_tract_audit = TRUE it also keeps the coverage
weight of each tract in the attribute cacs_tract_audit:
# Drive-time areas and ACS data bundled with the package
iso <- readRDS(system.file(
"extdata", "legacy_2025_isochrones.rds", package = "catchmentACS"
))
acs <- readRDS(system.file(
"extdata", "sample_alabama_subset.rds", package = "catchmentACS"
))
site_id <- "AL_SITE_17"
drive_time <- 10L
iso_one <- iso[
iso$site_id == site_id & iso$drive_time_min == drive_time, ,
drop = FALSE
]
weighted <- cacs_intersect_weight(
iso_sf = iso_one,
acs_sf = acs,
weight_method = "area",
keep_tract_audit = TRUE,
verbose = FALSE
)
propagated <- cacs_propagate_moe(weighted, verbose = FALSE)The number of people below the poverty level
(B17001_002) is a count, so its margin of error comes from
the formula for counts (Family A). The package forms that formula by
applying the handbook’s formula for a sum to the weighted tract
estimates (U.S. Census Bureau 2020, chap.
8). The code below takes the coverage weight w^{cov}_{sj} of each tract
(area_wt) from cacs_tract_audit and the
published margin of error M_j from the
ACS data:
VAR <- "B17001_002"
audit <- attr(weighted, "cacs_tract_audit") |>
transmute(GEOID, w_cov = area_wt)
published <- acs |>
sf::st_drop_geometry() |>
filter(GEOID %in% audit$GEOID, variable == VAR) |>
select(GEOID, M_j = moe)
hand_A <- audit |>
left_join(published, by = "GEOID") |>
arrange(desc(w_cov))
hand_A#> # A tibble: 3 × 3
#> GEOID w_cov M_j
#> <chr> <dbl> <dbl>
#> 1 01125001001 1.000 200
#> 2 01127001001 0.113 218
#> 3 01123001001 0.113 152One tract lies wholly inside the area, and the other two have coverage weights of about 0.11. The code below computes the formula for counts term by term at the 90 percent level, where the divisor and the multiplier z are both 1.645:
z <- 1.645 # at another level, qnorm(1 - (1 - level) / 2)
se_term <- hand_A$w_cov * hand_A$M_j / 1.645 # standard error of each term
hand_moe_A <- z * sqrt(sum(se_term^2)) # z * combined standard error
list(
per_tract_var = round(se_term^2, 4),
summed_var = sum(se_term^2),
hand_moe_A = hand_moe_A
)#> $per_tract_var
#> [1] 14781.8294 222.3397 108.0903
#>
#> $summed_var
#> [1] 15112.26
#>
#> $hand_moe_A
#> [1] 202.223The tract inside the area contributes most of the sum. The other two tracts have margins of error of 218 and 152, but their coverage weights reduce their squared terms to about 1.3 percent of what they would be if the tracts lay wholly inside the area. The package’s row for the same variable:
pkg_A <- propagated |>
filter(variable == VAR) |>
select(variable, estimate, moe, moe_formula_effective, moe_fallback)
pkg_A#> # A tibble: 1 × 5
#> variable estimate moe moe_formula_effective moe_fallback
#> <chr> <dbl> <dbl> <chr> <lgl>
#> 1 B17001_002 1237. 202. weighted_sum FALSE
all.equal(hand_moe_A, pkg_A$moe)#> [1] TRUEThe package records "weighted_sum" (Family A) in
moe_formula_effective, and its margin of error, 202.2, is
the value computed by hand.
A rate that falls back to the ratio formula
With formula_dispatch = "auto", the proportion formula
is chosen for poverty_rate and
labor_force_participation. In the area above it can be used
for both rates, but in the 10-minute area of site
AL_SITE_08 it cannot be used for
poverty_rate:
iso_08 <- iso[iso$site_id == "AL_SITE_08" & iso$drive_time_min == 10L, ,
drop = FALSE]
propagated_08 <- cacs_propagate_moe(
cacs_intersect_weight(iso_sf = iso_08, acs_sf = acs, verbose = FALSE),
verbose = FALSE
)
rates_08 <- cacs_derive_rates(propagated_08, formula_dispatch = "auto",
verbose = FALSE)#> Warning in cacs_derive_rates(propagated_08, formula_dispatch = "auto", verbose = FALSE): Some rates are `NA` or use a replacement formula:
#> • `moe_fallback_reason = "negative_variance"`: 1 rate row used the ratio
#> formula (C2), because the proportion formula (C1) gave a negative variance.
#> ℹ The failure_origin and moe_fallback_reason columns give the reason on each
#> rate row.The warning counts one row in which the ratio formula replaced the proportion formula. The rows of the two rates for which the proportion formula was chosen show which one:
rates_08 |>
filter(variable %in% c("poverty_rate", "labor_force_participation")) |>
select(variable, estimate, moe, moe_formula_requested,
moe_formula_effective, moe_fallback, moe_fallback_reason) |>
glimpse()#> Rows: 2
#> Columns: 7
#> $ variable <chr> "poverty_rate", "labor_force_participation"
#> $ estimate <dbl> 0.4276799, 0.6136355
#> $ moe <dbl> 0.12541387, 0.02354859
#> $ moe_formula_requested <chr> "proportion_subset", "proportion_subset"
#> $ moe_formula_effective <chr> "general_ratio_conservative", "proportion_subset"
#> $ moe_fallback <lgl> TRUE, FALSE
#> $ moe_fallback_reason <chr> "negative_variance", "n/a"For poverty_rate, the proportion formula was chosen and
the ratio formula used, with moe_fallback = TRUE and
moe_fallback_reason = "negative_variance";
labor_force_participation kept the proportion formula. The
estimates and margins of error of the two counts of the poverty rate
show why:
counts <- propagated_08 |>
filter(variable %in% c("B17001_002", "B17001_001")) |>
select(variable, estimate, moe)
counts#> # A tibble: 2 × 3
#> variable estimate moe
#> <chr> <dbl> <dbl>
#> 1 B17001_001 2563. 557.
#> 2 B17001_002 1096. 216.
A <- counts$estimate[counts$variable == "B17001_002"] # numerator
B <- counts$estimate[counts$variable == "B17001_001"] # denominator
M_A <- counts$moe[counts$variable == "B17001_002"]
M_B <- counts$moe[counts$variable == "B17001_001"]
R <- A / B
list(
relative_moe = c(numerator = M_A / A, denominator = M_B / B),
under_root_proportion = M_A^2 - R^2 * M_B^2,
ratio_formula = sqrt(M_A^2 + R^2 * M_B^2) / B,
package = rates_08$moe[rates_08$variable == "poverty_rate"]
)#> $relative_moe
#> numerator denominator
#> 0.1970576 0.2171622
#>
#> $under_root_proportion
#> [1] -10006.85
#>
#> $ratio_formula
#> [1] 0.1254139
#>
#> $package
#> [1] 0.1254139The relative margin of error of the numerator, 0.197, is smaller than
that of the denominator, 0.217, so r =
1.21 is greater than 1 and the value under the square root of the
proportion formula is negative. The ratio formula gives 0.1254, the
margin of error that the package reports for
poverty_rate.