Most published BRFSS prevalence estimates are age-adjusted, and most first analyses of the microdata are not, which is the usual reason the two disagree. Health outcomes vary strongly with age, and places differ in how old their residents are, so a crude comparison of two states, or of one state over twenty years, mixes the thing being measured with the age structure it sits on. Direct standardization removes that mix-up by re-weighting every estimate to one fixed population, and the fixed population everyone uses is the year-2000 projected U.S. population (Klein and Schoenborn, 2001).
This package ships that standard as data and leaves the arithmetic to
survey::svystandardize(), which already does it correctly.
The contribution here is the right table and the worked pattern, not
another implementation.
The standard population
brfss_std_pop_2000 carries the 2000 projected U.S.
population in two groupings. set = "age19" is NCHS’s 19
five-year groups covering all ages. set = "adult6"
collapses the adult population to the exact groups of BRFSS’s calculated
variable _AGE_G (18-24, 25-34, 35-44, 45-54, 55-64, 65 and
over), which is what adult BRFSS prevalence work needs.
library(brfssdata)
std <- subset(brfss_std_pop_2000, set == "adult6")
std
#> # A tibble: 6 × 6
#> set age_group age_min age_max std_pop std_weight
#> <chr> <chr> <int> <int> <dbl> <dbl>
#> 1 adult6 18-24 18 24 26258428 0.129
#> 2 adult6 25-34 25 34 37233437 0.183
#> 3 adult6 35-44 35 44 44659185 0.219
#> 4 adult6 45-54 45 54 37030152 0.182
#> 5 adult6 55-64 55 64 23961506 0.118
#> 6 adult6 65+ 65 NA 34709480 0.170The weights reproduce the age-18-and-over distribution in Klein and Schoenborn’s Statistical Note 20.
CDC’s own BRFSS adjustment is coarser than six groups. The BRFSS
Direct Age Adjustment users guide standardizes most BRFSS questions
to three categories, 18-44, 45-64, and 65 and over, weighted
0.530534557, 0.299194019, and 0.170271424. Those are the same projected
population collapsed further, so the adult6 rows give them
by addition: rows 1 through 3 for the first, rows 4 and 5 for the
second, row 6 for the third.
p3 <- c(sum(std$std_pop[1:3]), sum(std$std_pop[4:5]), std$std_pop[6])
p3 / sum(p3)
#> [1] 0.5305366 0.2991955 0.1702679The two groupings serve different purposes. Three categories is what reproduces CDC’s published age-adjusted tables. Six controls the age variation inside 18-44 and 45-64 that three leaves alone, and returns a slightly different number for it. Both are defensible; the write-up has to say which one produced the figure.
Crude and adjusted, one year
_AGE_G codes 1 through 6 correspond, in order, to the
six rows of the adult6 set, so the design carries
everything standardization needs. Two preparation steps go in one
mutate(). survey’s post-stratification
machinery wants syntactic, factor-typed grouping variables, and CDC’s
names are neither, so make a plain copy of the age variable first (the
same reason brfss_design() builds itself on
brfss_wt and friends). The outcome gets its own column at
the same time, because svystandardize() has to be told
which variable the estimate will use (see Cautions).
Take the factor labels from the standard-population table rather than
writing them out. svystandardize() matches
population to the levels of by by position and
never checks names, so a table sorted some other way is silently wrong
rather than an error. What the whole pattern rests on is that
std runs in ascending age order, which makes its row
i the row for _AGE_G code i. The
stopifnot() below costs nothing and says so out loud, which
matters if you filter or re-sort the table on the way in.
library(srvyr)
library(survey)
des <- brfss_design(2023, vars = c("GENHLTH", "_AGE_G"), quiet = TRUE) |>
mutate(
age_group = factor(`_AGE_G`, levels = 1:6, labels = std$age_group),
fair_poor = GENHLTH >= 4
)
stopifnot(!is.unsorted(std$age_min))The crude estimate first. Codes CDC uses for don’t know and refused
are already NA here (na = TRUE is the design
default), so fair_poor is fair-or-poor health over
substantive answers.
des |>
summarize(fair_poor = survey_mean(fair_poor, na.rm = TRUE))
#> # A tibble: 1 × 2
#> fair_poor fair_poor_se
#> <dbl> <dbl>
#> 1 0.194 0.00139svystandardize() re-weights the design so its age
distribution matches the standard, and everything downstream works as
before:
std_des <- svystandardize(
des,
by = ~age_group,
over = ~1,
population = std$std_pop,
excluding.missing = ~ fair_poor + age_group
)
svymean(~fair_poor, std_des, na.rm = TRUE)
#> mean SE
#> fair_poorFALSE 0.81538 0.0014
#> fair_poorTRUE 0.18462 0.0014Over the whole 2023 file the two differ modestly, because the sample’s weighted age distribution sits near the 2000 standard. The adjustment earns its keep in comparisons.
That call standardizes on the six groups. Passing p3
instead, with a matching three-level by variable, is the
version that lines up with CDC’s published tables; the next section runs
both.
Comparing states with different age structures
Florida’s adult population is much older than Utah’s, and
fair-or-poor health rises steeply with age, so the crude comparison
mixes health with demography. The states argument keeps the
extract small, and filtering by state before the design is built is
variance-exact because BRFSS strata nest within state.
two <- brfss_design(
2023,
vars = c("GENHLTH", "_AGE_G"),
states = c("FL", "UT"),
quiet = TRUE
) |>
mutate(
age_group = factor(`_AGE_G`, levels = 1:6, labels = std$age_group),
age3 = cut(`_AGE_G`, c(0, 3, 5, 6), labels = c("18-44", "45-64", "65+")),
fair_poor = as.numeric(GENHLTH >= 4),
state = factor(`_STATE`)
)
stopifnot(!is.unsorted(std$age_min))
two |>
group_by(state) |>
summarize(fair_poor = survey_mean(fair_poor, na.rm = TRUE))
#> # A tibble: 2 × 3
#> state fair_poor fair_poor_se
#> <fct> <dbl> <dbl>
#> 1 12 0.191 0.00739
#> 2 49 0.143 0.00451Crudely, Florida (FIPS 12) sits at 19.1 percent and Utah (FIPS 49) at
14.3. Standardizing within each state (over = ~state) puts
both on the 2000 age distribution. CDC’s three categories first:
two_std3 <- svystandardize(
two,
by = ~age3,
over = ~state,
population = p3,
excluding.missing = ~ fair_poor + age3
)
svyby(~fair_poor, ~state, two_std3, svymean, na.rm = TRUE, vartype = "ci")
#> state fair_poor ci_l ci_u
#> 12 12 0.1796502 0.1643326 0.1949677
#> 49 49 0.1436004 0.1348084 0.1523924CDC’s Prevalence & Trends tool publishes 18.0 percent (16.4 to 19.5) for Florida and 14.4 (13.5 to 15.2) for Utah on this measure in 2023. Both estimates and both intervals agree to rounding, which is the check to run before a published figure and your own are set side by side.
The six-group standard is the same call with age_group
and std$std_pop in place of age3 and
p3:
two_std6 <- svystandardize(
two,
by = ~age_group,
over = ~state,
population = std$std_pop,
excluding.missing = ~ fair_poor + age_group
)
svyby(~fair_poor, ~state, two_std6, svymean, na.rm = TRUE, vartype = "ci")
#> state fair_poor ci_l ci_u
#> 12 12 0.1773874 0.1623208 0.1924540
#> 49 49 0.1430210 0.1344115 0.1516305Florida comes back at 17.7 against 18.0, a quarter of a point below the three-category figure and outside publication rounding; Utah stays at 14.3. The finer grouping is not wrong; it answers a slightly different question. What it will not do is reproduce CDC’s table, and it should not be offered as if it had. Choose by what the number is for: three categories to reproduce or sit beside a CDC estimate, six when the comparison is internal and the extra control is worth departing from the published series.
Under either grouping Florida drops by more than a point once its older age structure no longer counts against it, while Utah barely moves; about a quarter of the crude gap between the two states was age, not health.
Cautions
excluding.missing has to name every variable the
estimate will use, the outcome included. svystandardize()
drops the rows that are missing on any variable in that formula and then
computes the age weights over the survivors, so a variable left out of
the formula is still carrying age weight at the moment
na.rm = TRUE throws its row away. What comes back is then
an average of the age-group means weighted by each group’s standard
share times its weighted response rate, where direct standardization
weights by the standard share alone. The two coincide only when item
nonresponse is flat across the age groups, and it generally is not.
Naming the age variable as well, as ~ fair_poor + age_group
does above, is also what makes the requirement of a complete age
variable actually bite; _AGE_G is imputed by CDC and is
complete in modern files, so nothing is dropped on that account
here.
The consequence is that a standardized design belongs to one outcome.
Build a fresh std_des for each measure instead of reusing
one across several, since the calibration is specific to the rows that
measure answers. The discrepancy is small when an item is nearly
complete and grows with age-patterned nonresponse, and it is largest for
questions that sit behind a skip pattern and are put to only part of the
sample. Standardizing such an outcome puts the eligible subpopulation on
the standard age distribution, which is the comparison that was
wanted.
Standardize before estimating subgroups, and note that a subgroup’s
standardized estimate uses the standard’s age distribution, not the
subgroup’s own. CDC applies the three categories above to most BRFSS
questions, and a few measures carry a grouping of their own, so check
the measure’s documentation before claiming an exact match; the
age19 set covers the finer groupings.
The whole-file estimates above are not fifty-state figures. The 2023
public-use file covers 52 reporting areas, 48 states plus the District
of Columbia, Guam, Puerto Rico, and the U.S. Virgin Islands, with
Kentucky and Pennsylvania absent that year.
brfss_year_info() reports the jurisdiction count for every
year in its states column, and the datasets
article traces how coverage moved across the series; pass
states = when the write-up claims a specific universe.
Age adjustment is for comparison, not description. The crude estimate remains the actual burden in the population, and reports usually show both.
Sources
Klein RJ, Schoenborn CA. Age adjustment using the 2000 projected U.S. population. Healthy People 2010 Statistical Notes, No. 20. Hyattsville, Maryland: National Center for Health Statistics, 2001.
Centers for Disease Control and Prevention. BRFSS Direct Age Adjustment, published with the annual BRFSS documentation (https://www.cdc.gov/brfss/annual_data/annual_data.htm), which gives the three age categories and their weights.
The table itself is aggregated from SEER’s single-age rendering of
the Census P25-1130 projection (https://seer.cancer.gov/stdpopulations/), with the
anchors verified against the published group tables; see
?brfss_std_pop_2000.