Calibrating a tax-benefit function¶

A tax-benefit system has hundreds of instruments. Macroeconomic models need it in a compact form: a function from a household’s gross income to the taxes it pays, net of the benefits it receives. euromod.eu_euromod.calibrate estimates such a function from a EUROMOD simulation, for a baseline or for a reform.

from euromod.eu_euromod.calibrate import prepare_sample, fit_hsv

MODEL = r"C:\EUROMOD_RELEASES_J2.0+"      # folder of the model
DATA = r"C:\EUROMOD_DATA"                 # folder with the input datasets

baseline = prepare_sample(MODEL, DATA, "PT_2022")
fit_hsv(baseline)

{'lambda': 0.1925, 'tau': 0.1631, 'se_intercept': 0.0015, 'se_lambda': 0.0012,
 'se_tau': 0.0024, 'r2': 0.9565, 'n': 5751}

The numbers on this page are from one release of the model and its input data. Yours will differ with another release.

How it works¶

The function¶

The function is the one of Heathcote, Storesletten and Violante (2017). For a household with gross income $y$, measured in multiples of mean income, the average net tax rate is

\[T(y) = 1 - (1 - \lambda)\, y^{-\tau}\]

$\lambda$ is the level: the net tax rate of a household at mean income.
$\tau$ is the progressivity: how fast the rate rises with income. With $\tau = 0$ every household pays the same rate.

Benefits are counted as negative taxes, so the rate is below zero for households that receive more than they pay.

The household sample¶

prepare_sample runs a system of the model and turns EUROMOD’s individual-level output into one row per household:

The head of a household is its main earner: the member with the highest earnings (ils_earns), the oldest one if earnings are equal, and the one listed first in the output if age is equal too.
Households are kept when the head is 25 to 60 years old and no member has self-employment income (yse). The self-employed face different rules, and their incomes are measured less reliably.
Incomes, taxes and benefits are summed over the members of the household. Gross income is original income plus replacement income (ils_origrepy); disposable income is ils_dispy.
Households with non-positive gross income, original income, disposable income or income after taxes are dropped. Households whose taxes are negative in total are kept: refundable tax credits are support like a benefit, and dropping their recipients would make the sample depend on whether a country pays its support through the tax system. Version 0.4 dropped them, which removed up to a fifth of the households in some countries; drop_negative_tax=True does so again.
Households below the 2.5th or above the 99th weighted percentile of gross income are dropped.
Incomes are divided by the weighted mean gross income of the households that remain, so that 1.0 is mean income.

Every choice in steps 2, 3, 4 and 5 is an argument of household_sample, e.g. prepare_sample(MODEL, DATA, "PT_2022", head_age=(25, 64), trim=(1, 99)).

The estimator¶

In logs the function is a straight line,

\[\log(\text{disposable income}) = \log(1 - \lambda) + (1 - \tau) \log(y),\]

which fit_hsv estimates by weighted least squares with the survey weights (dwt). The solution is in closed form: nothing is optimised numerically, so the result does not depend on starting values.

Using it¶

Household types¶

Household types are subsets of one sample. They keep its income scale, so their estimates can be compared with each other.

from euromod.eu_euromod.calibrate import HOUSEHOLD_TYPES, select_households

fit_hsv(select_households(baseline, "married"))["tau"]            # 0.1918
fit_hsv(select_households(baseline, "single_no_children"))["tau"] # 0.1309
fit_hsv(select_households(baseline, "hh_size == 1"))              # any condition

HOUSEHOLD_TYPES holds seven types: all, married, single, and the last two with and without children. A household is married when its head is (dms == 2); a child is a member aged 18 or under with a parent in the data.

Several countries, and the euro area¶

calibrate_hsv does everything in one call and returns a table with one row per system and household type:

from euromod.eu_euromod import euro_area, eu27_systems
from euromod.eu_euromod.calibrate import calibrate_hsv

table = calibrate_hsv(MODEL, DATA, eu27_systems(2022), pooled={"EA_2022": euro_area(2022)})
table.loc["PT_2022", ["lambda", "tau", "se_tau", "r2", "n"]]

                       lambda     tau  se_tau      r2     n
household_type
all                    0.1925  0.1631  0.0024  0.9565  5751
married                0.1783  0.1918  0.0030  0.9663  2616
single                 0.2012  0.1530  0.0036  0.9478  3135
...

pooled adds an estimate on the pooled sample of a group of countries. Each country keeps its own income scale and its own survey weights, so it counts in proportion to its population. euro_area(year) gives the members of a year: 19 up to 2022, 20 from 2023, when Croatia joined.

How well the function fits¶

from euromod.eu_euromod.calibrate import fit_diagnostics

fit = fit_hsv(baseline)
fit_diagnostics(baseline, fit)

Key	Meaning
`net_tax_rate_at_mean`	Observed net tax rate of households at mean income: the data counterpart of `lambda`.
`rmse`, `mae`	Errors of the net tax rate, household by household. They mostly reflect the spread of households around the function.
`rmse_binned`, `mae_binned`	The same against mean rates in 80 income bins. They measure the shape of the function.
`gini_gross`, `gini_disposable`	Observed Gini coefficients.
`gini_fitted`	Gini of the disposable income that the fitted function assigns to the same gross incomes.

Describing the sample¶

from euromod.eu_euromod.calibrate import (sample_characteristics, income_shares,
    gini_coefficients, tax_rates_by_quantile, benefit_ratio_by_quantile)

sample_characteristics(baseline)     # households, age of head, share married, ...
income_shares(baseline)              # quintile shares of gross and disposable income
tax_rates_by_quantile(baseline)      # mean net tax rate by decile and household type
benefit_ratio_by_quantile(baseline)  # means-tested benefits over gross income, by decile

Using the estimates¶

The published estimates¶

The estimates for the 27 EU member states and the euro area, policy years 2022 to 2025, ship with the package. They need neither the model nor any data:

from euromod.eu_euromod.calibrate import load_estimates, schedule

baseline = load_estimates()          # one row per system and household type
baseline.loc[("PT_2022", "all"), ["lambda", "tau"]]      # 0.1925, 0.1631
baseline.xs("married", level="household_type")["tau"]    # every system, married households
baseline.loc[[("PT_2022", "all"), ("PT_2025", "all")], "tau"]   # 0.1631, 0.1730
baseline.attrs                       # model release, policy years, sample definition

portugal = schedule(baseline, "PT_2022")                 # the fitted function, ready to use
portugal.net_tax_rate([0.5, 1.0, 2.0])                   # array([0.0958, 0.1925, 0.2788])
portugal.net_taxes(2.0)                                  # 0.5576: taxes net of benefits
portugal.disposable_income(2.0)                          # 1.4424
portugal.marginal_rate(2.0)                              # 0.3964
schedule(baseline, "PT_2022", "married")                 # the same for one household type

Incomes are in multiples of mean gross income: 2.0 is twice the mean, and a disposable income of 1.4424 is 1.4424 times mean gross income.

The systems are AT_2022 … SK_2025, every member state in every year, and EA_2022 … EA_2025, the pooled euro area of each year (19 members in 2022, 20 from 2023). Each system runs on the dataset the model marks as its best match. The columns are those of fit_hsv and fit_diagnostics. schedule works on any table of estimates, including your own.

Saving and loading your own¶

calibrate_hsv returns the same kind of table. save_estimates writes it to a JSON file, with any notes you give it, and load_estimates reads it back:

from euromod.eu_euromod.calibrate import calibrate_hsv, save_estimates, load_estimates, schedule

reform = {"$tin00_rate9": 0.58}
estimates = calibrate_hsv(MODEL, DATA, ["PT_2022"], constants=reform)
save_estimates(estimates, "pt_reform.json", constants=reform, note="top rate raised to 58%")

estimates = load_estimates("pt_reform.json")    # later, without running EUROMOD
estimates.attrs["note"]                         # 'top rate raised to 58%'
schedule(estimates, "PT_2022").disposable_income(2.0)    # used at once, like the published ones

Evaluating the function¶

hsv evaluates the function with numpy:

import numpy as np
from euromod.eu_euromod.calibrate import hsv, marginal_hsv

incomes = np.array([0.5, 1.0, 2.0])             # multiples of mean income
hsv(incomes, fit["lambda"], fit["tau"])           # array([0.0958, 0.1925, 0.2788])
marginal_hsv(incomes, fit["lambda"], fit["tau"])  # array([0.2433, 0.3242, 0.3964])

to_calibration turns a fit into a Calibration of the connector’s general calibration engine. It can be saved as JSON and gives a numba-compiled function, which another compiled function can call:

import numba
from euromod.eu_euromod.calibrate import to_calibration

cal = to_calibration(fit, system="PT_2022")
cal(incomes)                                    # array([0.0958, 0.1925, 0.2788])
cal.to_json("PT_2022_hsv.json")

schedule = cal.numba_function                   # schedule(income, lambda_, tau)

@numba.njit
def tax_bill(incomes, lambda_, tau):
    total = 0.0
    for y in incomes:
        total += y * schedule(y, lambda_, tau)
    return total

On arrays numpy and numba are equally fast. The compiled function matters when the schedule is called inside a model’s own compiled loop.

Reproducing a reform¶

A reform is either a counterfactual system defined in the model, or a baseline system with some constants overridden. Name a different system, or pass the constants:

reform = prepare_sample(MODEL, DATA, "PT_2022_ref")

reform = prepare_sample(MODEL, DATA, "PT_2022",
                        constants={"$tin00_rate6": 0.47, "$tin00_rate7": 0.535,
                                   "$tin00_rate8": 0.55, "$tin00_rate9": 0.58})

To see what a reform does to progressivity, hold the level fixed and let only tau move:

fit_hsv(baseline, lambda_fixed=0.192)["tau"]    # 0.1629
fit_hsv(reform, lambda_fixed=0.192)["tau"]      # 0.1733

This example raises the four highest income-tax rates of Portugal by ten percentage points. With lambda_fixed, tau comes from a regression through the origin of log(disposable income) - log(1 - lambda) on log(y).

Three things to keep in mind:

A baseline and a reform are only comparable on the same input dataset. A counterfactual system without a best-match dataset of its own needs dataset="...".
A constant that is defined per group takes the pair as key: constants={("$f_cpi", "2023"): 1.02}.
The reform sample is normalised by its own mean income. A reform that changes gross incomes therefore also moves the income scale.

What it needs¶

Paths. The folder of the model and the folder with the input datasets. A dataset is looked up as <data>/<name>.txt and as <data>/<CC>/<name>.txt. Paths are always arguments: nothing is read from the environment.

From the model.


Variables	`idhh`, `dag`, `dms`, `dwt`, `yse`, `idmother`, `idfather`
Income lists	`ils_earns`, `ils_origy`, `ils_origrepy`, `ils_dispy`, `ils_taxsim`, `ils_benmt`, `ils_bennt`
Conventions	`dms == 2` is married; a system name starts with its two-letter country code

If the output of a run lacks one of these, household_sample stops and names it.

Packages. numpy and pandas. Only to_calibration needs more: install euromod[calibrate] for sympy and numba.

Keeping simulation output¶

A run takes from seconds to minutes. With cache_dir the output is kept and reused on the next call:

baseline = prepare_sample(MODEL, DATA, "PT_2022", cache_dir="runs/baseline")
reform = prepare_sample(MODEL, DATA, "PT_2022", constants=REFORM, cache_dir="runs/reform")

Use one folder per scenario. Beside each output a small .json file records the dataset, constants, add-ons and switches that produced it. Asking for the same system with other settings is an error, so a stale file cannot pass for a reform. The files hold household microdata: keep them where that data is allowed to be.

Reference¶

Heathcote, J., Storesletten, K. and Violante, G. L. (2017). Optimal tax progressivity: an analytical framework. Quarterly Journal of Economics, 132(4), 1693-1754.