Marriage tax 2022

marriage_tax_2022.R

library(data.table)
library(rvest)
library(xml2)

brackets = "https://taxfoundation.org/data/all/federal/2022-tax-brackets/" |>
  read_html() |>
  xml_find_all("//table[contains(caption/text(), '2022 Federal Income')]") |>
  html_table()

brackets = as.data.table(brackets[[1L]])[-.N]

brackets[, marginal_rate := as.numeric(gsub("%", "", `Tax Rate`)) / 100]
lower_bound <- \(x) as.numeric(gsub("[$,]|\\s.*", "", x))
brackets[, lower_bound_single := lower_bound(`For Single Filers`)]
brackets[, lower_bound_married := lower_bound(`For Married Individuals Filing Joint Returns`)]

brackets = melt(
  brackets,
  id.vars = "marginal_rate",
  measure.vars = patterns("^lower_bound"),
  variable.name = "status",
  value.name = "lower_bound"
)
brackets[, status := gsub("lower_bound_", "", status)]
setkey(brackets, status)

brackets[, upper_bound := shift(lower_bound, -1L, fill = Inf), by = status]
brackets[, cum_tax_init := shift(cumsum(marginal_rate * (upper_bound - lower_bound)), fill = 0), by = status]

deductible <- function(status) switch(status, single = 12950, married = 25900)

tax <- function(income, STATUS) {
  taxable_income <- pmax(income - deductible(STATUS), 0)
  brackets[.(STATUS), {
    .SD[
      findInterval(taxable_income, lower_bound),
      cum_tax_init + marginal_rate * (taxable_income - lower_bound)
    ]
  }]
}

income_grid <- seq(0, 1e6, length.out = 250)

income_pairs <- CJ(income1 = income_grid, income2 = income_grid)

income_pairs[, tax_burden_separate := tax(income1, "single") + tax(income2, "single")]
income_pairs[, tax_burden_joint := tax(income1 + income2, "married")]

income_pairs[, plot(
  income1, income2,
  pch = 22L, col = fifelse(tax_burden_separate < tax_burden_joint, "red", "blue")
 )]