birdie: Bayesian Instrumental Regression for Disparity Estimation

Description

Bayesian models for accurately estimating conditional distributions by race, using Bayesian Improved Surname Geocoding (BISG) probability estimates of individual race. Implements the methods described in McCartan, Fisher, Goldin, Ho and Imai (2025) doi:10.1080/01621459.2025.2526695.

Author(s)

Maintainer: Cory McCartan mccartan@psu.edu

Other contributors:

• Kosuke Imai imai@harvard.edu [contributor]

• Daniel Ho dho@law.stanford.edu [contributor]

• Jacob Goldin jsgoldin@uchicago.edu [contributor]

• Robin Fisher robin.fisher@treasury.gov [contributor]

• The Stan Development Team (include/rstan) [copyright holder]

References

McCartan, C., Fisher, R., Goldin, J., Ho, D.E., & Imai, K. (2025). Estimating Racial Disparities when Race is Not Observed. Journal of the American Statistical Association. Available at doi:10.1080/01621459.2025.2526695.

See Also

Useful links:

• https://github.com/CoryMcCartan/birdie

• https://corymccartan.com/birdie/

• Report bugs at https://github.com/CoryMcCartan/birdie/issues

Fit BIRDiE Models

Description

Fits one of three possible Bayesian Instrumental Regression for Disparity Estimation (BIRDiE) models to BISG probabilities and covariates. The simplest Categorical-Dirichlet model (cat_dir()) is appropriate when there are no covariates or when all covariates are discrete and fully interacted with another. The more general Categorical mixed-effects model (cat_mixed()) is a supports any number of fixed effects and up to one random intercept. For continuous outcomes a Normal linear model is available (gaussian()).

Usage

birdie(
  r_probs,
  formula,
  data,
  family = cat_dir(),
  prior = NULL,
  weights = NULL,
  algorithm = c("em", "gibbs", "em_boot"),
  iter = 400,
  warmup = 50,
  prefix = "pr_",
  ctrl = birdie.ctrl()
)

Arguments

Details

By default, birdie() uses an expectation-maximization (EM) routine to find the maximum a posteriori (MAP) estimate for the specified model. Asymptotic variance-covariance matrices for the MAP estimate are available for the Categorical-Dirichlet and Normal linear models via bootstrapping. Full Bayesian inference is supported via Gibbs sampling for the Categorical-Dirichlet and Normal linear models as well.

Whatever model or method is used, a finite-population estimate of the outcome-given-race distribution for the entire observed sample is always calculated and stored as ⁠$est⁠ in the returned object, which can be accessed with coef.birdie() as well.

The Categorical-Dirichlet model is specified as follows:

Y_i \mid R_i, X_i, \Theta \sim \text{Categorical}(\theta_{R_iX_i}) \\ \theta_{rx} \sim \text{Dirichlet}(\alpha_r),

where Y is the outcome variable, R is race, X are covariates (fixed effects), and \theta_{rx} and \alpha_r are vectors with length matching the number of levels of the outcome variable. There is one vector \theta_{rx} for every combination of race and covariates, hence the need for formula to either have no covariates or a fully interacted structure.

The Categorical mixed-effects model is specified as follows:

Y_i \mid R_i, X_i, \Theta \sim \text{Categorical}(g^{-1}(\mu_{R_iX_i})) \\ \mu_{rxy} = W\beta_{ry} + Zu_{ry} \\ u_{r} \mid \vec\sigma_{r}, L_r \sim \mathcal{N}(0, \text{diag}(\vec\sigma_{r})C_r\text{diag}(\vec\sigma_{r})) \\ \beta_{ry} \sim \mathcal{N}(0, s^2_{r\beta}) \\ \sigma_{ry} \sim \text{Inv-Gamma}(4, 3s_{r\sigma}) \\ C_r \sim \text{LKJ}(2),

where \beta_{ry} are the fixed effects, u_{ry} is the random intercept, and g is a softmax link function. Estimates for \beta_{ry} and \sigma_{ry} are stored in the ⁠$beta⁠ and ⁠$sigma⁠ elements of the fitted model object.

The Normal linear model is specified as follows:

Y_i \mid R_i, \vec X_i, \Theta \sim \mathcal{N}(\vec X_i^\top\vec\theta, \sigma^2) \\ \sigma^2 \sim \text{Inv-Gamma}(n_\sigma/2, l_\sigma^2 n_\sigma/2) \\ \beta_{\text{intercept}} \sim \mathcal{N}(0, s^2_\text{int}) \\ \beta_k \sim \mathcal{N}(0, s^2_\beta), \\

where \vec\theta is a vector of linear model coefficients. Estimates for \theta and \sigma are stored in the ⁠$beta⁠ and ⁠$sigma⁠ elements of the fitted model object.

More details on the models and their properties may be found in the paper referenced below.

Value

An object of class birdie, for which many methods are available. The model estimates may be accessed with coef.birdie(), and updated BISG probabilities (conditioning on the outcome) may be accessed with fitted.birdie(). Uncertainty estimates, if available, can be accessed with ⁠$se⁠ and vcov.birdie().

References

McCartan, C., Fisher, R., Goldin, J., Ho, D.E., & Imai, K. (2025). Estimating Racial Disparities when Race is Not Observed. Journal of the American Statistical Association. Available at doi:10.1080/01621459.2025.2526695.

Examples

data(pseudo_vf)

r_probs = bisg(~ nm(last_name) + zip(zip), data=pseudo_vf)

# Process zip codes to remove missing values
pseudo_vf$zip = proc_zip(pseudo_vf$zip)

fit = birdie(r_probs, turnout ~ 1, data=pseudo_vf)
print(fit)
fit$se # uncertainty quantification

fit = birdie(r_probs, turnout ~ zip, data=pseudo_vf, algorithm="gibbs")

fit = birdie(r_probs, turnout ~ (1 | zip), data=pseudo_vf,
             family=cat_mixed(), ctrl=birdie.ctrl(abstol=1e-3))

summary(fit)
coef(fit)
fitted(fit)

Class "birdie" of BIRDiE Models

Description

The output of birdie() is an object of class birdie, which supports many generic functions. Notably coef.birdie() returns the main model estimates of outcome given race, and fitted.birdie() returns a table analogous to the output of bisg() with updated race probabilities.

Usage

## S3 method for class 'birdie'
coef(object, subgroup = FALSE, ...)

## S3 method for class 'birdie'
fitted(object, ...)

## S3 method for class 'birdie'
residuals(object, x_only = FALSE, ...)

## S3 method for class 'birdie'
predict(object, adj = NULL, ...)

## S3 method for class 'birdie'
simulate(object, nsim = 1, seed = NULL, ...)

## S3 method for class 'birdie'
plot(x, log = FALSE, ...)

## S3 method for class 'birdie'
tidy(x, subgroup = FALSE, ...)

## S3 method for class 'birdie'
glance(x, ...)

## S3 method for class 'birdie'
augment(x, data, ...)

## S3 method for class 'birdie'
formula(x, ...)

## S3 method for class 'birdie'
family(object, ...)

## S3 method for class 'birdie'
nobs(object, ...)

## S3 method for class 'birdie'
vcov(object, ...)

## S3 method for class 'birdie'
print(x, ...)

## S3 method for class 'birdie'
summary(object, ...)

Arguments

Details

The internal structure of birdie objects is not designed to be accessed directly. The generics listed here should be used instead.

Value

Varies, depending on the method. See generic functions' documentation for details.

Functions

• coef(birdie): Return estimated outcome-given-race distributions. When subgroup=FALSE this always returns a finite-population estimate of the outcome-given-race distribution for the observed sample.

• fitted(birdie): Return an updated race probability table. bisg() estimates Pr(R | G, X, S); this table is Pr(R | Y, G, X, S, Theta-hat).

• residuals(birdie): Return the residuals for the outcome variable as a matrix. Useful in sensitivity analyses and to get an idea of how well race, location, names, etc. predict the outcome.

• predict(birdie): Create point predictions of individual race. Returns factor vector of individual race labels. Strongly not recommended for any kind of inferential purpose, as biases may be extreme and in unpredictable directions.

• simulate(birdie): Simulate race from the posterior distribution Pr(R | Y, G, X, S, Theta-hat). Does not account for uncertainty in model parameters.

• plot(birdie): Visualize the estimated conditional distributions for a BIRDiE model. If available, marginal standard error estimates (⁠$se⁠) will be visualized with 95% confidence-level error bars.

• tidy(birdie): Put BIRDiE model coefficients in a tidy format.

• glance(birdie): Glance at a BIRDiE model.

• augment(birdie): Augment data with individual race predictions from a BIRDiE model.

• formula(birdie): Extract the formula used to specify a BIRDiE model.

• family(birdie): Return the BIRDiE complete-data model family.

• nobs(birdie): Return the number of observations used to fit a BIRDiE model.

• vcov(birdie): Return the estimated variance-covariance matrix for the BIRDiE model estimates, if available.

• print(birdie): Print a summary of the model fit.

• summary(birdie): Print a more detailed summary of the model fit.

Examples

methods(class="birdie")

BIRDiE Complete-Data Model Families

Description

BIRDiE supports a number of complete-data outcome models, including categorical regression models. Models specific to BIRDiE are listed here. See the Details section of birdie() for more information about each model.

Usage

cat_dir(link = "identity")

cat_mixed(link = "softmax")

Arguments

Value

A list of class family containing the specification.

Examples

cat_dir()
cat_mixed()

Control of BIRDiE Model Fitting

Description

Constructs control parameters for BIRDiE model fitting. All arguments have defaults.

Usage

birdie.ctrl(
  abstol = 1e-06,
  reltol = 1e-06,
  max_iter = 1000,
  fix_sigma = FALSE,
  accel = c("squarem", "anderson", "daarem", "none"),
  order = switch(match.arg(accel), none = 0L, anderson = -1L, daarem = -1L, squarem = 1L),
  anderson_restart = TRUE
)

Arguments

Value

A list containing the control parameters.

References

Varadhan, R., & Roland, C. (2004). Squared extrapolation methods (SQUAREM): A new class of simple and efficient numerical schemes for accelerating the convergence of the EM algorithm.

Walker, H. F., & Ni, P. (2011). Anderson acceleration for fixed-point iterations. SIAM Journal on Numerical Analysis, 49(4), 1715-1735.

Henderson, N. C., & Varadhan, R. (2019). Damped Anderson acceleration with restarts and monotonicity control for accelerating EM and EM-like algorithms. Journal of Computational and Graphical Statistics, 28(4), 834-846.

Examples

str(birdie.ctrl(max_iter=100))

Bayesian Improved Surname Geocoding (BISG)

Description

Calculates individual probabilities of belonging to racial groups given last name, location, and other covariates (optional). The standard function bisg() treats the input tables as fixed. An alternative function bisg_me(), assumes that the input tables are subject to measurement error, and uses a Gibbs sampler to impute the individual race probabilities, using the model of Imai et al. (2022).

Usage

bisg(
  formula,
  data = NULL,
  p_r = p_r_natl(),
  p_rgx = NULL,
  p_rs = NULL,
  save_rgx = TRUE
)

bisg_me(
  formula,
  data = NULL,
  p_r = p_r_natl(),
  p_rgx = NULL,
  p_rs = NULL,
  iter = 1000,
  warmup = 100,
  cores = 1L
)

## S3 method for class 'bisg'
summary(object, p_r = NULL, ...)

## S3 method for class 'bisg'
predict(object, adj = NULL, ...)

## S3 method for class 'bisg'
simulate(object, nsim = 1, seed = NULL, ...)

Arguments

Value

An object of class bisg, which is just a data frame with some additional attributes. The data frame has rows matching the input data and columns for the race probabilities.

Methods (by generic)

• summary(bisg): Summarize predicted race probabilities. Returns vector of individual entropies.

• predict(bisg): Create point predictions of individual race. Returns factor vector of individual race labels. Strongly not recommended for any kind of inferential purpose, as biases may be extreme and in unpredictable directions.

• simulate(bisg): Simulate race from the Pr(R | G, X, S) distribution.

Functions

• bisg(): The standard BISG model.

• bisg_me(): The measurement error BISG model.

References

Elliott, M. N., Fremont, A., Morrison, P. A., Pantoja, P., and Lurie, N. (2008). A new method for estimating race/ethnicity and associated disparities where administrative records lack self-reported race/ethnicity. Health Services Research, 43(5p1):1722–1736.

Fiscella, K. and Fremont, A. M. (2006). Use of geocoding and surname analysis to estimate race and ethnicity. Health Services Research, 41(4p1):1482–1500.

Imai, K., Olivella, S., & Rosenman, E. T. (2022). Addressing census data problems in race imputation via fully Bayesian Improved Surname Geocoding and name supplements. Science Advances, 8(49), eadc9824.

Examples

data(pseudo_vf)
bisg(~ nm(last_name), data=pseudo_vf)

r_probs = bisg(~ nm(last_name) + zip(zip), data=pseudo_vf)
summary(r_probs)
head(predict(r_probs))

data(pseudo_vf)
bisg_me(~ nm(last_name) + zip(zip), data=pseudo_vf)

Download Census Race Data

Description

Downloads and prepares race-by-geography tables from U.S. census data, using the easycensus package. Requires that an api key be set up through easycensus::cens_auth() in that package, usually by storing it in the CENSUS_API_KEY environment variable. Supports data from the decennial census and the American Community Survey at a variety of levels of geographic detail. The output of this function can be used directly in bisg().

Usage

census_race_geo_table(
  geo = c("us", "state", "county", "zcta", "tract"),
  ...,
  year = 2010,
  survey = c("dec", "acs1", "acs5"),
  GEOIDs = TRUE,
  counts = TRUE
)

Arguments

Value

A data frame with geographic identifier column(s) and six columns white, black, etc. containing the counts or proportion of residents in each racial group.

Examples

census_race_geo_table("zcta", year=2010)
## Not run: 
# Census API key required
census_race_geo_table("us", year=2010)
census_race_geo_table("state", year=2021, survey="acs1")
census_race_geo_table("county", state="NH", year=2020, GEOIDs=FALSE)

## End(Not run)

Compute Racial Disparities from Model Estimates

Description

This function lets you easily compute differences in conditional expectations between all pairs of specified racial groups.

Usage

disparities(x, subgroup = FALSE, races = TRUE)

Arguments

Value

A data frame containing a row with every possible disparity for the specified races, which are identified by columns race_1 and race_2. The reported disparity is estimate_1 - estimate_2.

Examples

data(pseudo_vf)
r_probs = bisg(~ nm(last_name) + zip(zip), data=pseudo_vf)
fit = birdie(r_probs, turnout ~ 1, data=pseudo_vf)

disparities(fit)
disparities(fit, races=c("white", "black"))

Calculate Weighted Estimate of (Discrete) Outcomes By Race

Description

Calculates the "standard" weighted estimator of conditional distributions of an outcome variable Y by race R, using BISG probabilities. This estimator, while commonly used, is only appropriate if Y \perp R \mid X, S, where S and X are the last names and covariates (possibly including geography) used in making the BISG probabilities. In most cases this assumption is not plausible and birdie() should be used instead. See the references below for more discussion as to selecting the right estimator.

Up to Monte Carlo error, the weighted estimate is equivalent to performing multiple imputations of the race vector from the BISG probabilities and then using them inside a weighted average or linear regression.

Usage

est_weighted(
  r_probs,
  formula,
  data = NULL,
  weights = NULL,
  prefix = "pr_",
  se_boot = 0
)

## S3 method for class 'est_weighted'
print(x, ...)

## S3 method for class 'est_weighted'
summary(object, ...)

Arguments

Value

An object of class est_weighted, inheriting from birdie, for which many methods are available. The model estimates may be accessed with coef(). Uncertainty estimates, if available, can be accessed with ⁠$se⁠ and vcov.birdie().

Methods (by generic)

• print(est_weighted): Print a summary of the model fit.

• summary(est_weighted): Print a more detailed summary of the model fit.

References

McCartan, C., Fisher, R., Goldin, J., Ho, D.E., & Imai, K. (2025). Estimating Racial Disparities when Race is Not Observed. Journal of the American Statistical Association. Available at doi:10.1080/01621459.2025.2526695.

Examples

data(pseudo_vf)

r_probs = bisg(~ nm(last_name) + zip(zip), data=pseudo_vf)

# Process zip codes to remove missing values
pseudo_vf$zip = proc_zip(pseudo_vf$zip)

est_weighted(r_probs, turnout ~ 1, data=pseudo_vf)

est = est_weighted(r_probs, turnout ~ zip, data=pseudo_vf)
tidy(est, subgroup=TRUE)

National Racial Demographics

Description

Returns the proportion of the U.S. population in six racial groups in a given year. Group definitions necessarily follow those used by the Census Bureau in its surname tables:

• white: Non-Hispanic White alone

• black: Non-Hispanic Black alone

• hisp: Hispanic, any race

• asian: Non-Hispanic Asian, Native Hawaiian, or Pacific Islander alone

• aian: Non-Hispanic American Indian/Alaska Native

• other: Non-Hispanic, two or more races, or other race

Usage

p_r_natl(year = 2021, vap = FALSE)

Arguments

Value

A named numeric vector of length 6.

Examples

p_r_natl(year=2010)

Preprocess Last Names and Geographic Identifiers

Description

These functions are called automatically by bisg() but may be useful, especially when geographic variables are included in a birdie() model. proc_zip() and proc_state() preprocess their corresponding geographic identifiers. States are partially matched to state names and abbreviations and are returned as FIPS codes. ZIP codes are crosswalked to Census ZCTAs. Missing identifiers are replaced with "<none>". proc_name() processes last names in accordance with Census processing rules (https://www2.census.gov/topics/genealogy/2010surnames/surnames.pdf). Names are converted to Latin characters, capitalized, stripped of prefixes and suffixes, and otherwise standardized.

Usage

proc_zip(x)

proc_state(x)

proc_name(x, to_latin = TRUE)

Arguments

Value

A processed character vector

Functions

• proc_zip(): Match ZIP codes to ZCTAs and fill in missing values.

• proc_state(): Match state names and abbreviations and fill in missing values.

• proc_name(): Process names to a Census-standardized format.

Examples

proc_name("Smith Jr.")
proc_zip("00501")
proc_state("Washington")

A pseudo-voterfile

Description

A dataset containing 5,000 fake voter records. Created by randomizing a subset of the North Carolina voter file. Turnout records are completely randomly generated.

Usage

pseudo_vf

Format

A data frame with 5,000 rows and 4 records:

last_name

Voter's last name

zip

5-digit ZIP code. May be NA

race

One of "white", "black", "hisp", "asian", "aian", or "other"

turnout

1 if the voter voted in the most recent election, 0 otherwise

Source

https://www.ncsbe.gov/results-data/voter-registration-data

Examples

data(pseudo_vf)
print(pseudo_vf)

Objects exported from other packages

Description

These objects are imported from other packages. Follow the links below to see their documentation.

generics

augment, glance, tidy